% AddZerosDFT


% x(n)=2*cos(2*pi*f1*n)+cos(2*pi*f2*n)

clc;
clear;
Repeat = 1;
disp('高密度譜VS高分辨譜')
disp('采用序列x(n)=2*cos(2*pi*f1*n)+cos(2*pi*f2*n)')
while Repeat == 1
    clc;
    clear;
    f1 = input('f1=');
    f2 = input('f2=');
    disp('為滿足Nyquist定理，采樣點數N>=(k1/f1,k2/f2)的最小公倍數');
    M = input('請輸入小于該值    :M=');
    N = input('請輸入大于等于該值:N=');
    
    n1 = 0:M-1;
    x1 = 2*cos(2*pi*f1*n1)+cos(2*pi*f2*n1);
    X1 = dft(x1,M);
    k1 = 0:M-1;
    w1 = 2*pi/M*k1;
    
    figure;
    subplot(2,3,1),stem(n1,x1),title('采樣點不足');
    subplot(2,3,4),stem(w1/pi,abs(X1)),title('頻譜');
    
    n2 = 0:N-1;
    x2(1:M) = x1;
    x2(M+1:N) = zeros(1,N-M);
    X2 = dft(x2,N);
    k2 = 0:N-1;
    w2 = 2*pi/N*k2;
    subplot(2,3,2),stem(n2,x2),title('補零');
    subplot(2,3,5),stem(w2/pi,abs(X2)),title('高密度頻譜');
    
    n3 = 0:N-1;
    x3 = 2*cos(2*pi*f1*n3)+cos(2*pi*f2*n3);
    X3 = dft(x3,N);
    k3 = 0:N-1;
    w3 = 2*pi/N*k3;
    
    subplot(2,3,3),stem(n3,x3),title('采樣點足夠');
    subplot(2,3,6),stem(w3/pi,abs(X3)),title('高分辨率頻譜');
    
    Repeat = input('Still want to repeat,Yes=1,No=0? Ans=');
end