离散傅立叶变换算法
实验项目名称 实验一:离散傅立叶变换的性质及应用
实验目的 1.了解DFT性质及应用。
2.熟悉MATLAB编程的特点。
实验要求 能够用MATLAB来进行数字信号的处理。
实验原理 1、DFT变换
正变换:
反变换:
2、序列卷积
设序列 的长度为N,序列 的长度为M。则分别对两个序列作 点的DFT得到 和 ,则两序列的线性卷积 等于 。即时域卷积频域为相乘关系。
实验仪器 计算机一台;
MATLAB软件
实验步骤 1、 用三种不同的DFT程序计算 的傅立叶变换 ,并比较三种程序计算机的运行时间。
(1)编制用for循环语句的M函数文件dft1.m,用循环变量逐点计算 ;
(2)编写用MATLAB矩阵运算的M函数文件dft2.m,完成下列矩阵运算:
(3)调用FFT库函数,直接计算 ;
(4) 分别利用上述三种不同方式编写的DFT程序计算序列 的傅立叶变换 ,并画出相应的幅频和相频特性,再比较各个程序的计算机运行时间。
2、利用DFT实现两序列的卷积运算,并研究DFT点数与混叠的关系。
给定 , 。用FFT和IFFT分别求线性卷积和混叠结果输出,并用函数stem(n,y)画出相应图形。选择不同的DFT点数进行对比,观察其混叠效应。
3、研究高密度频谱与高分辨率频谱
设有连续信号
以采样频率 对该信号采样,分析下列三种情况的幅频特性。
(1)采集数据长度 点,做 点的DFT,并画出幅频特性。
(2)采集数据长度 点,补零到256点的DFT,并画出幅频特性。
(3)采集数据长度 点,做 点的DFT,并画出幅频特性。
观察三幅不同频率特性图,分析和比较它们的特点以及形成的原因。
4、实现序列的内插和抽取所对应的傅立叶变换。
给定序列 ,做128点的傅立叶变换,并求
和 为整数
对应的傅立叶变换(128点)。比较这三个计算结果得到的幅频特性图,分析其差别产生的原因。选择不同的插值倍数和抽样倍数对比其幅频的变化
实验内容 1.(1)M函数文件dft1.m
function y=dft1(x)
N=length(x);
for k=1:N A=0;
for n=1:N
A=A+x(n)*exp((-j*2*pi/N)*(k-1)*(n-1));
y(k)=A;
end
end
end
(2)M函数文件dft2.m
function y=dft2(x)
N=length(x);
x1=x';
for r=1:N
for c=1:N
g(r,c)=exp((-j*2*pi/N)*(r-1)*(c-1));
end
end
y=g*x1;
end
2.
取23点的DFT
>> x=0:15;h=ones(1,8); y1=fft(x,23); y2=fft(h,23); y3=y1.*y2; y4=ifft(y3,23)
取24点的DFT
>> y1=fft(x,24); y2=fft(h,24); y3=y1.*y2; y4=ifft(y3,24)
取18点的DFT
>> y1=fft(x,18); y2=fft(h,18); y3=y1.*y2; y4=ifft(y3,18)
3
>> fs=32000;Ndata=16;N=16;
>> n=0:Ndata-1;t=n/fs;
>> x=cos(2*pi*6500*t)+cos(2*pi*9000*t);
>> y=fft(x,N);
>> A=abs(y);
>> f=(0:N-1)*fs/N;%真实的频率
>> subplot(3,1,1),plot(f(1:N/2),A(1:N/2)*2/N);
>> xlabel('频率/HZ'),ylabel('振幅'),grid on
>> title('数据16点 DFT 16点');
其他点数只需将Ndata与N修改即可。
4.>> n=0:127;
>> x=cos((pi/36)*n)+cos((1.5*pi/36)*n);
>> x1=chouqu(x,3)
>> x2=neicha(x,3);
>> y=fft(x,128);y1=fft(x1,128);y2=fft(x2,128);
>> subplot(3,1,1),stem(n,abs(y)),title('x(n)幅频特性')
>> subplot(3,1,2),stem(n,abs(y1)),title('x(3n)幅频特性')
>> subplot(3,1,3),stem(n,abs(y2)),title('x(n/3)幅频特性')
>> grid on
>> subplot(3,1,1),grid on
>> subplot(3,1,2),grid on
序列抽取函数chouqu.m
function y=chouqu(x,k)
n=length(x);
m=((n-1)-mod((n-1),k))/k;
for j=0:m
y(j+1)=x(k*j+1);
end
end
序列内插函数neicha.m
function y=neicha(x,k)
n=length(x);
751com.cn
((m/3)+1);
else y(m+1)=0;
end
end
end
为序列抽取和内插函数。174