%***********************************************************************************
%
%                              Jim Merriam
%                       University of Saskatchewan
%                          jim.merriam@usask.ca
%                              (June, 1995)  
%
%***********************************************************************************
%
clf
clear s;
subplot(211),axis([0 10 0 10]),axis('off')
s(1)={'\bf DISCRETE FOURIER TRANSFORMS II\rm'};
text(0,8,s,'FontName','times','FontSize',16),clear s;
s(1)={'This is the second demo on Discrete Fourier transforms. '};
s(2)={'It demonstrates:'};
s(3)={' '};
s(4)={'SPECTRAL LEAKAGE'};
s(5)={'ALIASING'};
s(6)={'GIBBS PHENOMENON'};
s(7)={'SPECTRUM OF RANDOM DATA'};
s(8)={'THE EFFECT OF LONG PERIODS'};
text(1,1,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton', 'units','normalized','position',[.8 .01 .1 .05], ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
axis([0 10 0 10]),axis('off')
s(1)={'In the last demo, the spectra contained isolated peaks at'};
s(2)={'exactly the frequency of the sine wave. This is not generally'};
s(3)={'the case, and was contrived by picking the sampling interval '};
s(4)={'and the number of samples so that one of the discrete Fourier'};
s(5)={'frequencies would be exactly equal to the frequency of the '};
s(6)={'sampled signal.'};
text(1,5,s,'FontName','times','FontSize',10),clear s


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end


clf
x=[1:128];
y=sin(x*pi/50);
Y=fft(y);
subplot(211),plot(imag(Y),'*'),title('IMAGINARY PART'),xlabel('FREQUENCY INDEX')
subplot(212),axis([0 10 0 10]),axis('off')
s(1)={'This is the also spectrum of a pure sinusoid. It looks different from'};
s(2)={'the spectra in the last examples because the "peaks" are composed'};
s(3)={'of not just one point, but several. In the previous examples, one'};
s(4)={'of the Fourier frequencies was exactly equal to the frequency of'};
s(5)={'the wave. In this case that is not true - none of the discrete'};
s(6)={'Fourier frequencies is exactly equal to the frequency of the signal,'};
s(7)={'so the spectrum peaks at the discrete frequency closest to the  actual '};
s(8)={'frequency and information is \bf LEAKED \rm to nearby frequencies.'};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
i=0:5;
j=[1:512];		% an apparent error fixed by inserting this line - I.M., Oct 6 2014
f2=sin(j*2*pi/32);
f1=f2(1:30:512);
F1=fft(f1);
subplot(211)
plot([1:30:512],f1,'*',j,f2,'--');
subplot(212);axis([0 10 0 10]);axis('off');
s(1)={'This is an example of \bf ALIASED \rm data. Notice that the samples'};
s(2)={'obviously do not reflect all of the intricacies in the data.'};
s(3)={'The signal has frequency 16 cy/18 samples =0.88cy/sample '};
s(4)={'but the sampled version looks like it has frequency of '};
s(5)={'about 1 cy/18 samples, or 0.055cy/sample.'};
s(6)={' '};
s(7)={'Which frequency will the transform indicate? '};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
hold off;
F1=fft(f1);
subplot(211)
plot(imag(F1),'*');
title('THIS IS THE IMAGINARY PART OF THE TRANSFORM'),ylabel('AMPLITUDE')
subplot(212);axis([0 10 0 10]);axis('off');
xlabel('THIS IS THE FREQUENCY INDEX m')
s(1)={'Here is the imaginary part of the transform.'};
s(2)={'What is the frequency indicated by the peak in the'};
s(3)={'transform? Is the 0.88 cy/sample signal indicated?'};
s(4)={' '};
s(5)={'The signal at 0.88 cy/sample has been aliased to '};
s(6)={'      fn-(0.88-fn)=0.12cy/sample'};
s(7)={'by undersampling. That is, our discrete samples were ' };
s(8)={'not close enough together to properly describe the signal.'};
text(1,5,s,'FontName','times','FontSize',10),clear s


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
x=ones(1,512);
stp=[x,-x,-x,x];
subplot(211),axis([1,length(stp),-2,2]),hold on,plot(stp)
subplot(212);axis([0 10 0 10]);axis('off');
s(1)={'Here is a square wave with a period of 2048 samples.'};
s(2)={'It will be resynthesized below by summing over a different'};
s(3)={'number of frequencies each time.'};
text(1,5,s,'FontName','times','FontSize',10),clear s


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
F1=fft(stp);
F2=F1;
F2(1,10:2040)=zeros(1,2031);
rs10=ifft(F2);
subplot(414),plot(stp,'-'),hold on,plot(real(rs10),':')
title('WITH 10 TERMS IN THE SUM')
F2(1,5:2045)=zeros(1,2041);
rs3=ifft(F2);
subplot(413),plot(stp,'-'),hold on,plot(real(rs3),':')
title('WITH 3 TERMS IN THE SUM')
F2(1,4:2046)=zeros(1,2043);
rs2=ifft(F2);
subplot(412),plot(stp,'-'),hold on,plot(real(rs2),':')
title('WITH 2 TERMS IN THE SUM')
F2(1,3:2047)=zeros(1,2045);
rs1=ifft(F2);
subplot(411),plot(stp,'-'),hold on, plot(real(rs1),':')
title('WITH 1 TERM IN THE SUM')


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
subplot(211)
axis([0 10 0 10]);axis('off');
s(1)={'Notice that as the number of terms used to resynthesize the signal is'};
s(2)={'increased the resynthesized signal looks more and more like the original,'};
s(3)={'and that by using only ten out of over a thousand terms, the resynthesized'}; 
s(4)={'curve looks pretty similar to the original. There is not much difference'}; 
s(5)={'between using one term and two terms, because term two is small '};
s(6)={'(as are all the even ones).'};  
s(7)={'The Gibbs phenomenon can be seen, but it will be even more'};
s(8)={'apparent if we plot the difference between the original and '};
s(9)={'the resynthesized data.'};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
subplot(411),plot(stp-real(rs1),'-')
title('ORIGINAL-RESYNTHESIZED WITH 1 TERM')
subplot(412),plot(stp-real(rs2),'-')
title('ORIGINAL-RESYNTHESIZED WITH 2 TERMS')
subplot(413),plot(stp-real(rs3),'-')
title('ORIGINAL-RESYNTHESIZED WITH 3 TERMS')
subplot(414),plot(stp-real(rs10),'-')
title('ORIGINAL-RESYNTHESIZED WITH 10 TERMS')


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end


clf;
x=randn(1,2046);
subplot(211),plot(x)
subplot(212),axis([0 10 0 10]);axis('off');
s(1)={'These are  random data, with a normal distribution, a mean of'};
s(2)={'zero, and a standard deviation of 1.'};
s(3)={'The Fourier transform of random data will be seen to be flat,'};
s(4)={'will have the same appearence, that is the (average) amplitude '};
s(5)={'eill be independent of frequency.'};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
X=fft(x);
subplot(211),plot(real(X));
subplot(212),axis([0 10 0 10]);axis('off');
s(1)={'Here is the transform of the random data. Notice that the mean'};
s(2)={'amplitude does not change across the spectrum (it is flat).'};
s(3)={'This is a simple FFT. From a Power Spectral Density, you would'};
s(4)={'be able to measure the standard deviation in the original data.'};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf;
x=[1:1024];
X=fft(x);
subplot(311),plot(x)
subplot(312),axis([0 10 0 10]);axis('off');
s(1)={'This is a linear ramp function, and below is the transform of the'};
s(2)={'ramp. Notice that the ramp has power mostly at low frequencies.'};
text(1,5,s,'FontName','times','FontSize',10),clear s
subplot(313),plot(imag(X));


h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end

clf; 
y=sin(x*1.2*pi/1024);
Y=fft(y);
subplot(311),plot(y)
subplot(312),axis([0 10 0 10]);axis('off');
s(1)={'This is less than a cycle of a sine wave, and below is the'};
s(2)={'transform of this function. Notice that the spectrum looks '}; 
s(3)={'somewhat like the spectrum of the ramp. '}; 
s(4)={'This is a perfect sinusoid, but less than a cycle is  observed,'};
s(5)={'so instead of getting a peak at the frequency of that wave,'};
s(6)={'the information is smeared out over a broad range of frequencies.'};
text(1,5,s,'FontName','times','FontSize',10),clear s
subplot(313),plot(real(Y));

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end


clf;


subplot(211),axis([0 10 0 10]),axis('off')
s(1)={'The lesson from these last two examples is to avoid having anything'};
s(2)={'with a "wavelength" longer than the length of the data set, because'};
s(3)={'the longest period in the transform (the lowest frequency)'};
s(4)={'will be only as long as the data set. Anything with a longer'};
s(5)={'period (ramps etc) will be smeared out over the frequencies'};
s(6)={'available to the transform. This is somewhat like a long '};
s(7)={'period equivalent of aliasing, but it is even deadlier than aliasing,'};
s(8)={'because with aliasing a high frequency will alias to (predictable)'};
s(9)={'lower frequencies without being smeared out over a frequency band.'};
text(1,5,s,'FontName','times','FontSize',10),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end


clf
axis([0 10 0 10]);axis('off');
s(1)={'\bf DISCRETE FOURIER TRANSFORMS II\rm'};
s(2)={''};
s(3)={'                 THE END'};
s(4)={' '};
text(1,6,s,'FontName','times','FontSize',16),clear s

h=uicontrol('style','pushbutton','position',[.8 .01 .1 .05], 'units','normalized', ...
'string','OK','Callback','OK=''continue '';');
OK='PAUSEaBIT';
while OK=='PAUSEaBIT'
   pause(0.01)
end


clf