数字信号处理教案（五篇）

第一篇：数字信号处理教案

“数字信号处理”教案

Digital Signal Processing —

Teaching Project

第一讲：信号的采集、基本DSP系统

Lecture 1 Conceptual introduction of DSP

 了解技术背景、各种信号的特征、A/D转换、采样与量化、Nyquist 定理

一、连续信号的采样与量化

信号的分类与特点、模拟信号到离散信号的转换、Nyquist采样定理以及量化。

二、采样前后频谱的变化

模拟信号以及相应离散信号频谱之间的关系。

三、从采样信号恢复连续信号

如何从采样后的离散信号恢复模拟信号。

Questions:(1)What is the advantage of DSP ?(2)Why generally put a LPF and a amplifier before the A/D conversion ?

第二讲：离散信号的描述与基本运算、线性卷积

Lecture 2 Discrete signal: its description and computations

 掌握离散信号的描述方法、典型信号的特征、信号之间的基本运算以及线性卷积 信号与系统分类

一、信号的分类

模拟信号、离散信号、数字信号

二、系统分类

模拟系统、离散系统、数字系统 连续时间信号的采样与量化 1 离散时间信号—序列

一、典型的序列

离散信号的时域描述；冲击信号、单位阶跃信号、指数信号、正弦信号等的描述。

二、序列的运算

信号序列之间的基本运算，能量的计算以及分解等。线性卷积

序列的线性卷积运算、具体步骤。

Questions:(1)What is absolute time for a time index n of x(n)?(2)In practical application, is determined signals such as sine need to be processed ? If not, what type of signal is we mostly faced ?

第三讲：系统的分类与描述

Lecture 3 Linear shift-invariant system and its description

 掌握LSI、因果、稳定、FIR、IIR系统的特征；LSI的I/O描述；线性常系数差分方程；系统结构描述 离散时间系统一、离散时间系统的类型

线性系统、移不变系统、因果系统、稳定系统、IIR与FIR系统。

二、离散时间系统的描述

LSI系统的I/O关系（线性卷积形式）、差分方程描述。

Questions:(1)Which system description is mostly used in practical application, why ?(2)Can a IIR system be replaced by a FIR system ?

第四讲：Z变换与系统函数

Lecture 4 Z transform

 掌握Z变换；系统函数以及零极点分析；系统函数与差分方程之间的转换 Z变换

一、Z变换的定义及其收敛域

双边Z变换、收敛域的概念、典型信号的Z变换；不同分布序列的收敛域特征。

二、逆变换

基本逆Z变换的定义、留数法以及幂级数法计算。

三、Z变换的性质

导数与极值等特性。离散时间系统的Z变换分析法

一、系统函数

系统函数定义；不同系统的系统函数特点；极点与零点的特性、与差分方程的关系等。

二、离散时间系统的信号流图描述

系统的结构框图、流图描述方法。

Questions:(1)why we need study Z transform, how important converge region is ?(2)why the condition for a causal stable LSI is that its converge region includes the unit circle ?

第五讲：离散信号的傅立叶变换

Lecture 5 Discrete time Fourier transform

 掌握离散信号的傅立叶变换DTFT；频谱、幅度谱与相位谱；离散信号DTFT的特征 离散信号的傅立叶变换

一、离散信号傅立叶变换的定义

离散信号DTFT与IDTFT的定义，典型信号的DTFT计算。

二、离散信号的傅立叶变换与Z变换的关系

单位圆上的Z变换。离散信号傅立叶变换的特点

Questions:(1)What a point on magnitude spectrum states ?(2)What is relation between frequency components of a signal and the points of its spectrum curve ?

第六讲：系统频率响应与频谱关系

Lecture 6 System frequency response and spectrum relations

 掌握LSI系统频率响应概念；零极点对频谱的影响；模拟信号频谱与对应离散信号频谱的关系。线性移不变系统的频率响应系统函数零极点与频率响应的关系离散信号频谱与模拟信号频谱之间的关系

一、离散时间傅立叶变换的导出

Questions: 从模拟信号以及频谱推导到离散信号的频谱。模拟信号频谱与对应离散信号频谱之间的关系。

二、DTFT与FT的关系 系统函数与频率响应的关系，零点和极点对系统频率响应的影响。由线性移不变系统对复指数信号的作用推导出系统的频率响应。对称、周期、卷积等特性，帕斯维尔（Parseval）定理。(1)What a point on magnitude frequency response states ?(2)What is response of a system to the points of spectrum of input signal ?

第七讲：频谱分析与应用

Lecture 7 Spectrum analysis and application

 掌握频谱的基本信息特征；频谱分析的典型应用；短时谱分析的概念 频谱分析与应用

一、频谱的基本特征

通过复正弦信号的频谱说明DTFT的意义以及频谱分析的意义。

二、信号调制与语音合成

通信中AM调制与语音合成中频谱的应用。

二、短时频谱分析

Questions:(1)propose some examples of spectrum analysis in application(2)what is the influence of short time processing for spectrums ?

第八讲：周期信号的傅立叶级数表示

Lecture 8 Fourier series of periodical discrete signal

 了解周期信号的DFS描述； DFS的频谱特征； 周期卷积 周期信号的离散傅立叶级数表示

一、离散傅立叶级数

周期信号的DFS定义及频谱分析。

二、周期卷积

从一个周期求和的线性卷积导出周期卷积。

第九讲：离散傅立叶变换 阐述实际应用中的频谱分析方法。Lecture 9 Discrete Fourier transform

 掌握DFT；DFT的基本前提与特征；频率取样定理；DFT与DFS和DTFT的关系 离散傅立叶变换离散傅立叶变换特性

一、有限长特性与频域采样定理

描述DFT的时频有限长特性；DFT作为DTFT采样的频域采样定理。

二、循环卷积特性

Questions:(1)Why we need DFT ?(2)What is the difference between DFT and spectrum sampling ?

第十讲：短时离散傅立叶变换

Lecture 10 Short-time DFT

 掌握循环卷积；STDFT的概念和实用意义；时间分辨率与频率分辨率 短时离散傅立叶变换分析

一、短时离散傅立叶变换的定义

非有限长信号的STDFT定义；STDFT与原始频谱之间的关系。

二、频率分辨率与时间分辨率

Questions:(1)why it is said, for non-stationary signal, short-time DFT is a unique selection ?(2)Is zero-padding enough for high frequency resolution ? 短时频谱的时间分辨率与频率分辨率，及其短时窗长的影响。有限长信号的循环卷积。DFT与IDFT的定义；DFT与短时谱；从DFT的信号完备恢复。

第十一讲：快速傅立叶变换与应用

Lecture 11 Fast Fourier transform ant application

 掌握基2运算的FFT算法；了解FFT在信号处理中的应用 快速傅立叶变换

一、基于时选的快速傅立叶变换

时域实行奇偶分解的FFT算法。

二、基于频选的快速傅立叶变换快速傅立叶变换的应用

一、信号去噪与语音识别

谱相减方法的去噪处理；应用频谱特征的语音识别应用。

二、利用FFT计算线性卷积

线性卷积与循环卷积的关系；通过循环卷积与DFT的对应关系得到FFT计算线性卷积的方法。

Questions:(1)Is there any difference between DFT and FFT ?(2)Can you propose a new fast algorithm of DFT ?

第十二讲：数字滤波器类型与技术指标

Lecture 12 Digital filter type and technical parameters

 了解IIR、FIR数字滤波器的结构特点；滤波器的设计技术指标；IIR数字滤波器的一般设计方法 数字滤波器的技术指标

频域实行奇偶分解的FFT算法。IFFT快速算法与FFT的关系。

三、傅立叶反变换的快速计算 通带、阻带、截止频率（3dB下降）、通带与阻带边界频率、阻带衰减。无限脉冲响应数字滤波器的结构模拟滤波器到数字滤波器的转换

一、脉冲响应不变法

从时域脉冲响应保持不变原理分析导出模拟滤波器到数字滤波器的转换。

二、双线性变换法

Questions:(1)how many technical parameters must be set for design of filter ?(2)what is advantages of bilinear transform ?

第十三讲：IIR数字滤波器的设计

Lecture 13 Design of IIR filter

 掌握Butterworth、Chebyshev和椭圆滤波器的设计方法；脉冲响应设计法与双线性设计法； LPF与HPF、BPF、BSF的转换 IIR滤波器的特性

一、巴特沃兹滤波器

Butterworth滤波器的特点；相应滤波器的设计方法。

二、切比雪夫滤波器 Chebyshev滤波器的特点；相应滤波器的设计方法。

三、椭圆滤波器

椭圆滤波器的特点以及设计方法。IIR滤波器设计的频率变换方法 从克服模拟滤波器到数字滤波器转换过程中频率畸变的问题，导出双线性频率变换方法。直接Ⅰ与Ⅱ型结构；级联与并联结构；全通滤波器。

一、模拟低通滤波器到其它滤波器的变换

模拟低通滤波器转换到高通、带通、带阻滤波器的方法。

二、数字低通滤波器到其它滤波器的变换

Questions:(1)do you think Butterworth is much easier than others ?(2)what is a general steps for design of IIR filters ?

第十四讲：IIR滤波器的应用与系数量化效应

Lecture 14 Application and coefficient effects of IIR filter

 了解IIR滤波器设计中的系数量化效应和实际应用 IIR滤波器实现与系数量化效应

一、IIR滤波器的实现

IIR滤波器的硬件与软件实现方法。

二、系数量化效应IIR滤波器应用

一、小循环阻抗容积信号处理

说明滤波器的具体应用与效果。

二、DTMF双音频信号的合成Questions:(1)Is it OK for use of IIR filter in image processing ?(2)Propose other IIR filter applications.介绍用一个IIR滤波器如何完成输出一个单频率信号。滤波器系数量化效应对性能的影响分析。数字低通滤波器转换到数字高通、带通、带阻滤波器的方法。第十五讲： 线性相位FIR滤波器及窗函数设计原理

Lecture 15 Linear phase FIR filter and principle of window method

 掌握FIR滤波器的特点；线性相位概念、意义及其实现条件；FIR滤波器窗函数设计法原理。FIR数字滤波器的特点

一、基本特点

脉冲响应、差分方程、系统函数以及系统结构等方面的特点。

二、线性相位特点

线性相位概念、系统设计中的意义，举例说明。

三、线性相位FIR滤波器的实现条件

如何实现线性相位，不同奇偶点数的区别。窗函数设计法原理

一、窗函数设计法原理

从时域逼近角度分析导出窗函数设计法，说明失真的情况。

二、理想低通滤波器

Questions:(1)What is the importance of linear phase for a filter ?(2)Can IIR be realized as a linear phase filter, why ?

第十六讲：窗函数设计分析与实例

Lecture 16 Design analysis and examples of window method

 掌握Hamming窗等5种基本窗函数的具体设计方法；特别是Kaiser窗设计实例 窗函数设计法分析

一、各种窗函数设计法 描述一个理想LPF的特点，特别是幅度特性。矩形窗、汉宁窗、哈明窗、布莱克曼窗、凯泽窗设计FIR的方法、特点。

二、窗函数设计法的进一步分析与总结

对窗长、窗的类型在设计中的影响做总结分析。利用凯泽窗设计FIR滤波器

一、低通滤波器设计

凯泽窗设计LPF的具体举例分析。

二、高通通滤波器设计

凯泽窗设计HPF的具体举例分析。

三、带通滤波器设计

凯泽窗设计BPF的具体举例分析。

四、带阻滤波器设计

凯泽窗设计BSF的具体举例分析。

Questions:(1)are you confident for design of FIR filter now ? why ?(2)If you are assigned to design a untypical filter, how can you do ?

第十七讲：频率取样设计与等波纹优化设计

Lecture 17 Frequency design and equal-ripple method of FIR filter design

 掌握频率取样设计方法；等波纹优化设计方法 频率取样设计法

一、频率取样设计法原理

从频率抽样形成DFT频谱，并进一步得到有限长脉冲响应的思路介绍，说明其实际失真。

二、设计实例分析等波纹逼近优化设计方法

举例说明频率取样设计法的具体过程、从不成功设计到成功设计的转变思路与方法。

一、最小均方误差优化设计 LMS准则下的优化设计介绍。

二、等波纹逼近优化设计法

Questions:(1)which is more excellent as a method ?(2)why some points must be set in transition band ?

第十八讲：系数量化效应与FIR滤波器应用

Lecture 18 Application and coefficient effects of FIR filter

 了解 FIR滤波器的系数量化效应以及实际应用 系数量化效应与溢出控制

一、系数量化效应

有限字长条件下滤波器系数的量化对频谱的影响，引起失真的情况。

二、溢出控制

怎样处理滤波器输出数据对D/A转换器或其他接收器的输入溢出问题。FIR滤波器应用

一、信号去噪

举例说明运用FIR实现限带噪声消除的实际应用。

二、信号的高频提升

Questions:(1)If to implement a FIR in a MCU, what should you consider ?(2)Propose some other application examples.最小误差意义下的频域的等波纹逼近，具体设计方法，MATLAB仿真设计举例。

一个简单的一阶高频FIR滤波器如何提升信号的高频部分。

第二篇：数字信号处理——第八讲(教案)

第八讲（3.6节 离散时间LTI系统的Z域分析 3.7节 梳状滤波器、全通滤波器和最小相位系统）

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

1、回顾第七讲内容：

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

2、本节内容概要：

<1> 系统函数的极点分布与系统因果性和稳定性的关系

<2> 系统函数的零极点分布对系统频率响应特性的影响

<3> 梳状滤波器 <4> 全通滤波器 <5> 最小相位系统 <6> 本章小结 <7> 本章作业

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <1> 系统函数的极点分布与系统因果性和稳定性的关系

在离散时间LTI系统的时域分析中（参见2.4节）已经证明：因果系统的充分必要条件是其单位脉冲响应hn满足

hn0, n0

(3.6.10)而稳定系统的充分必要条件是其单位脉冲响应hn绝对可和，即

nhn

(3.6.11)当系统的单位脉冲响应hn0（n0）时，则其系统函数Hz的收敛域一定包含点，即收敛域是某个圆的外部，表示为Rxz。因为收敛域中不含极点，所以因果系统其系统函数Hz的极点分布在某个圆内，收敛域是这个圆的外部。

《数字信号处理》教案——臧博 当系统的单位脉冲响应hn满足绝对可和的条件时，根据(3.2.3)式其DTFT存在，即

HejDTFThn

存在。由序列的Z变换与DTFT的关系，即

HejHz当序列hn的Hejzej

存在时，其Hz的收敛域一定包含单位圆。由此得出结论：稳定系统的系统函数Hz的收敛域包含单位圆，而收敛域中没有极点。

如果系统是因果稳定系统，则要求同时满足因果性和稳定性的条件，即系统的系统函数Hz的收敛域是包含单位圆的外部，表示为

rz，0r1

所以，因果稳定系统的系统函数Hz的所有极点一定分布在单位圆内。

例3.6.1 已知离散时间LTI系统的系统函数为

1a2，a1，a为实常数 Hz211aazaz试分析该系统的因果性和稳定性。

解：系统函数Hz的极点为p1a和p2a1,如图3.4.4所示（参见例题3.4.11）。针对三种可能的收敛域，分别讨论系统的因果性和稳定性。

(1)当收敛域为a1z时，Hz的极点在半径为a的圆的内部，对应的系统是

1因果系统，但由于收敛域不包含单位圆，因此系统是不稳定系统。系统的单位脉冲响应为

hnaannun

它是一个因果序列，但不收敛。

(2)当收敛域为0za时，Hz的极点在半径为a的圆的外部，且收敛域不包含单 位圆，对应的系统是非因果且不稳定的系统。系统的单位脉冲响应为

hnananun1

这是一个非因果且不收敛的序列。

(3)当收敛域aza1时，极点p1a在半径为a的圆内，而极点p2a1在半径为a1的圆外，即该圆的外部不是收敛域，对应的系统是一个非因果系统，但由于收敛域包含单位圆，因此它是稳定系统。系统的单位脉冲响应为

hnna

《数字信号处理》教案——臧博 这是一个收敛的双边序列。

最后说明，因果稳定的系统是物理可实现的稳定系统。如果系统是非因果的但是稳定的系统，如例3.6.1中的（3），严格讲它是物理不可实现的系统。但是在数字信号处理中，利用数字系统的存储功能，是可以近似实现的。基本方法是将非因果的单位脉冲响应hn，从nN到nN截取一段并存储，表示为h2Nn；把h2Nn作为具体实现系统的单位脉冲响应参与运算；运算的结果即系统的输出在时间上有所延迟。N愈大，h2Nn所表示的系统愈接近原hn表示的系统。这种非因果但稳定系统的数字实现，是数字信号处理技术优于模拟信号处理技术的又一体现。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <2> 系统函数的零极点分布对系统频率响应特性的影响

有理系统函数Hz经因式分解后，一般可表示为

HzA1zz1iM1pz1kk1i1N

(3.6.12)式中，Ab0a0,当a01时，Ab0；若Hz的零点zicii1,2,M，极点pkdkk1,2,N，则

HzA1cz1iM1dz1kk1i1N

(3.6.13)该式说明，除了反映系统幅度增益的比例常数A以外，整个系统函数可以由它的全部零点、极点唯一确定。

用零点和极点来表示系统函数的优点之一是它引导出一种获得系统频率响应的实用的几何方法。

设系统是稳定的，令zeMjj代入(3.6.13)式，得系统的频率响应。

HejA'eciedjkk1i1N

(3.6.14)式中，AAe'jNM，它只包含常系数的幅度增益和线性的相移，不影响Hejj的频率

j特性。在z平面上，eci可以用由零点ci指向单位圆上e点B的向量ciB来表示；同

《数字信号处理》教案——臧博

样，edk可以用由极点dk指向单位圆上点B的向量dkB来表示，如图3.6.1所示，即 jjImzdkdkBkBRezciBici1dk*

图3.6.1 系统频率响应的几何表示

ciBejci dkBejdk

ciB和dkB分别称为零点矢量和极点矢量。将它们用极坐标表示为

jciBciBei jdkBdkBek

将它们代入(3.6.14)式，得到

HejA'jHeedkBi1Nk1ciBMj

(3.6.15)式中

HejA'ciBMk1Ni1N

(3.6.16)dkBik

(3.6.17)i1k1M这样，系统的频率响应由(3.6.16)式和(3.6.17)式确定。其中，系统的幅频响应按(3.6.16)式由零点矢量长度之积与极点矢量长度之积的比值决定；而相频响应则按(3.6.17)式由零点矢量幅角之和与极点矢量幅角之和的差值决定。当频率ω从零变化到2时，这些向量的终

《数字信号处理》教案——臧博 点B沿单位圆逆时针旋转一周。从而可分别得到系统的幅频响应和相频响应。例如，图3.6.2(a)所示的具有一个零点和一个极点的系统，其幅频响应He和相频响应如

j图3.6.2(b)所示。

按照(3.6.16)式，知道系统函数的零极点分布后，可以分析零极点位置对系统频率特性的影响。当变化使B点转到极点附近时，该极点矢量长度短，因而幅频响应将出现峰值，且极点愈靠近单位圆，极点矢量长度愈短，峰值愈高愈尖锐。如果极点在单位圆上，该极点对应的幅频响应无穷大，系统是不稳定的，这与稳定系统的极点要求分布在单位圆的条件是一致的。对于零点，结果相反，当变化使B点转到零点附近时，该零点矢量长度最短，幅频响应将出现谷值，零点愈靠近单位圆，谷值愈接近零。当零点处在单位圆上时，谷值为零。

总结以上结论：系统函数的极点位置主要影响系统幅频响应的峰值位置及尖锐程度，零点位置主要影响系统幅频响应的谷值及谷深。

例3.6.2 设一阶系统的差分方程为

ynayn1xn, 0a1

试分析该系统的频率响应特性。

解：由系统的差分方程得到系统函数为

Hz1z，za

1az1za系统的单位脉冲响应为

hnanun

系统函数的零点z0，极点pa，如图3.6.2(a)所示。幅频响应特性和相频响应特性如图3.6.2(b)所示。

2|H(ej)|jImz1.510.500.20.40.60.81/1.21.41.61.820.20a1Rez()/0.10-0.1-0.200.20.40.60.81/1.21.41.61.82

(a)零极点分布

(b)频率响应特性曲线

图3.6.2 例3.6.2 系统函数的零极点分布和频率响应特性

例3.6.3 例3.6.2的单位脉冲响应hn是无限长因果序列，对应的系统是IIR数字滤波

《数字信号处理》教案——臧博 器。如果截取hn的一段，得到一个有限长单位脉冲响应

na,0nN1hn0,其他n0a

1对应的是FIR数字滤波器，试分析其频率响应特性。

解：系统的系统函数为

Hzazn0N1nn1aNzN，za 11az将Hz该写为

zNaN, za HzN1zza其N个零点为

ziaej2Ni1,N , i1,2它们等间隔分布在半径为a0a1的圆上。在za处有一个极点，其余N1个极点是处在原点的N1阶重极点。这样，在za处的极点被一个相同位置的零点所抵消。设N8，则其零、极点分布图和对应的频率响应特性分别如图3.6.3(a)、(b)所示。可以看出，系统幅频响应的峰值出现在0处，因为那里的零点被极点抵消了；另外，在每一零点附近的幅频响应均有陷落。随着截取长度N的增加，幅频响应特性的曲线将逐渐光滑而接近hnanun的结果。

64jImz|H(ej)|2000.20.40.60.81/1.21.41.61.820.4N1阶极点()/0aRez10.20-0.2-0.400.20.40.60.811.21.41.61.82/

(a)零极点分布(N=8)

(b)频率响应特性(N=8)

图3.6.3 例3.6.3 系统函数的零极点分布和频率响应特性

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <3> 梳状滤波器

《数字信号处理》教案——臧博 梳状滤波器(combfilter)、全通滤波器(whole pass filter)和最小相位系统(minimum phase system)是各具特色频率响应特性的离散时间LTI系统。本节分别讨论它们的基本特性。

滤波器通常只有一个通带或阻带，但在实际应用中，往往要求滤波器具有多个通（阻）带。梳状滤波器就是一个多通（阻）带的数字滤波器。梳状滤波器的频率响应是的周期函数，其周期为2N（N是一个正整数）。如果Hz是一个单通（阻）带滤波器，可用zN取代Hz中的z，从而获得梳状滤波器的系统函数HczHzN。称Hz为构成梳状滤波器的原型滤波器(prototype filter)。原型滤波器可以是FIR数字滤波器，也可以是IIR数字滤波器。

由系统函数为

Hz1z1，z0 的一阶高通FIR原型数字滤波器产生的梳状滤波器的系统函数为

Hcz1zNzN1N

(3.7.1)

zHcz的零点有N个，由分子多项式的决定，即

zN10

ziej2i1N,N，i1,2这N个零点等间隔分布在单位圆上。Hcz在z0处有N阶极点。

令zej，则梳状滤波器的频率响应Hce为

jHcej1zjN1cosnjsinn

(3.7.2)其幅频响应Hce为 jjHce21cosn1

2(3.7.3)相频响应c为

cargtansinn

(3.7.4)

1cosn设N8，Hcz的零极点分布如图3.7.1(a)所示。当从零变化到2时，每遇到一个零点，幅频响应为零，在两个零点中间幅频响应最大，形成峰值。幅频响应谷值点的频率为i2Ni1i1,2,N。N8时的幅频响应特性如图3.7.1(b)所示。通常将具有如图3.7.1(b)所示幅频响应特性的滤波器称为梳状滤波器。

《数字信号处理》教案——臧博 如果滤波器的系统函数Hcz是由一阶低通FIR原型数字滤波器Hz1z1产生的，即

Hcz1zN

(3.7.5)则

Hcej2ejN2cosN

(3.7.6)它也是梳状滤波器，但频率响应特性与(3.7.2)式所示的频率响应特性有所差别，请读者自己分析。

21.5|H(ej)|10.5000.20.40.60.81/1.21.41.61.82jImz0.5Rez()/10-0.500.20.40.60.81/1.21.41.61.82

(a)零极点分布(N=8)

(b)频率响应特性曲线(N=8)

图3.7.1 梳状滤波器的零极点分布和频率响应特性

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <4> 全通滤波器

全通滤波器是指系统的幅频响应恒为常数的数字滤波器。全通滤波器在滤波器结构、多采样速率信号处理、滤波器组设计等领域有着广泛的应用。

(1)一阶全通滤波器

一阶的因果稳定的全通滤波器的系统函数定义为

z1cHwp1z，c(3.7.7)1cz1系统的频率响应为

Hwp1ejej1cej

(3.7.8)j1ce由于(3.7.8)式中的分子分母互为共轭，故

Hwp1ej1

(3.7.9)所以将该系统称为全通滤波器。由(3.7.7)式可得全通滤波器Hwp1z的极点为

《数字信号处理》教案——臧博 p1crej

（3.7.10）系统的零点为

z111je

(3.7.11)cr称由(3.7.10)式和(3.7.11)式给出的一阶全通滤波器零点和极点的位置关于单位圆镜像对称，如图3.7.2(a)所示。

将一阶全通滤波器的系数c用幅度r和幅角的极坐标表示，则由(3.7.8)式可得一阶全通滤波器的相频响应为

rsinwp12argtan

(3.7.12)

1rcos将wp1对求导可得

dwp1d1r21rcosrsin220

(3.7.13)这说明一阶全通滤波器的相频响应是单调递减的。

由(3.7.12)式可知，当0，从0变到2时，一阶全通滤波器的相位将从0递减到2；当0时，记

12arctanrsin

1rcos一阶全通滤波的相位将从1递减到12。这说明一阶全通滤波器当从从0变到2，0变到2时，相位将减小2。

一阶全通滤波器的频率响应特性如图3.7.2(b)所示。

《数字信号处理》教案——臧博

21.5|H(ej)|10.5000.10.20.30.40.50.60.70.80.91/jImz1c1c0.5()/Rez001-0.5-100.20.40.60.811.21.41.61.82/

(a)零极点分布

(b)频率响应特性曲线

图3.7.2 一阶全通滤波器的零极点分布和频率响应特性

(2)N阶实系数全通滤波器

N阶实系数全通滤波器的系统函数为

HwpNzazii0Nii0NNiazNzNa1zN1a2zN2aN1a1z1a2z2aNzN

式中 zNANz1ANz，a01,a1,a2,aN为实数

（3.7..14）

ANz1a1z1a2z2aNzN

(3.7.15)为使系统稳定，实系数多项式ANz的根，即系统函数的极点pkk1,2,N必须全都在单位圆内。

由于系统函数中的系数a1,a2,aN是实系数，所以

ANz1zejANejAN(ej)

(3.7.16)j式中AN(ej)是ANe的共轭，它们两者的模是相等的。因此有

《数字信号处理》教案——臧博 HwpNejAN(ej)1

(3.7.17)jANe这就证明了(3.7.14)式系统函数所代表的系统具有全通滤波特性。

全通滤波器具有特殊的零点和极点分布规律。设zi是HwpNz的零点，按照(3.7.14)式，zi1必然是HwpNz的极点，记为pizi1，则pizi1，即全通滤波器的零点和极点互为共轭倒数关系。如果再考虑到ANz和ANz的系数均为实数，其零点、极点或者为实

1数，或者呈共轭复数对。其中，复数零点和复数极点必然以两对一组出现。例如，zi为HwpNz的复数零点，则必有复数零点zi，而复数极点为pizi1和pizi1。而实

数零极点则以两个一组出现，且零点与极点互为倒数关系。以一个实数极点和一对共轭复数极点为例，零极点位置示意图如图3.7.3所示。

jImzzipiRez0p1ppizi图3.7.3 全通滤波器的零极点分布

观察图3.7.3，如果将零点zi和极点pi组成一对，将零点zi和pi组成一对，那么全通滤波器的零点与极点便以共轭倒数关系出现。因此，实数全系通滤波器的系统函数也可以写成如下形式：

z1pi

(3.7.18)HwpNz11pzi1iN其中，极点pi和零点1pi构成一个一阶全通滤波器。所以N阶实系数全通滤波器可分解为N个一阶全通滤波器的级联。

由于N阶全通滤波器的相频响应是N个一阶全通滤波器相频响应的和，所以N阶全通

《数字信号处理》教案——臧博

12 滤波器的相频响应也是单调递减的。由于HwpNe点的相位为

1，N阶实系数全通滤波器在0j000

当从0变到2时，N阶全通滤波器的相频响应将从0递减到2N。

全通滤波器是一种纯相位滤波器，经常用于相位均衡。对幅频响应特性满足要求而相频响应特性有缺陷的滤波器，可以级联全通滤波器进行相位校正。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <5> 最小相位系统

系统零点、极点都在单位圆内的因果系统，称为最小相位系统，系统函数记为Hminz。任一实系数因果稳定系统的Hz都可表示为一个最小相位系统和一个全通系统的级联，即

HzHminzHwpNz

(3.7.19)为了证明上述结论，设系统Hz只有一个零点z1a在单位圆外，a1，而其余

零点都在单位圆内，那么Hz就能表示成

HzH1zz1a

(3.7.20)按定义H1z是一个最小相位系统。Hz也可等效地表示为

11az1a1zHzH1zzaH1z1az

(3.7.21)111az1az1因为a1，所以H1z1az也是最小相位系统，记为Hminz。由于

1z1a1az1是一阶全通系统（参见（3.7.7）式），故得

HzHminzHwp1z

(3.7.22)对单位圆外的每一个零点使用(3.7.21)式，即得(3.7.19)式。由于HwpNz是全通系统，因此(3.7.19)式中的Hz和Hminz具有相同的幅频响应。

由于全通系统的相位是非正的，在幅频响应与Hz相同的所有因果稳定系统中，最小相位系统的相位值最大，称之为最小相位滞后。因此，最小相位系统更为确切的术语应为最小相位滞后系统。由于最小相位系统是历史上早已惯用的名称，所以至今仍延用这一术语。

例3.7.1 已知实系数因果稳定系统的系统函数Hz为

《数字信号处理》教案——臧博 bz1，a1,b1 Hz1az1试判决该系统是否为最小相位系统。

解：系统函数Hz的零点z11b。由于b1，所以z11，零点在单位圆外，所以这不是一个最小相位系统。

由(3.7.21)式，与Hz具有相同幅频响应的最小相位系统为

1bz1，a1,b1 Hminz1az1图3.7.4示出了Hz和Hminz的相频响应特性曲线。

10.80.60.40.2()/0.30.20.10-0.2-0.4-0.6()/0-0.1-0.2-0.3-0.8-100.20.40.60.81/1.21.41.61.82-0.400.20.40.60.81/1.21.41.61.82

(a)Hz的相频响应

(b)Hminz的相频响应

图3.7.4

a0.9,b0.4时，例3.7.1的相频响应特性曲线

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <6> 本章小结

本章内容很多，下面我们将本章的内容做一个小结，一方面是帮助大家对过去两周学习的内容做一个回顾，加深一下印象，另一方面，通过这个小结，也帮助大家今后的考前复习。

《数字信号处理》教案——臧博

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <7> 本章作业

《数字信号处理》教案——臧博

第三篇：数字信号处理实验教学-电子教案

数字信号处理实验 电子教案

讲义1 Matlab简介及其安装使用说明..................................................2 讲义2 Matlab基本语句..........................................................................8 讲义3 Matlab基本数值运算................................................................13 讲义4 Matlab函数、及其调用方法....................................................16 实验1 常见离散信号产生和实现.........................................................20 实验2 离散系统的时域分析.................................................................22 实验3 FFT算法的应用..........................................................................24 实验4 离散系统的变换域分析.............................................................27 实验5 有限冲激响应数字滤波器设计.................................................32 实验6 无限冲激响应数字滤波器设计.................................................36 实验7 设计性和研究性实验.................................................................41讲义1 Matlab简介及其安装使用说明

一.MATLAB程序设计语言简介

MATLAB，Matrix Laboratory的缩写，是由Mathworks公司开发的一套用于科学工程计算的可视化高性能语言，具有强大的矩阵运算能力。与大家常用的Fortran和C等高级语言相比，MATLAB的语法规则更简单，更贴近人的思维方式，被称之为“草稿纸式的语言”。截至目前，MATLAB已经发展到7.x版，适用于所有32位的Windows操作系统，按NTFS(NT文件系统)格式下完全安装约需 850 MB。MATLAB软件主要由主包、仿真系统和工具箱三大部分组成。

二.MATLAB应用入门

1.MATLAB的安装与卸载

MATLAB软件在用户接口设计上具有较强的亲和力，其安装过程比较典型，直接运行光盘中的安装向导支撑程序SETUP.exe，按其提示一步步选择即可。MATLAB自身带有卸载程序，在其安装目录下有uninstall子目录，运行该目录下的uninstall.exe即可； 也可以通过Windows系统的安装卸载程序进行卸载。2.MATLAB的启动与退出

MATLAB安装完成后，会自动在Windows桌面上生成一个快捷方式，它是指向安装目录下binwin32matlab.exe的链接，双击它即可来到MATLAB集成环境的基本窗口，通常称之为命令窗口。MATLAB的退出与普通WIN32的程序一样，值得一提的是它有一个自身专有的快捷键Ctrl+Q。3.MATLAB界面简介

图 1 MATLAB基本界面——命令窗口

图2

图3

图4

图5

图6 指令历史

图7 1)菜单栏

菜单栏中包括File、Edit、View、Web、Window和Help六个菜单项。这里着重介绍File、help项。

File项：数据输入/输出的接口，包括10个子项，这里重点介绍其中的5个子项： New：新建文件项。有四个选择：M File(*.M，文本格式的MATLAB程序文件，可以直接通过文件名的方式在MATLAB环境下解释运行；Figure(图形);Model(仿真模型文件)和GUI(可视化界面文件)。

Open：打开所有MATLAB支持的文件格式，系统将自动识别并采用相应的程序对文件进行处理。例如, 打开一个.m文件，系统将自动打开M文件编辑器对它进行编辑。

Set Path...：设置工作路径。可以打开路径设置(Set Path)对话框(图2)，将用户自己建立的目录加入MATLAB的目录系统中，以便所编制的文件能够在MATLAB环境中直接调用。

图8 路径设置对话框

单击Add Folder...按钮可以将你的一个文件夹加入到系统路径中；Add with Subfolders...允许把一个文件夹包括其所有的子文件夹加入到系统路径中。这两种操作均可以直观地在右侧的路径栏内看到结果。

选中一个加入的文件夹，你可以利用Move to Top(移至所有路径的最前面),Move Up(上移一个),Move Down(下移一个),Move to Bottom(移至所有路径的最后面)等四个按钮将改变文件在系统路径中的排列位置以利于对文件的搜索使用，也可以利用Remove按钮将其删除。

对路径操作完毕后，按Save按钮予以保存;Help 项：

Matlab Help：打开以html超文本形式存储的的帮助文件主页； Demos：打开matlab演示窗主页；

About Matlab：Matlab注册图标、版本、制造商和用户信息；

图9 Help选项

图10 Help窗口

2)命令行区

进行通用操作，数值计算，编程和数据类型，输入输出，绘图，用户界面等命令，例如，命令：help函数名（*.m文件）； edit编辑函数、文件

对输入命令的解释MATLAB按以下顺序进行：

① 检查它是否是工作空间中的变量，是则显示变量内容。② 检查它是否是嵌入函数，是则运行之。③ 检查它是否是子函数。④ 检查它是否是私有函数。

⑤ 检查它是否是位于MATLAB搜索路径范围内的函数文件或脚本文件。

请注意，如果有两个以上的方案与输入的命令相匹配，MATLAB将只执行第一个匹配。

讲义2 Matlab基本语句

一.程序控制语句

1.循环语句

MATLAB的循环语句包括for循环和while循环两种类型。1)for循环 语法格式：

for 循环变量 = 起始值：步长：终止值 循环体 end 起始值和终止值为一整形数，步长可以为整数或小数，省略步长时，默认步长为1。执行for循环时，判定循环变量的值是否大于(步长为负时则判定是否小于)终止值，不大于(步长为负时则小于)则执行循环体，执行完毕后加上步长，大于(步长为负时则小于)终止值后退出循环。例1 给矩阵A、B赋值。MATLAB 语句及运行结果如下： k=5;a=zeros(k, k)%矩阵赋零初值 for m=1 : k for n=1: k a(m,n)=1/(m+n-1);end end for i=m :-1 : 1 b(i)=i;end 运行结果：  a= 1.0000 0.5000 0.3333 0.2500?0.2000 0.5000 0.3333 0.2500 0.2000 0.1667 0.3333 0.2500 0.2000 0.1667 0.1429 0.2500 0.2000 0.1667 0.1429 0.1250 0.2000 0.1667 0.1429 0.1250 0.1111 b= 1 2?3 4 5 2)while循环 语法格式：

while 表达式 循环体 end 其执行方式为：若表达式为真(运算值非0)，则执行循环体；结果为0)，则退出循环体，执行end后的语句。

例2 a=3;while a a=a-1 end 输出: a=2 a=1 a=0 2.条件转移语句

条件转移语句有if和switch两种。

若表达式为假(运算 1)if语句

MATLAB中if语句的用法与其他高级语言相类似，其基本语法格式有以下几种： 格式一： if 逻辑表达式 执行语句 end 格式二： if 逻辑表达式 执行语句1 else 执行语句2 end 格式三：

if 逻辑表达式1 执行语句1

else? if 逻辑表达式2 执行语句2 end 2)switch语句

switch语句的用法与其他高级语言相类似，其基本语法格式为：

switch表达式(标量或字符串)case 值1 语句1 case 值2 语句2 … otherwise 语句n end

二.绘图语句 常用的MATLAB绘图语句有figure、plot、subplot、stem等，图形修饰语句有title、axis、text等。

1．figure figure有两种用法，只用一句figure命令，会创建一个新的图形窗口，并返回一个整数型的窗口编号。figure(n)表示将第n号图形窗口作为当前的图形窗口，并将其显示在所有窗口的最前面;如果该图形窗口不存在，则新建一个窗口，并赋以编号n。2．plot 线型绘图函数。用法为plot(x,y,'s')。参数x为横轴变量，y为纵轴变量，s用以控制图形的基本特征如颜色、粗细等，通常可以省略，常用方法如表1所示。

表1

3．Stem 绘制离散序列图，常用格式stem(y)和stem(x,y)分别和相应的plot函数的绘图规则相同，只是用stem命令绘制的是离散序列图。4．subplot subplot(m,n,i)图形显示时分割窗口命令，把一个图形窗口分为m行，n列，m×n个小窗口，并指定第i个小窗口为当前窗口。5．绘图修饰命令

在绘制图形时，我们通常需要为图形添加各种注记以增加可读性。在plot语句后使用title('标题')可以在图形上方添加标题，使用xlabel('标记')或ylabel('标记')为X轴或Y轴添加说明，使用text(X值、Y值、'想加的标示')可以在图形中任意位置添加标示。例3 画图基本语句如图1所示。

MATLAB 语句及运行结果如下： x=0:0.1*pi:2*pi;%定义x向量

figure(1);%创建一个新的图形窗口，编号为1 subplot(2,2,1);%将窗口划分为2行，2列，在第1个窗口中作图 plot(x,sin(x));%画图

title('正弦线');%给图形加标题 subplot(2,2,2);%在第2个窗口中作图 plot(x,sin(x),'r');%画一正弦波，红色 xlabel('X');%给x轴加说明 ylabel('SIN(X)');%给y轴加说明 subplot(2,2,3);%在第3个窗口中作图 plot(x,sin(x),'--');%画一正弦波，破折线 subplot(2,2,4);%在第4个窗口中作图 plot(x,sin(x),'r+');%画一正弦波，红色+线 text(4,0,'注记');

图1 讲义3 Matlab基本数值运算

一.MATLAB内部特殊变量和常数

 MATLAB内部有很多变量和常数，用以表达特殊含义。常用的有：

1.变量ans： 指示当前未定义变量名的答案；

2.常数eps：表示浮点相对精度，其值是从1.0到下一个最大浮点数之间的差值。该变量值作为一些MATLAB函数计算的相对浮点精度，按IEEE标准，如：，近似为2.2204e-016；

3.常数Inf：表示无穷大。当输入或计算中有除以0时产生Inf；

4.虚数单位i,j：表示复数虚部单位，相当于5.NaN：表示不定型值，是由0/0运算产生的。

6.常数pi：表示圆周率π，其值为3.14***…；

；

 MATLAB中可表示的数字的近似范围从，1.有效数字表示的典型例子如下：1234.56789,123456.789E-2,1.23456789e3（format指令可以控制显示格式）

2.复数形式：3.5+4*j,-2.1-7.4*j（i 也可以）

取绝对值：abs()语法格式：abs(x)。当x为实数时计算x的绝对值；x为复数时得到的是复数的模值；x为字符串时得到各字符的ASCII码。取相角：angle()语法格式：angle(z)。求复矢量或复矩阵的相角，结果为一个以弧度为单位介于-π和+π之间的值。

二.变量

1.变量命名规则

MATLAB中对变量的命名应遵循以下规则：

（1）变量名可以由字母、数字和下划线混合组成，但必须以字母开头；（2）字符长度不能大于31；（3）变量命名区分大小写。2.局部变量和全局变量 局部变量是指那些每个函数体内自己定义的，不能从其他函数和MATLAB工作空间访问的变量。全局变量是指用关键字“global”声明的变量。全局变量名应尽量大写，并能反映它本身的含义。如果需要在工作空间和几个函数中都能访问一个全局变量，必须在工作空间和这几个函数中都声明该变量是全局的。

三.矩阵及其运算

MATLAB具有强大的矩阵运算和数据处理功能，对矩阵的处理必须遵从代数规则。1．矩阵生成1）一般矩阵的生成

对于一般的矩阵MATLAB的生成方法有多种。最简单的方法是从键盘直接输入矩阵元素。直接输入矩阵元素时应注意：各元素之间用空格或逗号隔开，用分号或回车结束矩阵行，用中括号把矩阵所有元素括起来。

例1 在工作空间产生一个3×3矩阵A可用MATLAB语言描述如下： A=［1 2 3;4 5 6;7 8 9］ 或 A=［1 2 3 4 5 6 7 8 9］ 运行结果： A= 1 2 3 4 5 6 7 8 9 Size(A)得到矩阵的大小，ans = 3 3 2）特殊矩阵的生成

对于特殊的矩阵可直接调用MATLAB的函数生成。

用函数zeros生成全0矩阵：格式 B=zeros(m,n)生成m×n的全0阵。用函数ones生成全1矩阵：格式 B=ones(m,n)生成m×n的全1阵。

用函数eye生成单位阵：格式 B=eye(m,n)生成m×n矩阵，其中对角线元素全为1，其他元素为0。2．矩阵的运算

矩阵的运算有基本运算和函数运算两种类型。基本运算包括矩阵的加、减、乘、除、乘方、求转置、求逆等，其主要特点是通过MATLAB提供的基本运算符＋、－、*、/（）、^等即可完成。函数运算主要是通过调用MATLAB系统内置的运算函数来求取矩阵的行列式（det(A)）,求秩（rank(A)）,求逆(inv(A))，求特征值和特征向量（［V，D］＝eig(A)）, 求Jordan标准形（jordan(A)）和矩阵分解等。需要用时可以参阅联机帮助和相关参考书。（.*,.,表示逐个元素的乘积和相除；矩阵X/Y 相当于X*inv(Y), XY相当于inv(Y)*X）

例2 矩阵的基本运算。A=［1, 2, 3;4, 5, 6］;B =［6, 5, 4;3, 2, 1］;C =A+B %计算两个矩阵的和 D =B' %计算矩阵B的转置

E=A*D %做矩阵乘法，必须要满足矩阵乘法的基本要求E应该是2阶方阵 F=det(E)%求E的行列式值 G=E^(-1)%求E的逆

输出结果： C= 7 7 7 7 7 7 D= 6 3 5 2 4 1 E= 28 10 73 28 F=54 G= 0.5185-0.1852-1.3519 0.5185 讲义4 Matlab函数、及其调用方法

在MATLAB语言中，M文件有两种形式：脚本和函数。

脚本没有输入/输出参数，只是一些函数和命令的组合。它可以在MATLAB环境下直接执行，也可以访问存在于整个工作空间内的数据。由脚本建立的变量在脚本执行完后仍将保留在工作空间中可以继续对其进行操作，直到使用clear命令对其清除为止。

函数是MATLAB语言的重要组成部分。MATLAB提供的各种工具箱中的M文件几乎都是以函数的形式给出的。函数接收输入参数，返回输出参数，且只能访问该函数本身工作空间中的变量，从命令窗或其他函数中不能对其工作空间的变量进行访问。

1．函数结构

MATLAB语言中提供的函数通常由以下五个部分组成：

（1）函数定义行: 以function开头，函数名(必须与文件名相同)及函数输入输出参数在此定义；（2）H1行：第一注释行，供lookfor和help在线帮助使用；

（3）函数帮助文件；通常包括函数输入输出参数的含义，调用格式说明；

（4）函数体：它包括进行运算和赋值的所有MATLAB程序代码。函数体中可以包括流程控制、输入/输出、计算、赋值、注释以及函数调用和脚本文件调用等。在函数体中完成对输出参数的计算；

（5）注释。

这五个部分中最重要的是函数定义行和函数体。

函数定义行是一个MATLAB函数所必需的，其他各部分的内容可以没有，这种函数称为空函数。

例1 function sa=circle(r,s)% Circle plot a circle of radii r in the line specified by s % r raddi % s line color % sa area of the circle % % circle(r)use blue line to draw a circle of radii r % circle(r,s)use 's' to draw circle % sa=circle(r)compute the area of the circle and draw it in blue % sa=circle(r,s)compute the area of the circle and draw it in color 's' if nargin>2 error('to many input ');end clf;t=0:pi/100:2*pi;x=r*exp(i*t);if nargout==0 plot(x,s);else sa=pi*r*r;fill(real(x),imag(x),s)end axis('square');%makes the current axis box square in size.例2 function y = imp_fun(n,n0)% IMP_FUN Unit impulse function.% % IMP_FUN(N,N0), where N is a vector of sequential integers, % returns a vector the same length as N with zeros everywhere except % N = N0.[m1,n1] = size(n);if ~(m1 == 1 | n1 == 1)error('The sample vector must be one-dimensional.');end y = zeros(m1,n1);i = find(n >= n0);if isempty(i)return end y(i(1))= 1;% where n = n0, set output to 1 2．函数调用 函数调用的过程实际上就是参数传递的过程。例如，在一个脚本文件里调用函数“max”可采用如下方式： n=1:20;

a=sin(2*pi*n/20); [Y,I]=max(a);该调用过程把变量“a”传给了函数中的输入参数“x”,然后把函数运算的返回值传给输出参数“Y”和“I”。其中，Y是a序列的最大值，I是最大值Y对应的坐标值。

例3 构造：y[n] = δ[n-3]: 调用函数： n = 0:6;y = imp_fun(n,3);stem(n,y)

图1 例4 构造：y[n] = 5δ[n]imp_fun(n,2);stem(n,y)；

图2 实验1 常见离散信号产生和实现

一、实验目的

1、加深对常用离散信号的理解；

2、熟悉使用MATLAB在时域中产生一些基本的离散时间信号。

二、实验原理

1、单位抽样序列

在MATLAB中可以利用

函数实现。

2、单位阶越序列

在MATLAB中可以利用

函数实现:

3、正弦序列

在MATLAB中实现过程如下：

4、复指数序列

在MATLAB中实现过程如下：

5、指数序列

在MATLAB中实现过程如下：

三、预习要求

1、预先阅读实验讲义（MATLAB基础介绍）；

2、讨论正弦序列、复指数序列的性质。

A．绘出信号，当、时、、时的信号实部和虚部图；当B．绘出信号

时呢？此时信号周期为多少？ 的频率是多少？周期是多少？产生一个数字频率为0.9的正弦序列，并显示该信号，说明其周期。

3、使用帮助功能学习square(方波),sawtooth(锯齿波)和sinc函数，并绘图。

四、实验内容

编制程序产生上述5种信号，长度可输入确定，函数需要的参数可输入确定，并绘出其图形。

实验2 离散系统的时域分析

一、实验目的

1、熟悉并掌握离散系统的差分方程表示法；

2、加深对冲激响应和卷积分析方法的理解。

二、实验原理

在时域中，离散时间系统对输入信号或者延迟信号进行运算处理，生成具有所需特性的输出信号，具体框图如下:

其输入、输出关系可用以下差分方程描述：

输入信号分解为冲激信号，记系统单位冲激响应，则系统响应为如下的卷积计算式：

当称系统为IIR系统。

时，h[n]是有限长度的（），称系统为FIR系统；反之，三、预习要求

1、在MATLAB中，熟悉利用函数

2、在MATLAB中，熟悉用函数 响应的过程。

实现差分方程的仿真；

计算卷积，用求系统冲激

四、实验内容

1、以下程序中分别使用conv和filter函数计算h和x的卷积y和y1，运行程序，并分析y和y1是否有差别，为什么要使用x[n]补零后的x1来产生y1；具体分析当h[n]有i个值，x[n]有j个值，使用filter完成卷积功能，需要如何补零？ % Program P2_7 clf;h = [3 2 1-2 1 0-4 0 3];%impulse response x = [1-2 3-4 3 2 1];%input sequence y = conv(h,x);n = 0:14;subplot(2,1,1);stem(n,y);xlabel('Time index n');ylabel('Amplitude');title('Output Obtained by Convolution');grid;x1 = [x zeros(1,8)];y1 = filter(h,1,x1);subplot(2,1,2);stem(n,y1);xlabel('Time index n');ylabel('Amplitude');title('Output Generated by Filtering');grid;

2、编制程序求解下列两个系统的单位冲激响应和阶跃响应，并绘出其图形。要求分别用 filter、conv、impz三种函数完成。,，给出理论计算结果和程序计算结果并讨论。实验3 FFT算法的应用

一、实验目的

1、加深对离散信号的DFT的理解；

2、在MATLAB中实现FFT算法。

二、实验原理

N点序列的DFT和IDFT变换定义式如下： , , 利用旋转因子具有周期性，可以得到快速算法（FFT）。

在MATLAB中，可以用函数反变换。

和计算N点序列的DFT正、三、预习要求

1、在MATLAB中，熟悉函数fft、ifft的使用；

2、阅读扩展练习中的实例，学习在MATLAB中的实现FFT算法的实现；

3、利用MATLAB编程完成计算，绘出相应图形。并与理论计算相比较，说明实验结果的原因。

四、实验内容 1、2N点实数序列

N=64。用一个64点的复数FFT程序，一次算出，形。

并绘出 的图

2、已知某序列在单位圆上的N=64等分样点的Z变换为:。

用N点IFFT程序计算出和。

五、扩展练习

例1：对连续的单一频率周期信号 按采样频率和N =16，观察其DFT结果的幅度谱。

采样，截取长度N分别选N =20解：此时离散序列算并作图，函数fft用于计算离散傅里叶变换DFT，程序如下： k=8;n1=[0:1:19];xa1=sin(2*pi*n1/k);subplot(2,2,1)plot(n1,xa1)xlabel('t/T');ylabel('x(n)');xk1=fft(xa1);xk1=abs(xk1);subplot(2,2,2)stem(n1,xk1)xlabel('k');ylabel('X(k)');n2=[0:1:15];，即k=8。用MATLAB计xa2=sin(2*pi*n2/k);subplot(2,2,3)plot(n2,xa2)xlabel('t/T');ylabel('x(n)');xk2=fft(xa2);xk2=abs(xk2);subplot(2,2,4)stem(n2,xk2)xlabel('k');ylabel('X(k)');

图 不同的截取长度的正弦信号及其DFT结果

计算结果示于图，(a)和(b)分别是N=20时的截取信号和DFT结果，由于截取了两个半周期，频谱出现泄漏；(c)和(d)分别是N=16时的截取信号和DFT结果，由于截取了两个整周期，得到单一谱线的频谱。上述频谱的误差主要是由于时域中对信号的非整周期截断产生的频谱泄漏。

实验4 离散系统的变换域分析

一、实验目的

1、熟悉对离散系统的频率响应分析方法；

2、加深对零、极点分布的概念理解。

二、实验原理

离散系统的时域方程为

其变换域分析方法如下： 频域:

系统的频率响应为:

Z域:

系统的转移函数为:

分解因式:，其中和称为零、极点。

三、预习要求

1.在MATLAB中，熟悉函数tf2zp、zplane、freqz、residuez、zp2sos的使用，其中：[z，p，K]=tf2zp（num，den）求得有理分式形式的系统转移函数的零、极点；zplane（z，p）绘制零、极点分布图；h=freqz(num,den,w)求系统的单位频率响应；[r，p，k]=residuez（num，den）完成部分分式展开计算；sos=zp2sos（z，p，K）完成将高阶系统分解为2阶系统的串联。

2.阅读扩展练习中的实例，学习频率分析法在MATLAB中的实现；

3.编程实现系统参数输入，绘出幅度频率响应和相位响应曲线和零、极点分布图。

四、实验内容

求系统 的零、极点和幅度频率响应和相位响应。

五、扩展练习

例1: 求下列直接型系统函数的零、极点，并将它转换成二阶节形式

解:用MATLAB计算程序如下：

num=[1-0.1-0.3-0.3-0.2];den=[1 0.1 0.2 0.2 0.5];[z,p,k]=tf2zp(num,den);m=abs(p);disp('零点');disp(z);disp('极点');disp(p);disp('增益系数');disp(k);sos=zp2sos(z,p,k);disp('二阶节');disp(real(sos));zplane(num,den)输入到“num”和“den”的分别为分子和分母多项式的系数。计算求得零、极点增益系数和二阶节的系数： 零点: 0.9615-0.5730-0.1443 + 0.5850i-0.1443-0.5850i 极点: 0.5276+0.6997i 0.5276-0.6997i-0.5776+0.5635i-0.5776-0.5635i 增益系数: 1 二阶节: 1.0000-0.3885-0.5509 1.0000 1.15520 0.6511 1.0000 0.28850 0.36300 1.0000-1.0552 0.7679 系统函数的二阶节形式为：

极点图见图:

图 系统函数的零、极点图

例2: 差分方程

所对应的系统的频率响应。

解:差分方程所对应的系统函数为：

用MATLAB计算的程序如下： k=256;num=[0.8-0.44 0.36 0.02];den=[1 0.7-0.45-0.6];w=0:pi/k:pi;h=freqz(num,den,w);subplot(2,2,1);plot(w/pi,real(h));grid title('实部')xlabel('omega/pi');ylabel('幅度')subplot(2,2,2);plot(w/pi,imag(h));grid title('虚部')xlabel('omega/pi');ylabel('Amplitude')subplot(2,2,3);plot(w/pi,abs(h));grid title('幅度谱')xlabel('omega/pi');ylabel('幅值')subplot(2,2,4);plot(w/pi,angle(h));grid title('相位谱')xlabel('omega/pi');ylabel('弧度')

图

实验5 有限冲激响应数字滤波器设计

一、实验目的：

1、加深对数字滤波器的常用指标理解。

2、学习数字滤波器的设计方法。

二、实验原理：

低通滤波器的常用指标：

（1）通带边缘频率；

（2）阻带边缘频率；

（3）通带起伏；

（4）通带峰值起伏，（5）阻带起伏，最小阻带衰减。

三、预习要求

1.在MATLAB中，熟悉函数fir1、kaiserord、remezord、remez的使用；

B = fir1(n,Wn,'high','noscale')设计滤波器；

[n,Wn,beta,ftype] = kaiserord(f,a,dev)估计滤波器阶数；

[n,fo,ao,w] = remezord(f,a,dev,fs)计算等波纹滤波器阶数n和加权函数w(ω)；

B=remez(n,f,a)进行等波纹滤波器的设计。

2.阅读扩展练习中的实例，学习FIR滤波器的设计方法及其在MATLAB中的实现；

3.给出FIR数字滤波器的冲激响应，绘出它们的幅度和相位频响曲线，讨论它们各自的实现形式和特点。

数字滤波器有IIR和FIR两种类型，它们的特点和设计方法不同。

四、实验内容：

利用MATLAB编程,分别用窗函数法和等波纹滤波器法设计两种FIR数字滤波器，指标要求如下：

通带边缘频率：，通带峰值起伏：。

阻带边缘频率：，最小阻带衰减：。

五、扩展练习

例1: 用凯塞窗设计一FIR低通滤波器，通带边界频率，阻带衰减

不小于50dB。，阻带边界频率解: 首先由过渡带宽和阻带衰减

来决定凯塞窗的N和

上图给出了以上设计的频率特性，(a)为N=30直接截取的频率特性(b)为凯塞窗设计的频率特性。凯塞窗设计对应的MATLAB程序为： wn=kaiser(30,4.55);nn=[0:1:29];alfa=(30-1)/2;hd=sin(0.4*pi*(nn-alfa))./(pi*(nn-alfa));h=hd.*wn';[h1,w1]=freqz(h,1);或者：b = fir1(29,0.4,kaiser(30,4.55));[h1,w1]=freqz(b,1);plot(w1/pi,20*log10(abs(h1)));axis([0,1,-80,10]);grid;xlabel('归一化频率/p');ylabel('幅度/dB');还可以使用[n,Wn,beta,ftype]=kaiserord(f,a,dev)函数来估计滤波器阶数等，得到凯塞窗滤波器：

fcuts = [0.3 0.5];%归一化频率omega/pi mags = [1 0];devs = [0.05 10^(-2.5)];[n,Wn,beta,ftype] = kaiserord(fcuts,mags,devs);%计算出凯塞窗N，beta的值 hh = fir1(n,Wn,ftype,kaiser(n+1,beta),'noscale');freqz(hh);实际中，一般调用MATLAB信号处理工具箱函数remezord来计算等波纹滤波器阶数N和加权函数W(ω)，调用函数remez可进行等波纹滤波器的设计，直接求出滤波器系数。函数remezord中的数组fedge为通带和阻带边界 例2：利用雷米兹交替算法设计等波纹滤波器，设计一个线性相位低通FIR数字滤波器，其指标为：通带边界频率fc=800Hz，阻带边界fr=1000Hz，通带波动At=40dB，采样频率fs=4000Hz。

阻带最小衰减解:

在MATLAB中可以用remezord 和remez两个函数设计，其结果如图2，MATLAB程序如下： fedge=[800 1000];mval=[1 0];dev=[0.0559 0.01];fs=4000;[N,fpts,mag,wt]=remezord(fedge,mval,dev,fs);b=remez(N,fpts,mag,wt);[h,w]=freqz(b,1,256);plot(w*2000/pi,20*log10(abs(h)));grid;xlabel('频率/Hz');ylabel('幅度/dB');所得图像如下所示：

实验6 无限冲激响应数字滤波器设计

一、实验目的

1、掌握双线性变换法及脉冲相应不变法设计IIR数字滤波器的具体设计方法;

2、熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。

二、实验原理

在MATLAB中，可以用下列函数辅助设计IIR数字滤波器：

1)利用buttord和cheb1ord可以确定低通原型巴特沃斯和切比雪夫滤波器的阶数和截止频率； 2)[num,den]=butter（N,Wn）(巴特沃斯)和[num,den]=cheby1（N,Wn）,[num,den]=cheby2（N,Wn）(切比雪夫1型和2型)可以进行滤波器的设计；

3）lp2hp，lp2bp，lp2bs可以完成低通滤波器到高通、带通、带阻滤波器的转换； 4)使用bilinear可以对模拟滤波器进行双线性变换，求得数字滤波器的传输函数系数； 5)利用impinvar可以完成脉冲响应不变法的模拟滤波器到数字滤波器的转换。

三、预习要求

1.在MATLAB中，熟悉函数butter、cheby1、cheby2的使用，其中：

[num,den]=butter（N,Wn）巴特沃斯滤波器设计； [num,den]=cheby1（N,Wn）切比雪夫1型滤波器设计； [num,den]=cheby2（N,Wn）切比雪夫2型滤波器设计。

2.阅读扩展练习中的实例，学习在MATLAB中进行数字滤波器的设计；

3.给出IIR数字滤波器参数和滤波器的冲激响应，绘出它们的幅度和相位频响曲线，讨论它们各自的实现形式和特点。

四、实验内容

利用MATLAB编程，用脉冲响应不变法和双线性变换法设计一个数字带通滤波器，指标要求如下： 通带边缘频率：阻带边缘频率：，，通带峰值起伏：，最小阻带衰减：

。

五、扩展练习例1：设采样周期T=250μs（采样频率fs =4kHz），用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器，其3dB边界频率为fc =1kHz。[B,A]=butter(3,2*pi*1000,'s');[num1,den1]=impinvar(B,A,4000);[h1,w]=freqz(num1,den1);[B,A]=butter(3,2/0.00025,'s');[num2,den2]=bilinear(B,A,4000);[h2,w]=freqz(num2,den2);f=w/pi*2000;plot(f,abs(h1),'-.',f,abs(h2),'-');grid;xlabel('频率/Hz ')ylabel('幅值/dB')程序中第一个butter的边界频率2π×1000，为脉冲响应不变法原型低通滤波器的边界频率；第二个butter的边界频率2/T=2/0.00025，为双线性变换法原型低通滤波器的边界频率.图1给出了这两种设计方法所得到的频响，虚线为脉冲响应不变法的结果；实线为双线性变换法的结果。脉冲响应不变法由于混叠效应，使得过渡带和阻带的衰减特性变差，并且不存在传输零点。同时，也看到双线性变换法，在z=-1即Ω=π或f=2000Hz处有一个三阶传输零点，这个三阶零点正是模拟滤波器在ω=∞处的三阶传输零点通过映射形成的。

下图给出了MATLAB计算的结果。

例2: 设计一数字高通滤波器，它的通带为400～500Hz，通带内容许有0.5dB的波动，阻带内衰减在小于317Hz的频带内至少为19dB，采样频率为1,000Hz。wc=2*1000*tan(2*pi*400/(2*1000));wt=2*1000*tan(2*pi*317/(2*1000));[N,wn]=cheb1ord(wc,wt,0.5,19,'s');[B,A]=cheby1(N,0.5,wn,'high','s');[num,den]=bilinear(B,A,1000);[h,w]=freqz(num,den);f=w/pi*500;plot(f,20*log10(abs(h)));axis([0,500,-80,10]);grid;xlabel('')ylabel('幅度/dB')下图给出了MATLAB计算的结果。

例3： 设计一巴特沃兹带通滤波器，其３dB边界频率分别为f2=110kHz和f1=90kHz，在阻带f3 = 120kHz处的最小衰减大于１０dB，采样频率fs=400kHz。w1=2*400*tan(2*pi*90/(2*400));w2=2*400*tan(2*pi*110/(2*400));wr=2*400*tan(2*pi*120/(2*400));[N,wn]=buttord([w1 w2],[0 wr],3,10,'s');[B,A]=butter(N,wn,'s');[num,den]=bilinear(B,A,400);[h,w]=freqz(num,den);f=w/pi*200;plot(f,20*log10(abs(h)));axis([40,160,-30,10]);grid;xlabel('频率/kHz')ylabel('幅度/dB')下图给出了MATLAB计算的结果。

例4: 一数字滤波器采样频率fs=1kHz，要求滤除100Hz的干扰，其３dB的边界频率为95Hz和105Hz，原型归一化低通滤波器为： w1=95/500;w2=105/500;[B,A]=butter(1,[w1, w2],'stop');[h,w]=freqz(B,A);f=w/pi*500;plot(f,20*log10(abs(h)));axis([50,150,-30,10]);grid;xlabel('频率/Hz')ylabel('幅度/dB')下图为MATLAB的计算结果。

实验7 设计性和研究性实验

设计性实验1 图像信号的抽取与插值

实验目的

1、熟悉图像处理常用函数和方法；

2、培养通过查阅文献解决问题的能力。实验要求

给出一个二维灰度图像，3、编程实现对该图像的任意比例的放大及缩小；

4、编程实现对该图像的任意角度旋转；

5、解决缩放及旋转时产生的锯齿等不图像不平滑问题。实验提示

6、利用上采样、下采样等方法对信号进行缩放变换；

7、观察对图像进行缩放或旋转时，图像是否会出现锯齿等不平滑现象？

8、分析产生锯齿现象的原因；

9、查阅文献了解解决锯齿现象的方法。（例如平滑滤波、双线性插值、双立方插值等处理）

设计性实验2 语音及音乐信号的采样、滤波

实验目的

1、理解采样率和量化级数对语音信号的影响；

2、设计滤波器解决实际问题。实验要求

利用电脑的声卡录一段语音信号及音乐信号，（1）观察使用不同采样率及量化级数所得到的信号的听觉效果，从而确定对不同信号的最佳的采样率；

（2）分析音乐信号的采样率为什么要比语音的采样率高才能得到较好的听觉效果；（3）注意观察信号中的噪声（特别是50hz交流电信号对录音的干扰，设计一个滤波器去除该噪声。实验提示

（1）推荐录音及播放软件：CoolEdit；

（2）分析语音及音乐信号的频谱，根据信号的频率特性理解采样定律对信号数字化的工程指导意义；

（3）可用带阻滤波器对50Hz交流电噪声进行去噪处理；

（4）也可研究设计自适应滤波器对50Hz噪声及其它随机环境噪声进行滤波处理。设计性实验3 双音多频（DTMF）信号的合成和识别

二、实验目的

1、了解电话按键音形成的原理，理解DTMF音频产生软件和DTMF解码算法；

2、利用FFT算法识别按键音；

三、实验要求

（1）设计音频产生函数，音频信号见下图，每个数据信号持续半秒；（2）实现解码函数：接受（1）产生的DTMF信号，识别信号的频率，并生成包含拨号数字的序列；

四、实验提示

（1）DTFT音频可以用两个正弦波按比例叠加产生；

（2）检测算法可以用FFT算法的DFT，或是用一组滤波器实现；

（3）Goertzel算法可以实现调谐滤波器；

设计性实验4 音乐信号处理

五、实验目的

1、了解回声的产生和梳妆滤波器；

2、混音效果的原理和均衡器的设计；

六、实验要求

（1）设计函数实现一段语音或音乐的回声产生；

（2）设计均衡器，使得得不同频率的混合音频信号，通过一个均衡器后，增强或削减某些频率区域，以便修正低频和高频信号之间的关系；

七、实验提示

（1）回声产生可以使用梳妆滤波器，y(n)=x(n)+ax(n-R), a<1(回声

zRH(z),1R1z衰减系数)；或者传输函数为的全通滤波器实现；比较这两种实现方式的区别，分析为什么会有这样的区别；

（2）可以用许多一阶和二阶参数可调的滤波器级联来实现均衡器的功能，滤波器的结构选择结构要求是调整方便，最好调一个参数只影响一个应用指标，且可调参数少；

第四篇：数字信号处理实验教案

数字信号处理实验教案

信息工程学院-通信工程教研室

数字信号处理是一门理论和实际密切结合的课程，为深入掌握课程内容，最好在学习理论的同时，做习题和上机实验。上机实验不仅可以帮助读者深入的理解和消化基本理论，而且能锻炼同学们的独立解决问题的能力。本讲义在第三版的基础上编写了五个实验，前2个实验属基础性的验证性实验，第3、4、5个实验属基本应用综合性实验。

实验一

离散时间信号的MATLAB实现

实验二

线性卷积与循环卷积的原理及应用

实验三

频率采样定理

实验四

离散系统的因果性和稳定性及频率响应特性

实验五

基于MATLAB的快速傅里叶变换

根据教学进度，理论课结束后进行相关实验。

实验一

时域离散信号的产生

一

实验目的

(1)了解常用的时域离散信号及其特点。

(2)掌握MATLAB产生常用时域离散信号的方法。二 实验内容

(1)编写程序，产生下列离散序列：

A.f(n)=δ(n)

(-3

B．f(n)=e(0.1+j1.6π)n

(0

（2）一个连续的周期性三角波信号频率为50Hz,信号幅度在0~+2V之间，在窗口上显示2个周期信号波形，对信号的一个周期进行16点采样来获取离散信号。试显示原连续信号和采样获得的离散信号波形。

(3)一个连续的周期性方波信号频率为200Hz,信号幅度在-1~+1V之间，在窗口上显示2个周期信号波形，用Fs=4kHz的频率对连续信号进行采样，试显示原连续信号和采样获得的离散信号波形。三 实验步骤

(1)在matlab命令窗口中逐行输入下列语句

>> n1=-3;n2=4;n0=0;%在起点n1、终点n2的范围内，于n0处产生冲激 >> n=n1:n2;

%生成离散信号的时间序列 >> x=[n==n0];

%生成离散信号x(n)>> stem(n,x,'filled');%绘制杆状图，且圆心处用实心圆表示 >> title('单位脉冲序列');>> xlabel('时间（n）');ylabel('幅度x（n）');

在上述语句输入完成之后，敲击回车键，弹出图形窗口，显示出如下图形，即已经满足题干所述条件，产生了 f(n)=δ(n)，(-3

>> n1=16;a=0.1;w=1.6*pi;>> n=0:n1;>> x=exp((a+j*w)*n);>>subplot(2,1,1),stem(n,real(x));%在指定位置描绘图像 >> title('复指数序列的实部');>> subplot(2,1,2),stem(n,imag(x));>> title('复指数序列的虚部');

在上述语句输入完成之后，敲击回车键，弹出图形窗口，显示出如下图形，即已经满足题干所述条件，产生了f(n)=e(0.1+j1.6π)n，(0

>> f=50;Um=1;nt=2;

%输入信号频率、振幅、显示周期 >> N=16;T=1/f;

%N为信号一个采样周期的采样点数，T为信号周期 >> dt=T/N;

%采样时间间隔 >> n=0:nt*N-1;

%建立离散时间的时间序列 >> tn=n*dt;

%确定时间序列样点在时间轴上的位置 >> f=Um*sawtooth(2*f*pi*tn)+1;>> subplot(2,1,1),stem(tn,f);%显示经采样的信号 >> title('离散信号');>> subplot(2,1,2),plot(tn,f);%显示原连续信号 >> title('连续信号');

在上述语句输入完成之后，敲击回车键，弹出图形窗口，显示出如下图形，即已经满足题干所述条件，显示了原连续信号和采样获得的离散信号波形(4)在matlab命令窗口中逐行输入下列语句

>> f=200;Um=1;nt=2;

%输入信号频率、振幅、显示周期 >> Fs=4000;N=Fs/f;T=1/f;%输入采样频率、求采样点数N、T为信号周期 >> dt=T/N;

%采样时间间隔 >> n=0:nt*N-1;

%建立离散时间的时间序列 >> tn=n*dt;

%确定时间序列样点在时间轴上的位置 >> f=Um*sin(2*f*pi*tn);>> subplot(2,1,2),plot(tn,f);%显示原连续信号 >> title('连续信号');>> subplot(2,1,1),stem(tn,f);%显示经采样的信号 >> title('离散信号');

在上述语句输入完成之后，敲击回车键，弹出图形窗口，显示出如下图形，即已经满足题干所述条件，显示了原连续信号和采样获得的离散信号波形 四 思考题

(1)如何在matlab下生产f(n)=3sin(nπ/4)(0

（2）改变实验步骤中最后两个实验的频率参数，分别重新生成相关的信号？

实验二

线性卷积与循环卷积的原理及应用

一、实验目的

（1）掌握两种卷积的原理和两者的异同。

（2）掌握MATLAB实现两种卷积的计算和比较。

二、实验内容

（1）用MATLAB设计线性卷积；（2）调试写出线性卷积和源代码；（3）用MATLAB设计循环卷积； 三 实验步骤 1 线性卷积定理

1）线性卷积的引入

在实际应用中，为了分析时域离散线性非移变系统或者对序列进行滤波处理等，需要计算两个序列的线性卷积。线性卷积既可以在时域中直接计算，也可以通过变换在频域中计算得到。

2）线性卷积的时域计算方法

计算卷积的基本运算是翻转、移位、相乘和相加，这类卷积称为序列的线性卷积。如果两个序列的长度为N和M，那么卷积结果的长度为N+M-1。线性卷积有四步运算：①卷积运算时，y(n)要先反折得到y(-n);②m>0表示y(-n)序列右移，m<0表示左移，不同的m表示不同的值。

假设h(n)和x(n)都是有限长序列，长度分别为N和M，它们的线性卷积可以表示如下：

y2 循环卷积定理

lh(n)x(n)m0h(m)x(nm)N1

MATLAB信号处理工具箱提供了conv函数，该函数用于计算两个有限序列的卷积。

1）循环卷积的引入

为了提高线性卷积的速度，希望用DFT(FFT)计算线性卷积。从而引入循环卷积来运用DFT快速计算线性卷积。循环卷积运用到离散傅立叶变换的循环移位性质，即时域循环移位定理。2）循环卷积的时域计算方法

假设h(n)和x(n)都是有限长序列，长度分别为N和M，它们的L点循环卷积可以表示如下：

[M, ]

LmaxNL称为循环卷积区间长度。n和m的变化区间均是[0,L-1]，直接计算该式比较麻烦。计cm0LLyh(n)x(n)h(m)L1x((nm))R(n)算机中采用矩阵相乘或快速傅里叶变换（FFT）的方法计算循环卷积。用矩阵相乘的方法计算两个序列的循环卷积，这里关键是先形成循环卷积矩阵。如果h(n)的长度N

3）线性卷积与循环卷积的关系

y上式说明，cc(n)qy(nqL)RlL(n)

ly(n)等于y(n)以L为周期的周期延拓序列的主值序列。y(n)的长度

l为NM1，因此只有当循环卷积长度LNM1时，cly(n)以L为周期进行周期延

l(n)y(n)y拓才无混叠现象。因此取其主值序列就满足=。即满足条件LNM1。

四 思考题

(1)比较线性卷积和循环卷积在序列长度不同时两者的联系？(2)试着写出循环卷积的源代码？

实验三 时域采样理论与频域采样定理验证

一、实验目的

1时域采样理论与频域采样理论是数字信号处理中的重要理论。要求掌握模拟信号采样前后频谱的变化，以及如何选择采样频率才能使采样后的信号不丢失信息；要求掌握频率域采样会引起时域周期化的概念，以及频率域采样定理及其对频域采样点数选择的指导作用。

二、实验原理及方法 时域采样定理的要点是：(a)对模拟信号

ˆ(j)xa(t)以间隔T进行时域等间隔理想采样，形成的采样信号的频谱X是原模拟信号频谱为： Xa(j)以采样角频率s（s2/T）为周期进行周期延拓。公式

1Xa(jjns)ˆˆX(j)FT[xa(t)]Tn

a

（b）采样频率s必须大于等于模拟信号最高频率的两倍以上，才能使采样信号的

频谱不产生频谱混叠。

利用计算机计算上式并不方便，下面我们导出另外一个公式，以便用计算机上进行实验。

理想采样信号ˆa(t)xx(t)之间的关系为： 和模拟信号a

ˆa(t)xa(t)(tnT)xn

对上式进行傅立叶变换，得到：

ˆXa(j)[xa(t)(tnT)]ejtdtn

＝nxa(t)(tnT)ejtdt

在上式的积分号内只有当tnT时，才有非零值，因此：

ˆ(j)Xa

上式中，在数值上

nxa(nT)ejnT

xa(nT)＝x(n)，再将T代入，得到：

ˆ(j)Xa

nx(n)ejn

jX(e)，即 上式的右边就是序列的傅立叶变换

ˆ(j)X(ej)XaT

上式说明理想采样信号的傅立叶变换可用相应的采样序列的傅立叶变换得到，只要将自变量ω用T代替即可。

频域采样定理的要点是：

a)对信号x(n)的频谱函数X(ejω)在[0，2π]上等间隔采样N点，得到

XN(k)X(ej)

则N点IDFT[列，公式为： 2kN , k0,1,2,,N1

XN(k)]得到的序列就是原序列x(n)以N为周期进行周期延拓后的主值区序



xN(n)IDFT[XN(k)]N[x(niN)]RN(n)i

(b)由上式可知，频域采样点数N必须大于等于时域离散信号的长度M(即N≥M)，才能使时域不产生混叠，则N点IDFT[果N>M，XN(k)]得到的序列xN(n)就是原序列x(n),即xN(n)=x(n)。如xN(n)比原序列尾部多N-M个零点；如果N

了时域混叠失真，而且

在数字信号处理的应用中，只要涉及时域或者频域采样，都必须服从这两个采样理论的要点。

对比上面叙述的时域采样原理和频域采样原理，得到一个有用的结论，这两个采样理论具有对偶性：“时域采样频谱周期延拓，频域采样时域信号周期延拓”。因此放在一起进行实验。

三 实验步骤

（1）时域采样理论的验证。

tx(t)Aesin(0t)u(t)

a给定模拟信号，式中A=444.128，=502π，0=502πrad/s，它的幅频特性曲线如图10.2.1

图10.2.1

xa(t)的幅频特性曲线

现用DFT(FFT)求该模拟信号的幅频特性，以验证时域采样理论。安照xa(t)的幅频特性曲线，选取三种采样频率，即Fs=1kHz，300Hz，200Hz。观测时间选Tp50ms。为使用DFT，首先用下面公式产生时域离散信号，对三种采样频率，采样序列按顺序用x1(n)，x2(n)，x3(n)表示。

nTx(n)x(nT)Aesin(0nT)u(nT)a

因为采样频率不同，得到的x1(n)，x2(n)，公式

x3(n)的长度不同，长度（点数）用NTpFs计算。选FFT的变换点数为M=64，序列长度不够64的尾部加零。

X(k)=FFT[x(n)]，k=0,1,2,3,-----,M-1

式中k代表的频率为

k2kM。

要求： 编写实验程序，计算x1(n)、x2(n)和析频谱混叠失真。

（2）频域采样理论的验证 给定信号如下：

x3(n)的幅度特性，并绘图显示。观察分

n10n13x(n)27n14n260其它

jX(e)FT[x(n)]在区间[0,2]上等间隔采样32 编写程序分别对频谱函数和16点，得到X32(k)和X16(k)：

X32(k)X(ej)

2k32 , k0,1,2,31

X16(k)X(ej)

再分别对

216k , k0,1,2,15

X32(k)和X16(k)进行32点和16点IFFT，得到x32(n)和x16(n)：

x32(n)IFFTX[32k(3)2] n , 0, 1,2,x16(n)IFFTX[16k(1)6] n , 0, 1,2，31,15jX(e)、X32(k)和X16(k)的幅度谱，并绘图显示x(n)、x32(n)和x16(n)的波形，分别画出进行对比和分析，验证总结频域采样理论。

提示：频域采样用以下方法容易变程序实现。

jX(k)FFT[x(n)]X(e)在[0,2]的3232① 直接调用MATLAB函数fft计算就得到32点频率域采样

② 抽取X32(k)的偶数点即可得到X(ej)在[0,2]的16点频率域采样X16(k)，即X16(k)X32(2k), k0,1,2,,15。

○3 当然也可以按照频域采样理论，先将信号x(n)以16为周期进行周期延拓，取其主值区（16

jX(e)在[0,2]的16点频率域采样X16(k)。点），再对其进行16点DFT(FFT),得到的就是

四 思考题

如果序列x(n)的长度为M，希望得到其频谱()jXeω在]2,0[π上的N点等间隔采样，当N

实验四 离散系统的零极点分析

一、实验目的

1． 熟悉MATLAB的仿真及应用环境

2． 在MATLAB的环境下研究控制系统稳定性

二、实验内容和要求

1.了解离散系统的零极点与系统因果性和稳定性的关系。2.观察离散系统零极点对系统冲激响应的影响。

3.熟悉MATLAB中进行离散系统零极点分析的常用子函数。三 实验步骤

一）MATLAB子函数

1.zplane

功能：显示离散系统的零极点分布图。

调用格式：

zplane(z，p)；绘制由列向量z确定的零点、列向量p确定的极点构成的零极点分布图。

zplane(b，a)；绘制由行向量b和a构成的系统函数确定的零极点分布图。

［hz，hp，ht］＝zplane(z，p)；执行后可得到3个句柄向量：hz为零点线句柄，hp为极点线句柄，ht为坐标轴、单位圆及文本对象的句柄。

2.roots 功能：求多项式的根。

调用格式：

r＝roots(a)；由多项式的分子或分母系数向量求根向量。其中，多项式的分子或分母系数按降幂排列，得到的根向量为列向量。

二）实验原理

1.离散系统的因果性和稳定性

1)因果系统

由理论分析可知，一个离散系统的因果性在时域中必须满足的充分必要条件是：

h(n)＝0

n<0

即系统的冲激响应必须是右序列。

在变换域，极点只能在z平面上一个有界的以原点为中心的圆内。如果系统函数是一个多项式，则分母上z的最高次数应大于分子上z的最高次数。

2)稳定系统

在时域中，离散系统稳定的充分必要条件是：它的冲激响应绝对可加，即



h(n)

n0

在变换域，则要求所有极点必须在z平面上以原点为中心的单位圆内。

3)因果稳定系统

综合系统的因果性和稳定性两方面的要求可知，一个因果稳定系统的充分必要条件是：系统函数的全部极点必须在z平面上以原点为中心的单位圆内。

2.系统极点的位置对系统响应的影响

系统极点的位置对系统响应有着非常明显的影响。下面举例说明系统的极点分别是实数

和复数时的情况，使用MATLAB提供的zplane子函数制作零极点分布图进行分析。

3.系统的因果稳定性实例分析

在MATLAB中提供了roots子函数，用于求多项式的根。配合使用zplane子函数制作零极点分布图，可以帮助我们进行系统因果稳定性的分析。

4.实验任务

(z0.3)H3(z)(z1j)(z1j)

41.6z11.6z24z3H2(z) 10.4z10.35z20.4z3

求该系统的零极点及零极点分布图，并判断系统的因果稳定性。

四 思考题

1结合本次实验与书本上相关原理，对书本后面的习题进行相关的matlab软件仿真？

2因果稳定的离散系统必须满足的充分必要条件是什么？MATLAB提供了哪些进行零极点求解的子函数？如何使用？

实验五 基于MATLAB的快速傅里叶变换

一 实验目的

学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的误差及其原因，以便正确应用FFT。

二 实验原理

用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是N/2π，因此要求DN≤/2π。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。

三 实验步骤及内容

（1）对以下序列进行谱分析。

x1(n)R4(4)

n+1 0≤n≤3 x2(n)={ 8-n

4≤n≤7

0

其它n

4-n

0≤n≤3 X3(n)={ n-3

4≤n≤7

0

其它n

选择FFT的变换区间N为8和16 两种情况进行频谱分析。分别打印其幅频特性曲线。并进行对比、分析和讨论。

（2）对以下周期序列进行谱分析。

x3(n)cosn4

选择FFT的变换区间N为8和16 两种情况分别对以上序列进行频谱分析。分别打印其幅频特性曲线。并进行对比、分析和讨论。

（3）对模拟周期信号进行谱分析

x4(n)cos8tcos16tcos20t

选择 采样频率HzFs64=，变换区间N=16,32,64 三种情况进行谱分析。分别打印其幅频特性，并进行分析和讨论。

四 思考题

（1）对于周期序列，如果周期不知道，如何用FFT进行谱分析？

（2）如何选择FFT的变换区间？（包括非周期信号和周期信号）

第五篇：数字信号处理课程设计..

课程设计报告

课程名称： 数字信号处理 课题名称： 语音信号的处理与滤波

姓 名： 学 号： 院 系： 专业班级： 指导教师： 完成日期： 2013年7月2日

目录

第1部分 课程设计报告………………………………………3 一．设计目的……………………………………………3 二．设计内容……………………………………………3 三．设计原理……………………………………………3 四．具体实现……………………………………………5 1.录制一段声音…………………………………5 2.巴特沃斯滤波器的设计………………………8 3.将声音信号送入滤波器滤波…………………13 4.语音信号的回放………………………………19 5.男女语音信号的频谱分析……………………19 6.噪声的叠加和滤除……………………………22 五． 结果分析……………………………………………27 第2部分 课程设计总结………………………………28 一． 参考文献……………………………………………28

第1部分 课程设计报告

一．设计目的

综合运用本课程的理论知识进行频谱分析以及滤波器设计，通过理论推导得出相应结论，并利用MATLAB作为工具进行实现，从而复习巩固课堂所学的理论知识，提高对所学知识的综合应用能力，并从实践上初步实现对数字信号的处理。

二．设计内容

录制一段个人自己的语音信号，并对录制的信号进行采样；画出采样后语音信号的时域波形和频谱图；给定滤波器的性能指标，采用窗函数法和双线性变换法设计滤波器，并画出滤波器的频率响应；然后用自己设计的滤波器对采集的信号进行滤波，画出滤波后信号的时域波形和频谱，并对滤波前后的信号进行对比，分析信号的变化；回放语音信号；换一个与你性别相异的人录制同样一段语音内容，分析两段内容相同的语音信号频谱之间有什么特点；再录制一段同样长时间的背景噪声叠加到你的语音信号中，分析叠加前后信号频谱的变化，设计一个合适的滤波器，能够把该噪声滤除；

三．设计原理

1.在Matlab软件平台下，利用函数wavrecord()，wavwrite(),wavread(),wavplay()对语音信号进行录制，存储，读取，回放。

2.用y=fft(x)对采集的信号做快速傅立叶变换，并用[h1,w]=freqz(h)进行DTFT变换。

3.掌握FIR DF线性相位的概念，即线性相位对h(n)、H()及零点的约束，了解四种FIR DF的频响特点。

4.在Matlab中，FIR滤波器利用函数fftfilt对信号进行滤波。

5.抽样定理

连续信号经理想抽样后时域、频域发生的变化（理想抽样信号与连续信号频谱之间的关系）

理想抽样信号能否代表原始信号、如何不失真地还原信号即由离散信号恢复连续信号的条件（抽样定理）

理想采样过程描述： 时域描述：

ˆa(t)xa(t)T(t)xa(t)(tnT)xa(nT)(tnT)xnnT(t)频域描述：利用傅氏变换的性质，时域相乘频域卷积，若

n(tnT)ˆa(t)Xa(j)xXa(j)xa(t)T(j)T(t)

则有

ˆ(j)1X(j)(j)XaaT2121ˆXa(j)Xa(jjk)Xa(jjks)TkTTkˆ(j)与X(j)的关系：理想抽样信号的频谱是连续信号频谱的Xaa

周期延拓,重复周期为s(采样角频率)。如果：

X(j)Xa(j)a0s/2s/2即连续信号是带限的，且信号最高频率不超过抽样频率的二分之一，则可不失真恢复。

奈奎斯特采样定理：要使实信号采样后能够不失真还原，采样频率必须大于信号最高频率的两倍：s2h 或 fs2fh

四．具体实现

1.录制一段声音

1.1录制并分析

在MATLAB中用wavrecord、wavread、wavplay、wavwrite对声音进行录制、读取、回放、存储。

程序如下：

Fs=8000;%抽样频率 time=3;%录音时间 fprintf('按Enter键录音%ds',time);%文字提示 pause;%暂停命令 fprintf('录音中......');x=wavrecord(time*Fs,Fs,'double');%录制语音信号 fprintf('录音结束');%文字提示 fprintf('按Enter键回放录音');pause;%暂停命令

wavplay(x,Fs);%按任意键播放语音信号

wavwrite(x,Fs,'C:UsersacerDesktop数字信号sound.wav');%存储语音信号

N=length(x);%返回采样点数 df=fs/N;%采样间隔 n1=1:N/2;f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 figure(2);subplot(2,1,1);plot(x);%录制信号的时域波形 title('原始信号的时域波形');%加标题 ylabel('幅值/A');%显示纵坐标的表示意义 grid;%加网格

y0=fft(x);%快速傅立叶变换 figure(2);subplot(2,1,2);plot(f,abs(y0(n1)));%原始信号的频谱图 title('原始信号的频谱图');%加标题 xlabel('频率w/pi');%显示横坐标表示的意义 ylabel('幅值 ');%显示纵坐标表示的意义 title('原始信号的频谱图');%加标题

grid;%加网格

图1.1 原始信号的时域与频谱图

1.2滤除无效点

针对实际发出声音落后录制动作半拍的现象，如何拔除对无效点的采样的问题: 出现这种现象的原因主要是录音开始时，人的反应慢了半拍，导致出现了一些无效点，而后而出现的无效的点，主要是已经没有声音的动作，先读取声音出来，将原始语音信号时域波形图画出来，根据己得到的信号，可以在第二次读取声音的后面设定采样点，取好有效点，画出滤除无效点后的语音信号时域波形图，对比可以看出。这样就可以解决这个问题。

x=wavread('C:UsersacerDesktop数字信号sound.wav', 7

[4000,24000]);%从4000点截取到24000结束 plot(x);%画出截取后的时域图形 title('截取后的声音时域图形');%标题 xlabel('频率');ylabel('振幅');grid;%画网格

图1.2 去除无效点

2.巴特沃斯滤波器的设计

2.1设计巴特沃思低通滤波器

MATLAB程序如下。滤波器图如图3.3所示。

%低通滤波

fp=1000;fs=1200;Fs=22050;rp=1;rs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Fs1=1;wap=2*tan(wp/2);was=2*tan(ws/2);[N,wc]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs1);figure(1);[h,w]=freqz(Bz,Az,512,Fs1*22050);plot(w,abs(h));title('巴特沃斯低通滤波器');xlabel('频率（HZ）');ylabel('耗损（dB）');gridon;9

图2.1 巴特沃思低通滤波器

2.2设计巴特沃思高通滤波器

MATLAB程序如下。滤波器图如图3.5所示。%高通滤波

fp=4800;fs=5000;Fs=22050;rp=1;rs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;T=1;Fs1=1;wap=2*tan(wp/2);was=2*tan(ws/2);10

[N,wc]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'high','s');[Bz,Az]=bilinear(B,A,Fs1);figure(1);[h,w]=freqz(Bz,Az,512,Fs1*22050);plot(w,abs(h));title('巴特沃斯高通滤波器');xlabel('频率（HZ）');ylabel('耗损（dB）');grid on;

图2.2巴特沃思高通滤波器

2.3设计巴特沃思带通滤波器

MATLAB程序如下。滤波器图如图3.7所示。%带通滤波

fp=[1200,3000];fs=[1000,3200];Fs=8000;rp=1;rs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;T=1;Fs1=1;wap=2*tan(wp/2);was=2*tan(ws/2);[N,wc]=buttord(wap,was,rp,rs,'s');[B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs1);figure(4);[h,w]=freqz(Bz,Az,512,Fs1*1000);plot(w,abs(h));title('巴特沃斯带通滤波器');xlabel('频率（HZ）');ylabel('耗损（dB）');grid on;12

图2.3巴特沃思带通滤波器

3．将声音信号送入滤波器滤波

x=wavread('C:UsersacerDesktop数字信号sound.wav');%播放原始信号

wavplay(x,fs);%播放原始信号 N=length(x);%返回采样点数 df=fs/N;%采样间隔 n1=1:N/2;f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 figure(4);subplot(4,2,1);plot(x);%录制信号的时域波形

title('原始信号的时域波形');%加标题 ylabel('幅值/A');%显示纵坐标的表示意义 grid;%加网格

y0=fft(x);%快速傅立叶变换 subplot(4,2,3);plot(f,abs(y0(n1)));%原始信号的频谱图 title('原始信号的频谱图');%加标题 xlabel('频率w/pi');%显示横坐标表示的意义 ylabel('幅值 ');%显示纵坐标表示的意义 title('原始信号的频谱图');%加标题 grid;%加网格

3．1低通滤波器滤波 fs=8000;beta=10.056;wc=2*pi*1000/fs;ws=2*pi*1200/fs;width=ws-wc;wn=(ws+wc)/2;n=ceil(12.8*pi /width);h=fir1(n,wn/pi,'band',kaiser(n+1,beta));[h1,w]=freqz(h);

ys=fftfilt(h,x);%信号送入滤波器滤波，ys为输出 fftwave=fft(ys);%将滤波后的语音信号进行快速傅立叶变换 figure(4);subplot(4,2,2);%在四行两列的第二个窗口显示图形 plot(ys);%信号的时域波形

title('低通滤波后信号的时域波形');%加标题 xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义 grid;%网格

subplot(4,2,4);%在四行两列的第四个窗口显示图形 plot(f, abs(fftwave(n1)));%绘制模值 xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义

title('低通滤波器滤波后信号的频谱图');%标题 grid;%加网格

wavplay(ys,8000);%播放滤波后信号

3．2高通滤波器滤波 fs=8000;beta=10.056;ws=2*5000/fs;wc=2*4800/fs;

width=ws-wc;wn=(ws+wc)/2;n=ceil(12.8*pi/width);h=fir1(n,wn/pi, 'high',kaiser(n+2,beta));[h1,w]=freqz(h);ys=fftfilt(h,x);%将信号送入高通滤波器滤波 subplot(4,2,5);%在四行两列的第五个窗口显示图形 plot(ys);%信号的时域波形 xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义 title('高通滤波后信号的时域波形');%标题 ylabel('幅值/A');%显示纵坐标的表示意义 grid;%网格

fftwave=fft(ys);%将滤波后的语音信号进行快速傅立叶变换 subplot(4,2,7);%在四行两列的第七个窗口显示图形 plot(f,abs(fftwave(n1)));%绘制模值 axis([0 1 0 50]);xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义

title('高通滤波器滤波后信号的频谱图');%标题 grid;%加网格

wavplay(ys,8000);%播放滤波后信号

3．3带通滤波器 fs=8000;beta=10.056;wc1=2*pi*1000/fs;wc2=2*pi*3200/fs;ws1=2*pi*1200/fs;ws2=2*pi*3000/fs;width=ws1-wc1;wn1=(ws1+wc1)/2;wn2=(ws2+wc2)/2;wn=[wn1 wn2];n=ceil(12.8/width*pi);h=fir1(n,wn/pi,'band',kaiser(n+1,beta));[h1,w]=freqz(h);ys1= fftfilt(h,x);%将信号送入高通滤波器滤波 figure(4);subplot(4,2,6);%在四行两列的第六个窗口显示图形 plot(ys1);%绘制后信号的时域的图形 title('带通滤波后信号的时域波形');%加标题 xlabel('频率w/pi');ylabel('幅值/A');%显示纵坐标表示的意义 grid;%网格

fftwave=fft(ys1);%对滤波后的信号进行快速傅立叶变换 subplot(4,2,8);%在四行两列的第八个窗口显示图形

plot(f, abs(fftwave(n1)));%绘制模值 axis([0 1 0 50]);xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义 title('带通滤波器滤波后信号的频谱图');%加标题 grid;%网格

wavplay(ys1,8000);%播放滤波后信号 图形如下：

原始信号的时域波形幅值/A0-1012x 10原始信号的频谱图34幅值/A1低通滤波后信号的时域波形0.50-0.5012频率w/pi3400.51频率w/pi高通滤波后信号的时域波形幅值/A0幅值/A0幅值/Ax 10高通滤波器滤波后信号的频谱图5012频率w/pi34幅值/A0.20-0.2幅值/A2001000x 10低通滤波器滤波后信号的频谱图200100000.51频率w/pi带通滤波后信号的时域波形0.50-0.501234频率w/pix 10带通滤波器滤波后信号的频谱图50幅值 00.5频率w/pi1000.5频率w/pi1

分析：三个滤波器滤波后的声音与原来的声音都发生了变化。其中低

通的滤波后与原来声音没有很大的变化，其它两个都又明显的变化

4.语音信号的回放

sound(xlow,Fs,bits);%在Matlab中，函数sound可以对声音进行回放，其调用格式: sound(xhigh, Fs,bits);%sound(x, Fs, bits);sound(xdaitong, Fs,bits);5.男女语音信号的频谱分析

5.1 录制一段异性的声音进行频谱分析

Fs=8000;%抽样频率 time=3;%录音时间 fprintf('按Enter键录音%ds',time);%文字提示 pause;%暂停命令 fprintf('录音中......');x=wavrecord(time*Fs,Fs,'double');%录制语音信号 fprintf('录音结束');%文字提示 fprintf('按Enter键回放录音');pause;%暂停命令 wavplay(x,Fs);%按任意键播放语音信号

wavwrite(x,Fs,'C:UsersacerDesktop数字信号sound2.wav');%存储语音信号

5.2 分析男女声音的频谱

x=wavread(' C:UsersacerDesktop数字信号sound2.wav ');%播放原始信号,解决落后半拍

wavplay(x,fs);%播放原始信号 N=length(x);%返回采样点数 df=fs/N;%采样间隔 n1=1:N/2;

f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 figure(1);subplot(2,2,1);plot(x);%录制信号的时域波形

title('原始女生信号的时域波形');%加标题 ylabel('幅值/A');%显示纵坐标的表示意义 grid;%加网格

y0=fft(x);%快速傅立叶变换 subplot(2,2,2);plot(f,abs(y0(n1)));%原始信号的频谱图 title('原始女生信号的频谱图');%加标题 xlabel('频率w/pi');%显示横坐标表示的意义 ylabel('幅值 ');%显示纵坐标表示的意义 grid;%加网格

[y,fs,bits]=wavread(' C:UsersacerDesktop数字信号sound.wav ');% 对语音信号进行采样

wavplay(y,fs);%播放原始信号 N=length(y);%返回采样点数 df=fs/N;%采样间隔 n1=1:N/2;f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 subplot(2,2,3);plot(y);%录制信号的时域波形

title('原始男生信号的时域波形');%加标题 ylabel('幅值/A');%显示纵坐标的表示意义 grid;%加网格

y0=fft(y);%快速傅立叶变换

subplot(2,2,4);%在四行两列的第三个窗口显示图形 plot(f,abs(y0(n1)));%原始信号的频谱图 title('原始男生信号的频谱图');%加标题 xlabel('频率w/pi');%显示横坐标表示的意义 ylabel('幅值 ');%显示纵坐标表示的意义 grid;%加网格

5.3男女声音的频谱图

原始女生信号的时域波形0.50-0.5-1150100原始女生信号的频谱图幅值/A幅值 012345000x 10原始男生信号的时域波形0.50.5频率w/pi原始男生信号的频谱图1300200幅值/A0幅值 012x 1034100-0.5000.5频率w/pi1

图5.3男女声音信号波形与频谱对比

分析：就时域图看，男生的时域图中振幅比女生的高，对于频谱图女生的高频成分比较多

6.噪声的叠加和滤除

6.1录制一段背景噪声

Fs=8000;%抽样频率 time=3;%录音时间 fprintf('按Enter键录音%ds',time);%文字提示 pause;%暂停命令 fprintf('录音中......');x=wavrecord(time*Fs,Fs,'double');%录制语音信号

fprintf('录音结束');%文字提示 fprintf('按Enter键回放录音');pause;%暂停命令 wavplay(x,Fs);%按任意键播放语音信号 wavwrite(x,Fs,'C:UsersacerDesktop数字信号噪音.wav');%存储语音信号

6.2 对噪声进行频谱的分析

[x1,fs,bits]=wavread(' C:UsersacerDesktop数字信号噪音.wav ');%对语音信号进行采样

wavplay(x1,fs);%播放噪声信号 N=length(x1);%返回采样点数 df=fs/N;%采样间隔

n1=1:N/2;f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 figure(5);subplot(3,2,1);plot(x1);%信号的时域波形 title('噪声信号的时域波形');grid;ylabel('幅值/A');y0=fft(x1);%快速傅立叶变换

subplot(3,2,2);plot(f,abs(y0(n1)));%噪声信号的频谱图 ylabel('幅值');title('噪声信号的频谱图');

6.3原始信号与噪音的叠加

fs=8000;[x,fs,bits]=wavread(' C:UsersacerDesktop数字信号sound.wav ');%对录入信号进行采样

[x1,fs,bits]=wavread(' C:UsersacerDesktop数字信号噪音.wav ');%对噪声信号进行采样

yy=x+x1;%将两个声音叠加

6.4叠加信号的频谱分析：

wavplay(yy,fs);%播放叠加后信号 N=length(yy);%返回采样点数 df=fs/N;%采样间隔 n1=1:N/2;f=[(n1-1)*(2*pi/N)]/pi;%频带宽度 figure(5);subplot(3,2,3);plot(yy,'LineWidth',2);%信号的时域波形

title('叠加信号的时域波形');xlabel('时间/t');ylabel('幅值/A');grid;y0=fft(yy);%快速傅立叶变换 subplot(3,2,4);plot(f,abs(y0(n1)));%叠加信号的频谱图 title('叠加信号的频谱图');xlabel('频率w/pi');ylabel('幅值/db');grid;

6.5 设计一个合适的滤波器将噪声滤除 fs=18000;%采样频率 Wp=2*1000/fs;%通带截至频率 Ws=2*2000/fs;%阻带截至频率 Rp=1;%最大衰减 Rs=100;%最小衰减

[N,Wn]=buttord(Wp,Ws,Rp,Rs);%buttord函数（n为阶数，Wn为截至频率）

[num,den]=butter(N,Wn);%butter函数（num为分子系数den为分母系数）

[h,w]=freqz(num,den);%DTFT变换

ys=filter(num,den,yy);%信号送入滤波器滤波，ys为输出 fftwave=fft(ys);%将滤波后的语音信号进行快速傅立叶变换 figure(5);subplot(3,2,5);plot(ys);%信号的时域波形

title('低通滤波后信号的时域波形');%加标题 ylabel('幅值/A');%显示标表示的意义 grid;%网格 subplot(3,2,6);plot(f, abs(fftwave(n1)));%绘制模值 title('低通滤波器滤波后信号的频谱图');%标题 xlabel('频率w/pi');ylabel('幅值/A');%显示标表示的意义 grid;%加网格

wavplay(ys,8000);%播放滤波后信号 grid;图形如下：

噪声信号的时域波形1100噪声信号的频谱图幅值/A0-1幅值0123450000.5叠加信号的频谱图1x 10叠加信号的时域波形10-101时间/t2200幅值/db34幅值/A100000.5频率w/pi1x 10低通滤波后信号的时域波形0.5低通滤波器滤波后信号的频谱图200幅值/A0-0.5幅值/A012x 1034100000.5频率w/pi1

图6.1噪音的叠加与滤除前后频谱对比

7.结果分析

1.录制刚开始时，常会出现实际发出声音落后录制动作半拍，可在[x,fs,bits]=wavread('d:matlavworkwomamaaiwo.wav')加 窗[x,fs,bits]=wavread('d:matlavworkwomamaaiwo.wav'，[100 10000])，窗的长度可根据需要定义。

2.语音信号通过低通滤波器后，把高频滤除，声音变得比较低沉。当通过高通滤波器后，把低频滤除，声音变得比较就尖锐。通过带通滤波器后，声音比较适中。

3.通过观察男生和女生图像知：时域图的振幅大小与性别无关，只与说话人音量大小有关，音量越大，振幅越大。频率图中，女生高 27

频成分较多。

4.叠加噪声后，噪声与原信号明显区分,但通过低通滤波器后，噪声没有滤除，信号产生失真。原因可能为噪声与信号频率相近无法滤除。

第2部分 课程设计总结

通过本次课程设计，使我们对数字信号处理相关知识有了更深刻的理解，尤其是对各种滤波器的设计。在设计的过程中遇到了很多问题，刚刚开始时曾天真的认为只要把以前的程序改了参数就可以用了，可是问题没有我想象中的那么简单，单纯的搬程序是不能解决问题的。通过查阅资料和请教同学收获了很多以前不懂的理论知识。再利用所学的操作，发现所写的程序还是没有能够运行，通过不断地调试，运行，最终得出了需要的结果。整个过程中学到了很多新的知识，特别是对Matlab的使用终于有些了解。在以后的学习中还需要深入了解这方面的内容。在这次的课程设计中让我体会最深的是：知识来不得半点的马虎。也认识到自己的不足，以后要进一步学习。

八．参考文献

[1]数字信号处理教程(第三版)程佩青 清华大学出版社 [2]MATLAB信号处理 刘波 文忠 电子工业出版社 [3]MATLAB7.1及其在信号处理中的应用 王宏 清华大学出版社

[4]MATLAB基础与编程入门 张威 西安电子科技大学出版社

[5] 数字信号处理及其MATLAB实验 赵红怡 张常 化学工业出版社

[6]MATLAB信号处理详解 陈亚勇等 人民邮电出版社 [7] 数字信号处理

钱同惠 机械工业出版社 29