数字信号实验报告(全)

数字信号处理实验报告实验一：用FFT做谱分析一、实验目的1、进一步加深DFT算法原理和基本性质的理解。2、熟悉FFT算法原理和FFT子程序的应用。3、学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。二、实验原理用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是2π/N≤D。可以根据此时选择FFT的变换区间N。误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。对模拟信号的频谱时，首先要按照采样定理将其变成时域离散信号。如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。三、实验内容和步骤对以下典型信号进行谱分析：其它nnnnnnx其它nnnnnnxnRnx,074,330,4)(,074,830,1)()()(32414()cos4xnn5()cos(/4)cos(/8)xnnn6()cos8cos16cos20xtttt对于以上信号，x1(n)~x5(n)选择FFT的变换区间N为8和16两种情况进行频谱分析。分别打印其幅频特性曲线。并进行对比、分析和讨论;；x6(t)为模拟周期信号，选择采样频率HzFs64，变换区间N=16,32,64三种情况进行谱分析。分别打印其幅频特性，并进行分析和讨论。四、实验代码如下所示：clear;%清除变量closeall;%关闭全部绘图窗口b=menu('请选择信号x1(n)--x6(n)','x1(n)=R4(n)','x2(n)=[12344321]',...'x3(n)=[43211234]','x4(n)=cos(npi/4)','x5(n)=sin(npi/8)',...'x6(n)=cos(8pit)+cos(16pit)+cos(20pit)','Exit');i=0;A=[8,16,32,64];while(b~=7)%当选择EXIT时，返回值7，则退出循环ifb==6m=menu('请选择FFT变换区间长度N','N=16','N=32','N=64');N=A(m+1);fs=64;n=0:(N-1);x=cos(8*pi*n/fs)+cos(16*pi*n/fs)+cos(20*pi*n/fs);elsem=menu('请选择FFT变换区间长度N','N=8','N=16','N=32');N=A(m);n=0:(N-1);ifb==1x=[1,1,1,1,0,0,0,0,zeros(1,N-8)];elseifb==2x=[1,2,3,4,4,3,2,1,zeros(1,N-8)];elseifb==3x=[4,3,2,1,1,2,3,4,zeros(1,N-8)];elseifb==4x=cos(n*pi/4);elseifb==5x=sin(n*pi/8);endend%先画出信号源图i=i+1;figure(i);%创建绘图窗口subplot(2,2,1);%指定1号子图xlabel('n');%标记X坐标stem(n,x,'.r');ylabel('x(n)');title(['x',num2str(b),'(n)的波形']);%进行FFTf=fft(x,N);%再画出FFT波形subplot(2,2,3);stem(n,abs(f),'.b');xlabel('k');ylabel('|X(k)|');title(['x',num2str(b),'(n)的N=',num2str(N),'点FFT']);b=menu('请选择信号x1(n)--x6(n)','x1(n)=R4(n)','x2(n)=[12344321]',...'x3(n)=[43211234]','x4(n)=cos(npi/4)','x5(n)=sin(npi/8)',...'x6(n)=cos(8pit)+cos(16pit)+cos(20pit)','Exit');endcloseall;五、实验结果图及分析1、实验结果图分析：（1）x1(n)的波形如图1-1、图1-2和图1-3所示，由3张图可知道,N值越大，频率分辨率越高。（2）x2(n)的波形如图1-4、图1-5和图1-6所示，由3张图可知道,N值越大，频率分辨率越高。（3）x3(n)的波形如图1-7、图1-8和图1-9所示，由3张图可知道,N值越大，频率分辨率越高。（4）x4(n)的波形如图1-10、图1-11和图1-12所示。根据参数可得出X4(t)的频率f=8Hz，当N=8、16、32时，频率分辨率为F0=fs/N=8Hz、4Hz、2Hz，因此在FFT图里分别在N=1、2、4有高幅值,因为截取的为周期序列的整数倍，所以所得出的谱正确。（5）x5(n)的波形如图1-13、图1-14和图1-15所示。根据参数可得出X5(t)的频率f=4Hz。当N=8时，频率分辨率F0=fs/N=8Hz，因为截取的不是为周期序列的整数倍，而且频率分辨率不够，所得出的谱有较大的误差，所以FFT图包含一些频率分量，不能清楚看清原信号的频率f。当N=16及32时，频率分辨率F0=fs/N=4Hz、2Hz，因此在FFT图里在N=1、2有高幅值,因为截取的为周期序列的整数倍，所以所得出的谱正确。（6）x6(n)的波形如图1-16、图1-17和图1-18所示。根据参数可得出X6(t)里包含3个频率，分别为f1=4，f2=8，f3=10。当N=16,频率分辨率F0=fs/N=4Hz，因为截取的不是x6里各周期序列的整数倍，所得出的谱有频谱泄漏，FFT图里可以看出信号cos(8pit)和cos(16pit)的频率f1=4,f2=8（在点N=1，2处有较大的幅值），而且频率分辨率不够高，不能分辨开第三个信号cos(20pit)的频率f3。当N=32，频率分辨率F0=fs/N=2Hz，因此在FFT图里的点N=2有高幅值，在N=4有高幅值N=5也有高幅值。因为截取的为周期序列的整数倍，所以所得出的谱正确。当N=64，频率分辨率F0=fs/N=1Hz，因此在FFT图里的点N=4有高幅值，在N=8有高幅值N=10也有高幅值。因为截取的为周期序列的整数倍，所以所得出的谱正确。变换区间N=64时频谱幅度是变换区间N=32时2倍，这种结果正好验证了用FFT对中期序列谱分析的理论。2、误差分析误差产生的原因：（1）对周期序列的截取不当，造成频谱泄漏（2）抽样点数N太少，频率分辨率不够用FFT做谱分析时参数的选择:（1）抽样频率要满足奈奎斯特准则，不小于信号最高频率的2倍（2）在抽样频率一定的情况下，抽样点数N要适当。太小会造成频率分辨力不够，太大会造成数据冗余。对周期序列，最好截取周期的整数倍进行谱分析图1-1：x1(n)的波形:N=8图1-2：x1(n)的波形N=16图1-3：x1(n)的波形N=32图1-4：x2(n)的波形N=8图1-5：x2(n)的波形:N=16图1-6：x2(n)的波形N=32图1-7：x3(n)的波形N=8图1-8：x3(n)的波形N=16图1-9：x3(n)的波形N=32图1-10：x4(n)的波形N=8图1-11：x4(n)的波形N=16图1-12：x4(n)的波形N=32图1-13：x5(n)的波形N=8图1-14：x5(n)的波形N=16图1-15：x5(n)的波形N=32图1-16：x6(n)的波形N=16图1-17：x6(n)的波形N=32图1-18：x6(n)的波形N=64实验二：用双线性变换法设计IIR数字滤波器一、实验目的1、熟悉用双线性变换法设计IIR数字滤波器的原理与方法。2、掌握数字滤波器的计算机仿真方法。二、实验内容及步骤1、用双线性变换法设计一个butterworth低通IIR数字滤波器。设计指标参数为：在通带内频率低于0.2时，最大衰减小于1dB，在阻带内[0.3,]频率区间上，最小衰减大于15dB。2、打印出数字滤波器在频率区间[0,]上的幅频响应特性曲线。3、运用MATLAB产生两个正弦信号，信号频率为50Hz和400Hz，采样频率为1000Hz。两个正弦信号相叠加为输入信号y(t)。设计一滤波器，保留源信号中50Hz的低频信号，对y(t)信号进行滤波。观察滤波前后信号的频谱特性，评价滤波器效果。三、实验步骤1、双线性变换法设计butterworth低通IIR数字滤波器复习有关butterworth模拟滤波器设计和用双线性变换法设计IIR数字滤波器的内容，用双线性变换法设计数字滤波器系统函数zH。其中满足本实验要求的数字滤波器系统函数为：212121612155.09044.013583.00106.117051.02686.1110007378.0zzzzzzzzHzHkk31（2.1）式中：3211212121，，，kzCzBzzAzHkkk（2.2）2155.09044.03583.00106.17051.02686.109036.0332211CBCBCBA，，，根据设计指标，调用MATLAB信号处理工具箱buttord和butter，也可以得到zH。由公式（2.1）和（2.2）可见，滤波器zH由三个二阶滤波器zH1、zH2和zH3级联而成，如图2-1所示。zH1zH2zH3nxny1ny2nyny3图2-1滤波器zH的组成2、编写滤波器仿真程序，计算H(z)对心电信号采样序列x(n)的响应序列y(n)。设yk(n)为第k级二阶滤波器Hk(z)的输出序列，yk－1(n)为输入序列，如图1所示。由式2.2可以得到差分方程：当k=1时，()()。所以H(z)对x(n)的总响应序列y(n)可以用顺序迭代算法得到。即以此对k=1,2,3，求解差分方程式2.3，最后得到()()。仿真程序就是实现上述求解差分方程和顺序迭代算法的通用程序。也可以直接调用Matlab的filter库函数实现仿真。四、实验代码。1、双线性变换法设计butterworth低通IIR数字滤波器的代码如下所示：T=1;Fs=1/T;wpz=0.2;wsz=0.3;wp=2*tan(wpz*pi/2);ws=2*tan(wsz*pi/2);rp=1;rs=15;[N,wc]=buttord(wp,ws,rp,rs,'s');[B,A]=butter(N,wc,'s');fk=0:1/512:1;wk=2*pi*fk;Hk=freqs(B,A,wk);subplot(2,2,1);plot(fk,20*log10(abs(Hk)));gridon;xlabel('f/Hz');ylabel('幅度(dB)');axis([0,1,-100,5]);title('(模拟滤波器的幅频响应曲线)');[N,wc]=buttord(wpz,wsz,rp,rs);[Bz,Az]=butter(N,wc);wk=0:pi/512:pi;Hz=freqz(Bz,Az,wk);subplot(2,2,4);plot(wk/pi,20*log10(abs(Hz)));gridon;xlabel('\omega/\pi');ylabel('幅度(dB)');axis([0,1,-100,5]);title('(数字滤波器的幅频响应曲线)');2、根据叠加输入信号y(t)设计的滤波器代