DSP语音信号处理课程设计最新版
目录
第一章绪论 (1)
1.1 MATLAB简介 (1)
1.2 数字信号处理简介 (1)
1.3语音信号处理简介 (2)
1.4 GUI简介 (2)
第二章方案论证 (4)
2.1 论证方案 (4)
2.2语音的录入与打开 (4)
第三章系统的总体设计 (5)
3.1 FFT的MATLAB实现 (5)
3.2设计原理 (6)
3.2.1 运用自相关方法估计语音信号的声道参数原理 (6)
3.2.2解决噪声污染的原理 (6)
3.2.3 自相关检测原理 (8)
3.2.4 中心消波法检测的原理 (8)
3.2.5 三电平中心消波法原理 (8)
第四章模块实现过程 (10)
4.1 DFT和DTFT设计 (10)
4.1.1 DFT定点分析 (10)
4.1.2 DTFT设计 (11)
4.2 滤波器设计 (13)
4.2.1原理 (13)
4.2.2简介IIR数字滤波器和FIR数字滤波器 (13)
4.2.3设计内容 (14)
4.3 分离观察 (19)
第五章总结与展望 (21)
参考文献 (22)
附录 (23)
第一章绪论
1.1 MATLAB简介
MATLAB的名称源自Matrix Laboratory，它是一种科学计算软件，专门以矩阵的形式处理数据。MATLAB 将高性能的数值计算和可视化集成在一起，并提供了大量的内置函数，从而被广泛地应用于科学计算、控制系统、信息处理等领域的分析、仿真和设计工作，而且利用MA TLAB 产品的开放式结构，可以非常容易地对MATLAB 的功能进行扩充，从而在不断深化对问题认识的同时，不断完善MATLAB 产品以提高产品自身的竞争能力Matlab的数据分析和处理功能十分强大，运用它来进行语音信号的分析、处理和可视化相当便捷。MATLAB是一种科学计算软件,主要适用于矩阵运算和信息处理领域的分析设计,它使用方便,输入简捷,运算高效,内容丰富,并且很容易由用户自行扩展。MATLAB当前已成为美国和其他发达国家在大学教学和教学研究中最常用而必不可少的工具。其主要有信号处理(signal processing)、系统控制(control system)、神经网络(neural network)、图形处理(image processing)、鲁棒控制(robust control)、非线形系统控制设计(nonlinear control)、系统辨识(system identification)、最优化(optimisation)、µ分析与综合（µ analysis and synthesis）、模糊逻辑(fuzzy logic)、小波(wavlet)、样条(spline)等工具箱,而且工具箱还在不断增加。这些工具箱给各个领域的研究和工具应用提了有利的工具。借助于这些“巨人肩上的工具”，各个层次研究人员可直观、方便地进行分析计算及设计工作,从而大大地节省时间。
1.2 数字信号处理简介
数字信号处理的主要研究对象是数字信号，且是采用运算的方法达到处理的目的。因此，其实现方法，基本上分成两种实现方法，即软件和硬件实现方法。软件实现方法指的是按照原理和算法，自己编写程序或者采用现成的程序在通用计算机上实现，硬件实现指的是按照具体的要求和算法，设计硬件结构图，用乘法器加法器延时器、控制器、存储器以及输入输出接口部件实现的一种方法。显然前者灵活，只要改变程序中的有关参数，但是运算速度慢，一般达不到实时处理，因此，这种方法适合于科研和教学。后者运算速度快，可以达到实时处理要求，但是不灵活。目前DSP芯片已进入市场，且正在高速发展，速度高，体积小，性能优良，价格也在不断下降。可以说，用DSP芯片实现数字信号处理，正在变成工程技术领域的主要方法。用合适的DSP芯片，配有合适的芯片语言及任务要求的软件，来实现信号处理功能是一种最佳的数字信号处理系统。
1.3语音信号处理简介
语音信号处理是研究用数字信号处理技术和语音学知识对语音信号进行处理的新兴的学科，是目前发展最为迅速的信息科学研究领域的核心技术之一。通过语音传递信息是人类最重要、最有效、最常用和最方便的交换信息形式。同时，语言也是人与机器之间进行通信的重要工具，它是一种理想的人机通信方式，因而可为计算机、自动化系统等建立良好的人机交互环境，进一步推动计算机和其他智能机器的应用，提高社会的信息化和自动化程度。
语音信号处理是一门新兴的学科，同时又是综合性的多学科领域和涉及面很广的交叉学科。虽然从事这一领域研究的人员主要来自信号与信息处理及计算机应用等学科，但是它与语音学、语言学、声学、认知科学、生理学、心理学等许多学科也有非常密切的联系。
20世纪60年代中期形成的一系列数字信号处理的理论和算法，如数字滤波器、快速傅立叶变换（FFT）等是语音信号数字处理的理论和技术基础。随着信息科学技术的飞速发展，语音信号处理取得了重大的进展：进入70年代之后，提出了用于语音信号的信息压缩和特征提取的线性预测技术（LPC），并已成为语音信号处理最强有力的工具，广泛应用于语音信号的分析、合成及各个应用领域，以及用于输入语音与参考样本之间时间匹配的动态规划方法；80年代初一种新的基于聚类分析的高效数据压缩技术—矢量量化（VQ）应用于语音信号处理中；而用隐马尔可夫模型（HMM）描述语音信号过程的产生是80年代语音信号处理技术的重大发展，目前HMM已构成了现代语音识别研究的重要基石。近年来人工神经网络(ANN)的研究取得了迅速发展，语音信号处理的各项课题是促进其发展的重要动力之一，同时，它的许多成果也体现在有关语音信号处理的各项技术之中。
语音信号处理技术的应用极其广泛，涉及工业、军事、交通、医学、民用等各个领域，其中最重要的包括语音编码、语音合成、语音识别以及语音增强等。1.4 GUI简介
图形用户界面（Graphical User Interface，简称GUI，又称图形用户接口）是指采用图形方式显示的计算机操作用户介面。与早期计算机使用的命令行界面相比，图形界面对于用户来说在视觉上更易于接受。
如Windows是以图形界面方式操作的，因为你可以用鼠标来点击按钮来进行操作，很直观。而DOS就不具备GUI，所以他只能输入命令。DOS 的这种界面叫CUI (Command line User Interface ) 命令行模式的人机接口。
GUI 是Graphical User Interface 的简称，即图形用户界面，通常人
机交互图形化用户界面设计经常读做“goo-ee”，准确来说GUI 就是屏幕产品的视觉体验和互动操作部分。
GUI 是一种结合计算机科学、美学、心理学、行为学，及各商业领域需求分析的人机系统工程，强调人—机—环境三者作为一个系统进行总体设计。
这种面向客户的系统工程设计其目的是优化产品的性能，使操作更人性化，减轻使用者的认知负担，使其更适合用户的操作需求，直接提升产品的市场竞争力。
GUI 即人机交互图形化用户界面设计。纵观国际相关产业在图形化用户界面设计方面的发展现状，许多国际知名公司早已意识到GUI 在产品方面产生的强大增值功能，以及带动的巨大市场价值，因此在公司内部设立了相关部门专门从事GUI 的研究与设计，同业间也成立了若干机构，以互相交流GUI 设计理论与经验为目的。随着中国IT 产业，移动通讯产业，家电产业的迅猛发展，在产品的人机交互界面设计水平发展上日显滞后，这对于提高产业综合素质，提升与国际同等业者的竞争能力等等方面无疑起了制约的作用。
GUI的广泛应用是当今计算机发展的重大成就之一，他极大地方便了非专业用户的使用人们从此不再需要死记硬背大量的命令，取而代之的是可用用通过窗口、菜单、按键等方式来方便地进行操作。而嵌入式GUI具有下面几个方面的基本要求：轻型、占用资源少、高性能、高可靠性、便于移植、可配置等特点。
GUI 准则：
（1）减少用户的认知负担
（2）保持界面的一致性。
（3）满足不同目标用户的创意需求
（4）用户界面友好性
（5）图标识别平衡性
（6）图标功能的一致性
（7）建立界面与用户的互动交流
第二章方案论证
2.1 论证方案
利用MATLAB中的wavread命令来读入（采集）语音信号，将它赋值给某一向量。再将该向量看作一个普通的信号，对其进行FFT变换实现频谱分析，再依据实际情况对它进行滤波。对于波形图与频谱图（包括滤波前后的对比图）都可以用MATLAB画出。我们还可以通过sound命令来对语音信进行回放，以
便在听觉上来感受声音的变化。
选择设计此方案，是对数字信号处理的一次实践。在数字信号处理的课程学习过程中，我们过多的是理论学习，几乎没有进行实践方面的运用。这个课题正好是对数字语音处理的一次有利实践，而且语音处理也可以说是信号处理在实际应用中很大众化的一方面。
这个方案用到的软件也是在数字信号处理中非常通用的一个软件——MATLAB软件。所以这个课题的设计过程也是一次数字信号处理在MATLAB中应用的学习过程。课题用到了较多的MATLAB语句，而由于课题研究范围所限，真正与数字信号有关的命令函数却并不多。
2.2语音的录入与打开
[y,fs,bits]=wavread('Blip',[N1 N2]);用于读取语音，采样值放在向量y 中，fs表示采样频率(Hz)，bits表示采样位数。[N1 N2]表示读取从N1点到N2点的值（若只有一个N的点则表示读取前N点的采样值）。
sound(x,fs,bits); 用于对声音的回放。向量y则就代表了一个信号（也即一个复杂的“函数表达式”）也就是说可以像处理一个信号表达式一样处理这个声音信号。
第三章系统的总体设计
3.1 FFT的MATLAB实现
在MATLAB的信号处理工具箱中函数FFT和IFFT用于快速傅立叶变换和逆变换。下面介绍这些函数。
函数FFT用于序列快速傅立叶变换。
函数的一种调用格式为y=fft(x)其中，x是序列，y是序列的FFT，x可以为一向量或矩阵，若x为一向量，y是x的FFT。且和x相同长度。若x为一矩阵，则y是对矩阵的每一列向量进行FFT。
如果x长度是2的幂次方，函数fft执行高速基－2FFT算法；否则fft执行一种混合基的离散傅立叶变换算法，计算速度较慢。
函数FFT的另一种调用格式为 y=fft(x,N)式中，x，y意义同前，N为正整数。函数执行N点的FFT。若x为向量且长度小于N，则函数将x补零至长度N。若向量x的长度大于N，则函数截短x使之长度为N。若x 为矩阵，按相同方法对x进行处理。经函数fft求得的序列y一般是复序列，通常要求其幅值和相位。MATLAB提供求复数的幅值和相位函数：abs，angle，这些函数一般和FFT 同时使用。
函数abs(x)用于计算复向量x的幅值，函数angle(x)用于计算复向量的相角，介于和之间，以弧度表示。
函数unwrap(p)用于展开弧度相位角p ,当相位角绝对变化超过时，函数把它扩展至。
（2）频率计算
若N点序列x(n)(n=0,1,…,N-1)是在采样频率下获得的。它的FFT也是N 点序列，即X(k)(k=0,1,2,…,N-1)，则第k点所对应实际频率值为f=k*f /N.
（3）作FFT分析时，幅值大小与FFT选择点数有关，但不影响分析结果。
3.2设计原理
3.2.1 运用自相关方法估计语音信号的声道参数原理
由均方预测误差最小的得到正则方程
∑==Φ=Φp k n
n k p i i k i a 1,...3,2,1],0,[],[ （3-1）
其中，
∑∞-∞=≤≤--=
Φm n
n n p k i k m s i m s k i ,1],[][],[ （3-2） 在最佳解时的误差为
∑=Φ-Φ=p k n k n n k a E 1]
,0[]0,0[ （3-3）
在自相关法中式2-1，式2-3变为
∑==-=-p k n
n k p i i r k i r a 1,...3,2,1],0[][ （3-4）
∑=-=p k n k n n k R a R E 1][]0[
（3-5）
由 式2-4 可列出方程组式2-6
???????
?????????=????????????????????????????------][]3[]2[]1[]0[]3[]2[]1[...]3[...]0[]1[]2[]2[...]1[]0[]1[]1[...]2[]1[]0[421p R R R R R p R p R p R p R R R R p R R R R p R R R R n n n n n n n n n n n n n n n n n n n n M M K M M M M ααα（3-6） 解方程组 式2-6 求出线性预测系数，通过误差式2-5可求出增益G
∑=-==p k n k n n k R R E G 12
][]0[α
（3-7）
3.2.2解决噪声污染的原理
噪声来源于实际的应用环境，因而其特性变化很大。噪声可以是加性的，也
可以是非加性的。对于非加性噪声，有些可以通过变换转变为加性噪声。例如，
乘积性噪声或卷积性噪声可以通道同态变换而成为加性噪声。
（1）通过直接变换频谱可以消除噪声的周期性成分，这可以用数字信号处理的方法来实现。变换滤波器如图2-1(a)所示。信号要经过离散傅里叶变换(DFT)，变换到频域，在领域进行处理，然后用反变换IDFT来重建语音信号。
图3-1通过频谱整形削去周期噪声
图3-1（b）所示为频谱整形器，它可以是简单的一系列选通门。它可将噪声成分变换到零值，则反变换后的信号周期性干扰将被滤除。
（2）同态滤波法的关键部分具有非线性处理性质，它应用于语音识别中，着眼于将语音信息(基音、频谱)中的乘性噪声或干扰分离，或者将已减少了噪声、干扰的信息重新合成得到降噪时域信号再进行识别。
同态滤波法的原理框图如图2-2所示，在傅里叶变换所得到的复倒谱中基音信息明显地出现，故可将其经适当时窗再作平滑、滤波及博里叶变换得到降噪处理后的共振峰等。它们可作为识别用的特征。也可以令其与得到的音调信号合成以还原成降噪的语音信号，然后进入识别系统作其他的特址提取相识别。
图3-2 非线性处理中的同态滤波法框图
3.2.3 自相关检测原理
对于离散的数字语音信号序列X(n),自相关函数定义如下：
-∞
R(k)=∑ x(m)x(m+k) （3-8）M=+∞
式中，k为信号的延迟点数。对于随机性信号序列或周期性信号序列，自相关函数定义为：
R(k)=lim 1／2N+1∑X(m)X(m+k) （3-9）自相关函数具有以下性质：如果序列X(n)具有周期Np，其自相关函数也是同周期的周期函数。即X(n)=X(n+Np)则R(k)=R(k+Np)清音信号没有周期性，它的自相关函数也没有周期，R(k)会随着k的增大迅速衰减。浊音信号具有准周期性，它的自相关函数R(k)具有与｛X(m)｝相同周期。自相关法基音检测正是利用R(k)的这一性质对语音信号进行基音检测的。
3.2.4 中心消波法检测的原理
中心消波处理是使用如图所示的中心消波函数进行处理：
x(n)-L 当x(n)>CL时
Y(n)=C(n)= 0 当|x(n)|≤CL时（3-10） x(n)+L 当x(n)<-CL时
其中削波电平CL一般取最大信号幅度的60%-70%图2-3所给出中心削波处理后的结果。中心削波后，再用自相关法检测基音频率，错判为倍频或分频的情况就可大大减少了。
图3-3
3.2.5 三电平中心消波法原理
为了减少自相关计算中的乘法运算，可以把上述中心削波以后的信号{y(n)}的自相关用两个信号的互相关代替，其中一个信号是{y(n)}另一个信号是对{y(n)}进行三电平量化产生的结果{y’(n)}即：
+1 当y(n)>0时
y’(n)=C’(n)=0 当y(n)=0时（3-11） -1 当y(n)<0时
显然y’(n)只有-1，0，+1三种肯能的取值，因而这里的互相关计算至需做加减法，而这个互相关序列的周期性与{y(n)}的自相关序列是近似相同的。
第四章模块实现过程
4.1 DFT和DTFT设计
4.1.1 DFT定点分析
定点分析：已知一个语音信号，数据采样频率为100Hz，试分别绘制N＝128点DFT的幅频图和N＝1024点DFT幅频图。
编程如下：
x=wavread('ah02.wav');
sound(x);
fs=100;N=128;
y=fft(x,N);
magy=abs(y);
f=(0:length(y)-1)'*fs/length(y);
subplot(221);plot(f,magy);
xlabel('频率（Hz）');ylabel('幅值');
title('N=128(a)');grid
subplot(222);plot(f(1:N/2),magy(1:N/2));
xlabel('频率（Hz）');ylabel('幅值');
title('N=128(b)');grid
fs=100;N=1024;
y=fft(x,N);
magy=abs(y);
f=(0:length(y)-1)'*fs/length(y);
subplot(223);plot(f,magy);
xlabel('频率（Hz）');ylabel('幅值');
title('N=1024(c)');grid
subplot(224);plot(f(1:N/2),magy(1:N/2));
xlabel('频率（Hz）');ylabel('幅值');
title('N=1024(d)');grid
运行结果如图：
图4-1DFT
上图(a)、(b)为N=128点幅频谱图，(c)、(d)为N=1024点幅频谱图。由于采样频率f =100Hz，故Nyquist频率为50Hz。(a)、(c)是0～100Hz频谱图，(b)、(d)是0～50Hz频谱图。由(a)或(c)可见，整个频谱图是以Nyquist频率为轴对称的。因此利用fft对信号作频谱分析，只要考察0～Nyquist频率（采样频率一半）范围的幅频特性。比较(a)和(c)或(b)和(d)可见，幅值大小与fft选用点数N 有关，但只要点数N足够不影响研究结果。从上图幅频谱可见，信号中包括15Hz 和40Hz的正弦分量。
4.1.2 DTFT设计
若信号长度T=25.6s，即抽样后x(n)点数为T/Ts=256，所得频率分辨率为Hz，以此观察数据长度N的变化对DTFT分辨率的影响：
编程如下：
[x,fs,bits]=wavread('ah02.wav');
N=256;
f=0:fs/N:fs/2-1/N;
X=fft(x);
X=abs(X);
subplot(211)
plot(f(45:60),X(45:60));grid
xlabel('Hz'),ylabel('|H(ejw)|')
%数据长度N扩大4倍后观察信号频谱
N=N*4;
f=0:fs/N:fs/2-1/N;
X=fft(x);
X=abs(X);
subplot(212)
plot(f(45*4:4*60),X(4*45:4*60));grid xlabel('Hz'),ylabel('|H(ejw)|')
结果如图：
图4-2 DTFT
4.2 滤波器设计
4.2.1原理
设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统，并使系统函数H(z)具有指定的频率特性。
数字滤波器从实现的网络结构或者从单位冲激响应分类，可以分成无限长单位冲激响应（IIR）数字滤波器和有限长单位冲激响应（FIR）数字滤波器。
数字滤波器频率响应的三个参数：
（1）幅度平方响应：
（2）相位响应
（3）群时延响应
4.2.2简介IIR数字滤波器和FIR数字滤波器
（1）IIR数字滤波器：
IIR数字滤波器的系统函数为的有理分数，即IIR数字滤波器的逼近问题就是求解滤波器的系数和，使得在规定的物理意义上逼近所要求的特性的问题。如果是在s平面上逼近，就得到模拟滤波器，如果
是在z平面上逼近，则得到数字滤波器。
（2）FIR数字滤波器：
设FIR的单位脉冲响应h(n)为实数，长度为N，则其z变换和频率响应分别为按频域采样定理FIR数字滤波器的传输函数H(z)和单位脉冲响应h(n)可由它的N个频域采样值H(k)唯一确定。
MATLAB中提供了几个函数，分别用于实现IIR滤波器和FIR滤波器。（1）卷积函数conv
卷积函数conv的调用格式为 c=conv(a,b)
该格式可以计算两向量a和b的卷积，可以直接用于对有限长信号采用FIR滤波器的滤波。
（2）函数filter
函数filter的调用格式为y=filter(b,a,x)该格式采用数字滤波器对数据进行滤波，既可以用于IIR滤波器，也可以用于FIR滤波器。其中向量b和a分别表示系统函数的分子、分母多项式的系数，若a＝1，此时表示FIR滤波器，否则就是IIR 滤波器。该函数是利用给出的向量b和a，对x中的数据进行滤波，结果放入向量y。
（3）函数fftfilt函数fftfilt的调用格式为y=fftfilt(b,x)
该格式是利用基于FFT的重叠相加法对数据进行滤波，这种频域滤波技术只对FIR滤波器有效。该函数是通过向量b描述的滤波器对x数据进行滤波。
关于用butter函数求系统函数分子与分母系数的几种形式。
[b,a]=butter(N,wc,'high'):设计N阶高通滤波器，wc为它的3dB边缘频率，以单位，故[b,a]=butter(N,wc):当wc为具有两个元素的矢量wc=[w1,w2]时，它设计2N阶带通滤波器，3dB通带为，w的单位为。[b,a]=butter(N,wc,'stop'):若wc=[w1,w2]，则它设计2N阶带阻滤波器，3dB通带为，w的单位为。如果在这个函数输入变元的最后，加一个变元“s”，表示设计的是模拟滤波器。这里不作讨论。
为了设计任意的选项巴特沃斯滤波器，必须知道阶数N和3dB边缘频率矢量wc。这可以直接利用信号处理工具箱中的buttord函数来计算。如果已知滤波器指标，，和，则调用格式为[N,wc]=buttord(wp,ws,Rp,As)
对于不同类型的滤波器，参数wp和ws有一些限制：对于低通滤波器，wp 4.2.3设计内容 （1）滤波器示例： 在这里为了说明如何用MATLAB来实现滤波，特举出一个简单的函数信号滤波实例（对信号x(n)=sin( n/4)+5cos( n/2)进行滤波，信号长度为500点），从中了解滤波的实现过程。 程序如下： Wn=0.2*pi; N=5; [b,a]=butter(N,Wn/pi); n=0:499; x=sin(pi*n/4)+5*cos(pi*n/2); X=fft(x,4096); subplot(221);plot(x);title('滤波前信号的波形'); subplot(222);plot(X);title('滤波前信号的频谱'); y=filter(b,a,x); Y=fft(y,4096); subplot(223);plot(y);title('滤波后信号的波形'); subplot(224);plot(Y);title('滤波后信号的频谱'); 结果如图： 图4-3 函数信号滤波器 （2）N阶高通滤波器的设计（在这里，以5阶为例，其中wc为其3dB边缘频率，以为单位） 程序设计如下： x=wavread('ah02.wav'); sound(x); N=5;wc=0.3; [b,a]=butter(N,wc,'high'); X=fft(x); subplot(322);plot(X);title('滤波前信号的频谱'); y=filter(b,a,x); Y=fft(y); subplot(323);plot(y);title('IIR滤波后信号的波形'); subplot(324);plot(Y);title('IIR滤波后信号的频谱'); z=fftfilt(b,x); Z=fft(z); subplot(325);plot(z);title('FIR滤波后信号的波形'); subplot(326);plot(Z);title('FIR滤波后信号的频谱'); 等到结果如图： 图4-4 N阶高通滤波器 （3）N阶低通滤波器的设计（在这里，同样以5阶为例，其中wc为其3dB边缘频率，以为单位）， 程序设计如下： x=wavread('ah02.wav'); sound(x); N=5;wc=0.3; [b,a]=butter(N,wc); X=fft(x); subplot(321);plot(x);title('滤波前信号的波形'); y=filter(b,a,x); Y=fft(y); subplot(323);plot(y);title('IIR滤波后信号的波形'); subplot(324);plot(Y);title('IIR滤波后信号的频谱'); z=fftfilt(b,x); Z=fft(z); subplot(325);plot(z);title('FIR滤波后信号的波形'); subplot(326);plot(Z);title('FIR滤波后信号的频谱'); 得到结果如图： 图4-5 N阶低通滤波器 （4）2N阶带通滤波器的设计（在这里，以10阶为例，其中wc为其3dB边缘频率，以为单位，wc=[w1,w2],w1 wc w2） 程序设计如下： x=wavread('ah02.wav'); sound(x); N=5;wc=[0.3,0.6]; [b,a]=butter(N,wc); X=fft(x); subplot(321);plot(x);title('滤波前信号的波形'); y=filter(b,a,x); Y=fft(y); subplot(323);plot(y);title('IIR滤波后信号的波形'); subplot(324);plot(Y);title('IIR滤波后信号的频谱'); z=fftfilt(b,x); Z=fft(z); subplot(325);plot(z);title('FIR滤波后信号的波形'); subplot(326);plot(Z);title('FIR滤波后信号的频谱'); 得到结果如图： 图4-6 2N阶带通滤波器 （5）2N阶带阻滤波器的设计（在这里，以10阶为例，其中wc为其3dB边缘频率，以为单位，wc=[w1,w2],w1 wc w2） 程序设计如下： x=wavread('ah02.wav'); sound(x); N=5;wc=[0.2,0.7]; [b,a]=butter(N,wc,'stop'); X=fft(x); subplot(321);plot(x);title('滤波前信号的波形'); subplot(322);plot(X);title('滤波前信号的频谱'); y=filter(b,a,x); Y=fft(y); subplot(323);plot(y);title('IIR滤波后信号的波形'); subplot(324);plot(Y);title('IIR滤波后信号的频谱'); z=fftfilt(b,x); Z=fft(z); subplot(325);plot(z);title('FIR滤波后信号的波形'); subplot(326);plot(Z);title('FIR滤波后信号的频谱'); 得到结果如图： 图4-7 2N阶带阻滤波器 4.3 分离观察 通过函数的调用，也可以将信号的频谱进行“分离观察”，如显出信号的幅值或相位。下面，通过改变系统函数的分子与分母系数比，来观察信号滤波前后的幅值与相位。并且使结果更加明显，使人耳得以很容易的辨听。 x=wavread('ah02.wav'); sound(x); b=100;a=5; y=filter(b,a,x); X=fft(x,4096); subplot(221);plot(x);title('滤波前信号的波形');