set(0,'defaultaxesfontsize',22);
clear all
clf
a=0.8;
N=41; off=(N-1)/2;
D=off;
Nfft=1024;
Xa=zeros(1,Nfft);
%wd=linspace(0,2*pi-2*pi/Nfft,Nfft);
wd=linspace(-pi,pi,Nfft);
wda=linspace(0,2*pi-2*pi/N,N);
%
Xw=ones(1,N)./(1-.8*exp(-j*wda))+ones(1,N)./(1-.9*exp(-j*wda));
Xwd=ones(1,Nfft)./(1-.8*exp(-j*wd))+ones(1,Nfft)./(1-.9*exp(-j*wd));
%Xw=((1-a^2)./((1-2*a*cos(wda)+a^2))).*exp(-j*wda*off);
%Xwd=((1-a^2)./((1-2*a*cos(wd)+a^2))).*exp(-j*wd*off);
xa=real((ifft(Xw,N)));
%
Xrec=abs(fft(xa,Nfft));
for k=0:N-1,
wk=2*pi*k/N;
wvec(1,k+1)=wk;
Xk=1/(1-.8*exp(-j*wk))+1/(1-.9*exp(-j*wk));
%Xk=(1-a^2)./((1-2*a*cos(wk)+a^2)).*exp(-j*wk*off);
Xvec(1,k+1)=Xk;
Xa=Xa+Xk*(sin((wd-wk)*N/2).*exp(-j*(wd-wk)*(N-1)/2))./(N*sin((wd-wk)/2));
end
%
plot(wd,abs(Xa),'b','Linewidth',3)
%axis([0 2*pi 0 max(abs(Xwd))])
axis([-pi pi 0 max(abs(Xwd))])
hold on
plot(wd,abs(Xwd),'r','Linewidth',3)
plot(wvec-(pi-pi/N),abs(fftshift(Xvec)),'bx','MarkerSize',16,'Linewidth',3)
%plot(wd,Xrec,'k','Linewidth',3)
%legend('Reconstructed Spectrum from N samples','True Spectrum','Zero-Padded FFT of IFFT of N samples')
legend('Reconstructed Spectrum from N samples','True Spectrum')
hold off
pause
n=-off+D:off+D;
n=0:(N-1);
stem(n,xa,'b','Linewidth',3);
%xalias=a.^abs(n-D)+a.^abs(n-D-N)+ a.^abs(n-D+N);
xalias=0.8.^n+0.8.^(n+N)+ 0.8.^(n+2*N);
xalias=xalias+0.9.^n+0.9.^(n+N)+ 0.9.^(n+2*N);
xalias=xalias+0.9.^(n+3*N)+ 0.8.^(n+3*N);
%
hold on
stem(n,xalias,'k','Linewidth',3)
stem(n,0.8.^n+0.9.^n,'r','Linewidth',3)
%stem(n,0.8.^abs(n-D),'r','Linewidth',3)
legend('N-pt IFFT of sampled X(w)','Time-Domain Aliasing Formula','True x[n]')
hold off