数字信号处理DFTMATLAB程序
DFT和DTFTmatlab程序

已知序列xn=[1,1,1,1],试用MATLAB编写程序,计算该序列的离散付里叶变换及逆离散付里叶变换。
(1)MATLAB程序:function xk=dft(xn,N) %dftn=[0:1:N-1];k=n;WN=exp(-j*2*pi/N); %旋转因子nk=n'*k;WNnk=WN.^nk;xk=xn*WNnk;function xn=idft(xk,N) %idftn=[0:1:N-1];k=n;WN=exp(-j*2*pi/N);nk=n'*k;WNnk=WN.^(-nk);xn=xk*WNnk/N;xn=[1,1,1,1]; %计算dftN=4;xk=dft(xn,N)'xk=[4,0,0,0]; %计算idftN=4;xn=idft(xk,N)'仿真结果:DFT:(2)MATLAB程序:xn=[1,1,1,1];N=length(xn);n=0:N-1;k=0:N-1;Xk=xn*exp(-j*2*pi/N).^(n'*k); %计算DFTx=(Xk*exp(j*2*pi/N).^(n'*k))/N; %计算IDFTsubplot(1,2,2);stem(k,abs(Xk));title('|X(k)|'); %画图axis([-1,N,1.1*min(abs(Xk)),1.1*max(abs(Xk))]); %加坐标subplot(1,2,1);stem(n,xn);title('x(n)');axis([-1,N,1.1*min(xn),1.1*max(xn)]);仿真结果:图3-1 序列x(n)及其DFT变换(3)MATLAB程序:xn=[1,1,1,1];N=length(xn);n=0:N-1;subplot(2,2,1);stem(n,xn);title('x(n)');k=0:N-1;Xk=fft(xn,N); %计算Xksubplot(2,1,2);stem(k,abs(Xk));title('Xk=DFT(xn)');xn1=ifft(Xk,N); %计算xnsubplot(2,2,2);stem(n,xn1);title('x(n)=IDFT(Xk)');仿真结果:图3-2 序列的正逆离散傅里叶变换取一个周期的正弦信号,作8点采样,求它的连续频谱。
数字信号处理: MATLABdft对称性验证以及应用

数字信号处理: MATLABdft对称性验证以及应用武汉理工大学《数字信号处理》课程设计说明书目录1 MATLAB基本操作及常用命令介绍 (1)1(1 MATLAB的启动 .....................................................................11(2桌面平台 ..................................................................... . (1)1.3 基本平面图形绘制命令plot (2)2 理论分析 ..................................................................... (3)2.1实验内容 ..................................................................... . (3)2.2序列对称性的理论验证 (3)3 程序验证 ............................................................................................. 4 4 结果分析 ..................................................................... ........................ 7 5 对称性的应用 ..................................................................... .. (10)5.1 FFT算法的基本思想 (10)5.2 对称性应用的程序实现 (11)6 心得体会 ..................................................................... ...................... 15 参考文献 ..................................................................... .. (16)武汉理工大学《数字信号处理》课程设计说明书1 MATLAB基本操作及常用命令介绍1(1 MATLAB的启动启动MATLAB有多种方式,最常用的方法就是双击系统桌面的MATLAB图标,也可以在开始菜单的程序选项中选择MATLAB快捷方式。
matlab实现DFT

DFT 基于Matlab 的实现一、实验目的1.掌握DFT 函数的用法。
2. 利用DFT 进行信号检测及谱分析。
3.了解信号截取长度对谱分析的影响。
二、实验内容1.利用DFT 计算信号功率谱。
实验程序:t=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t)+randn(1,length(t));Y=dft(x,512);P=Y.*conj(Y)/512;f=1000*(0:255)/512;plot(f,P(1:256))2. 进行信号检测。
分析信号频谱所对应频率轴的数字频率和频率之间的关系。
模拟信号)8cos(5)4sin(*2)(t t t x ππ+=,以n t 01.0= 10-≤≤N n 进行取样,求N 点DFT 的幅值谱。
实验程序:subplot(2,2,1)N=45;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=dft(x,N);plot(q,abs(y));title('DFT N=45')subplot(2,2,2)N=50;n=0:N-1;t=0.01*n; q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=dft(x,N);plot(q,abs(y));title('DFT N=50')subplot(2,2,3)N=55;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=dft(x,N);plot(q,abs(y));title('DFT N=55')subplot(2,2,4)N=60;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=dft(x,N);plot(q,abs(y));title('DFT N=60')3. 对2,进一步增加截取长度和DFT点数,如N加大到256,观察信号频谱的变化,分析产生这一变化的原因。
用FFT对信号作频谱分析Matlab程序

对以下序列进行FFT 分析x 1(n)=R 4(n)x 2(n)=x 3(n)=x1n=[ones(1,4)]; %产生R4(n)序列向量X1k8=fft(x1n,8); %计算x1n 的8点DFTX1k16=fft(x1n,16); %计算x1n 的16点DFT%以下绘制幅频特性曲线N=8;f=2/N*(0:N-1); (不懂)figure(1);subplot(1,2,1);stem(f,abs(X1k8),'r','、'); %绘制8点DFT 的幅频特性图,abs 求得Fourier 变换后的振幅title('(1a) 8点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');N=16;f=2/N*(0:N-1);subplot(1,2,2);stem(f,abs(X1k16),'、'); %绘制8点DFT 的幅频特性图title('(1b) 16点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');%x2n 与 x3nM=8;xa=1:(M/2); xb=(M/2):-1:1; %从M/2到1每次递减1x2n=[xa,xb]; %产生长度为8的三角波序列x2(n)x3n=[xb,xa];X2k8=fft(x2n,8);X2k16=fft(x2n,16);X3k8=fft(x3n,8);X3k16=fft(x3n,16);figure(2);N=8;f=2/N*(0:N-1);subplot(2,2,1);stem(f,abs(X2k8),'r','、'); %绘制8点DFT 的幅频特性图n+1 0≤n ≤3 8-n 4≤n ≤7 0 其它n 4-n 0≤n ≤3 n-3 4≤n ≤70 其它ntitle('(2a) 8点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,3);stem(f,abs(X3k8),'r','、'); %绘制8点DFT的幅频特性图title('(3a) 8点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');N=16;f=2/N*(0:N-1);subplot(2,2,2);stem(f,abs(X2k16),'、'); %绘制8点DFT的幅频特性图title('(2b) 16点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,4);stem(f,abs(X3k16),'、'); %绘制8点DFT的幅频特性图title('(3b) 16点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');%x4n 与 x5nN=8;n=0:N-1;x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k8=fft(x4n,8);X4k16=fft(x4n,16);X5k8=fft(x5n,8);X5k16=fft(x5n,16);figure(3);N=8;f=2/N*(0:N-1);subplot(2,2,1);stem(f,abs(X4k8),'r','、'); %绘制8点DFT的幅频特性图title('(4a) 8点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,3);stem(f,abs(X5k8),'r','、'); %绘制8点DFT的幅频特性图title('(5a) 8点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');N=16;f=2/N*(0:N-1);subplot(2,2,2);stem(f,abs(X4k16),'、'); %绘制8点DFT的幅频特性图title('(4b) 16点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,4);stem(f,abs(X5k16),'、'); %绘制8点DFT的幅频特性图title('(5b) 16点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');%x8nFs=64; T=1/Fs;N=16;n=0:N-1; %对于N=16的情况nT = n*T;x8n=cos(8*pi*nT)+cos(16*pi*nT)+cos(20*pi*nT)X8k16=fft(x8n,16);N=16;f=2/N*(0:N-1);figure(4);subplot(2,2,1);stem(f,abs(X8k16),'、'); %绘制8点DFT的幅频特性图title('(6a) 16点DFT[x_8(n)]');xlabel('ω/π');ylabel('幅度');N=32;n=0:N-1; %对于N=16的情况nT = n*T;x8n=cos(8*pi*nT)+cos(16*pi*nT)+cos(20*pi*nT)X8k32=fft(x8n,32);N=32;f=2/N*(0:N-1);subplot(2,2,2);stem(f,abs(X8k32),'、'); %绘制8点DFT的幅频特性图title('(6b) 32点DFT[x_8(n)]');xlabel('ω/π');ylabel('幅度'); N=64;n=0:N-1; %对于N=16的情况nT = n*T;x8n=cos(8*pi*nT)+cos(16*pi*nT)+cos(20*pi*nT)X8k64=fft(x8n,64);N=64;f=2/N*(0:N-1);subplot(2,2,3);stem(f,abs(X8k64),'、'); %绘制8点DFT的幅频特性图title('(6c) 64点DFT[x_8(n)]');xlabel('ω/π');ylabel('幅度');。
matlab中实现dft算法代码

matlab中实现dft算法代码
在MATLAB中实现DFT(离散傅里叶变换)算法的代码如下:
```matlab
function X = myDFT(x)
N = length(x); % 输入信号的长度
X = zeros(1, N); % 存储DFT结果的数组
for k = 0:N-1
for n = 0:N-1
X(k+1) = X(k+1) + x(n+1) * exp(-1i*2*pi*k*n/N);
end
end
end
```
在这段代码中,`x`是输入信号的数组,`N`是输入信号的长度,`X`是用于存储DFT结果的数组。
通过双重循环计算每个频率点的复数值,然后将其存储在数组`X`中。
最后,函数返回计算得到的DFT结果数组`X`。
要使用这个函数进行DFT计算,可以按照以下步骤:
```matlab
x = [1, 2, 3, 4]; % 输入信号
X = myDFT(x); % 调用自定义的DFT函数进行计算
disp(X); % 显示DFT结果
```
在这个例子中,输入信号`x`是一个包含了[1, 2, 3, 4]的数组。
然后,通过调用`myDFT`函数计算DFT结果,并将结果存储在`X`中。
最后,通过使用`disp`函数来显示计算得到的DFT结果`X`。
需要注意的是,这只是一个简单的DFT算法实现代码,可能没有考虑到性能优化和其他复杂情况。
在实际应用中,可以使用MATLAB内置的`fft`函数来进行更高效和准确的DFT计算。
数字信号处理DFTMATLAB程序

实验三 频域信号处理1. 实验目的(1) 学习信号DFT 变换的matlab 实现; (2) 学习fft 的matlab 实现; (3) 验证DFT 的相关性质。
2. 思考题(1) 若()()()sin sin 4x n n n ππ=+是一个128点的有限长序列,求其128点DFT 结果; 程序如下:求DFT 变换矩阵A : clc; clear; N=128;A=dftmtx(N)Ai=conj(dftmtx(N)); n=0:(N-1); k=0:(N-1); nk=n'*k;Wn=(sin(pi/8)+sin(pi/4)).^nk Wk=conj(Wn)/N;求128点的DFT (分别用FFT 函数和dftmtx 函数)clc; clear; N=128; n=0:N-1;x=sin(pi/8*n)+sin(pi/4*n); subplot(3,1,1) plot(n,x); grid on title('原图') y1=fft(x,N); A=dftmtx(N); y2=x(1:N)*A; subplot(3,1,2) plot(n,y1) grid on title('FFT')subplot(3,1,3) plot(n,y2) grid ontitle('dftmtx')程序运行结果如图1所示:图 1(2) 对模拟信号()()()2sin 45sin 8x t t t ππ=+,以0.01t n =,()0:1n N =-进行采样,求a ) N=40点的FFT 幅度谱,从图中能否观察出两个频谱分量;b ) 提高采样点数值N=128,再求该信号的幅度频谱,此时幅度频谱发生了什么变化?信号的两个模拟频率和数字频率分别为多少?FFT 频谱分析结果和理论上是否一致? 程序如下:clc; clear; N=40; n=0:N-1; t=0.01*nx=2*sin(4*pi*t)+5*sin(8*pi*t); subplot(2,1,1)plot(x(1:N)) grid on title('原图') y1=fft(x,N); subplot(2,1,2) plot(n,y1(1:N))原图-13FFT-13dftmtxgrid ontitle('40点采样')结果如图2:原图40点采样图2 clc;clear;N=128;n=0:N-1;t=0.01*nx=2*sin(4*pi*t)+5*sin(8*pi*t);subplot(2,1,1)plot(x(1:N))title('原图')y1=fft(x,N);subplot(2,1,2)plot(n,y1(1:N))title('128点采样')结果如图3:图 3从图2,3中可以看出:N=40和128点的FFT 幅度谱中均出现了两个频谱分量.对于128点采样,11128,128f T f ===,模拟频率分别为124,8ππΩ=Ω=,那么,数字频率分别为1122114,81283212816T T ππωπωπ=Ω=⨯==Ω=⨯=。
(完整word版)用matlab实现DFTFFT

用matlab实现DFT FFT目录实验目的 (2)实验内容 (2)1.用MATLAB实现DFT (2)2.用MATLAB实现FFT,分析有限离散序列的FFT (3)3.通过分别计算时间,得出DFT与FFT的算法差异 (7)实验原理 (8)1. 离散傅里叶变换的快速算法FFT (8)2. FFT提高运算速度的原理 (9)3. 理论分析DFT与FFT算法差异 (11)实验步骤 (12)实验结果 (13)实验分析 (27)实验结论 (33)实验体会 (33)实验目的1.通过研究DFT,FFT性质,用语言实现DFT, FFT。
不使用MATLAB现有的FFT函数,自己编写具体算法。
2.掌握FFT基2时间抽选法,理解其提高减少乘法运算次数提高运算速度的原理。
3.设计实验,得出DFT和FFT算法差异的证明,如复杂度等(精度、不同长度的序列等)。
实验内容1. 用MATLAB实现DFTN点序列x(n) 的DFT为:DFT的矩阵为:根据DFT公式与矩阵展开,通过MATLAB实现DFT:2.用Matlab实现FFT编程思想及程序框图:●原位计算因为DIT-FFT与DIF-FFT的算法类似,这里我们以DIT-FFT为例。
N=2M点的FFT共进行M级运算,且每一级都由N/2个蝶形运算组成,后一级的节点数据由前一级同处一条水平线位置的节点数据产生,所以我们同样可以将后一级的节点数据储存到前一级的节点中,这样的方法叫做原位计算,它大大节省了内存资源,降低了成本,简化了运算。
●序列的倒序无论是进行DIT-FFT还是DIF-FFT都需要进行倒序,包括输入倒序与输出倒序,以一定的方式将数组进行重新排列。
倒序的方法:首先由于N=2M,我们就可以用M位二进制数来表示节点的顺序,并且按照奇偶时域抽取。
然后,如图1所示,第一次按最低位n0的0、1值分解为奇偶组,第二次按次低位n1的0、1值分解为奇偶组,以此类推。
最后,所得二进制数所对应的十进制数即为序列倒序后产生的序列。
数字信号处理基于matlab(用DFT作谱分析,窗函数的设计)

实验一:用DFT作谱分析x1=[1 1 1 1];x2=[1 2 3 4 4 3 2 1];n1=0:8;x3=cos(n1*pi/4);n2=0:8;x4=sin(n2*pi/8);figuresubplot(2,2,1);stem(x1);subplot(2,2,2);stem(x2);subplot(2,2,3);stem(x3);subplot(2,2,4);stem(x4);N1=8;F1x1=fft(x1,N1);F1x2=fft(x2,N1);F1x3=fft(x3,N1);F1x4=fft(x4,N1);figuresubplot(2,2,1);stem(abs(F1x1));subplot(2,2,2);stem(abs(F1x2));subplot(2,2,3);stem(abs(F1x3));subplot(2,2,4);stem(abs(F1x4));N3=256;F3x1=fft(x1,N3);F3x2=fft(x2,N3);F3x3=fft(x3,N3);F3x4=fft(x4,N3);w=2/N3*[0:N3-1];figuresubplot(2,2,1);plot(w,abs(F3x1));subplot(2,2,2);plot(w,abs(F3x2));subplot(2,2,3);plot(w,abs(F3x3));subplot(2,2,4);plot(w,abs(F3x4));N6=64;t=0:1/64:1/64*(N6-1);x6=cos(8*pi*t)+cos(16*pi*t)+cos(20*pi*t);figure;plot(x6)Fx6=fft(x6,N6);N=64w6=2/N6*[0:N6-1];figure;plot(w6,Fx6)实验二:用双线性法设计IIR数字滤波器clear all;close all;clcT=1;Wp=2/T*tan(0.2*pi/2);Ws=2/T*tan(0.3*pi/2);Rp=1;Rs=15;[N,Wc]=buttord(Wp,Ws,Rp,Rs,'s');[B,A]=butter(N,Wc,'s');[C,D]=bilinear(B,A,1/T);W=0:0.001*pi:0.5*pi;H=freqs(B,A,W);Hd=freqz(C,D,W);figuresubplot(3,1,1);plot(W/pi,abs(H))grid ontitle('模拟巴特沃斯滤波器')xlabel('Frequency/Hz')ylabel('Magnitude')subplot(3,1,2);plot(W/pi,abs(Hd))grid ontitle('数字巴特沃斯滤波器')xlabel('Didital Frequency/pi')ylabel('Magnitude')Hd_db=-20*log(abs(Hd(1)./(abs(Hd)+eps)));subplot(3,1,3);plot(W/pi,Hd_db)grid ontitle('数字巴特沃斯滤波器波特图')xlabel('Didital Frequency/pi')ylabel('bd_Magnitude')x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4.8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];M=512;fx=fft(x,M);y=filter(C,D,x);fy=fft(y,M);f=1/512*[0:511]*250figuresubplot(2,1,1);plot(x);subplot(2,1,2);plot(f,abs(fx));figuresubplot(2,1,1);plot(y);subplot(2,1,2);plot(f,abs(fy));心电图谱分析1、滤波前波形和频谱2、滤波后波形和频谱实验三:用窗函数法设计FIR数字滤波器一、N=33,用四种窗设计滤波器function hd=ideal(N,wc)for n=0:N-1if n==(N-1)/2hd(n+1)=wc/pi;else hd(n+1)=sin(wc*(n-(N-1)/2))/(pi*(n-(N-1)/2));endendclear all;close all;clcN=33;wc=pi/4;hd=ideal(N,wc);w1=boxcar(N);w2=hamming(N);w3=hann(N);w4=blackman(N);h1=hd.*w1';h2=hd.*w2';h3=hd.*w3';h4=hd.*w4';M=512;fh1=fft(h1,M);db1=-20*log10(abs(fh1(1)./(abs(fh1)+eps)));fh2=fft(h2,M);db2=-20*log10(abs(fh2(1)./(abs(fh2)+eps)));fh3=fft(h3,M);db3=-20*log10(abs(fh3(1)./(abs(fh3)+eps)));fh4=fft(h4,M);db4=-20*log10(abs(fh4(1)./(abs(fh4)+eps)));w=2/M*[0:M-1];figure;subplot(2,2,1);stem(h1);grid on;title('矩形窗');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh1));grid on;title('矩形窗');xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db1);grid on;title('矩形窗');xlabel('w');ylabel('20lg(H(k)/H(0))'); subplot(2,2,4);plot(w,angle(fh1)); grid on;title('矩形窗');xlabel('w');ylabel('相位特性');figure;subplot(2,2,1);stem(h2);grid on;title('汉宁窗');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh2)); grid on;title('汉宁窗');xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db2);grid on;title('汉宁窗');xlabel('w');ylabel('20lg(H(k)/H(0))'); subplot(2,2,4);plot(w,angle(fh2)); grid on;title('汉宁窗');xlabel('w');ylabel('相位特性');figure;subplot(2,2,1);stem(h3);grid on;title('海明窗');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh3)); grid on;title('海明窗');xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db3);grid on;title('海明窗');xlabel('w');ylabel('20lg(H(k)/H(0))'); subplot(2,2,4);plot(w,angle(fh3)); grid on;title('海明窗');xlabel('w');ylabel('相位特性');figure;subplot(2,2,1);stem(h4);grid on;title('布莱克曼窗');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh4)); grid on;title('布莱克曼窗');xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db4);grid on;title('布莱克曼窗');xlabel('w');ylabel('20lg(H(k)/H(0))'); subplot(2,2,4);plot(w,angle(fh4)); grid on;title('布莱克曼窗');xlabel('w');ylabel('相位特性');二,用一种窗采用N=15或N=33、w=pi/4设计滤波器clear all;close all;clcN=15;wc=pi/4;hd=ideal(N,wc);w1=blackman(N);h1=hd.*w1';M=512;fh1=fft(h1,M);db1=-20*log10(abs(fh1(1)./(abs(fh1)+eps)));w=2/M*[0:M-1];figure;subplot(2,2,1);stem(h1);grid on;title('布莱克曼窗(N=15)');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh1));grid on;title('布莱克曼窗(N=15)');xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db1);grid on;title('布莱克曼窗(N=15)');xlabel('w');ylabel('20lg(H(k)/H(0))');subplot(2,2,4);plot(w,angle(fh1));grid on;title('布莱克曼窗(N=15)');xlabel('w');ylabel('相位特性');N=33;wc=pi/4;hd=ideal(N,wc);w2=blackman(N);h2=hd.*w2';M=512;fh2=fft(h2,M);db2=-20*log10(abs(fh2(1)./(abs(fh2)+eps)));w=2/M*[0:M-1];figure;subplot(2,2,1);stem(h2);grid on;title('布莱克曼窗(N=33)');xlabel('n');ylabel('h(n)');subplot(2,2,2);plot(w,abs(fh2)); grid on;title('布莱克曼窗(N=33)'); xlabel('w');ylabel('H(k)');subplot(2,2,3);plot(w,db2);grid on;title('布莱克曼窗(N=33)'); xlabel('w');ylabel('20lg(H(k)/H(0))'); subplot(2,2,4);plot(w,angle(fh2)); grid on;title('布莱克曼窗(N=33)'); xlabel('w');ylabel('相位特性');三、提取50HZ基频信号N=512;t=0:1/512:1/512*(N-1);x=sin(100*pi*t)+sin(200*pi*t)+sin(300*pi*t); Fx=fft(x,N);w3=2/N*[0:N-1];figure;subplot(2,1,1);plot(x);grid on;title('抽样信号');xlabel('n');ylabel('y(n)');subplot(2,1,2);plot(w3,abs(Fx));grid on;title('抽样信号频谱');xlabel('pi');ylabel('Magnitude');N=33;wc=0.4;hd=ideal(N,wc);w1=boxcar(N);w2=blackman(N);h1=hd.*w1';y1=conv(x,h1);h2=hd.*w2';y2=conv(x,h2); figure;subplot(2,1,1);plot(y1(1:64)); grid on;title('矩形窗');xlabel('n');ylabel('y(n)');subplot(2,1,2);plot(y2(1:64)); grid on;title('布莱克曼窗');xlabel('n');ylabel('y(n)');Fy1=fft(y1,N);Fy2=fft(y2,N);figure;subplot(2,1,1);plot(abs(Fy1)); subplot(2,1,2);plot(abs(Fy2));。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
页脚内容1
实验三 频域信号处理
1. 实验目的
(1) 学习信号DFT 变换的matlab 实现;
(2) 学习fft 的matlab 实现;
(3) 验证DFT 的相关性质。
2. 思考题
(1) 若()()()sin sin 4x n n n ππ=+是一个128点的有限长序列,求其128点DFT 结果; 程序如下:
求DFT 变换矩阵A :
clc;
clear;
N=128;
A=dftmtx(N)
Ai=conj(dftmtx(N));
n=0:(N-1);
k=0:(N-1);
nk=n'*k;
Wn=(sin(pi/8)+sin(pi/4)).^nk
Wk=conj(Wn)/N;
求128点的DFT(分别用FFT函数和dftmtx函数)
clc;
clear;
N=128;
n=0:N-1;
x=sin(pi/8*n)+sin(pi/4*n);
subplot(3,1,1)
plot(n,x);
grid on
title('原图')
y1=fft(x,N);
A=dftmtx(N);
页脚内容2
y2=x(1:N)*A;
subplot(3,1,2)
plot(n,y1)
grid on
title('FFT')
subplot(3,1,3)
plot(n,y2)
grid on
title('dftmtx')
程序运行结果如图1所示:
原图
-13FFT
-13dftmtx
页脚内容3
页脚内容4
图 1
(2) 对模拟信号()()()2sin 45sin 8x t t t ππ=+,以0.01t n =,()0:1n N =-进行采样,求
a ) N =40点的FFT 幅度谱,从图中能否观察出两个频谱分量;
b ) 提高采样点数值N=128,再求该信号的幅度频谱,此时幅度频谱发生了什么变化?信号的两个模拟频率和数字频率分别为多少?FFT 频谱分析结果和理论上是否一致?
程序如下:
clc;
clear;
N=40;
n=0:N-1;
t=0.01*n
x=2*sin(4*pi*t)+5*sin(8*pi*t);
subplot(2,1,1)
plot(x(1:N))
grid on
title('原图')
y1=fft(x,N);
subplot(2,1,2)
plot(n,y1(1:N))
grid on
title('40点采样')
结果如图2:
原图
40点采样
图2 clc;
clear;
N=128;
n=0:N-1;
页脚内容5
页脚内容6 t=0.01*n
x=2*sin(4*pi*t)+5*sin(8*pi*t);
subplot(2,1,1)
plot(x(1:N))
title('原图')
y1=fft(x,N);
subplot(2,1,2)
plot(n,y1(1:N))
title('128点采样')
结果如图3:
原
图
128点采样
页脚内容7
图 3
从图2,3中可以看出:N=40和128点的FFT 幅度谱中均出现了两个频谱分量.对于128点采样,11128,128
f T f ===,模拟频率分别为124,8ππΩ=Ω=,那么,数字频率分别为1122114,81283212816
T T ππωπωπ=Ω=⨯==Ω=⨯=。
(3) 输入为()2sin(/4)()x n n u n π=,020n ≤≤
系统采样响应为h(n)=δ(n)+2.5δ(n -1)+2.5δ(n -2)+δ(n -3) 010n ≤≤
验证 请用时域卷积及频域两种方法计算系统输出响应,并验证20点的循环卷积定理 程序如下:
clc;
clear;
N=20;
n1=0:19;
x=2*sin(n1*pi/4);
X1=fft(x,N)
n2=0:9;
h=zeros(1,length(n2));
h(1)=1;
h(2)=2.5;
h(3)=2.5;
h(4)=1;
X2=fft(h,N)
y1=conv(x,h);
y2=circonv(x,h,20)
subplot(2,1,1)
X=fft(y2,N)
stem(X);
title('X=fft(y2,N)')
xlabel('n')
ylabel('X')
subplot(2,1,2)
XK=X1.*X2
stem(XK);
title('XK=X1.*X2')
页脚内容8
页脚内容9 xlabel('n')
ylabel('XK')
suptitle('时域卷积')
X =fft(y2,N)n X 时域卷积
X K=X 1.*X 2n X K 图 4 clc;
clear;
N=20;
n1=0:19;
x=2*sin(n1*pi/4);
X1=fft(x,N)
n2=0:19;
h=zeros(1,length(n2));
h(1)=1;
h(2)=2.5;
h(3)=2.5;
h(4)=1;
X2=fft(h,N)
subplot(2,1,1)
X=fft(x.*h,N)
stem(X)
title('X=fft(x.*h,N)')
xlabel('n')
ylabel('X')
subplot(2,1,2)
XK=1/N*circonv(X1,X2,20)
stem(XK)
title('XK=1/N*circonv(X1,X2,20)')
xlabel('n')
页脚内容10
页脚内容11
ylabel('XK')
suptitle('频域卷积')
X =fft(x.*h,N)
n X
频域卷积
X K=1/N*circonv(X 1,X 2,20)n X K
图 5
(4) 若()
1944j n x n e π⎛⎫+ ⎪⎝⎭=,请计算实部与虚部对应的频域,并验证共轭对称性。
程序如下:
clc;
clear;
N=512;
e=2.73;
n=0:N-1;
x=4.*e.^(n/9.0).*(cos(pi/4*n)+i*sin(pi/4*n));
yx=fft(x,N);
xr=4.*e.^(n/9.0).*cos(pi/4*n);
xi=4.*e.^(n/9.0).*sin(pi/4*n);
yr=fft(xr,N);
yi=fft(i*xi,N);
subplot(2,1,1)
plot(n,yx)
grid on
title('x(n)的FFT')
subplot(2,1,2)
plot(n,yr+yi)
grid on
title('x(n)实部的FFT与虚部的FFT之和')
页脚内容12
26x(n)的FFT
26x(n)实部的FFT与虚部的FFT之和
图6
页脚内容13。