使用频率采样法设计地频率采样型FIR滤波器课程设计

1、长沙理工大学数字信号处理课程设计报告学院计算机与通信工程专业通信工程班级_学号学生姓名 指导教师 课程成绩 完成日期 课程设计成绩评定学 院计算机与通信工程专 业通信工程班 级学 号学生姓名指导教师完成日期 2015年3月13日指导教师对学生在课程设计中的评价评分项目优良中及格不及格课程设计中的创造性成果学生掌握课程内容的程度课程设计完成情况课程设计动手能力文字表达学习态度规范要求课程设计论文的质量指导教师对课程设计的评定意见综合成绩指导教师签字2015年3月15 日课程设计任务书计算机与通信工程 系通信工程 专业课程名称数字信号处理课程设计时间 20142015学年第二学期12周学生姓名指导

2、老师题目语音信号滤波去噪一使用频率米样法设计的频率米样型FIR滤波器主要内容：用麦克风米集一段语音信号，绘制波形并观察其频谱特点，加入一个带 外单频噪声，使用频率米样法设计 个满足指标的 FIR滤波器，对该含噪语音 信号进行滤波去噪处理，画出频率采样型滤波器结构图。比较滤波前后的波形 和频谱并进行分析，根据结果和学过的理论得出合理的结论。 与不冋滤波方法 的同学比较各种滤波方法性能的优劣。要求：（1）滤波器指标必须符合工程实际。（2）设计完后应检查其频率响应曲线是否满足指标。（3 ）处理结果和分析结论应该一致，而且应符合理论（4）独立完成课程设计并按要求编写课程设计报告书 应当提交的文件：（1

3、 ）课程设计学年论文。（2）课程设计附件（主要是源程序）。语音信号滤波去噪一一使用频率采样法设计的频率采样型FIR滤波器学生姓名：指导老师：摘要此次课程设计是使用频率采样法设计一个频率采样型滤波器，对一段加噪的语音信号进行滤波去噪处理，并且对滤波前后的时域波形和频谱分析滤波性能。此次课程设计仿真平台为MATLAB7.0，开发工具是M语言编程。首先利用Windows录音机录 制一段语音信号，并人为加入单频噪声，然后对信号进行频谱分析以确定所加噪声频率， 再设计滤波器进行滤波去噪处理，最后比较滤波前后的波形和频谱并进行分析，并用 VISIO绘制滤波器的结构图。由分析结果可知，滤波后的语音信号与原始

4、语音信号基本 一致，即设计的FIR滤波器能够去除信号中所加单频噪声，即达到了设计目的。关键词 MATLAB ;滤波去噪；频率采样法；频率采样型滤波器；VISIO1引言滤波去噪是信号处理中一种最基本但十分重要的技术。利用滤波可以从复杂的信 号中提取所需的信号，抑制不需要的信号。滤波器就是这样一种可以在时域和频域对信 号进行滤波处理的系统。通常情况下，有用信号和干扰信号是在不同频段上的，于是通 过对滤波器的频率特性精心设计就能达到滤波的目的。本课程设计是采用频率采样法设计频率抽样型滤波器，从而对语音信号滤波去噪。通过对比滤波前后的波形图及回放滤波前后的语音信号，来判断滤波器对噪声信号确实有滤除作用

5、。1.1课程设计目的(1) 熟悉使用MATLAB ；(2) 了解FIR滤波器原理及结构；(3) 利用所学数字信号处理相关知识用 MATLAB设计一个FIR滤波器；提高自己动手的能力。(5)对加噪的语音信号进行滤波去噪处理，比较滤波前后的时域波形和频谱并进行分 析；1.2课程设计的要求(1) 滤波器指标必须符合工程设计。(2) 设计完后应检查其频率响应曲线是否满足指标。(3) 处理结果和分析结论应该一致，而且应符合理论。(4) 独立完成课程设计并按要求编写课程设计报告。1.3设计平台本课程设计仿真平台为 MATLAB7.0。MATLAB的名称源自Matrix Laboratory ，1984年由

6、美国Mathworks 公司推向市场。它是一种科学计算软件，专门以矩阵的形 式处理数据。MATLAB将高性能的数值计算和可视化集成在一起，并提供了大量的内 置函数，从而被广泛地应用于科学计算、控制系统、信息处理等领域的分析、仿真和设 计工作。1993年MathWorks 公司从加拿大滑铁卢大学购得 MAPLE软件的使用权，从而以MAPLE为“引擎”开发了符号数学工具箱(Symbolic Math Toolbox ) 22设计原理用windows自带的录音机录制一段语音信号，绘制波形并观察其频谱，给定相应 技术指标，用频率采样法设计的一个满足指标的频率采样型FIR滤波器，对该语音信号进行滤波去噪

7、处理，比较滤波前后的波形和频谱并进行分析。2.1 FIR滤波器的设计FIR(Finite Impulse Response)滤波器3:有限长单位冲激响应滤波器，又称为非递归型滤波器，是数字信号处理系统中最基本的元件，它可以在保证任意幅频特性的同 时具有严格的线性相频特性，同时其单位抽样响应是有限长的，因而滤波器是稳定的系 统。因此，FIR滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。有限长单位冲激响应(FIR)滤波器有以下特点：(1) 系统的单位冲激响应h (n)在有限个n值处不为零(2) 系统函数H(z)在|z|>0处收敛，极点全部在z = 0处(因果系统)(3) 结构上主要

8、是非递归结构，没有输出到输入的反馈，但有些结构中(例如频率 抽样结构)也包含有反馈的递归部分。2.2频率采样型结构把一个有限长序列(长度为N点)的z变换H (z)在单位圆上作N等分抽样，就得 到H (k)，其主值序列就等于h (n)的离散傅里叶变换H (k)。那里也说到用H (k)表示的 H (z)的内插公式为H(z)(1N、1 N 1 H (k)N k o1 Wn z 1(2.1)这个公式就为FIR滤波器提供了另外一种结构，这种结构由两部分级联组成。1 N 1H(z) Hc(z)-Hk(z)(2.2)N k o其中级联的第一部分为梳状滤波器，其结构如下图所示：Hc(z)(1 z N)(2.3

9、)X(n)仁Y(n)m-Z -N图2-1梳状滤波器结构图第二部分由N个谐振器组成的谐振柜。它是由N个一阶网络并联组成，而这每一个一阶网络都是一个谐振器Hk(z)1 WNkz1(2.4)其结构如下图所示：H(k)Hk(z>图2-2 一阶谐振器频率抽样型结构特点：2(1) 它的系数H(k)直接就是滤波器在Wk k处的频率响应。因此，控制滤波器的N频率响应是很直接的。(2) 结构有两个主要缺点：(a) 所有的相乘系数及 H(k)都是复数，应将它们先化成二阶的实数，这样乘起来较复杂，增加乘法次数，存储量。(b) 所有谐振器的极点都是在单位圆上，由WNk决定考虑到系数量化的影响，当系数量化时，极点

10、会移动，有些极点就不能被梳状滤波器的零点所抵消。(零点由延时单元决 定，不受量化的影响)系统就不稳定了。(3) 将一阶网络合并为二阶网络(a)第k和第N-k个谐振器合并为一个实系数的二阶网络，因为h(n)是实数，它的DFT也是圆周共轭对称的(2.5)H(k) H*(N k) k 1,2,3丄 N 1因此，可以将第k和第N-k个谐振器合并为一个二阶网络Hk(z)H(k)1 rWNkz 1H (N k)H (k)(NT)1k11 Wnz 1 Wn z*H (k)1 rWNk*z 1H(k)1 z 1Wn k* WNkr r2WNkWNk*z210k 1kz1 z 12r cos(2 k)N其中：0

11、k 2ReH( k), 1k2rReH( k)W,(2.6)(b)第k和第N-k个谐振器合并为一个二阶网络的极点在单位圆内，而不是在单位圆上，因而从频率响应的几何解释可知，它相当于一个有限Q的谐振器。其谐振频率Wk 为：丨 .f1 z1丁 图2-3二阶网络结构图除了共轭复根外，还有实根。当N=偶数时，有一对实根,H°(z)0它们分别为严和H1 rz 1(2.7)当N=奇数时，只有一个实根z=r（k=O ）,即只有HO（z）.（c）修正频率抽样结构流图（N=偶数）x(n)图2-4修正频率抽样结构流图（N=偶数）H(z)N N、(1 r z )1 H(0)N1 rz 1H(j)N/211

12、 rz 1Ok11kzz12rco畤k)(2.8)修正频率抽样结构流图（N=奇数）x(n)7y(n),v -y-rz2rcos()"nJ1-1-rz .Ar2rcos().V/(0J图2-5修正频率抽样结构流图（N=奇数）N N、1H(z) (1 r z)1H(0)1 rz 1N/2 11Ok 1kz2k 1 1 z 12r cos(k)N(2.9)2.3频率采样法设计思路:这种设计方法是从频域进行设计的一种方法，首先给定一个希望逼近的频率响应。H(k)Hd(ej )| 2k 0,1,2,., N 1(2.10)kN知道H(k)后,由IDFT定义，可以用这N个采样值H(k)来惟一确定

13、有限长序列h(n)， 即N 11nkh(n)H(k)WNn 0,1,2,., N 1(2.11)N k oN 1H (z)h(n)z n(2.12)n 0N 1H (ejw ) h(n )e jwn(2.13)n 0内插公式H(ej )N 12H (k)kk 0N(2.14)式中,0(3 )是内插函数()sin( N / 2) j (n 1)/2eN sin(/2)(2.15)四种线性相位的FIR滤波器如下表2-1所示表2-1四种线性相位的FIR滤波器类型甲®)1型及关于疋 亦偶对称®厘型刃山庆于=Ck Qjt偶对茹*关于tu=ic奇对称2型仙戶也Mid). M旗迫做血关于

14、=0、JL亦奇时称V- 1J74型辺)关于-0、加苛对称关于厲=討禺对称2 2逼近误差及其改进措施这样设计所得到的频响H(ej 3)与要求的理想频响Hd(ej 3)会有一定的差别(1)采样点上滤波器的实际频率响应是严格地和理想频率响应数值相等的。但是在采样点之间的频响则是由各采样点的加权内插函数的延伸叠加而成的,因而有一定的逼近误差，误差大小取决于理想频率响应曲线形状。(3) 理想频率响应特性变化越平缓，则内插值越接近理想值，逼近误差越小。(4) 如果采样点之间的理想频率特性变化越陡，则内插值与理想值的误差就越大，因而在理想频率特性的不连续点附近，就会产生肩峰和起伏。滤波器性能的改善：(1)

15、增加过渡带采样点，它可以大大减少振荡，阻带衰减也可以得到进一步改善。一般一点到二点的过渡带采样即可得到满意的结果。(2) 增加采样点密度，过渡带的宽度与采样点数N成反比。但N值意味着或长度的增加，滤波器运算量必然增大4。3设计步骤3.1设计流程图录制语音信号ht.wav。并用wavread函数采集加入2500Hz的噪声。进行快速傅里叶变换，绘制加噪前后的频谱对比图，和时域对比图设定滤波器的性能指标回放语音信号图3-1流程图3.2录制语音信号使用电脑录音机录制一段语音信号属性设置为8000Hz ,8位单声道的音频格式。如图3.1-2所示:用频率采样法设计频率采样型的FIR滤波器丿是否、符合标准/

16、VN用设计的滤波器进行滤波处理， 回放滤波后的信号，并保存1比较滤波前后语音信言号的波形和频谱图3-2语音信号设置然后在MATLAB软件平台下，利用函数 wavread对语音信号进行采样，源程序为：x,fs,bits=wavread('ht.wav')，记住采样频率和采样点数，matlab 实现得：fs=8000 ；bits=8。3.3语音信号的频谱分析在MATLAB中编写加噪之后时域，频域对比图的源程序如下：x,fs,bits=wavread('D:ht.wav');%输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样

17、率，bits是波形文件每样本的编码位数。sound(x,fs,bits); %按指定的采样率和每样本编码位数回放N=length(x); %计算信号x的长度fn=2500;%单频噪声频率，此参数可改t=0:1/fs:(N-1)/fs; %计算时间范围，样本数除以采样频率x=x' y=x+0.1*si n(fn*2*pi*t); % 加噪sound(y,fs,bits); %应该可以明显听出有尖锐的单频啸叫声X=abs(fft(x); Y=abs(fft(y); %对原始信号和加噪信号进行fft变换，取幅度谱X=X(1:N/2); Y=Y(1:N/2);%截取前半部分deltaf=fs/

18、N;%计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf;%计算频谱频率范围运行结果如下图所示:图3-3加噪前后时域频域对比图由图3.3-1可以看出，在频域为2500Hz处加入了一个单频噪声，而加入噪声之后, 时域的波形出现了明显失真，通过听取原声音信号 x和加噪声音信号y，可以明显听到 y语音信号中有一明显尖锐噪声。3.4滤波器设计设计一个带阻滤波器，将单频信号滤出去，源程序如下所示：fpd=2400;fsd=2450;fsu=2550;fpu=2600;Rp=0.5;As=45;% 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=mi n

19、( (fsd-fpd),(fpu-fsu);%计算上下边带中心频率fcd fcu，和频率间隔dfwp=fcd/fs*2*pi;ws=fcu/fs*2*pi;dw=df/fs*2*pi;% 将 Hz 为单位的模拟频率换算为rad为单位的数字频率M=41;alpha=(M-1)/2;%设置抽样点数l=0:M-1;wl=2*pi*l/M;T1=0.6025,T2=0.127% 过渡带的抽样值NP仁fix(M*fpd/8000);NS仁fix(M*fsd/8000);NS2=ceil(M*fsu/8000);NP2=ceil(M*fpu/8000);%各边界频率的抽样值 1Hrs=o nes(1,NP

20、1-2),0.6025,0.127,zeros(1,NP2-NP1+1),0.127,0.6025, on es(1,M-2 *NP2-5),0.6025,0.127,zeros(1,NP2-NP1+1),0.127,0.6025,o nes(1,NP1-3);Hdr=1,1,0,0,1,1;% 带阻滤波器的参数设置wdl 仁 fpd*2/fs;wdl2=fsd*2/fs;wdl3=fsu*2/fs;wdl4=fpu*2/fs;%rad wdl=0,wdl1,wdl2,wdl3,wdl4,1;k1=0:floor(M-1)/2);k2=floor(M-1)/2)+1:M-1;an gH=-al

21、pha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2);H=Hrs.*exp(j*a ngH);h=real(ifft(H,M);db,mag,pha,grd,w=freqz_m(h,1);Hr,ww,a,L=hr_type2(h);程序运行结果如下图3.3-2所示:M ig j re 1Fil Edit Vi曹训n榊rt 口Qk 总衣hfiop' WwHfow 屮*P."：' X固 i 旦遞液器覃怔冲击响N图3-4滤波器参数图由图3-4可以看出，滤波器的衰减大于设定值As=45,满足性能指标，滤波器的衰减可由增加过度带宽来得到。上图中的抽样值M二

22、41，因此阻带的带宽比较宽。可以通过调大M值来使得阻带的带宽变窄，但这样会使得滤波器的阶数相应的增加，使得 其结构图标的复杂，为了是滤波器能达到相应的性能要求，又使得其结构变得简单，最 终取M值为41 03.5信号滤波处理原程序如下所示：%进行滤波处理z=filter(h,1,y); %用设计好的滤波器对y进行滤波Z=fft(z);Z=Z(1:length(Z)/2);% 计算频谱取前一半程序运行结果如下图3.3-3所示:图3-5滤波前后时域，频域对比图频谱取前一半由上图3-5可以看出，在滤波之后时域图能得到恢复，频域图中的单频噪声信号也的到滤除，说明了设计的滤波器能滤除加入的噪声信号，因此说

23、达到了 设计的要求。3.6结果分析开始通过分析决定设计一个带阻的滤波器来滤除加入的单频噪声，根据噪声的频率来设计阻带的范围。在采集语音信号后，按照步骤用频率采样法设计频率采样型滤波器。由图3-4可知,设计的滤波器达到要求。我们观察到图3-5滤波前后语音信号的波形对比图，发现时域波形中加干扰噪声后有明显的变化，不过经过滤波后几乎没有变化，说明设计的滤波 器达到要求。再通过听取原始语音信号，加噪信号 y,和滤波之后的信号z。对比之后, 发现滤波器确实滤除了噪声。从理想的角度考虑，该带阻滤波器的阻带带宽应该可以变得更窄， 让滤波效果更好, 但是这样的采样值会变得非常大，考虑到实际的情况，通过权衡之后

24、，决定牺牲带宽来 使得滤波器的阶数降低，因此在上图3-5中我们可以看出，在噪声频谱左右两边的信号 也被滤除了。3.6滤波器结构设计通过计算得出滤波器的系数如下图所示:B -2. QQdQ1, 0000-E97661. OOOQ乱 gotof10.W11. 000 Q-1.S0681.00002.叽 9083-0,98831. 0000-1.79231. 00002.0000-0-97370. 9737昇 0000-L«36»1.00002. ODOQ0. 95347. 9534】* 0000-U44101,00002. 0000-0.»2750.A275L 000

25、0-1. 212&1,00002, 00000, 8»«27” B962.0000-0.95541. 00002. ODCO-0, 8&960,95961. 0000-0.4760l.ooao2- 00000. 8179-0.町"h 0D00-CL 38081. 0000-0*77151,77161” 00000,07661. 0000队 2540叭 72057+ 7205L 00000.22941. 00000 . 0000663Q. 6653.00000.53001.0000a. oata0. 3E62» »«h

26、0D000.81B11, 0000此 oodo(L 00000.»1351,0000Loan1.0000仇 26400. 117+ "酣1.00001.3307E_ 0000U2050-0.40110” 40911. 00001.54301.00002. 00000.-0.33M1.00<0L719I1.00002. 0000一0. 2«500. 2临益0I.OUOQL>5501.00002.00(10-0*10041. 00001.94741.00002-0000-0. 11470, IU71. 00001.99H1.00D01. 0000S 283

27、-0,0383L 0000-LOQOO0图3-6滤波器系数根据以上求出的系数，利用 VISIO软件来绘制滤波器的结构图如下图所示:Z1Z1-11-1Z1-11.6359Z1-10.9275 Az1z1才-0.89621.2125-1.1*z1-10.3808 z1.1-1z1z1-0.6653 0.25400.665307-1z1L-1.7191z1-1z111-11.0000z1-0.7205z1 -0.9534 I-11.9068Z1 -0.98830.98832.0000-0.7715 2.00000.95342.00001.44104-1-0.26502.0000-0.4091 0.2

28、540 b 0.9971 i0.89622.00000.2650z1>1.97661.7923Z1 0.97370.33801.2050-0.8596 2.0000Z10.8596-0.9737 2.0000-0.8179J-0.99712.0000>0.81792.0000Z1 0.4091 1z1-0.92752.00000.72051.2050-0.11472.0000/ 0.1147：0.9554-12.0000 0.0383z1Z-0.0383I -1.9474z10.19042.0000-0.1904 Ay（n）-0.2294-1Z1-1.5430 z1 -033801

29、-1-1.8550-1Z1F1图3-7滤波器的结构图4出现的问题及解决办法在这次课程设计当中，由于基础不扎实，出现了很多问题，即有MATLAB软件操作不当，也有知识掌握程度不够出现的各种问题。1，在一开始录制语音信号时，没有修改参数值，导致频率抽样过高，后在老师的 提醒之下将频率修改为8000hz。2，在调用音乐文件时没有将文件放在 MATLAB的工作文件夹下面，到时文件找不 到，还有程序中用到的各种函数都没有放入正确的位置。3，在绘制加噪前后频率对比图时，留白过多，对比不明显，通过axis函数对横纵坐标进行限定4，在使用visio软件绘制滤波器的结构图时，没有合理规划整个页面，导致重新绘 制

30、。5，在使用visio绘制结构图时，前面滤波器参数设置出现错误，所以滤波器的系数 发生了改变，只能重新修改参数。6，在设计滤波器的过程中，没有准确理解每一条指令代表的含义，导致程序前后 不对应，出现很多错误，通过 matlab中的错误提示，准确找到错误的那一行代码，进 行修改。7,在设计滤波器的过程当中，滤波器的衰减小于开始所设置的值，通过牺牲过渡 带和调节过渡带的采样值，即T1和T2来使得衰减大于所设定值As。8，频谱图没有进行去前半段的处理。9，最后听取滤波后声音，没有将其保存，在请教了同学之后，学会了如何保存滤 波后的声音文件。10,在设计滤波器的过程当中，抽样点数过高，是的滤波器的阶数

31、很高，导致绘图 困难。通过降低抽样的点数，来平衡滤波器阶数和阻带带宽之间的矛盾。5结束语这是第三次进行课程设计了，在前面的课程设计当中我们学习到了如何使用 matlab，所以对于matlab软件的使用并没有那么陌生了，尽管如此，在使用 matlab 的过程当中还是出现了很多错误，比如说忘记添加函数文件，参数前后不对应之类的低 级错误。在这次课程设计当中，还学习到了另外一个绘图软件，VISIO，相对而言，是比较容易操作的，在这次课程设计中主要用来绘制各类结构图，以及流程图。此次课程设计，让我更深入的了解到频率抽样法以及频率抽样型的滤波器，由开始 的无从下手，再翻阅书上的例子，首先决定使用低通的滤

32、波器将噪声滤除，但由于是加 入的是一个单频的信号噪声，就决定使用带阻滤波器，在网上查找相关资料时，大多是 低通滤波器，没有带阻滤波器的例子，所以就通过书籍信号与系统上高通的例题来 进行修改，当然这过程不是一蹴而就的。经过了一次又一次的错误，才得出来最后的模 型。还有这也是我和其他同学一起讨论出来的。这次课程设计让我对滤波器的类型有了一个更加完整了解，在设计中也使我对一些概念有了更深刻的认识。例如：在滤波器分类方面，我深刻的了解了低通，高通滤波器 与带通，带阻滤波器的特性区别。还有在课程设计中每一次的数据输入都有其重要意义, 用MATLAB编译程序时，可以根据滤波器指标的要求实时知道对滤波器的影

33、响。通过 一次次的调试和权衡使滤波器的性能达到最佳。课程设计不仅要求对滤波器理论的研究，更重要的是培养一种遇到问题解决问题 的思维。因为有了这次课程设计，我懂得了书本知识只是实际应用的理论指导。如果仅 仅只学习书本知识，不去在实践中运用，那只是停留在只知其表，不知其因的层面。比 如在数学计算上，可以将噪声完全滤除，而在这次设计中，若要完全滤除噪音，滤波器 的阶数就会增高，在现实生活中是很难实现的，所以噪声是不能完全滤除的。课程设计结束了，我相信这次课程设计对今后的学习是很有帮助的，它让我将理论 更好地和实践相结合，提高了动手的能力，也填补了自己学习上的一些不足。这次课程设计能取得成功，不仅仅是

34、我一个人的努力的结果，更离不开指导老师与 同学的帮助，在此向老师和同学们表示衷心的感谢。参考文献1 吴镇扬，数字信号处理M.高等教育出版社，20042 张圣勤.MATLAB7.0 实用教程M.北京：机械工程出版社，20063 程佩青数字信号处理教程M 北京：清华大学出版社，20024 高西全，丁玉美.数字信号处理M.第三版.西安：西安科大出版社，1994附录：源程序%程序名称：htl.m%程序功能：设计一频率采样法设计的频率采样型的滤波器，然后对一加单凭频噪声信 号去噪%程序作者：贺婷%最后修改日期：2015-3-10x,fs,bits二wavread('D:ht.wav');

35、%输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。 sound(x,fs,bits);%按指定的采样率和每样本编码位数回放N=length(x); %计算信号x的长度fn=2500;%单频噪声频率，此参数可改t=0:1/fs:(N-1)/fs; % 计算时间范围，样本数除以采样频率x=x' y=x+0.1*si n(fn*2*pi*t); % 添加噪声sound(y,fs,bits); %可以明显听出有尖锐的单频啸叫声X=abs(fft(x);Y二abs(fft(y);%对原始信号和加噪信号进行fft变

36、换,取幅度谱X=X(1:N/2); Y=Y (1:N/2);% 截取前半部分deltaf=fs/N;%计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf;%计算频谱频率范围%设计滤波器fpd=2400;fsd=2450;fsu=2550;fpu=2600;Rp=0.5;As=45;% 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df二min(fsd-fpd),(fpu-fsu);% 计算上下边带中心频率fcd fcu，和频率间隔dfwp=fcd/fs*2*pi;ws=fcu/fs*2*pi;dw=df/fs*2*pi;% 将 Hz 为单位的模拟

37、频率换算为rad为单位的数字频率M=41;alpha=(M-1)/2;%设置抽样点数l=0:M-1;T1=0.6025,T2=0.127% 过渡带的抽样值NP1= fix(M*fpd/8000);NS 1= fix(M*fsd/8000);NS2二ceil(M*fsu/8000);NP2二ceil(M*fpu/8000);%各边界频率的抽样值1Hrs=o nes(1,NP1-2),0.6025,0.127,zeros(1,NP2-NP1+1),0.127,0.6025,on es(1,M-2*NP2-5),0.6025,0.127,zeros(1,NP2-NP1+1),0.127,0.6025

38、,on es(1,NP1-3);k1=0:floor(M-1)/2);k2=floor(M-1)/2)+1:M-1;an gH=-alpha*(2*pi)/M*k1,alpha*(2*pi)/M*(M-k2);H=Hrs.*exp(j*a ngH);h=real(ifft(H,M);% 求 FIR 单位冲击响应db,mag,pha,grd,w=freqz_m(h,1);% 幅频特性曲线Hr,ww,a,L=hr_type2(h);%画图subplot(2,2,1);stem(l,h);title(' 滤波器 单位冲 击响应');xlabel(' n');ylabe

39、l('y( n)');axis(O 40 -1 1);grid on; subplot(2,2,2);plot(w/pi,pha);title('滤 波器相 位响 应');xlabel('w/pi');ylabel('相位 pha');axis(0 1-4 4);grid on;subplot(2,2,3);plot(w/pi,db);title('滤 波 器 幅 度 响 应');xlabel('w/pi');ylabel('dB');axis(O 1 -50 20);grid on

40、;%画线X_l二0,0,wp/pi,ws/pi;1,1,wp/pi,ws/pi;Y_l二-As,-Rp,-50,-50;-As,-Rp,20,20;%在wp , ws , Rp , As处画线以更直观判断设计是否达标，每列参数是 每个线条的端点坐标line(X_l,Y_l,'Color','r','LineWidth',2,'LineStyle','-')% 添加线宽为 2的红色线条subplot(2,2,4);plot(w/pi,mag);title('滤 波器幅 度响应');xlabel('w/pi');ylabel(' 幅度 mag');axis(0 1 0 1.2);grid on;%进行滤波处理z=filter(h,1,y); % 用设计好的滤波器对x1进行滤波Z=fft(z);Z=Z(1:length(Z)/2);% 计算频谱取前一半%比较滤波前后的性能%在同一张大图里分别绘制原始信号x，加噪信号y，滤波去噪信号z的时域波形和频谱，以便比较和分析。figure(2);subplot(321);plot(t,x);axis tight;grid on;title('原始语音信号时间x