%this demo shows the frequencies of the sampled
%Cos(w0*t)=Cos(w0*k*h) in the frequency range [0 w_s]
% as a function of the samping period h

clear
clf

w0=2;   %input frequency, integer
w=0:pi/100:2*pi;

ww=[6 4.5 3] %sampling frequencies
hh=2*pi./ww;  %sampling period

for i=1:3
  h=hh(i);  %sampling period
  ws=ww(i);  %sampling frequency
  num=[1 -cos(w0*h) 0];
  den=[1 -2*cos(w0*h) 1];  %num and den of Cos(w0*k*h)
  y=freqz(num,den,w);  %discrete time Fourier transform

  if i == 1
    subplot(311)
    elseif i == 2
    subplot(312)
    else
    subplot(313)
  end

  plot(w/h,abs(y))
  hold on
  plot([ws ws],[0 max(abs(y))],'-.')
  axis([0 1.1*2*pi/h 0 max(abs(y))])
  set(gca,'YTickLabel',[])

  if i == 1
    ss='Fourier transform of Cos(q*k*h)';
    ss(26)=num2str(w0);
    title(ss,'fontsize',12)
  end
end

xlabel('frequency \omega  (dash-dot line is at the sampling frequency \omega_s)','fontsize',12)