数字信号处理实验报告.doc

数字信号处理实验报告课程名称：数字信号处理学 院：信息科学与工程学院专业班级：通信1502班学生姓名：侯子强学 号：0905140322指导教师：李宏2017年5月28日实验一 离散时间信号和系统响应一. 实验目的1. 熟悉连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解2. 掌握时域离散系统的时域特性3. 利用卷积方法观察分析系统的时域特性4. 掌握序列傅里叶变换的计算机实现方法，利用序列的傅里叶变换对离散信号及系统响应进行频域分析二、实验原理 1. 采样是连续信号数字化处理的第一个关键环节。对采样过程的研究不仅可以了解采样前后信号时域和频域特性的变化以及信号信息不丢失的条件，而且可以加深对离散傅里叶变换、Z变换和序列傅里叶变换之间关系式的理解。对连续信号以T为采样间隔进行时域等间隔理想采样，形成采样信号：式中为周期冲激脉冲，为的理想采样。的傅里叶变换为： 上式表明将连续信号采样后其频谱将变为周期的，周期为s=2/T。也即采样信号的频谱是原连续信号xa(t)的频谱Xa(j)在频率轴上以s为周期，周期延拓而成的。因此，若对连续信号进行采样，要保证采样频率fs2fm，fm为信号的最高频率，才可能由采样信号无失真地恢复出原模拟信号计算机实现时，利用计算机计算上式并不方便，因此我们利用采样序列的傅里叶变换来实现，即而为采样序列的傅里叶变换2. 时域中，描述系统特性的方法是差分方程和单位脉冲响应，频域中可用系统函数描述系统特性。已知输入信号，可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应。本实验仅在时域求解，对于差分方程可用Matlab中的工具箱函数filter()函数求解 一个时域离散线性时不变系统的输出与输入间的关系为：可用Matlab中的工具箱函数conv()函数求解三、实验内容及步骤1. 时域采样定理的验证给定模拟信号：式中 。其幅频特性如图所示：选择三种采样频率Fs=1kHz, 300Hz, 200Hz, 生成采样序列分别用序列表示。编写程序计算三个序列的幅频特性曲线，并绘图显示。观察在折叠频率附近与连续信号频谱有无明显差别，分析频谱混叠现象。实验程序如下%时域采样定理的验证%Fs=1KHzTp=64/1000; %Tp=64ms Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); %MFFT yn=xa(nT);subplot(3,2,1); stem(xnt); %box on;title(a) Fs=1000Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk);title(a) T*FTxa(nT),Fs=1000Hz); xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk)%Fs=300HzTp=64/1000; Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); yn=xa(nT);subplot(3,2,1); stem(xnt); box on;title(a) Fs=300Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk),r);title(a)T*FTxa(nT),Fs=300Hz); xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk)%Fs=200HzTp=64/1000; %64ms Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1; A=444.128;alph=pi*50*20.5;omega=pi*50*20.5;xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); yn=xa(nT);subplot(3,2,1); stem(xnt,.); box on;title(a) Fs=200Hz); k=0:M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk);title(a) T*FTxa(nT),Fs=200Hz); xlabel(f(Hz);ylabel(幅度);axis(0,Fs,0,1.2*max(abs(Xk);2. 给定一个低通滤波器的差分方程为：输入序列 （1）分别求出和的系统响应，并画出其波形（2） 求出系统的单位脉冲响应，画出其波形A=1,-0.9;B=0.05,0.05;x1n=ones(1,8),zeros(1,50)x2n=ones(1,200);hn=impz(B,A,50);subplot(3,1,1);stem(hn);title(（1）系统单位脉冲响应h(n);y1n=filter(B,A,x1n);subplot(3,1,2);stem(y1n);title(2)系统对R（8）的响应 y1(n);y2n=filter(B,A,x2n);subplot(3,1,3);stem(y2n);title(3)系统对u(n)的响应)y2(n);3. 给定系统的单位脉冲响应为 用线性卷积法求分别对系统和的输出响应，并画出波形x1n=ones(1,8);h1n=ones(1,10) zeros(1,20);h2n=1,2.5,2.5,1,zeros(1,10);y11n=conv(h1n,x1n);y22n=conv(h2n,x1n);subplot(2,2,1);stem(h1n,.b);title(4)系统单位脉冲响应h1(n);subplot(2,2,2);stem(y11n,.b);title(5)h1(n)与R8(n)的卷积y11(n);subplot(2,2,3);stem(h2n,.b);title(6)系统单位脉冲响应h2(n);subplot(2,2,4);stem(y22n,.b);title(7)h2(n)与R8（n）的卷积y22(n);四、实验思考1. 在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？它们所对应的模拟频率是否相同？为什么？答：当采样频率不同时，数字度量不同，但是模拟频率相同。因为数字频率W是模拟角频率用采样频率FS归一化频率。数字频率和模拟角频率之间的关系是W=T,模拟信号的模拟角频率不变，当采样频率不同时，T不同，所以数字频率不同。因此，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量不相同，但是它们所对应的模拟频率相同。2. 如果输入信号为无线长序列，系统的单位脉冲响应是有限长序列，可否用线性卷积法求系统的响应？如何求？答：（1）对输入信号序列分段； （2）求单位脉冲响应与各段的卷积； （3）将各段卷积结果相加。3. 如果信号经过低通滤波器，把信号的高频分量滤掉，时域信号会有何变化？用前面第二个实验结果进行分析说明答：把信号经过低通滤波器，把信号的高频成分滤掉，时域信号的剧烈将变得平滑。五、 实验心得及体会通过本次实验我重新温习了MATLAB这个软件的使用方法，运行环境。通过这款软件使我们的学习更加便利。实验二 用FFT对信号作频谱分析一、实验目的1. 进一步加深DFT算法原理和基本性质的理解2. 掌握用FFT对连续信号和时域离散信号进行频谱分析的方法3. 了解用FFT进行频谱分析时可能出现的分析误差及其原因，以便在实际中正确应用FFT二、实验原理用FFT对信号作频谱分析是学习数字信号处理的重要内容，经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率F和分析误差。频谱分辨率直接和FFT的变换区间N有关，FFT能够实现的频率分辨率是2p/N，因此要求2p/NF。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时，离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。对模拟信号进行谱分析时，首先要按照采样定理将其变为时域离散信号。如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。三、实验步骤及内容1. 对以下给出的各序列进行谱分析： 选择FFT的变换区间N为8和16两种情况进行频谱分析。分别打印其幅频特性曲线，并进行对比、分析、讨论。x1n=ones(1,4); %产生R4(n)序列向量X1k8=fft(x1n,8); %计算x1n的8点DFTX1k16=fft(x1n,16); %计算x1n的16点DFTN=8;f=2/N*(0:N-1);figure(1);subplot(1,2,1);stem(f,abs(X1k8),.); %绘制8点DFT的幅频特性图title(1a) 8点DFTx_1(n);xlabel(/);ylabel(幅度);N=16;f=2/N*(0:N-1);subplot(1,2,2);stem(f,abs(X1k16),.); %绘制8点DFT的幅频特性图title(1a) 16点DFTx_1(n);xlabel(/);ylabel(幅度);%x2n 和 x3nM=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=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),.); %绘制8点DFT的幅频特性图title(2a) 8点DFTx_2(n);xlabel(/);ylabel(幅度);subplot(2,2,3);stem(f,abs(X3k8),.); %绘制8点DFT的幅频特性图title(3a) 8点DFTx_3(n);xlabel(/);ylabel(幅度);N=16;f=2/N*(0:N-1);subplot(2,2,2);stem(f,abs(X2k16),.); %绘制8点DFT的幅频特性图title(2a) 16点DFTx_2(n);xlabel(/);ylabel(幅度);subplot(2,2,4);stem(f,abs(X3k16),.); %绘制8点DFT的幅频特性图title(3a) 16点DFTx_3(n);xlabel(/);ylabel(幅度);2. 对以下各周期序列进行频谱分析 选FFT的变换区间N为8和16两种情况分别对以上序列进行频谱分析。分别打印其幅频特性曲线，并进行对比、分析、讨论。n=0:8; xn4=cos(pi.*n)/4); subplot(2,3,1);stem(n,xn4,.); X8k4 = fft(xn4,8); n21 = 0:length(X8k4)-1; subplot(2,3,2);stem(n21,X8k4,.); X16k4 = fft(xn4,16); n22 = 0:length(X16k4)-1; subplot(2,3,3);stem(n22,X16k4,.);%endn=0:16; xn5=cos(pi.*n)/4)+cos(pi.*n)/8); subplot(2,3,4);stem(n,xn5,.); X8k5 = fft(xn5,8); n21 = 0:length(X8k5)-1; subplot(2,3,5);stem(n21,X8k5,.); X16k5 = fft(xn5,16); n22 = 0:length(X16k5)-1; subplot(2,3,6);stem(n22,X16k5,.);3. 对模拟周期信号进行频谱分析选择样频率Fs=64Hz，对变换区间N=16，32，64三种情况进行谱分析。分别打印其幅频特性曲线，并进行对比、分析、讨论。程序如下：Fs=64;T=1/Fs; N=16;n=0:N-1; %FFT的变换区间N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)16点采样 X6k16=fft(x6nT); %计算x6nT的16点DFT X6k16=fftshift(X6k16); %将零频率移到频谱中心 Tp=N*T;F=1/Tp; %频率分辨率Fk=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率（以零频率为中心） subplot(3,1,1);stem(fk,abs(X6k16),.);box on %绘制8点DFT的幅频特性图 title(6a) 16点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度); axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16) N=32;n=0:N-1; %FFT的变换区间N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)32点采样 X6k32=fft(x6nT); %计算x6nT的32点DFT X6k32=fftshift(X6k32); %将零频率移到频谱中心 Tp=N*T;F=1/Tp; %频率分辨率F k=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率（以零频率为中心） subplot(3,1,2);stem(fk,abs(X6k32),.);box on %绘制8点DFT的幅频特性图 title(6b) 32点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度); axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32) N=64;n=0:N-1; %FFT的变换区间N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)64点采样 X6k64=fft(x6nT); %计算x6nT的64点DFT X6k64=fftshift(X6k64); %将零频率移到频谱中心 Tp=N*T;F=1/Tp; %频率分辨率F k=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率（以零频率为中心） subplot(3,1,3);stem(fk,abs(X6k64),.); box on%绘制8点DFT的幅频特性图 title(6a) 64点|DFTx_6(nT)|);xlabel(f(Hz);ylabel(幅度); axis(-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64)四、实验思考1. 在N=8时， 和的幅频特性会相同吗？为什么？N=16时呢？答：在N=8时， 和 的幅频特性相同，而N=16时不相同。因为 = ,所以 和 的8点DFT的模相等。但当N=16时， 和 不满足循环移位关系，所以两者幅频特性不相同。2. 对于周期序列，如果周期不知道，如何用FFT进行谱分析？答：周期信号的周期预先不知道时,可先截取M点进行DFT,再将截取长度扩大1倍截取,比较结果,如果二者的差别满足分析误差要求,则可以近似表示该信号的频谱,如果不满足误差要求就继续将截取长度加倍,重复比较,直到结果满足要求 。五、实验总结及心得体会通过实验进一步加深DFT算法原理和基本性质的理解，掌握用FFT对连续信号和时域离散信号进行频谱分析的方法，了解用FFT进行频谱分析时可能出现的分析误差及其原因，以便在实际中正确应用FFT。实验三 用双线性变换法设计IIR数字滤波器一、实验目的1. 熟悉用双线性变换法设计IIR数字滤波器的原理和方法2. 掌握IIR数字滤波器的Matlab实现方法3. 通过观察对实际心电图信号的滤波作用，获得数字滤波的感性认识二、实验原理设计IIR数字滤波器一般采用间接设计法脉冲响应不变法和双线性变换法，应用最广泛的是双线性变换法。脉冲响应不变法的基本思想是：使数字滤波器的单位脉冲响应h(n)近似于模拟滤波器的单位脉冲响应ha(t)，即使 其S平面和Z平面的映射关系为：双线性变换法的基本思想是：使描述数字滤波器的差分方程近似描述模拟滤波器的微分方程S平面和Z平面的映射关系为：双线性变换法中的频率变换是一种非线性变换，这种非线性引起的幅频特性畸变可通过预变形矫正法而得到校正。设计IIR 数字滤波器的一般步骤: （1）确定所需类型数字滤波器的技术指标：通带截止频率p、通带衰减p、阻带截止频率s、阻带衰减s。（2）将所需类型数字滤波器的技术指标转换成相应类型模拟滤波器的技术指标。（3）设计该类型模拟滤波器（4）通过复频率变换将模拟滤波器转换成所需类型的数字滤波器。三、实验内容1. 分别用脉冲响应不变法和双线性变换法设计一个巴特沃斯低通IIR数字滤波器，设计指标参数为：在通带内频率低于0.2p时，最大衰减小于1dB，在阻带内0.3p，p频率区间上，最小衰减大于15dB。观察并画出所设计数字滤波器的幅频特性曲线和相频特性曲线，记录带宽和衰减量，检查是否满足要求。比较这两种方法的优缺点。Matlab程序为：%脉冲响应法T=1;wp=0.2*pi;ws=0.3*pi;rp=1;as=15; %输入低通滤波器要求n,wpo=buttord(wp,ws,rp,as,s); %计算阶数B,A=butter(n,wpo,s); %计算表达式分子分母的系数矩阵B1,A1=impinvar(B,A); Hk,w=freqz(B1,A1);subplot(2,1,1);plot(w/pi,20*log10(abs(Hk); %画出滤波器损耗函数曲线grid on;title(1)脉冲响应不变法衰减曲线);xlabel(频率（w/pi）);ylabel(幅度（dB）);%双线性法wp=2*tan(0.2*pi/2);ws=2*tan(0.3*pi/2);rp=1;as=15;n,wpo=buttord(wp,ws,rp,as,s); %计算阶数B,A=butter(n,wpo,s); %计算表达式分子分母的系数矩阵B1,A1=bilinear(B,A,1);Hk,w=freqz(B1,A1);subplot(2,1,2);plot(w/pi,20*log10(abs(Hk); %画出滤波器损耗函数曲线title(1)双线性法衰减曲线);grid on;xlabel(频率（w/pi）);ylabel(幅度（dB）);优缺点比较：（1）脉冲响应不变法会产生频谱混叠，但具有很好的线性特性，其单位脉冲响应完全模仿模拟滤波器的单位冲激响应波形，时域逼近性好。适合于带通、低通滤波器的设计。（2）双线性法很好地消除了频谱混叠，但是其数字频率与模拟频率之间不具有线性关系。2. 用双线性变换法设计一个切比雪夫高通IIR数字滤波器，设计指标参数为：在通带内频率高于0.3KHz, 最大衰减小于1dB，在阻带内频率低于0.2 KHz，最小衰减大于20dB，T=1ms。 画出所设计数字滤波器的幅频特性曲线和相频特性曲线，观察其通带损耗和阻带衰减是否满足要求。Matlab程序如下：%切比雪夫高通滤波器的设计fp=2*pi*300*0.001;fs=2*pi*200*0.001;wp=2000*tan(fp/2);ws=2000*tan(fs/2); %进行频率变换rp=1;as=20; %输入高通滤波器要求n,wpo=buttord(wp,ws,rp,as,s); %计算阶数B,A=butter(n,wpo,high,s); %计算表达式分子分母的系数矩阵B1,A1=impinvar(B,A);Hk,w=freqz(B,A);plot(w/pi,20*log10(abs(Hk); %画滤波器损耗函数曲线grid on;title(3)切比雪夫高通IIR数字滤波器);xlabel(频率（HZ）);ylabel(幅度（dB）);3. 人体心电图信号在测量过程中往往受到工业高频干扰，所以必须经过低通滤波处理后，才能作为判断心脏功能的有用信息。下面给出一实际心电图信号采样序列样本x(n)，其中存在高频干扰。 用1所计的滤波器对心电图信号采样序列x(n)进行仿真滤波处理，画出处理前后的信号波形。xn=-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;subplot(2,1,1);plot(xn); title(xn);y1n=filter(B1,A1,xn);subplot(2,1,2);plot(y1n);title(xn);四、实验思考1. 用双线性变换法设计数字滤波器过程中，变换公式 中T取值，对设计结果有无影响？为什么？答：没有。因为在第一步和第二步中从数字角频率W变到模拟角频率，再从模拟角频率变到数字角频率W，两次变换是对称的，只要两次变换过程的T是相同的即可，T的取值是无关紧要的2. 双线性变换法中和之间的关系是非线性的，在实验中你注意到这种非线性关系了吗？从哪几种数字滤波器的幅频特性曲线中可以观察到这种非线性关系？答：注意到了。双线性变换是从S平面映射到S1平面，再从S1平面映射到Z平面，一个线性相位的模拟滤波器经过双线性法变换后，就变成了非线性的了。切比雪夫的幅频特性是非线性的。五、实验总结及心得体会 通过实验学会了用双线性变换法设计IIR数字滤波器的原理和方法，并且掌握IIR数字滤波器的Matlab实现方法，通过观察对实际心电图信号的滤波作用，获得数字滤波的感性认识。实验四 用窗函数法设计FIR数字滤波器一、实验目的1. 掌握用窗函数法设计FIR数字滤波器的原理和方法2. 熟悉线性相位FIR数字滤波器特性3. 了解各种窗函数对滤波特性的影响二、实验原理 窗函数法设计 FIR 滤波器的步骤为： （1）构建希望逼近的理想频率响应函数及技术指标（2）求滤波器的单位脉冲响应 如果复杂，可对从采样M个点，采样值为，则：（3）根据对过渡带及阻带衰减的要求，选择窗函数的形式，并估计窗口宽度N，设要求的过渡带宽为，则（4）计算滤波器的单位脉冲响应： （5）求H(ej)，分析其幅频特性，若不满足要求，可适当改变窗函数形式或长度N，重复上述设计过程，以得到满意的结果。窗函数傅里叶变换W(ej)的主瓣决定了H(ej)过渡带宽，W(ej)的旁瓣大小和多少决定了H(ej)在通带和阻带范围内波动幅度，常用的几种窗函数有：矩形窗；Hanning窗；Hamming窗；Blackmen窗；Kaiser窗三、实验内容及步骤1. 用升余弦窗设计一线性相位低通FIR 数字滤波器，截止频率 。窗口长度N=15，33。要求在两种窗口长度情况下，分别求出h(n)，打印出相应的幅频特性和相频特性曲线，观察3dB带宽和20dB 带宽，总结窗口长度N对滤波特性的影响%用升余弦窗设计一线性相位低通FIR 数字滤波器%N=15hn1=fir1(14,pi/4,hanning(15);figure(1);Hn1,w=freqz(hn1,1);subplot(2,1,1);plot(w/pi,abs(Hn1);title(N=15的h(n)的幅频曲线)xlabel(w);ylabel(幅度);w1=angle(Hn1);subplot(2,1,2);plot(w/pi,w1);title(N=15的h(n)的相频曲线)xlabel(w);ylabel(angle(Hn1)%N为33hn2=fir1(32,pi/4,hanning(33);figure(2);Hn2,w=freqz(hn2,1);subplot(2,1,1);plot(w/pi,abs(Hn2);title(N=33的h(n)的幅频曲线)title(N=33的h(n)的相频曲线)xlabel(w);ylabel(幅度);w2=angle(Hn2);subplot(2,1,2);plot(w/pi,w2);xlabel(w);ylabel(angle(Hn2)对窗口长度N对滤波特性的影响的总结： 调整窗口长度N只能有效地控制过渡带的宽度，当N增大，主瓣幅度加高，同时旁瓣也加高，WRg（）的主瓣和旁瓣幅度变窄。导致波动频率加快。因此，加大N，并不是解决吉布斯效应的有效方法。对窗口长度N对滤波特性的影响的总结： 调整窗口长度N只能有效地控制过渡带的宽度，当N增大，主瓣幅度加高，同时旁瓣也加高，WRg（）的主瓣和旁瓣幅度变窄。导致波动频率加快。因此，加大N，并不是解决吉布斯效应的有效方法。2. N=33, , 用四种窗函数设计线性相位低通滤波器，绘制相应的幅频特性曲线，观察3dB带宽和20dB 带宽以及阻带最小衰减，比较四种窗函数对滤波器特性的影响%矩形窗hn1=fir1(32,pi/4,boxcar(33);figure(1);Hn1,w=freqz(hn1,1);subplot(4,1,1);plot(w/pi,abs(Hn1);title(（1）（2）矩形窗的幅频曲线)xlabel(w);ylabel(幅度);w1=angle(Hn1);subplot(4,1,2);plot(w/pi,w1);xlabel(w);ylabel(angle(Hn)%三角窗hn2=fir1(32,pi/4,bartlett(33);Hn2,w=freqz(hn2,1);subplot(4,1,3);plot(w/pi,abs(Hn2);title(（3）(4）三角窗的幅频曲线)xlabel(w);ylabel(幅度);w2=angle(Hn2);subplot(4,1,4);plot(w/pi,w2);xlabel(w);ylabel(angle(Hn2)%汉宁窗hn3=fir1(32,pi/4,hanning(33);figure(2);Hn3,w=freqz(hn3,1);subplot(4,1,1);plot(w/pi,abs(Hn3);title(（1）（2）汉宁窗的幅频曲线)xlabel(w);ylabel(h(n）幅频曲线);w3=angle(Hn3);subplot(4,1,2);plot(w/pi,w3);xlabel(w);ylabel(angle(Hn3)%哈明窗hn4=fir1(32,pi/4,hamming(33);Hn4,w=freqz(hn4,1);subplot(4,1,3);plot(w/pi,abs(Hn4);title(（3）(4）哈明窗的幅频曲线)xlabel(w);ylabel(h(n）幅频曲线);w4=angle(Hn4);subplot(4,1,4);plot(w/pi,w4);xlabel(w);ylabel(angle(Hn4)3. 调用信号产生函数xtg（见教材P295）产生具有加性噪声的信号xt，显示xt波形及其频谱。设计一FIR 数字低通滤波器，从高频噪声中提取xt中的单频调幅信号，要求信号幅频失真小于0.1dB,噪声频谱衰减不小于60dB。（1）观察xt的频谱，确定滤波器指标参数。（2）根据滤波器指标选择合适的窗函数，计算窗函数的长度N，设计一个FIR低通滤波器，绘图显示滤波器的频响特性曲线（3）用设计的FIR低通滤波器对xt进行滤波，绘图显示滤波器输出信号的时域和频域波形图Matlab程序为：function xt=xtg%信号x(t)产生函数，并显示信号的时域波形和幅频特性曲线%xt=xtg产生一个长度为N，有加性高频噪声的单调调幅信号xt，n=1000，%采样频率fs=1000HZ%载波频率fc=fs/10=100HZ,调制正弦波频率f0=fc/10=10HZN=1000;Fs=1000;T=1/Fs;Tp=N*T;t=0:T:(N-1)*T;fc=Fs/10;f0=fc/10;mt=cos(2*pi*f0*t);ct=cos(2*pi*fc*t);xt=mt.*ct;nt=2*rand(1,N)-1;%设计高通滤波器hn,用于滤除噪声nt中的低频成分，生成高通噪声fp=150;fs=200;rp=0.1;as=70; %滤波器指标fb=fp,fs;m=0,1;dev=10(-as/20),(10(rp/20)-1)/(10(rp/20)+1);n,fo,mo,W=remezord(fb,m,dev,Fs); %确定remez函数所需参数hn=remez(n,fo,mo,W); %调用remez函数进行设计，用于滤除噪声nt中的低频成分yt=filter(hn,1,10*nt); %滤除随机噪声中的低频成分，生成高通噪声yt%以下为绘图成分xt=xt+yt;%噪声加信号fst=fft(xt,N);k=0:N-1;f=k/Tp;subplot(2,1,1);plot(t,xt);grid;xlabel(t/s);ylabel(x(t);axis(0,Tp/5,min(xt),max(xt);title(a)xt信号加噪声波形);subplot(2,1,2);plot(f,abs(fst)/max(abs(fst);grid;title(b)xt信号加噪声频谱);axis(0,Fs/2,0,1.2);xlabel(f/HZ);ylabel(幅度);%阻带衰减as不低于60，所以选用布莱克曼窗hn=fir1(119,pi/10,blackman(120);figure(2);Hn,w=freqz(hn,1);subplot(2,1,1);plot(w/pi,abs(Hn);title(滤波器的频响特性);xlabel(w);ylabel(h(n）频响特性曲线);%滤波过程yn=conv(xt,hn);figure(3);subplot(2,1,1);plot(yn);title(xt经过低通后的时域波形);Hn,w=freqz(yn,1);subplot(2,1,2);plot(w/pi,abs(Hn)/max(abs(Hn);xlabel(w);title(xt经过低通后的频域波形);四、实验思考1. 如果给定通带截止频率和阻带截止频率以及阻带最小衰减，如何用窗函数法设计线性相位低通滤波器？答：1）根据阻带衰减确定窗函数，然后根据通带与阻带截止频率确定窗的长度N2）写出低通滤波器的频域函数3）求低通滤波器的时域函数4）求=2. 如果要求用窗函数法设计带通滤波器，且给定上、下边带截止频率为和，试求理想带通的单位脉冲响应答：1）根据阻带衰减确定窗函数，然后根据通带与阻带截止频率确定窗的长度N2）写出带通滤波器的频域函数3）求带通滤波器的时域函数4）求=五. 实验总结IIR滤波器的优点是可利用模拟滤波器设计的结果，缺点是相位是非线性的。而FIR滤波器具有良好的线性相位。FIR滤波器线性相位的特点： 如果FIR滤波器的单位抽样响应H(N)为实数，而且满足以下任一条件： 偶对称H(N)H(N-1-N) 奇对称H(N)-H(N-1-N) 且其对称中心在N(N-1)/2处则滤波器具有准确的线性相位。以下对FIR及IIR进行对比：FIR滤波器IIR滤波器设计方法一般无解析的设计公式，要借助计算机程序完成利用AF的成果，可简单、有效地完成设计设计结果线性相位(最大优点)只能得到幅频特性，相频特性未知， 如需要线性相位，须用全通网络校 准，但增加滤波器阶数和复杂性稳定性极点全部在原点(永远稳定)有稳定性的问题 心得体会：1.更加熟悉掌握用窗函数法设计FIR数字滤波器的原理和方法；2. 熟悉了线性相位FIR数字滤波器特性；3. 了解了各种窗函数对滤波特性的影响；4.对FIR和IIR 有了更加深刻的理解.5.实验过程中也遇到了很多困难，在同学帮助下得到解决。