%  ex_04_10_1.m :  習題 4.10

clear all;
% Specifications: sample frequency and data length
fs = 40;
N = 64;
t_int = 0:1/fs:(1/fs * (N - 1));
 
% Frequency components of the continuous time signal
f1 = 10.8;
f2 = 11.75;
f3 = 12.55;
% Calculate the discrete time signal - x(n)
x = sin(2 * pi * f1 * t_int) + sin(2 * pi * f2 * t_int) + sin(2 * pi * f3 * t_int);
figure('Color', 'w');
plot(x);grid on;

X=fft(x);
figure('Color', 'w');
f=-20:40/N:20-20/N;
subplot(221);stem(f,fftshift(abs(X)),'.');grid on; 
x3N = [x, zeros(1, 3 * N)];
X=fft(x3N);
f=-20:40/(4*N):20-20/(4*N);
subplot(222);
plot(f,fftshift(abs(X)));grid on; 
x7N = [x, zeros(1, 7 * N)];
X=fft(x7N);
f=-20:40/(8*N):20-20/(8*N);
subplot(223);
plot(f,fftshift(abs(X)));grid on; 
x15N = [x, zeros(1, 15 * N)];
X=fft(x15N);
f=-20:40/(16*N):20-20/(16*N);
subplot(224);
plot(f,fftshift(abs(X)));grid on; 
 
% Specifications of CZT: length, starting point and step
M = 50;
W = exp(- i * 2 * pi * 0.2 / fs);
A = exp(i * 2 * pi * 9 / fs);
x_czt = czt(x, M, W, A);f=9:0.2:9+0.2*(M-1);
figure('Color', 'w');
subplot(221);
plot(f,abs(x_czt));grid on;


M = 60;
W = exp(- i * 2 * pi * 0.12 / fs);
A = exp(i * 2 * pi * 8/ fs);
x_czt = czt(x, M, W, A);f=8:0.12:8+0.12*(M-1);
subplot(222);
plot(f,abs(x_czt));grid on;