数字信号处理上机实习报告.pdf

仅供个人参考 不得用于商业用途 专题一 离散卷积的计算 一、实验内容 设线性时不变(LTI)系统的冲激响应为 h(n)，输入序列为 x(n)1、h(n)=(0.8)n，0 nW 4;x(n)=u(n)-u(n-4)For pers onal use only in study and research;not for commercial use 2、h(n)=(0.8)nu(n),x(n)=u(n)-u(n-4)n 3、h(n)=(0.8)u(n),x(n)=u(n)求以上三种情况下系统的输出 y(n)。For pers onal use only in study and research;not for commercial use 二、实验目的 1、掌握离散卷积计算机实现。2、进一步对离散信号卷积算法的理解 三、原理及算法概要 算法：把冲激响应 h(n)与输入序列 x(n)分别输入到程序中，然后调用离散卷积函数 y=conv(x.,h)即可得到所要求的结果。原理：离散卷积定义为 cd y(n)二 x(k)h(n-k)当序列为有限长时，则 n y(n)二 x(k)h(n-k)k=e 仅供个人参考 不得用于商业用途 四理论计算 1、h(n)=(0.8)n,0 W nW 4;x(n)=u(n)-u(n-4)仅供个人参考 不得用于商业用途 C30 y(n)=x(n)-h(n)二x(m)h(n-m)m=.：:(a)当 n7 时，y(n)=0 理论结果与上图实验结果图中所示吻合。2、h(n)=(0.8)nu(n),x(n)=u(n)-u(n-4)oC y(n)=x(n)-h(n)二、x(m)h(n _ m)m 二 3(a)当 n0 时，y(n)=0(b)n 当0乞n乞3时，y(n八(0.8)n mq0 n(c)当 4 23 时，y(n)=0 在卷积序列中，因为 n 趋于无穷，在实验中无法实现，而 0.8 的 20 次幕基本接近于 0,所以这里只取到 n 为 20 就可以了。由以上可以看出，理论结果与上图实验结果图中所 示吻合。3、h(n)=(0.8)nu(n),x(n)=u(n)QO y(n)=x(n)-h(n)=、x(m)h(n _ m)m=:(a)当 n0 时，y(n)=0(b)n 当 0 140 时，y(n)=0 在卷积序列中，因为 n 趋于无穷，在实验中无法实现，而 0.8 的 70 次幕的叠加基本接 近于最大值，所以这里只取到 n 为 70 就可以了。由以上可以看出，理论结果与上图实验结 果图中所示吻合。五程序 x1=1 1 1 1;n x 仁 0:3;h1=1 0.8 0.64 0.8A3 0.8A4;nh1=0:4;y1=c on v(x1,h1);subplot(3,3,1);stem(nx1,x1);title(序列 x1);xlabel(n);ylabel(x1(n);subplot(3,3,2);stem(nh1,h1);title(序列 h1);xlabel(n);ylabel(h1(n);subplot(3,3,3);stem(y1);title(序列 y1);xlabel(n);ylabel(y1(n);x2=1 1 1 1;nx2=0:3;nh2=0:1:20;h2=(0.8).A nh2;y2=c on v(x2,h2);subplot(3,3,4);stem(x2);title(序列 x2);xlabel(n);ylabel(x2(n);subplot(3,3,5);stem(h2);title(序列 h2);xlabel(n);ylabel(h2(n);subplot(3,3,6);stem(y2);title(序列 y2);xlabel(n);ylabel(y2(n)n x3=0:1:20;X3=1.A nx3;nh3=0:1:20;h3=(0.8).A nh3;y3=c on v(x3,h3);subplot(3,3,7);stem(nx3,x3);title(序列 x3);xlabel(n);ylabel(x3(n);subplot(3,3,8);stem(nh3,h3);title(序列 h3);xlabel(n);ylabel(h3(n);subplot(3,3,9);stem(y3);title(序列 y3);xlabel(n);ylabel(y3(n)六.程序运行结果 仅供个人参考 不得用于商业用途 七结果分析 有限长序列的离散卷积计算结果与理论值一致，而存在无限长序列做卷积时，由于在程序 处理时是用比较长有限长序列代替的，所以与理论值基本相同。八.专题实习心得 该专题主要进行的是有 matlab 实现离散卷积的计算，分为三个类型：冲激响应 h(n)与输 入序列 x(n)都为有限长，一个序列为有限长一个序列为无限长和两个序列均为无限 长。这部分实习内容不难，是原来做过的，主要是为了拾拣起原来所学的知识。第一 种类型是最简单，也是最基本的，直接调用函数 y=conv(x.,h)即可。在第二种和第三 种类型的计算中遇到了一些困难，在输入序列时，由于存在无限长的序列，不知道该 怎么输入，查过了相关的题目才弄明白。通过本专题的实习，让我对数字信号处理中 的一些基础知识有了一些回顾，对原来所学过的知识也熟悉了一些。仅供个人参考 不得用于商业用途 专题二 离散傅里叶变换及其应用 分析下列三种情况下的幅频特性。(1)采集数据长度 N=16,分析 16 点的频谱，并画出幅频特性。(2)采集数据长度 N=16,并补零到 64 点，分析其频谱，并画出幅频特性。(3)采集数据长度 N=64,分析 46 点的频谱，并画出幅频特性。观察三幅不同的幅频特性图，分析和比较它们的特点及形成原因。2、实现序列的内插和抽取所对应的离散傅里叶变换。(选做)4对应的128点傅里叶变换，比较三个结果所得到的幅 频特性，分析其差别产生的原因。二、实验目的 1、了解 DFT 及 FFT 的性质和特点 2、利用 FFT 算法计算信号的频谱。三、关键算法 算法 1:读入离散序列 x(n)=cos(0.48 n n)+cos(0.52 n n)，采集长度为 N=16 的数据，调 用 matlab中的函数fft(x,16)与 fft(x,64)对其作离散傅里叶变换得到 16 点、64 点 的频谱。采集数据长度为 N=64,，调用 matlab 中的函数fft(x,46)对其作离散傅里叶变 换得到46 点的频谱。实验内容 设有离散序列 x(n)=cos(0.48 n n)+cos(0.52 n n)x n I-cos 36 x,n ix 4n 和 x2 n i=x 仅供个人参考 不得用于商业用途 算法 2:读入离散序列信号 x n=cos n亠cos n R128 n，每隔 N=4采集-(36丿 J 36丿 个数值，得到新的抽取离散序列 x,n=x4n，每隔K=0.25采集一个数值得到新的内插 f n、离散序列X2(n)=x 1。分别对其做离散傅氏变换，得到对应的频谱图。原理概要：当抽样数”=労时，以下为算法蝶形图。X(0)W社X=X(2)曲.X(6)WN0 X0)X(5)X(3)0 WN#X(7)1、当 N=2M 时，则要进行 M 次分解，即进行 M 级蝶形单元的计算 2、按自然顺序输入，输出是码位倒置。3、每一级包含 N/2 个基本蝶形运算 4、第 L 级有 2L-1个蝶群，蝶群间隔为 N/2 如果是 Matlab 实现的话，可用以下两种方法计算信号频谱 1、调用库函数为：fft(),直接计算 X(k)-2、进行矩阵运算 X(0)=x(1)*1 WN0 WN Wr WN 9 9 WN TT X(0)-W畀 X(1)*-X(N-1)1 WN N-1 WN WNH_ x(N-1)四.程序 程序 1:n=0:1:15;x1=cos(0.48*3.14*n)+cos(0.52*3.14*n);g1=abs(fft(x1,16);subplot(3,2,1);stem(x1);title(x1);subplot(3,2,2);stem(g1);title(g1);x(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7)般规律如下:仅供个人参考 不得用于商业用途 n2=0:1:15;x2=cos(0.48*3.14*n2)+cos(0.52*3.14*n2);x2=x2 zeros(1,48);g2=abs(fft(x2,64);subplot(3,2,3);stem(x2);title(x2);subplot(3,2,4);stem(g2);title(g2);n3=0:1:64;x3=cos(0.48*3.14*n3)+cos(0.52*3.14*n3);g3=abs(fft(x3,46);subplot(3,2,5);stem(x3);title(x3);subplot(3,2,6);stem(g3);title(g3);程序 2：x=;for n=0:127 x=x(cos(1/36*pi*n)+cos(1.5/36*pi*n)end X=abs(fft(x);subplot(3,2,1);stem(x);title(序列 x);xlabel(n);ylabel(x(n);subplot(3,2,2);stem(X);title(序列 X);xlabel(omega);ylabel(X);x1=;for n=0:4:508 x1=x1(cos(1/36*pi*n)+cos(1.5/36*pi*n)end X1=abs(fft(x1);subplot(3,2,3);stem(x1);title(序列 x1);xlabel(n);ylabel(x1(n);subplot(3,2,4);stem(X1);title(序列 X1);xlabel(omega);ylabel(X1);x2=;for n=0:0.25:31.75 x2=x2(cos(1/36*pi*n)+cos(1.5/36*pi*n)end X2=abs(fft(x2);subplot(3,2,5);stem(x2);title(序列 x2);xlabel(n);ylabel(x2(n);subplot(3,2,6);stem(X2);title(序列 X2);仅供个人参考 不得用于商业用途 xlabel(omega);ylabel(X2);五程序运行结果 80 60 x 40 20 0 1 1 1 0 50 100 150 C9 结果 1：结果 2:序列爲 50 100 150 序列仅供个人参考 不得用于商业用途 六结果分析 N 点 DFT 的频谱分辨率是 2 n/N。一节指出可以通过补零观察到更多的频点，但是 这并不意味着补零能够提高真正的频谱分辨率。这是因为 xn实际上是 x(t)采样的 主值序列，而将 xn补零得到的 xn 周期延拓之后与原来的序列并不相同，也不是 x(t)的采样。因此是不同离散信号的频谱。对于补零至 M 点的 x的 DFT，只能说它的 分辨率 2 n/M 仅具有计算上的意义，并不是真正的、物理意义上的频谱。频谱分辨率的 提高只能通过提高采样频率实现。七.专题实习心得 离散傅里叶变换是一种快速算法，由于有限长序列在其频域也可离散化为有限长序列，因 此离散傅里叶变换在数字信号处理中是非常有用的。DFT 是重要的变换，在分析有限长序 列的有用工具、信号处理的理论上有重要意义、运算方法上起核心作用，谱分析、卷积、相关都可以通 DFT 在计算机上实现。DFT 解决了两个问题：一是离散与量化，二是快速运 算。通过编程实践体会到了时域、频域信号的对应关系，也对采样频率的含义有了深刻的 认识，同时也加深了对采样信号频谱周期性的理解。n 0 50 100 150 仅供个人参考 不得用于商业用途 八.思考题 1 直接计算 DFT 与用 FFT 计算，理论上的速度差别应是多少？整个 DFT 运算总共需要 4N 2 次实数乘法和 2NN2 N-1)次实数加法。整个 FFT 运算总共需要(N2/2)+(N2)=NN+1)/2衣N2/2 次复数乘法和 N(N/2-1)+N=N2/2 次复数加法。由此可见，理论上用 FFT 计算运算速度是 DFT 的二倍。2、什么是高密度频谱及高分辨率频谱？实验 1 中至少应采集几个样本才能区分两个频率 分量？DFT 的频谱分辨率一定时，在尾部补零可以得到 DFT 的高密度频谱，可以细化当前分辨率 下的频谱，但不能改变 DFT 的分辨率。高分辨率是指对信号中两个靠得较近的频谱分量的 识别能力比较高，在采样频率不变时，采样点数 N 越大，频谱分辨率越高。在实验 1 中至 少要采集 8 个点才能区分两个频率分量。专题三 IIR 滤波器的设计 一、实验内容 1、设计一个 Butterworth 数字低通滤波器，设计指标如下：通带截止频率：0.2 n，幅度衰减不大于 1 分贝 阻带截止频率：0.3 n，幅度衰减大于 15 分贝 2、让不同频率的正弦波通过滤波器，验证滤波器性能。3、分析不同滤波器的特点和结果。4、编程设计实现 IIR 滤波器。二、实验目的 掌握不同 IIR 滤波器的性质、特点。通过实验学习如何设计各种常用的 IIR 滤波器，以便在实际工作中能根据具体情况使用 IIR 滤波器。三、原理及算法概要 算法：输入通带截止频率 Wp 阻带截止频率 Ws 通带衰减 Rp,阻带衰减 Rs,通过这些数值 调用N Wn=buttord(Wp,Ws,Rp,Rs)函数计算巴特沃斯数字滤波器的阶数 N 和截止频率 wn,再根据阶数 N 通过函数b,a=butter(N,Wn)，即可得到所要的巴特沃斯滤波器。设计一个正 弦波信号，再调用函数 A=filter(b,a,I)让正弦波信号通过滤波器，得到滤波信号。原理概要：1、滤波器类型 Butterworth 滤波器 Butterworth 滤波器的特点是在通带内的频率特性是平坦的，并且随着频率的增加而 衰减。Butterworth 滤波器又是最简单的滤波器。仅供个人参考 不得用于商业用途 2 1|Ha 2N N 阶低通 Butterworth 滤波器的幅度平方函数为：1 Chebyshev 滤波器 Chebyshev 滤波器的在通带内的响应是等纹波的，而在阻带内是单调下降的，或在 通带内是单调下降的，在阻带内是等纹波的特性|Ha(j)J.2N 其幅度平方函数为 1.；PN(/C)其中 VN 是 Chebyshev 多项式。椭圆滤波器 椭圆滤波器在通带和阻带内都是等纹波振荡。2 1|Ha(jO)|2=“2 宀 c/c 2N 椭圆滤波器的特性函数为：1 U N八C)其中 UN 是 N 阶雅可比函数。模拟滤波器设计完成后，就可以进行典型滤波器设计的第二个步骤，通过冲激响应不变法和 双线性变换转变为相应的数字滤波器。2、变换方法(1)冲激响应不变法 冲激响应不变法的基本原理是从滤波器的冲激响应出发，对模拟滤波器冲激响应 h(t)进 行取样，所得到的离散序列 h(nT)作为数字滤波器的单位取样响应。H(z)是由 H(s)通过下式的对应关系得到。1 _ 1 s-Sk 1 (2)双线性变换是在所得到满足性能指标要求的模拟滤波器的基础上 ，通过变换 1 1-Z S=丄二二,从而得到相应的数字滤波器。T 1 z4 四.程序 Wp=0.2;Ws=0.3;Rp=1;Rs=15;N Wn=buttord(Wp,Ws,Rp,Rs)%用于计算巴特沃斯数字滤波器的阶数 N 和截止频率 wn 仅供个人参考 不得用于商业用途 b,a=butter(N,Wn);%计算 N 阶巴特沃斯数字滤波器系统函数分子、分母多项式的系数 向量 b、a,设计所需的低通滤波器 h,omega=freqz(b,a,512);%返回量 h 包含了离散系统频响，调用中若 N 默认，默认值为 512。plot(omega/pi,20*log10(abs(h);grid;xlabel(omega/pi);ylabel(Gain,dB);title(IIR Butterworth Lowpass Filter);Wp=0.2;Ws=0.3;Rp=1;Rs=15;N1,Wn1=buttord(Wp,Ws,Rp,Rs);%用于确定阶次 b,a=butter(N,Wn);%用于直接设计巴特沃兹数字滤波器，即为 IIR 滤波器%freqz(b,a);t=1:300 I=sin(0.1*pi*t)+sin(0.4*pi*t);%设计正弦波 plot(I);figure;A=filter(b,a,I);%正弦波通过滤波器 plot(A);五程序运行结果 N=Wn=0.2329 N1=6 Wn1=仅供个人参考 不得用于商业用途 0.2329 HR Butterworth Lowpass Filter 50 0-50-100-150-200-250-300-350 400 0 0.1 02 G.3 0.4 0.6 Q.0.7 08 0.9 1 仅供个人参考 不得用于商业用途 六结果分析 Butterworth 滤波器在通带内的频率特性是平坦的,并且随着频率的增加而衰减。正弦信号 在经过 仅供个人参考 不得用于商业用途 IIR 滤波器滤波后，高频信号被滤除，低频信号被保留了下来。七.专题实习心得 通过本专题的学习，熟悉和巩固了 Butterworth 数字低通滤波器的设计方法和原理，实现 滤波器设计的有关经典算法，更重要的是熟练掌握使用 MATLAB言设计各种要求的数字滤 波器。这一部分的内容相对来说是有些难度了，做起来花费的精力也多了一些，不过对数 字信号处理内容的掌握上又加深了一层。八思考题 1、比较双线性变换与 冲激响应不变法的优缺点。冲激响应不变法的特点：模仿模拟滤波器的单位抽样相应，时域模仿性好 产生频率响应的混叠失真。保持系 统的相位特性 换不改 变系统的因果性和 稳定性;因为存在频率响应的混 叠效应，只适用于限 带的模拟滤波器，高通和 带阻滤波器 不宜采用冲激 响应不变法，冲激 响应不变法是一 种时域模仿，缺点是 产生频率响应的 混叠失真。冲激 响应不变法产生混叠失真的原因是映射是多 对一的映射。双线性变换 的特点是一 对一的映射，消除了多 值变换 性，也就消除了 频谱混叠。1、和频率变换相比，数字变换有什么优点？数字滤波器的幅 频响应相对于模 拟滤波器的幅 频响应有畸变 3、实验中求取相应模拟滤波器原形时，为什么要对模拟频率 Qc 进行归一化？用归一化巴特沃思低通 滤波器 为 原型 滤波器，一旦 归一化低通 滤波器的系 统函数 确定后，其 它巴特沃思低通 滤波、高通、带通、带阻滤波器的传递函数都可以通 过变 换法从归一化低通原型的 传递函数 Ha(S)得到。归一化原型滤波器是指截止 频率 Qc 已 经归一化成 Qc=1 的低通滤波器。对于截止频率为某个 Qc 的低通滤波器，则令 S/Qc 代替归一化原型滤波器系统函数中的 S人，即 S人 S/Qc对于其他高通、带通、带阻滤波器，可 应用后面讨论到的频带变换 法，由其 变换得出。专题四 用窗函数设计 FIR 滤波器 一、实验内容 选取合适窗函数设计一个线性相位 FIR 低通滤波器，使它满足如下性能指标：通带截止频率：3 p-0.5 n,通带截止频率处的衰减不大于 3 分贝；阻带截止频率：3 s=0.66 n，阻带衰减不小于 40 分贝。二、实验目的 1 掌握用窗函数法设计 FIR 滤波器的原理和方法。仅供个人参考 不得用于商业用途 2、熟悉线性相位滤波器特性。3、了解各种窗函数对滤波器特性的影响。三、原理及算法概要 算法：通过其通带截止频率 3 P 与阻带截止频率3 s 算出其过渡带的宽度与滤波器的长度，从 而得到理想滤波器的截止频率，根据所要求的理想滤波器，得到 hd(n)。由于其通带截止 频率处的衰减不大于 3 分贝与阻带衰减不小于 40 分贝，我选择最接近的汉宁窗，最后调用函 数h=hd.*win 及 freqz(h,1,512)得到实际汉宁窗的响应和实际滤波器的幅度响应。原理概要：1 利用窗函数法设计 FIR 滤波器 FIR 滤波器的最大特点是其相位特性可以设计为严格的线性，而其幅值可以任意设 置，这样输出波形就不会相位失真。理想低通滤波器的单位取样响应 hd(n)是无限长的，所以要用一个有限长的因果序列 h(n)进行逼近，最有效的方法是截断 hd(n),即用有限长的窗函数 w(n)来截取 hd(n),表示为 h(n)=hd(n)w(n)为获得线性相位的 FIR 滤波器，h(n)必须满足中心对称条件，序列 h(n)应 有一定的延迟a,且a=(N-1)/2 频率响应逼近 hd(ejw)的 FIR 滤波器，最简单的窗函数为矩形窗 0 n(N-1)/2 加窗后的频谱 H(e)二 Hd(e)“W(ejJ 二 1 _W(ej)H(ej_)加窗后使 2 J!实际频响偏离理想频响，影响主要有两个方面：(1)通带和阻带之间存在过渡带，过渡带宽度取决于窗函数频响的主瓣宽度。(2)通带和阻带区间有纹波，这是由窗函数的旁瓣引起的，旁瓣越多，纹波越多。增加窗函数的宽度 N 其主瓣宽度减小，但不改变旁瓣的相对值。为了改善滤波器的性能，要求窗函数的主瓣宽度尽可能窄，以获得较窄的过渡带；旁瓣衰减尽可能大，数量尽可能大，从而改善纹波状况，使实际频响 H(ej3)更好地逼近 理想频响 Hd(ej3)仅供个人参考 50 不得用于商业用途 除了矩形窗外，一般还可以采用以下几种窗函数 五.程序运行结果 N=Battlet:w(n)=丄sin空厲 M sin 临/2)2 二 2 二 X1 汉宁窗：W(n)=0.5WR(J 0.25WR()WR()N-1 N-1 2 辽 2 海明窗 W(n)=0.54WR()0.23WR()WR()N-1 N-1 2 下 2TT 布来克曼窗 W(n)=0.42WR(，)0.23WR()WR()N-1 N-1 其中：WR 是矩形窗的频谱 四.程序 wp=0.5*pi;ws=0.66*pi;wdelta=ws-wp;%N=ceil(8*pi/wdelta)%if rem(N,2)=0 N=N+1;end Nw=N;wc=(wp+ws)/2;%n=0:N-1;alpha=(N-1)/2;m=n-alpha+0.00001;hd=sin(wc*m)./(pi*m);%win=(ha nnin g(Nw);%h=hd.*wi n;%freqz(h,1,512);%过渡带宽度 滤波器长度 理想低通滤波器的截止频率 一个响应 汉宁窗 实际汉宁窗的响应 实际滤波器的幅度响应 仅供个人参考 不得用于商业用途 418 Edi t View Insert Tools Desktop Window Help、裁凰謠民S 0 六结果分析 上图中第一个图为 FIR 滤波器的相频特性，第二个图为幅频特性；根据滤波器特性：通带截止频率处的衰减不大于 3 分贝；阻带衰减不小于 40 分贝。选用汉宁窗。FIR 滤波器 的最大特点是其相位特性可以设计为严格的线性，而其幅值可以任意设置，这样输出波形 就不会相位失真。理想低通滤波器的单位取样响应 hd(n)是无限长的，所以要用一个有 限长的因果序列 h(n)进行逼近，最有效的方法是截断 hd(n),即用有限长的窗函数 w(n)来截取 hd(n)，表示为 h(n)=hd(n)w(n)为获得线性相位的 FIR 滤波器，h(n)必须满足中心对称条件，序列 h(n)应 有一定的延迟a,且a=(N-1)/2 频率响应逼近 hd(ejw)的 FIR 滤波器，最简单的窗函数为矩形窗 1 n(N-1)/2-100 1 1 1 II 1 i i 4 1 1*1111 4 1 1 1 il 1 1 i liii*1 _ i|1 1 I V M V V v H k I 1 -a-1-T 1-厂 li ii 1 i II b I I i P P 1 4 1 1 1 P f H d ii 1 1 i i a i I 1 i I Y A y 八 MM 亠!_ 1 _ _ _ 1-50 0 0 1 0.2 0.3 0.4 0,5 0.6 0.7 03 09 Normalized Frequency rad/sample)1-100D-2000 3000 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 NormNhEEcI Frequency(XK rad/sample)仅供个人参考 不得用于商业用途 加窗后的频谱 H（e）二 Hd（e）“W（e，）二1 二 W（ejr）H（ej（）加窗后使 2-兀 实际频响偏离理想频响，影响主要有两个方面：（3）通带和阻带之间存在过渡带，过渡带宽度取决于窗函数频响的主瓣宽度。（4）通带和阻带区间有纹波，这是由窗函数的旁瓣引起的，旁瓣越多，纹波越多。增加窗函数的宽度 N 其主瓣宽度减小，但不改变旁瓣的相对值。为了改善滤波器的性能，要求窗函数的主瓣宽度尽可能窄，以获得较窄的过渡带；旁瓣衰减尽可能大，数量尽可能大，从而改善纹波状况，使实际频响 H（ejw）更好地逼近 理想频响 Hd（ejw）七.专题实习心得 FIR 有限冲击响应数字滤波器由于设计灵活，在滤波效果和响应时间之间能较好的折 衷，因此在数字信号处理领域应用广泛。用窗函数法设计 FIR 数字滤波器在截断点依旧取模 拟样本冲激响应的全值更利于理解和仿真。但是，窗函数法设计 FIR 数字滤波器的频率特性 与模拟样本频率特性不完全一样，且具有线性相移的特点。窗函数法设计 FIR 数字滤波器是 傅立叶变换的典型运用。通过学习有益于加深对数字滤波器的理解，并提高使用傅立叶变 换分析解决实际问题的能力。八：思考题 1 在这个实验中，你选择哪一种窗函数？试说明原因。选用汉宁窗。根据所给定滤波器特性：通带截止频率处的衰减不大于 减不小于 40 分贝。参照书本上的各滤波器的特性。2、用窗函数法设计 FIR 滤波器时，滤波器的过渡带宽度和阻带衰减和哪些因素有关?窗函数的主瓣宽度决定滤波器过渡带宽度，旁瓣幅值决定阻带衰减。选择窗函数 应使其频谱满足：主瓣宽度尽可能的窄，以使 过渡带尽可能陡；最大旁瓣相 对于主瓣 尽可能小，使能量 尽可能集中在主瓣 内，这样可以使肩峰和波 动减小。3 分贝；阻带衰 仅供个人参考 不得用于商业用途 综合 一、实验内容 录制一段自己的语音信号，时间为 10s 左右，并对录制的信号进行采样；画出采样后 语音信号的时域波形和频谱图；给定滤波器的性能指标，采用窗函数法或双线性变换设计 滤波器，并画出滤波器的频率响应；然后用自己设计的滤波器对采集的语音信号进行滤 波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变 化；回放语音信号；最后，用 MATLAB 设计一信号处理系统界面。二 算法 调用函数 function pushbutton1_Callback(hObject,eventdata,handles)实现一个信号 处理系统界面。选择左键时，用双线性变换法设计滤波器来对信号进行处理，选择右键 时，用窗函数法设计滤波器来对信号进行处理。读取语音信号，对语音信号进行 f=8000 的 频率进行采样，调用函数 y 仁 fft(x1,2048)对所采集的点做 2048 点 FFT 变换。先设计 butterworth 模拟滤波器，再用双线性变换法实现模拟滤波器到数字滤波器的转换。最后 调用函数 f1=filter(bz,az,x2)对加了噪声的语音信号进行滤波，得到滤波后的频谱图。三程序 function pushbutton1_Callback(hObject,eventdata,handles)%hObject handle to pushbutton1(see GCBO)%eventdata reserved-to be defined in a future version of MATLAB%handles structure with handles and user data(see GUIDATA)fs=8000;%语音信号采样频率为 8000 实习 信息工程实习 第五题 语音.wav);A1=0.05;A2=0.10;d=A1*cos(2*pi*3800*t)+A2*sin(2*pi*3600*t);x2=x1+d;wp=0.8*pi;ws=0.85*pi;t=(0:length(x1)-1)/8000;y1=fft(x1,2048);%f=fs*(0:1023)/2048;figure(1)subplot(2,2,1)plot(t,x1)%grid on;axis tight;title(原始语音信号);xlabel(time(s);ylabel(幅度);subplot(2,2,2)plot(f,abs(y1(1:1024)%grid on;axis tight;title(原始语音信号的 FFT 频谱)xlabel(Hz);ylabel(幅度);%双线性变换法设计的巴特沃斯滤波器 对信号做 2048 点 FFT 变换 做原始信号的时域波形 做原始信号的 FFT 频谱 仅供个人参考 不得用于商业用途 Rp=1;Rs=15;Fs=8000;Ts=1/Fs;wp1=2/Ts*tan(wp/2);%ws1=2/Ts*tan(ws/2);N,Wn=buttord(wp1,ws1,Rp,Rs,s);%Z,P,K=buttap(N);%Bap,Aap=zp2tf(Z,P,K);b,a=lp2lp(Bap,Aap,Wn);bz,az=bilinear(b,a,Fs);%H,W=freqz(bz,az);%subplot(2,2,3)plot(W*Fs/(2*pi),abs(H)grid on;axis tight;xlabel(频率(Hz)ylabel(频率响应)title(Butterworth)f1=filter(bz,az,x2);figure(2)subplot(2,2,1)plot(t,x2)%grid on;axis tight;title(滤波前的时域波形);subplot(2,2,2)plot(t,f1);%grid on;axis tight;title(滤波后的时域波形);y3=fft(f1,2048);y2=fft(x2,2048);subplot(2,2,3);plot(f,abs(y2(1:1024);%grid on;axis tight;title(滤波前的频谱)xlabel(Hz);ylabel(幅度);subplot(2,2,4)plot(f,abs(y3(1:1024);%grid on;axis tight;title(滤波后的频谱)xlabel(Hz);ylabel(幅度);将模拟指标转换为数字指标 选择滤波器最小阶数 创建 butterworth 模拟滤波器 用双线性法实现模拟到数字的转换 绘制频率响应曲线 画出滤波前的时域图 画出滤波后的时域图 画出滤波前的频谱图 画出滤波后的频谱图 function pushbutton2_Callback(hObject,eventdata,handles)%hObject handle to pushbutton2(see GCBO)%eventdata reserved-to be defined in a future version of MATLAB%handles structure with handles and user data(see GUIDATA)fs=8000;%语音信号采样频率为 8000 实习 信息工程实习 第五题 语音.wav);t=(0:length(x1)-1)/8000;y1=fft(x1,2048);f=fs*(0:1023)/2048;figure(1)subplot(2,1,1)plot(t,x1)%做原始信号的时域波形 grid on;axis tight;title(原始语音信号);xlabel(time(s);ylabel(幅度);subplot(2,1,2)plot(f,abs(y1(1:1024)%做原始信号的 FFT 频谱 grid on;axis tight;对信号做 2048 点 FFT 变换 仅供个人参考 不得用于商业用途 title(原始语音信号的 FFT 频谱)xlabel(Hz);ylabel(幅度);%窗函数设计滤波器 t=(0:length(x1)-1)/8000;f=fs*(0:2047)/4096;A1=0.05;A2=0.10;d=A1*cos(2*pi*3600*t)+A2*sin(2*pi*3800*t);x2=x1+d;wp=0.8*pi;ws=0.85*pi;wdelta=ws-wp;N=ceil(6.6*pi/wdelta);%取整 wn=(0.8+0.85)*pi/2;bz,az=fir1(N,wn/pi,hamming(N+1);%选择窗函数，并归一化截止频率 figure(2)freqz(bz,az);grid on;axis tight;f2=filter(bz,az,x2);figure(3)subplot(2,2,1)plot(t,x2);grid on;axis tight;title(滤波前的时域波形);subplot(2,2,2)plot(t,f2);grid on;axis tight;title(滤波后的时域波形);y3=fft(f2,4096);f=fs*(0:2047)/4096;figure(8)y2=fft(x2,4096);subplot(2,2,3);plot(f,abs(y2(1:2048);grid on;axis tight;title(滤波前的频谱)xlabel(Hz);ylabel(幅度);subplot(2,2,4)plot(f,abs(y3(1:2048);仅供个人参考 不得用于商业用途 grid on;axis tight;title(滤波后的频谱)xlabel(Hz);ylabel(幅度);四.实验结果 实验结果：Push Button 选择左边的按键:0 1000 2Q00 3000 频率(Hz)原始语音信号 原始语音信号FFT频谱 0 2 0 1 0-0 1-0 2 0 10 20 time(s)Butterworth 12 10 鰹8 1g 6 2 0 1000 2000 3000 Hz 1 8-6/2 QO-O-I Push Button 仅供个人参考 不得用于商业用途 滤波前的频谱 选择右边的按键:原始语音信号 timeja)原始语音信号FFT频谙 滤波前的时域波形 0 10 20 滤波后的时域波形 0 10 20 60 擲 40-g 20 一 虹 0 1000 2000 3000 O 00 1 O O O 2 0 6 10 15 20 25 仅供个人参考 不得用于商业用途 汉明窗：50 O 50-O O-150 0 9 0-5000-10Q0Q 0.1 0.2 0.3 04 0 5 0 6 0.7 Normalised Frequency(XT radsample 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.9 0.9 1 Normalized frequency(XT：radsample)-16000 滤波前的吋域波形 滤波后的时域波形 0 10 20 0 10 20 仅供个人参考 不得用于商业用途 五结果分析 因为录音是自己用电脑自己说话来录音的，所以从语音的频谱可以看出，几乎没有高频 信号，之前所设计的滤波器都是低通滤波，为了再使用低通滤波器，所以人为加了高频信 号，通过低通滤波器之后，滤波后的时域波形对滤波前的时域波形看起来相对光滑一些，滤波后的频谱图，滤除了添加的高频噪声信号，保留了低频信号，语音效果也改善了不 少。通过最后的频谱图可以看出，与之前语音信号的频谱图几乎一样。而通过比较双线性 变换法与窗函数设计的滤波器可以看出，它们的滤波