Fs = 515; Fstep = 7; Fmin = 10;
while Fs>=Fmin
%for Fs = 515:-7:5
%Fs = 250;
F0 = 100;
f0 = F0/Fs;
f0 = f0 - floor(f0+0.5); % alias
Nfactor = 100;
T = 1/Fs;
Nfreq = 1e4;
maxfreq = 3.5;
freq = (2*[0:Nfreq-1]/Nfreq - 1)*maxfreq*Fs;

timelen = 5*Fs/F0;
ctime = [0:timelen*Nfactor-1]/(Fs*Nfactor);
dtime = [0:timelen-1];
xc = cos(2*pi*F0*ctime);
xd = cos(2*pi*f0*dtime);
xdc = cos(2*pi*f0*Fs*ctime);

subplot(311);
plot(ctime,xc,'b'); hold on;
stem(dtime/Fs,xd,'.r'); hold on;
plot(ctime,xdc,'r'); hold off;
title(sprintf('F0 = %f, Fs = %f, f0 = %f',F0,Fs,f0),'FontSize',20);

subplot(312);
stem(F0*[-1 1],[1 1],'b.');
xlim(maxfreq*Fs*[-1 1]);
ylim([0 1.2]);

subplot(313);
stem(F0*[-1 1],[1 1],'b.','LineWidth',2); hold on;
ylim([0 1.2]);
num = 10;
map = winter(num);
for i=1:num
  stem(F0*[-1 1]-i*Fs,[1 1],'o','Color',map(end+1-i,:)); hold on;
  stem(F0*[-1 1]+i*Fs,[1 1],'x','Color',map(i,:)); hold on;
end
ax = axis;
plot(-Fs/2*[1 1],ax(3:4),'--k');
plot(+Fs/2*[1 1],ax(3:4),'--k');
hold off;
xlim(maxfreq*Fs*[-1 1]);
title(sprintf('Sample Rate = %d samples/second',Fs),'FontSize',20);

if(Fs<=2*F0)
    Fstep = 2;
end
Fs = Fs - Fstep;

pause;
end;
