要求

有一单频信号y(t)=sin(2πft)，其中f=100Hz，和f=500Hz，分别用DFT求y(t)的谱。
 抽样频率fs以不发生混叠
 记录长度T0取整数周期
 用stem语句绘出幅度谱，横坐标为模拟频率f，用plot画出时域波形

一、连续函数抽样

f=100
fs=1000;                
dt=1/fs;               
T=0:dt:0.05-dt;            
xn=sin(2*pi*f*T);
f=500
fs=8000;                
dt=1/fs;               
T=0:dt:0.01-dt;            
xn=sin(2*pi*f*T);

二、抽样长度确定

由奈奎斯特抽样定理，fs≥2f，
即f=100时，fs≥200，取fs=400；
f=500时，fs≥1000，取fs=2000

三、DFT

1、使用一次循环实现DFT

clear
f=100;                
t=(0:0.001:0.1);      
x=sin(2*pi*f*t);  

t1=clock;
fs=1000;                
dt=1/fs;               
T=0:dt:0.1-dt;            
xn=sin(2*pi*f*T);
N=length(xn);
WN=exp(-1i.*2.*pi./N);    
xk=zeros(1,N);
xk1=zeros(1,N);
k=[0:N-1];
 for n=0:1:N-1
   xk(n+1)=xn*WN.^(n.*k');  %此处下标一定得从1开始，因为matlab的下标是从1开始的
 end    
t2=clock;
etime(t2,t1)

subplot(2,1,1)
plot(T,xn)
subplot(2,1,2)
stem(k,abs(xk))

运行结果

2、使用矩阵实现DFT

function xk=dt_2(xn);
N=length(xn);

WN=exp(-j*2*pi/N);   
n=0:1:N-1;     %定义一个一维矩阵，即行向量,从0到N-1
k=0:1:N-1;
nk=k'*n;      %行向量k变换为列向量 乘上 行向量n ，得到一个N x N的矩阵
WNnk=WN.^(nk);      %做幂运算后的参数仍为一个 N x N的系数矩阵
xk=xn*WNnk;       %行向量 乘以  N x N的系数矩阵 即为DFT变换后的矩阵

总结

1、n的取值范围是[0，N-1]，有0时不应该包括N，否则会得到不一样的结果（如图）。可以看到因为取点不准确造成频谱泄露。

2、在连续信号进行傅里叶频域变换后，正弦信号的频域只有一个冲击，为w。
对离散信号而言，进行离散傅里叶变换以后，由离散傅里叶变换的共轭对称性可以推出时域实序列在频域为共轭对称序列，幅频函数为偶对称。同时可以将傅里叶频域转换到z域后用z反变换快速求出两个冲击信号和。