% file:      	chirp_spectrogram.m, (c) Matthew Roughan, Thu Oct  7 2004
%
N = 1000; K=2;
Fs = N;
t = (1:K*N)/N;
b = 15.0;
c = 0.0;
chirp = cos(2*pi*(b*t.^2 + c * t));
% b = -15.0;
% c = 60.0;
% chirp = cos(2*pi*(b*t.^2 + c * t));
z = fft(chirp);
w = abs(fftshift(z));
w = 20*log10(w/max(w));

figure(12)
hold off
plot(t, chirp, 'linewidth', 1);
set(gca, 'ytick', [], 'ylim', [-1.2 1.2], 'xlim', [0 2], 'fontsize', 18, 'linewidth', 3);
set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
print('-depsc', 'Plots/chirp_spect_1.eps');

figure(15)
plot(-Fs/2:Fs/(K*N):Fs/2-Fs/(K*N), w, 'linewidth', 3);
set(gca, 'fontsize', 18, 'linewidth', 3);
set(gca, 'xlim', [0 200]);
ylabel('dB');
set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
print('-depsc', 'Plots/chirp_spect_5.eps');

figure(16)
s = -Fs/2+Fs/(K*N):Fs/(K*N):Fs/2;
result = calc_sinc(s);
hold off
plot(s, w, 'linewidth', 3);
hold on
plot(s, 20*log10(abs(result)), 'g');
% plot(s, 20*log10(abs(cos(s.^2/(8*pi*b) + pi/4))), 'm');
plot(s, 20*log10(abs(conv2(cos(s.^2/(pi*b) + pi/4), result, 'same'))), 'r');
set(gca, 'fontsize', 18, 'linewidth', 3);
% set(gca, 'xlim', [0 200]);
ylabel('dB');


% figure(121)
% hold off
% Timefreq = myspecgram(real(chirp)', N/2, Fs, ones(N/2,1));
% %% say yes to plot, and then plot max freq 100
% colormap('default')
% brighten(-0.2);
% shading faceted;
% hold on
% plot(t, (2*b*t + c), 'linewidth', 3);
% set(gca, 'fontsize', 18, 'linewidth', 3, 'xlim', [0 2]);
% ylabel('frequency');
% set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
% print('-depsc', 'Plots/chirp_spect_2.eps');
% YY = get(gca,'ylim');
% XX = get(gca,'xlim');

M = 200;
[STFT, X] = stft(chirp, M);
figure(221)
hold off
I = 20*log10(abs(fftshift(STFT,1)));
I = I - max(max(I));
imagesc(t(1:M:N*K)+M/(2*N), -Fs/2:Fs/M:Fs/2-Fs/M, I);
colormap('default')
cc = colorbar('EastOutside');
set(get(cc,'YLabel'),'String','dB','fontsize', 18);
hold on
plot(t, (2*b*t + c), 'linewidth', 3);
set(gca,'ylim', [0 100]);
set(gca,'ydir', 'normal');
set(gca,'xlim', [0 K]);
set(gca, 'fontsize', 18, 'linewidth', 3, 'xlim', [0 2]);
xlabel('time');
ylabel('frequency (Hz)');
title(sprintf('sample freq = 1000Hz, block size = %d', M));
set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
print('-depsc', 'Plots/chirp_spect_2.eps');

figure(321);
hold off
for i=1:N*K/M
  subplot(N*K/M,1,i);
  hold off
  plot(-Fs/2:Fs/M:Fs/2-Fs/M, 20*log10(abs(fftshift(STFT(:,i),1))));
  set(gca,'xlim', [0 100]);
  YY = get(gca, 'ylim');
  hold on
  plot([2*b*t(M*(i-1)+1) + c, 2*b*t(M*(i-1)+1) + c], YY, 'r')
  plot([2*b*t(M*i) + c,       2*b*t(M*i) + c], YY, 'r')
end


% figure(122)
% hold off
% Fs = N;
% Timefreq = myspecgram(real(chirp)', N/5, Fs, ones(N/5,1));
% %% say yes to plot, and then plot max freq 100
% colormap('default')
% brighten(-0.2);
% shading faceted;
% hold on
% plot(t, (2*b*t + c), 'linewidth', 3);
% set(gca, 'fontsize', 18, 'linewidth', 3);
% ylabel('frequency');
% set(gca, 'xlim', XX, 'ylim', YY);
% set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
% print('-depsc', 'Plots/chirp_spect_3.eps');

M = 100;
[STFT, X] = stft(chirp, M);
figure(222)
hold off
I = 20*log10(abs(fftshift(STFT,1)));
I = I - max(max(I));
imagesc(t(1:M:N*K)+M/(2*N), -Fs/2:Fs/M:Fs/2-Fs/M, I);
colormap('default')
cc = colorbar('EastOutside');
set(get(cc,'YLabel'),'String','dB','fontsize', 18);
hold on
plot(t, (2*b*t + c), 'linewidth', 3);
set(gca,'ylim', [0 100]);
set(gca,'ydir', 'normal');
set(gca,'xlim', [0 K]);
set(gca, 'fontsize', 18, 'linewidth', 3, 'xlim', [0 2]);
xlabel('time');
ylabel('frequency (Hz)');
title(sprintf('sample freq = 1000Hz, block size = %d', M));
set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
print('-depsc', 'Plots/chirp_spect_3.eps');

% figure(123)
% hold off
% Fs = N;
% Timefreq = myspecgram(real(chirp)', N/10, Fs, ones(N/10,1));
% %% say yes to plot, and then plot max freq 100
% colormap('default')
% brighten(-0.2);
% shading faceted;
% hold on
% plot(t, (2*b*t + c), 'linewidth', 3);
% % set(gca, 'xtick', [], 'ytick', []);
% set(gca, 'fontsize', 18, 'linewidth', 3);
% ylabel('frequency');
% set(gca, 'xlim', XX, 'ylim', YY);
% set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
% print('-depsc', 'Plots/chirp_spect_4.eps');

M = 400;
[STFT, X] = stft(chirp, M);
figure(223)
hold off
I = 20*log10(abs(fftshift(STFT,1)));
I = I - max(max(I));
imagesc(t(1:M:N*K)+M/(2*N), -Fs/2:Fs/M:Fs/2-Fs/M, I);
colormap('default')
cc = colorbar('EastOutside');
set(get(cc,'YLabel'),'String','dB','fontsize', 18);
hold on
plot(t, (2*b*t + c), 'linewidth', 3);
set(gca,'ylim', [0 100]);
set(gca,'ydir', 'normal');
set(gca,'xlim', [0 K]);
set(gca, 'fontsize', 18, 'linewidth', 3, 'xlim', [0 2]);
xlabel('time');
ylabel('frequency (Hz)');
title(sprintf('sample freq = 1000Hz, block size = %d', M));
set(gcf, 'PaperUnits', 'centimeters', 'PaperPosition', [0 0 32 16])
print('-depsc', 'Plots/chirp_spect_4.eps');
