%用一個N點的FFT計算一個2N點的實序列的DFT
clear    %必要的清除變量
N=32;
n0=0:1:31;
n=0:1:63;  %2N=64
x=n.^2;
%抽偶序列
for i=0:1:(N-1)
    x1(i+1)=x(2*i+1);
end
%抽奇序列
for i=0:1:(N-1)
    x2(i+1)=x(2*i+2);
end
%FFT
Xk1=fft(x1,N);
Xk2=fft(x2,N);
%
Yk=Xk1+j*Xk2;
Ykr=real(Yk);
Yki=imag(Yk);
Yk=Ykr+j*Yki;
%%%建立X1(k)、Yr(k)和X2(k)、Yi(k)之間的關系
Xk1(0+1)=[Ykr(1)+Ykr(1)]/2+j*[Yki(1)-Yki(1)]/2;
Xk2(0+1)=[Yki(1)+Yki(1)]/2-j*[Ykr(1)-Ykr(1)]/2;
for k=1:1:N-1
    Xk1(k+1)=[Ykr(k+1)+Ykr((N-k)+1)]/2+j*[Yki(k+1)-Yki((N-k)+1)]/2;    %%%%??疑問：若按課件上的關系Ykr(N-k)和Yki(N-k)
    Xk2(k+1)=[Yki(k+1)+Yki((N-k)+1)]/2-j*[Ykr(k+1)-Ykr((N-k)+1)]/2;    %%%%??對應k+1的話,則必然不對？但是前面又必然要加。
end

%%%%
for k=0:1:(N-1)
    Xk(k+1)=Xk1(k+1)+Xk2(k+1)*exp(-j*pi*k/N);     %前半部分
    Xk(k+1+N)=Xk1(k+1)-Xk2(k+1)*exp(-j*pi*k/N);   %后半部分
end
%%%%將前后兩部分組合起來得到2N點的FFT
Xk;      %%%%???好象不用組合了？
%%%
subplot(2,1,1)
plot(n,abs(Xk))
title('幅頻特性')
%%%驗證
Xk=fft(x,2*N);
subplot(2,1,2)
plot(n,abs(Xk))
title('檢驗')    %%%???結果不吻合