华北电力大学数字信号处理实验六共29页
数字信号处理实验六报告
实验六 频域抽样定理和音频信号的处理实验报告 （一）频域抽样定理给定信号1, 013()27, 14260, n n x n n n +≤≤⎧⎪=-≤≤⎨⎪⎩其它 1．利用DTFT 计算信号的频谱()j X e ω，一个周期内角频率离散为M=1024点，画出频谱图，标明坐标轴。
n=0:100; %设定n 及其取值范围for n1=0:13 %对于n 处于不同的取值范围将n 代入不同的表达式xn(n1+1)=n1+1;endfor n2=14:26xn(n2+1)=27-n2;endfor n3=27:100xn(n3+1)=0;endM=1024; %设定抽样离散点的个数k=0:M-1; %设定k 的取值范围w=2*pi*k/M; %定义数字角频率[X,w] = dtft2( xn,n, M ) %调用dtft2子程序求频谱plot(w,abs(X)); %画出幅度值的连续图像xlabel('w/rad');ylabel('|X(exp(jw))|');title(' M=1024时的信号频谱图像'); %标明图像的横纵坐标和图像标题function [X,w] = dtft2(xn, n, M ) %定义x(n)的DTFT 函数w=0:2*pi/M:2*pi-2*pi/M; %将数字角频率w 离散化L=length(n); %设定L 为序列n 的长度 for (k=1:M) %外层循环，w 循环M 次sum=0; %每确定一个w 值，将sum 赋初值为零for (m=1:L) %内层循环,对n 求和，循环次数为n 的长度sum=sum+xn(m)*exp(-j*w(k)*n(m)); %求和X(k)=sum; %把每一次各x(n)的和的总值赋给X ，然后开始对下一个w 的求和过程end %内层循环结束end%外层循环结束M=1024时的信号频谱图像如图1-1所示：图1-1 M=1024时的信号频谱图像2．分别对信号的频谱()jX eω在区间π[0,2]上等间隔抽样16点和32点，得到32()X k和16()X k。
数字信号处理实验报告
数字信号处理实验报告引言数字信号处理（Digital Signal Processing，DSP）是一门研究数字信号的获取、分析、处理和控制的学科。
在现代科技发展中，数字信号处理在通信、图像处理、音频处理等领域起着重要的作用。
本次实验旨在通过实际操作，深入了解数字信号处理的基本原理和实践技巧。
实验一：离散时间信号的生成与显示在实验开始之前，我们首先需要了解信号的生成与显示方法。
通过数字信号处理器（Digital Signal Processor，DSP）可以轻松生成和显示各种类型的离散时间信号。
实验设置如下：1. 设置采样频率为8kHz。
2. 生成一个正弦信号：频率为1kHz，振幅为1。
3. 生成一个方波信号：频率为1kHz，振幅为1。
4. 将生成的信号通过DAC（Digital-to-Analog Converter）输出到示波器上进行显示。
实验结果如下图所示：（插入示波器显示的正弦信号和方波信号的图片）实验分析：通过示波器的显示结果可以看出，正弦信号在时域上呈现周期性的波形，而方波信号则具有稳定的上下跳变。
这体现了正弦信号和方波信号在时域上的不同特征。
实验二：信号的采样和重构在数字信号处理中，信号的采样是将连续时间信号转化为离散时间信号的过程，信号的重构则是将离散时间信号还原为连续时间信号的过程。
在实际应用中，信号的采样和重构对信号处理的准确性至关重要。
实验设置如下：1. 生成一个正弦信号：频率为1kHz，振幅为1。
2. 设置采样频率为8kHz。
3. 对正弦信号进行采样，得到离散时间信号。
4. 对离散时间信号进行重构，得到连续时间信号。
5. 将重构的信号通过DAC输出到示波器上进行显示。
实验结果如下图所示：（插入示波器显示的连续时间信号和重构信号的图片）实验分析：通过示波器的显示结果可以看出，重构的信号与原信号非常接近，并且能够还原出原信号的形状和特征。
这说明信号的采样和重构方法对于信号处理的准确性有着重要影响。
数字信号处理实验报告_完整版
实验1 利用DFT 分析信号频谱一、实验目的1.加深对DFT 原理的理解。
2.应用DFT 分析信号的频谱。
3.深刻理解利用DFT 分析信号频谱的原理，分析实现过程中出现的现象及解决方法。
二、实验设备与环境 计算机、MATLAB 软件环境 三、实验基础理论1.DFT 与DTFT 的关系有限长序列 的离散时间傅里叶变换 在频率区间 的N 个等间隔分布的点 上的N 个取样值可以由下式表示：212/0()|()()01N jkn j Nk N k X e x n eX k k N πωωπ--====≤≤-∑由上式可知，序列 的N 点DFT ,实际上就是 序列的DTFT 在N 个等间隔频率点 上样本 。
2.利用DFT 求DTFT方法1：由恢复出的方法如下：由图2.1所示流程可知：101()()()N j j nkn j nN n n k X e x n eX k W e N ωωω∞∞----=-∞=-∞=⎡⎤==⎢⎥⎣⎦∑∑∑ 由上式可以得到：IDFTDTFT( )12()()()Nj k kX e X k Nωπφω==-∑ 其中为内插函数12sin(/2)()sin(/2)N j N x eN ωωφω--= 方法2：实际在MATLAB 计算中，上述插值运算不见得是最好的办法。
由于DFT 是DTFT 的取样值，其相邻两个频率样本点的间距为2π/N ，所以如果我们增加数据的长度N ，使得到的DFT 谱线就更加精细，其包络就越接近DTFT 的结果，这样就可以利用DFT 计算DTFT 。
如果没有更多的数据，可以通过补零来增加数据长度。
3.利用DFT 分析连续信号的频谱采用计算机分析连续时间信号的频谱，第一步就是把连续信号离散化，这里需要进行两个操作：一是采样，二是截断。
对于连续时间非周期信号,按采样间隔T 进行采样，阶段长度M ，那么：1()()()M j tj nT a a a n X j x t edt T x nT e ∞--Ω-Ω=-∞Ω==∑⎰对进行N 点频域采样，得到2120()|()()M jkn Na a M kn NTX j T x nT eTX k ππ--Ω==Ω==∑因此，可以将利用DFT 分析连续非周期信号频谱的步骤归纳如下： （1）确定时域采样间隔T ，得到离散序列（2）确定截取长度M ，得到M 点离散序列,这里为窗函数。
华北电力大学 数字信号处理课程设计实验报告
（4） 计算 Y(k)=X1(k)*X2(k);
实
（5） 计算 Y（k）的反变换，即 y(n)=IFFT[X1(k)*X2(k)]。
直接计算 DFT 共需 N*N 次复数乘法和 N(N-1)次复数加法。而 FFT 仅需
验
计算 0.5M 次复数乘法和 M*N 次复数加法。由于在计算机上计算乘法所需的 原 时间比计算加法多得多，所以 FFT 的运算量比 DFT 要少的少。
设计滤波器，首先要对模拟频率进行数字与转换、归一化，采用双线性变 换法还要进行预畸变。接下来依次采用非归一化巴特沃斯模拟滤波器设计函 数、设计双线性变换法的函数、计算离散系统频率响应的函数，最后画出幅 频特性、相频特性、群延迟等图像即可完成图像。 2、源程序： fp=1000; fs=1500; Fs=10000; Rp=1; As=40; T=1/Fs; wp=2*pi*fp*T; ws=2*pi*fs*T; 实 omegap=(2/T)*tan(wp/2); 验 omegas=(2/T)*tan(ws/2); [cs,ds]=afd_buttap(omegap,omegas,Rp,As); 内 [b,a]=bilinear(cs,ds,Fs); [db,mag,pha,grd,w]=freqz_m(b,a); 容 subplot(2,2,1);plot(w/pi,mag);ylabel('幅度'); xlabel('以π为单位的频率');title('幅度响应');axis([0,0.8 0 1]); subplot(2,2,3);plot(w/pi,db);title('幅度响应（dB）');grid; xlabel(' 以 π 为 单位 的 频 率 ');ylabel(' 对 数幅 度 dB');axis([0,0.8 -60 0]); subplot(2,2,2);plot(w/pi,pha);title('相位响应');grid; xlabel('以π为单位的频率');ylabel('相位');axis([0,0.8 -4 4]); subplot(2,2,4);plot(w/pi,grd);title('群延迟');grid; xlabel('以π为单位的频率');ylabel('样本');axis([0,0.8 0 10]);
课程大作业——数字信号处理实验报告
实验一 信号、系统及系统响应一．实验目的1．熟悉理想采样的性质，了解信号采用前后的频谱变化，加深对采样定理的理解。
2．熟悉离散信号和系统的时域特性。
3．熟悉线性卷积的计算编程方法：利用卷积的方法，观察、分析系统响应的时域特性。
4．掌握序列傅氏变换的计算机实现方法，利用序列的傅氏变换对离散信号、系统及系统响应进行频域分析。
二．实验原理1．连续时间信号的采样采样是从连续时间信号到离散时间信号的过渡桥梁，对采样过程的研究不仅可以了解采样前后信号时域和频域特性发生的变化以及信号内容不丢失的条件，而且有助于加深对拉氏变换、傅氏变换、z 变换和序列傅氏变换之间关系的理解。
对一个连续时间信号进行理想采样的过程可以表示为该信号和个周期冲激脉冲的乘积，即)()()(ˆt M t x t xa a = （1－1） 其中)(ˆt xa 是连续信号)(t x a 的理想采样，)(t M 是周期冲激脉冲 ∑+∞-∞=-=n nT t t M )()(δ （1－2）它也可以用傅立叶级数表示为：∑+∞-∞=Ω=n tjm s e T t M 1)( （1－3）其中T 为采样周期，T s /2π=Ω是采样角频率。
设)(s X a 是连续时间信号)(t x a 的双边拉氏变换，即有：⎰+∞∞--=dt e t xs X st aa )()( （1－4）此时理想采样信号)(ˆt xa 的拉氏变换为 ∑⎰+∞-∞=+∞∞--Ω-===m s a sta a jm s X T dt e t x s X )(1)(ˆ)(ˆ （1－5）作为拉氏变换的一种特例，信号理想采样的傅立叶变换[]∑+∞-∞=Ω-Ω=Ωm s a a m j X T j X )(1)(ˆ （1－6）由式（1－5）和式（1－6）可知，信号理想采样后的频谱是原信号频谱的周期延拓，其延拓周期等于采样频率。
根据Shannon 采样定理，如果原信号是带限信号，且采样频率高于原信号最高频率分量的2倍，则采样以后不会发生频率混淆现象。
数字信号处理实验报告
一、实验目的1. 理解数字信号处理的基本概念和原理。
2. 掌握离散时间信号的基本运算和变换方法。
3. 熟悉数字滤波器的设计和实现。
4. 培养实验操作能力和数据分析能力。
二、实验原理数字信号处理（Digital Signal Processing，DSP）是利用计算机对信号进行采样、量化、处理和分析的一种技术。
本实验主要涉及以下内容：1. 离散时间信号：离散时间信号是指时间上离散的信号，通常用序列表示。
2. 离散时间系统的时域分析：分析离散时间系统的时域特性，如稳定性、因果性、线性等。
3. 离散时间信号的变换：包括离散时间傅里叶变换（DTFT）、离散傅里叶变换（DFT）和快速傅里叶变换（FFT）等。
4. 数字滤波器：设计、实现和分析数字滤波器，如低通、高通、带通、带阻滤波器等。
三、实验内容1. 离散时间信号的时域运算（1）实验目的：掌握离散时间信号的时域运算方法。
（2）实验步骤：a. 使用MATLAB生成两个离散时间信号；b. 进行时域运算，如加、减、乘、除等；c. 绘制运算结果的时域波形图。
2. 离散时间信号的变换（1）实验目的：掌握离散时间信号的变换方法。
（2）实验步骤：a. 使用MATLAB生成一个离散时间信号；b. 进行DTFT、DFT和FFT变换；c. 绘制变换结果的频域波形图。
3. 数字滤波器的设计和实现（1）实验目的：掌握数字滤波器的设计和实现方法。
（2）实验步骤：a. 设计一个低通滤波器，如巴特沃斯滤波器、切比雪夫滤波器等；b. 使用MATLAB实现滤波器；c. 使用MATLAB对滤波器进行时域和频域分析。
4. 数字滤波器的应用（1）实验目的：掌握数字滤波器的应用。
（2）实验步骤：a. 采集一段语音信号；b. 使用数字滤波器对语音信号进行降噪处理；c. 比较降噪前后的语音信号，分析滤波器的效果。
四、实验结果与分析1. 离散时间信号的时域运算实验结果显示，通过MATLAB可以方便地进行离散时间信号的时域运算，并绘制出运算结果的时域波形图。
数字信号处理实验报告
数字信号处理实验报告一、实验目的本次数字信号处理实验的主要目的是通过实际操作和观察，深入理解数字信号处理的基本概念和方法，掌握数字信号的采集、处理和分析技术，并能够运用所学知识解决实际问题。
二、实验设备与环境1、计算机一台，安装有 MATLAB 软件。
2、数据采集卡。
三、实验原理1、数字信号的表示与采样数字信号是在时间和幅度上都离散的信号，可以用数字序列来表示。
在采样过程中，根据奈奎斯特采样定理，为了能够准确地恢复原始信号，采样频率必须大于信号最高频率的两倍。
2、离散傅里叶变换（DFT）DFT 是将时域离散信号变换到频域的一种方法。
通过 DFT，可以得到信号的频谱特性，从而分析信号的频率成分。
3、数字滤波器数字滤波器是对数字信号进行滤波处理的系统，分为有限冲激响应（FIR）滤波器和无限冲激响应（IIR）滤波器。
FIR 滤波器具有线性相位特性，而 IIR 滤波器则在性能和实现复杂度上有一定的优势。
四、实验内容与步骤1、信号的采集与生成使用数据采集卡采集一段音频信号，或者在 MATLAB 中生成一个模拟信号，如正弦波、方波等。
2、信号的采样与重构对采集或生成的信号进行采样，然后通过插值算法重构原始信号，观察采样频率对重构信号质量的影响。
3、离散傅里叶变换对采样后的信号进行DFT 变换，得到其频谱，并分析频谱的特点。
4、数字滤波器的设计与实现（1）设计一个低通 FIR 滤波器，截止频率为给定值，观察滤波前后信号的频谱变化。
（2）设计一个高通 IIR 滤波器，截止频率为给定值，比较滤波前后信号的时域和频域特性。
五、实验结果与分析1、信号的采集与生成成功采集到一段音频信号，并在MATLAB 中生成了各种模拟信号，如正弦波、方波等。
通过观察这些信号的时域波形，对不同类型信号的特点有了直观的认识。
2、信号的采样与重构当采样频率足够高时，重构的信号能够较好地恢复原始信号的形状；当采样频率低于奈奎斯特频率时，重构信号出现了失真和混叠现象。
华北电力大学数字信号处理实验报告
华北电力大学实验报告实验环境MATLAB 6.5实验名称提高性实验实验三IIR数字低通滤波器的设计实验目的1.熟悉MATLAB在数字信号处理方面的应用2.深刻理解频域上的频谱特性3.学会数字滤波器的设计4.学会比较不同方法，不同参数设计的滤波器性能比较实验原理设计IIR数字滤波器一般采用间接法（冲激响应不变法和双线性变换法），该设计采用双线性不变法。
IIR的基本设计过程是：首先，将给定的数字滤波器的指标转换成过渡的模拟滤波器的指标；其次，设计过渡模拟滤波器；最后，将过渡模拟滤波器系统函数转换成数字滤波器的系统函数。
设计过程要事先预畸变。
由于从s平面到z平面的映射具有多值性，使得设计出来的数字滤波器不可避免的出现频谱混迭现象。
为了克服脉冲响应不变法可能产生的频谱混叠效应的缺点，使用一种新的变换——双线性变换。
双线性变换法可认为是基于对微分方程的积分，利用对积分的数值逼近的思想。
仿真滤波器的传递函数H(s)为将展开为部份分式的形式，并假设无重复极点，则那么，对于上述函数所表达的数字信号处理系统来讲，其仿真输入x（t）和模拟输出y(t)有如下关系实验环境MATLAB 6.5实验名称提高性实验实验四利用窗函数法实现线性相位FIR数字低通滤波器设计实验目的1、根据ALPF指标，构建物理可实现的线性相位FIR滤波器的冲激响应函数；2、采用多种窗函数，设计线性相位型FIR滤波器4、对比分析多种窗函数法设计的数字滤波器性能。
实验原理FIR滤波器通常采用窗函数方法来设计。
窗设计的基本思想是，首先选择一个适当的理想选频滤波器（它总是具有一个非因果，无限持续时间脉冲响应），然后截取（加窗）它的脉冲响应得到线性相位和因果FIR滤波器。
我们用Hd（e^jw）表示理想的选频滤波器，它在通带上具有单位增益和线性相位，在阻带上具有零响应。
一个带宽wc<pi的低通滤波器由下式给定：为了从hd(n)得到一个FIR滤波器，必须同时在两边截取hd（n）。
数字信号处理实验报告
数字信号处理实验报告实验⼀信号、系统及系统响应⼀、实验⽬的1、熟悉连续信号经理想采样前后的频谱变化关系，加深对时域采样定理的理解；2、熟悉时域离散系统的时域特性；3、利⽤卷积⽅法观察分析系统的时域特性；4、掌握序列傅⽴叶变换的计算机实现⽅法，利⽤序列的傅⽴叶变换对连续信号、离散信号及系统响应进⾏频域分析。
⼆、实验原理及⽅法采样是连续信号数字处理的第⼀个关键环节。
对采样过程的研究不仅可以了解采样前后信号时域和频域特性发⽣变化以及信号信息不丢失的条件，⽽且可以加深对傅⽴叶变换、Z 变换和序列傅⽴叶变换之间关系式的理解。
对⼀个连续信号进⾏理想采样的过程可⽤下式表⽰：，其中为的理想采样，p(t)为周期脉冲，即的傅⽴叶变换为上式表明为的周期延拓。
其延拓周期为采样⾓频率（）。
只有满⾜采样定理时，才不会发⽣频率混叠失真。
在实验时可以⽤序列的傅⽴叶变换来计算。
公式如下：离散信号和系统在时域均可⽤序列来表⽰。
为了在实验中观察分析各种序列的频域特性，通常对在[0，2]上进⾏M点采样来观察分析。
对长度为N的有限长序列x(n)，有：其中，，k=0,1,……M-1时域离散线性⾮移变系统的输⼊/输出关系为上述卷积运算也可在频域实现三、实验内容及步骤1、认真复习采样理论，离散信号与系统，线性卷积，序列的傅⽴叶变换及性质等有关内容，阅读本实验原理与⽅法。
2、编制实验⽤主程序及相应⼦程序。
①信号产⽣⼦程序，⽤于产⽣实验中要⽤到的下列信号序列：xa(t)=Ae-at sin(Ω0t)u(t)进⾏采样，可得到采样序列xa(n)=xa(nT)=Ae-anT sin(Ω0nT)u(n), 0≤n<50其中A为幅度因⼦，a为衰减因⼦，Ω0是模拟⾓频率，T为采样间隔。
这些参数都要在实验过程中由键盘输⼊，产⽣不同的xa(t)和xa(n)。
b. 单位脉冲序列：xb(n)=δ(n)c. 矩形序列：xc(n)=RN(n), N=10②系统单位脉冲响应序列产⽣⼦程序。
数字信号处理实验报告
数字信号处理实验报告实验一：混叠现象的时域与频域表现实验原理：当采样频率Fs不满足采样定理，会在0.5Fs附近引起频谱混叠，造成频谱分析误差。
实验过程：考虑频率分别为3Hz，7Hz，13Hz 的三个余弦信号，即：g1(t)=cos(6πt), g2(t)=cos(14πt), g3(t)=cos(26πt),当采样频率为10Hz 时，即采样间隔为0.1秒，则产生的序列分别为：g1[n]=cos(0.6πn), g2[n]=cos(1.4πn), g3[n]=cos(2.6πn)对g2[n]，g3[n] 稍加变换可得：g2[n]=cos(1.4πn)＝cos((2π-0.6π)n)= cos(0.6πn)g3[n]=cos(2.6πn)= cos((2π+0.6π)n)=cos(0.6πn)利用Matlab进行编程：n=1:300;t=(n-1)*1/300;g1=cos(6*pi*t);g2=cos(14*pi*t);g3=cos(26*pi*t);plot(t,g1,t,g2,t,g3);k=1:100;s=k*0.1;q1=cos(6*pi*s);q2=cos(14*pi*s);q3=cos(26*pi*s);hold on; plot(s(1:10),q1(1:10),'bd');figuresubplot(2,2,1);plot(k/10,abs(fft(q1)))subplot(2,2,2);plot(k/10,abs(fft(q2)))subplot(2,2,3);plot(k/10,abs(fft(q3)))通过Matlab软件的图像如图所示：如果将采样频率改为30Hz，则三信号采样后不会发生频率混叠，可运行以下的程序，观察序列的频谱。
程序编程改动如下：k=1:300;q=cos(6*pi*k/30);q1=cos(14*pi*k/30);q2=cos(26*pi*k/30);subplot(2,2,1);plot(k/10,abs(fft(q)))subplot(2,2,2);plot(k/10,abs(fft(q1)))subplot(2,2,3);plot(k/10,abs(fft(q2)))得图像：问题讨论：保证采样后的信号不发生混叠的条件是什么？若信号的最高频率为17Hz，采样频率为30Hz，问是否会发生频率混叠？混叠成频率为多少Hz的信号？编程验证你的想法。
华北电力大学科技学院数字信号处理课内实验任务书以及代码
《数字信号处理基础》课内实验任 务 书一、 目的与要求1． 掌握《数字信号处理基础》课程的基本理论；2． 掌握应用MATLAB 进行信号处理方面的分析和程序设计方法。
二、 主要内容1. 序列的产生和运算熟悉MATLAB 环境，掌握基本编程方法；熟悉MA TLAB 中序列产生和运算的基本函数。
完成实验指导书中P4中的“练习四”和“练习五”。
练习四：n=0:3;x=[2 3 4 5];y=(x.^3-2*x.^2+x-6.3)./(x.^2+0.05*x-3.14);plot(n,y);练习五：clc,clear;%清屏，清除命令窗口的作用N=64;n=0:N-1;w=randn(1,N);x=2*sin(4*pi*n)+5*sin(8*pi*n)+0.8*w;plot(n,x);2. 因果性数字系统的时域实现掌握实现差分方程()()()k ry n b x n k a x n r =-+-∑∑的程序编写方法。
编程完成下列信号滤波： 512()sin()()3x n n R n π=，滤波器的系统函数为1212132()10.50.6z z H z z z -----+=-+，求()y n 。
代码：clc;clear;n=0:511;x=sin(pi/3*n);b=[1,-3,2];a=[0.6,-0.5];y=dfilter(x,b,a,511);n1=0:length(y)-1;stem(n1,y);3. 离散傅里叶变换（DFT）及其快速算法（FFT）学习DFT和FFT的原理及Matlab编程实现方法。
完成实验指导书中P16中的“练习”。
代码：N=32; %抽取32点t=0:N-1;x2=cos(pi/9*t)+cos(1.5*pi/9*t);Xk2=fft(x2,128); %做128点傅里叶变换Am2=abs(Xk2);n=[0:127];w=(2*pi/128)*n;figure(2)plot(w,Am2,'b');title('Magnitude part');xlabel('frequency in radians');ylabel('|X(exp(jw))|');x3=zeros(1,128); %初始化数组为零N=32;for i=0:N-1 %n=4i时赋值x3(1,4*i+1)=cos(pi/36*i)+cos(1.5*pi/36*i);endXk3=fft(x3); %做128点fftAm3=abs(Xk3);n=[0:127];w=(2*pi/128)*n;figure(3)plot(w,Am3,'b');title('Magnitude part');xlabel('frequency in radians');ylabel('|X(exp(jw))|');%以下是练习改的要求x4=zeros(1,128);N=128;for i=0:N-1x4=cos(pi/36*n)+cos(1.5*pi/36*n);endXk4=fft(x4);Am4=abs(Xk4);n=[0:127];w=(2*pi/128)*n;figure(4)plot(w,Am4,'b');title('Magnitude part');xlabel('frequency in radians'); ylabel('|X(exp(jw))|');附图一张：需要用到的调用函数有（1）function[Am,pha]=dft1(x)N=length(x);w=exp(-j*2*pi/N);for k=1:Nsum=0;for n=1:Nsum=sum+x(n)*w^((k-1)*(n-1));endAm(k)=abs(sum);pha(k)=angle(sum);end（2）function [Am,pha]=dft2(x)N=length(x);n=[0:N-1];k=[0:N-1];w=exp(-j*2*pi/N);nk=n'*k;wnk=w.^(nk);Xk=x*wnk;Am=abs(Xk);pha=angle(Xk);（3）function [Am,pha]=dft3(x)Xk=fft(x);Am=abs(Xk);pha=angle(Xk);4. FFT 的典型应用学习FFT 应用于快速卷积和谱分析时的编程实现方法。
华北电力大学科技学院数字信号处理课程设计
课程设计(综合实验)报告( 2015-- 2016年度第一学期)名称：数字信号处理课程设计题目：MATLAB编程院系：信息工程系班级：13K2学号：31学生姓名：指导教师：孙老师设计周数： 2成绩：日期：2015年12月18日一、课程设计(综合实验)的目的与要求一、 目的与要求1． 掌握《数字信号处理基础》课程的基本理论；2． 掌握应用MATLAB 进行数字信号处理的程序设计方法。
二、 主要内容设计题目及设计要求：已知低通数字滤波器的性能指标如下：0.26p ωπ=，0.75dB p R =，0.41s ωπ=，50dB s A =要求：1. 选择合适的窗函数，设计满足上述指标的数字线性相位FIR 低通滤波器。
用一个图形窗口，包括四个子图，分析显示滤波器的单位冲激响应、相频响应、幅频响应和以dB 为纵坐标的幅频响应曲线。
2. 用双线性变换法，设计满足上述指标的数字Chebyshev I型低通滤波器。
用一个图形窗口，包括三个子图，分析显示滤波器的幅频响应、以dB 为纵坐标的幅频响应和相频响应。
3. 已知模拟信号1234()2sin(2)5sin(2)8cos(2)7.5cos(2)x t f t f t f t f t ππππ=+++其中10.12f kHz =，2 4.98f kHz =，3 3.25f kHz =，4 1.15f kHz =，取采样频率10s f kHz =。
要求：（1）以10s f kHz =对()x t 进行取样，得到()x n 。
用一个图形窗口，包括两个子图，分别显示()x t 以及()x n (0511n ≤≤)的波形；（2）用FFT 对()x n 进行谱分析，要求频率分辨率不超过5Hz 。
求出一个记录长度中的最少点数x N ，并用一个图形窗口，包括两个子图，分别显示()x n 以及()X k 的幅值；（3）用要求1中设计的线性相位低通数字滤波器对()x n 进行滤波，求出滤波器的输出1()y n ，并用FFT 对1()y n 进行谱分析，要求频率分辨率不超过5Hz 。
数字信号处理实验报告
数字信号处理实验报告数字信号处理实验报告一、实验目的本实验旨在通过数字信号处理的方法，对给定的信号进行滤波、频域分析和采样率转换等操作，深入理解数字信号处理的基本原理和技术。
二、实验原理数字信号处理（DSP）是一种利用计算机、数字电路或其他数字设备对信号进行各种处理的技术。
其主要内容包括采样、量化、滤波、变换分析、重建等。
其中，滤波器是数字信号处理中最重要的元件之一，它可以用来提取信号的特征，抑制噪声，增强信号的清晰度。
频域分析是指将时域信号转化为频域信号，从而更好地理解信号的频率特性。
采样率转换则是在不同采样率之间对信号进行转换，以满足不同应用的需求。
三、实验步骤1.信号采集：首先，我们使用实验室的信号采集设备对给定的信号进行采集。
采集的信号包括噪声信号、含有正弦波和方波的混合信号等。
2.数据量化：采集到的信号需要进行量化处理，即将连续的模拟信号转化为离散的数字信号。
这一步通常通过ADC（模数转换器）实现。
3.滤波处理：将量化后的数字信号输入到数字滤波器中。
我们使用不同的滤波器，如低通、高通、带通等，对信号进行滤波处理，以观察不同滤波器对信号的影响。
4.频域分析：将经过滤波处理的信号进行FFT（快速傅里叶变换）处理，将时域信号转化为频域信号，从而可以对其频率特性进行分析。
5.采样率转换：在进行上述处理后，我们还需要对信号进行采样率转换。
我们使用了不同的采样率对信号进行转换，并观察采样率对信号处理结果的影响。
四、实验结果及分析1.滤波处理：经过不同类型滤波器处理后，我们发现低通滤波器可以有效抑制噪声，高通滤波器可以突出高频信号的特征，带通滤波器则可以提取特定频率范围的信号。
这表明不同类型的滤波器在处理不同类型的信号时具有不同的效果。
2.频域分析：通过FFT处理，我们将时域信号转化为频域信号。
在频域分析中，我们可以更清楚地看到信号的频率特性。
例如，对于噪声信号，我们可以看到其频率分布较为均匀；对于含有正弦波和方波的混合信号，我们可以看到其包含了不同频率的分量。
数字信号处理实验文档
实验一：离散信号的MATLAB实现一、实验目的：1、掌握离散时间信号的一般表示方法。
2、熟悉连续信号经理想采样后的频谱变化关系，加深对时域采样定理的理解。
3、掌握离散信号序列的操作。
二、实验内容：M1-1 已知g1(t)=cos(6*pi*t), g1(t)=co 14*pi*t), g1(t)=cos(26*pi*t),以抽样频率fsam=10Hz对上述三个信号进行抽样。
在同一张图上画出g1(t),g2(t)和g3(t)及其抽样点，对所得结果进行讨论。
解：代码如下：100:100)*1/100;g1t=cos(6*pi*t);g2t=cos(14*pi*t);g3t=cos(26*pi*t);subplot(3,1,1);plot(g1t);subplot(3,1,2);plot(g2t);subplot(3,1,3);plot(g3t);绘出的图形如图1_1所示:图1_1采样频率为fsam=10Hz，采样时间为0.1s，而f1=3Hz，f2=7Hz，f3=13Hz，使得三个信号的采样图形相似，这样不能很好还原原来的信号图像。
所以对信号的采样频率应足够大，应满足fsam>=2fm.M1-2利用MATLAB的filter函数，求出下列系统的单位脉冲响应，并判断系统是否稳定。
讨论题所获得的结果。
代码1：k=1:300;x=zeros(1,300);x(1)=1;b1=[1];a1=[1,-1.845,0.850586];y1=filter(b1,a1,x);subplot(1,2,1);plot(k,y1);xlabel('k');ylabel('幅度y1');b2=[1];a2=[1,-1.85,0.85];y2=filter(b2,a2,x);subplot(1,2,2);plot(k,y2);xlabel('k');ylabel('幅度y2');图1_2_1代码2：x=zeros(1,500);x(1)=1;b1=[1];a1=[1,-1.845,0.850586];y1=filter(b1,a1,x);plot(k,y1);b2=[1];a2=[1,-1.85,0.85];y2=filter(b2,a2,x);plot(k,y1,k,y2,':');xlabel('k');ylabel('幅度');legend('y1''y2');图1_2_2结论：H1（z）的两个极点都在单位圆内，所以系统稳定，从图中可以看出响应曲线升高后有回落，系统最终趋向于0；H2（z）的一个极点在单位圆内，另一个在单位圆上，所以系统最终临界稳定，从图中可以看出响应曲线上升后没有回落，系统最终趋向于6.7左右。
数字信号处理课程设计实验报告
华北电力大学实验报告实验环境MATLAB 7.1实验名称实验一：FFT的应用实验目的1、熟悉MATLAB在数字信号处理中的应用。
2、掌握利用FFT计算序列线性卷积的基本原理及编程实现。
3、掌握对连续信号进行采样的基本原理和方法，并利用FFT对信号进行频谱分析。
实验原理1.线性卷积和圆周卷积设x(n)为L点序列，h(n)为M点序列，x(n)和h(n)的线性卷积为的长度为L+M-1。
x(n)和h(n)的N点圆周卷积为圆周卷积与线性卷积相等而不产生交叠的必要条件为圆周卷积定理：根据DFT的性质，x(n)和h(n)的N点圆周卷积的DFT等于它们DFT的乘积2.快速卷积快速卷积算法用圆周卷积实现线性卷积，根据圆周卷积定理利用FFT算法实现圆周卷积。
可以将快速卷积的步骤归纳如下：(1)为了使线性卷积可以用圆周卷积来计算，必须选择;同时为了能使用基-2FFT完成卷积运算，要求。
采用补零的办法使x(n)实验原理和h(n)的长度均为N。
(2)计算x(n)和h(n)的N点FFT(3)组成乘积(4)利用IFFT计算Y(K)的IDFT，得到线性卷积y(n)3.分段卷积我们考察单位取样响应为h(n)的线性系统，输入为x(n),输出为y(n)，则当输入序列x(n)极长时，如果要等x(n)全部集齐时再开始进行卷积，会使输出相对输入有较大的延时，再者如果序列太长，需要大量存贮单元。
为此我们把x(n)分段，分别求出每段的卷积，合在一起得到最后总的输出。
这种方法称为分段卷积分段卷积可细分为重叠保留法和重叠相加法。
重叠保留法：设x(n)的长度为，h(n)的长度为M。
我们把序列x(n)分成多段N点序列(n),每段与前一段重叠M-1个样本。
由于第一段没有前一段保留信号，为了修正，我们在第一个输入端前面填充M-1个零。
计算每一段与h(n)的圆周卷积，则其每段卷积结果的前M-1个样本不等于线性卷积值，不是正确的样本值。
所以我们将每段卷积结果的前M-1个样本社区，只保留后面的N-M+1个正确输出样本，把这些输出样本合起来，得到总输出。
数字信号处理实验报告(全)
实验一、离散时间系统及离散卷积1、单位脉冲响应源程序：function pr1() %定义函数pr1a=[1,-1,0.9]; %定义差分方程y(n)-y(n-1)+0.9y(n-2)=x(n) b=1;x=impseq(0,-20,120); %调用impseq函数n=[-40:140]; %定义n从-20 到120h=filter(b,a,x); %调用函数给纵座标赋值figure(1) %绘图figure 1 (冲激响应) stem(n,h); %在图中绘出冲激title('冲激响应'); %定义标题为:'冲激响应'xlabel('n'); %绘图横座标为nylabel('h(n)'); %绘图纵座标为h(n)figure(2) %绘图figure 2[z,p,g]=tf2zp(b,a); %绘出零极点图zplane(z,p)function [x,n]=impseq(n0,n1,n2) %声明impseq函数n=[n1:n2];x=[(n-n0)==0];结果：Figure 1:Figure 2:2、离散系统的幅频、相频的分析源程序：function pr2()b=[0.0181,0.0543,0.0543,0.0181];a=[1.000,-1.76,1.1829,-0.2781];m=0:length(b)-1; %m从0 到3l=0:length(a)-1; %l从0 到3K=5000;k=1:K;w=pi*k/K; %角频率wH=(b*exp(-j*m'*w))./(a*exp(-j*l'*w));%对系统函数的定义magH=abs(H); %magH为幅度angH=angle(H); %angH为相位figure(1)subplot(2,1,1); %在同一窗口的上半部分绘图plot(w/pi,magH); %绘制w(pi)-magH的图形grid;axis([0,1,0,1]); %限制横纵座标从0到1xlabel('w(pi)'); %x座标为 w(pi)ylabel('|H|'); %y座标为 angle(H)title('幅度，相位响应'); %图的标题为:'幅度，相位响应' subplot(2,1,2); %在同一窗口的下半部分绘图plot(w/pi,angH); %绘制w(pi)-angH的图形grid; %为座标添加名称xlabel('w(pi)'); %x座标为 w(pi)ylabel('angle(H)'); %y座标为 angle(H)结果：3、卷积计算源程序：function pr3()n=-5:50; %声明n从-5到50u1=stepseq(0,-5,50); %调用stepseq函数声用明u1=u(n)u2=stepseq(10,-5,50); %调用stepseq函数声用明u2=u(n-10) %输入x(n)和冲激响应h(n)x=u1-u2; %x(n)=u(n)-u(n-10)h=((0.9).^n).*u1; %h(n)=0.9^n*u(n)figure(1)subplot(3,1,1); %绘制第一个子图stem(n,x); %绘制图中的冲激axis([-5,50,0,2]); %限定横纵座标的范围title('输入序列'); %规定标题为:'输入序列'xlabel('n'); %横轴为nylabel('x(n)'); %纵轴为x(n)subplot(3,1,2); %绘制第二个子图stem(n,h); %绘制图中的冲激axis([-5,50,0,2]); %限定横纵座标的范围title('冲激响应序列'); %规定标题为:'冲激响应序列'xlabel('n'); %横轴为nylabel('h(n)'); %纵轴为h(n)%输出响应[y,ny]=conv_m(x,n,h,n); %调用conv_m函数subplot(3,1,3); %绘制第三个子图stem(ny,y);axis([-5,50,0,8]);title('输出响应'); %规定标题为:'输出响应'xlabel('n');ylabel('y(n)'); %纵轴为y(n)%stepseq.m子程序%实现当n>=n0时x(n)的值为1function [x,n]=stepseq(n0,n1,n2)n=n1:n2;x=[(n-n0)>=0];%con_m的子程序%实现卷积的计算function [y,ny]=conv_m(x,nx,h,nh)nyb=nx(1)+nh(1);nye=nx(length(x))+nh(length(h));ny=[nyb:nye];y=conv(x,h);结果：实验二、离散傅立叶变换与快速傅立叶变换1、离散傅立叶变换（DFT）源程序：function pr4()F=50;N=64;T=0.000625;n=1:N;x=cos(2*pi*F*n*T); %x(n)=cos(pi*n/16)subplot(2,1,1); %绘制第一个子图x(n)stem(n,x); %绘制冲激title('x(n)'); %标题为x(n)xlabel('n'); %横座标为nX=dft(x,N); %调用dft函数计算x(n)的傅里叶变换magX=abs(X); %取变换的幅值subplot(2,1,2); %绘制第二个子图DFT|X|stem(n,X);title('DFT|X|');xlabel('f(pi)'); %横座标为f(pi)%dft的子程序%实现离散傅里叶变换function [Xk]=dft(xn,N)n=0:N-1;k=0:N-1;WN=exp(-j*2*pi/N);nk=n'*k;WNnk=WN.^nk;Xk=xn*WNnk;结果：F=50,N=64,T=0.000625时的波形F=50,N=32,T=0.000625时的波形:2、快速傅立叶变换（FFT）源程序：%function pr5()F=50;N=64;T=0.000625;n=1:N;x=cos(2*pi*F*n*T); %x(n)=cos(pi*n/16) subplot(2,1,1);plot(n,x);title('x(n)');xlabel('n'); %在第一个子窗中绘图x(n)X=fft(x);magX=abs(X);subplot(2,1,2);plot(n,X);title('DTFT|X|');xlabel('f(pi)'); %在第二个子图中绘图x(n)的快速傅%里叶变换结果：3、卷积的快速算法源程序：function pr6()n=0:14;x=1.^n;h=(4/5).^n;x(15:32)=0;h(15:32)=0;%到此 x(n)=1, n=0~14; x(n)=0,n=15~32% h(n)=(4/5)^n, n=0~14; h(n)=0,n=15~32subplot(3,1,1);stem(x);title('x(n)');axis([1,32,0,1.5]); %在第一个子窗绘图x(n)横轴从1到32,纵轴从0到1.5 subplot(3,1,2);stem(h);title('h(n)');axis([1,32,0,1.5]); %在第二个子窗绘图h(n)横轴从1到32,纵轴从0到1.5 X=fft(x); %X(n)为x(n)的快速傅里叶变换H=fft(h); %H(n)为h(n)的快速傅里叶变换Y=X.*H; %Y(n)=X(n)*H(n)%Y=conv(x,h);y=ifft(Y); %y(n)为Y(n)的傅里叶反变换subplot(3,1,3) %在第三个子窗绘图y(n)横轴从1到32,纵轴从0到6 stem(abs(y));title('y(n=x(n)*h(n))');axis([1,32,0,6]);结果：实验三、IIR数字滤波器设计源程序：function pr7()wp=0.2*pi;ws=0.3*pi;Rp=1;As=25;T=1;Fs=1/T;OmegaP=(2/T)*tan(wp/2); %OmegaP(w)=2*tan(0.1*pi) OmegaS=(2/T)*tan(ws/2); %OmegaS(w)=2*tan(0.15*pi)ep=sqrt(10^(Rp/10)-1);Ripple=sqrt(1/(1+ep.^2));Attn=1/10^(As/20);N=ceil((log10((10^(Rp/10)-1)/(10^(As/10)-1)))/(2*log10(OmegaP/OmegaS) ));OmegaC=OmegaP/((10.^(Rp/10)-1).^(1/(2*N)));[cs,ds]=u_buttap(N,OmegaC);[b,a]=bilinear(cs,ds,Fs);[mag,db,pha,w]=freqz_m(b,a);subplot(3,1,1); %在第一个子窗绘制幅度响应的图形plot(w/pi,mag);title('幅度响应');xlabel('w(pi)');ylabel('H');axis([0,1,0,1.1]);set(gca,'XTickmode','manual','XTick',[0,0.2,0.35,1.1]);set(gca,'YTickmode','manual','YTick',[0,Attn,Ripple,1]);grid;subplot(3,1,2); %在第二个子窗以分贝为单位绘制幅度响应的图形plot(w/pi,db);title('幅度响应(dB)');xlabel('w(pi)');ylabel('H');axis([0,1,-40,5]);set(gca,'XTickmode','manual','XTick',[0,0.2,0.35,1.1]);set(gca,'YTickmode','manual','YTick',[-50,-15,-1,0]);grid;subplot(3,1,3); %在第三个子窗绘制相位响应的图形plot(w/pi,pha);title('相位响应');xlabel('w(pi)');ylabel('pi unit');%axis([0,1,0,1.1]);set(gca,'XTickmode','manual','XTick',[0,0.2,0.35,1.1]);set(gca,'YTickmode','manual','YTick',[-1,0,1]);grid;function [b,a]=u_buttap(N,OmegaC)[z,p,k]=buttap(N);p=p*OmegaC;k=k*OmegaC.^N;B=real(poly(z));b0=k;b=k*B;a=real(poly(p));function [mag,db,pha,w]=freqz_m(b,a)[H,w]=freqz(b,a,1000,'whole');H=(H(1:501))';w=(w(1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);结果：实验四、FIR数字滤波器的设计源程序：function pr8()wp=0.2*pi;ws=0.35*pi;tr_width=ws-wp;M=ceil(6.6*pi/tr_width)+1;n=0:M-1;wc=(ws+wp)/2;alpha=(M-1)/2;m=n-alpha+eps;hd=sin(wc*m)./(pi*m);w_ham=(hamming(M))';h=hd.*w_ham;[mag,db,pha,w]=freqz_m(h,[1]);delta_w=2*pi/1000;Rp=-(min(db(1:wp/delta_w+1)));As=-round(max(db(ws/delta_w+1:501)));subplot(2,2,1);stem(n,hd);title('理想冲激响应');axis([0,M-1,-0.1,0.3]);ylabel('hd(n)');subplot(2,2,2);stem(n,h);title('实际冲激响应');axis([0,M-1,-0.1,0.3]);ylabel('h(n)');subplot(2,2,3);plot(w/pi,pha);title('滤波器相位响应');axis([0,1,-pi,pi]);ylabel('pha');set(gca,'XTickmode','manual','XTick',[0,0.2,0.3,1.1]); set(gca,'YTickmode','manual','YTick',[-pi,0,pi]); grid;subplot(2,2,4);plot(w/pi,db);title('滤波器幅度响应');axis([0,1,-100,10]);ylabel('H(db)');set(gca,'XTickmode','manual','XTick',[0,0.2,0.3,1.1]); set(gca,'YTickmode','manual','YTick',[-50,-15,0]);function [mag,db,pha,w]=freqz_m(b,a)[H,w]=freqz(b,a,1000,'whole');H=(H(1:501))';w=(w(1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);结果：。
实验三-华北电力大学-数字信号处理实验
实验报告实验名称____________ ____课程名称____________ ____院系部:专业班级：学生姓名：学号:同组人:实验台号:指导教师：成绩：实验日期:华北电力大学1.实验目的分析常用窗函数的时域和频域特性，灵活运用窗函数分析信号频谱和设计FIR 数字滤波器。
2.实验原理在确定信号谱分析、随机信号功率谱估计以及FIR 数字滤波器设计中，窗函数的选择起着重要的作用。
在信号的频谱分析中，截短无穷长的序列会造成频率泄漏，影响频谱分析的精度和质量。
合理选取窗函数的类型，可以改善泄漏现象。
在FIR 数字滤波器设计中，截短无穷长的系统单位脉冲序列会造成FIR 滤波器幅度特性的波动，且出现过渡带。
3.实验内容及步骤(1) 1. 分析并绘出常用窗函数的时域特性波形。
2. 利用fft 函数分析常用窗函数的频域特性, 并从主瓣宽度和 旁瓣相对幅度两个角度进行比较分析。
3. 研究凯塞窗(Kaiser)的参数选择对其时域和频域的影响。
(1) 固定beta=4，分别取N =20, 60, 110；(2) 固定N =60，分别取beta=1,5,11。
4. 序列 ⎪⎪⎭⎫ ⎝⎛+⎪⎪⎭⎫ ⎝⎛=k k k x 20π9cos 20π11cos 5.0][ ，分析其频谱。
(1) 利用不同宽度N 的矩形窗截短该序列, N 分别为20,40,160，观察不同长度N 的窗对谱分析结果的影响；(2) 利用哈明窗重做 (1)；(3) 利用凯塞窗重做 (1)；(4) 比较和分析三种窗的结果；(5) 总结不同长度或类型的窗函数对谱分析结果的影响。
4.数据处理与总结1.分析并绘出常用窗函数的时域特性波形。
程序如下:clear;subplot(2,3,1);N=51;w=boxcar(N);stem(w)title('矩形窗')subplot(2,3,2);w=hanning(N);stem(w)title('Hanning窗')subplot(2,3,3);w=hamming(N);stem(w)title('Hamming窗')subplot(2,3,4);w=blackman(N);stem(w)title('blackman窗')subplot(2,3,5);w=bartlett(N);stem(w)title('三角形窗')subplot(2,3,6);w=kaiser(N);stem(w)title('kaiser窗')2,利用fft函数分析常用窗函数的频域特性clear;N=51;w=boxcar(N);y=fft(w,200);subplot(3,3,1);stem([0:N-1],w);title('时域波形');subplot(3,3,2);y0= abs(fftshift(y));plot([-100:99],y0);title('矩形窗频域'); subplot(3,3,3);w=hanning(N);y=fft(w,200);y0= abs(fftshift(y));plot([-100:99],y0);title('hanning窗频域'); subplot(3,3,4);w=hamming(N);y=fft(w,200);y0= abs(fftshift(y));plot([-100:99],y0);title('哈明窗频域');subplot(3,3,5);w=blackman(N);y=fft(w,200);y0= abs(fftshift(y));plot([-100:99],y0);title('布莱克曼窗频域');subplot(3,3,6);w=bartlett(N);y=fft(w,200);y0= abs(fftshift(y));plot([-100:99],y0);title('三角形窗频域');subplot(3,3,7);w=kaiser(N);y=fft(w,200);y0= abs(fftshift(y));plot([-100:99],y0);title('kaiser窗频域');3. 研究凯塞窗(Kaiser)的参数选择对其时域和频域的影响。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验六 IIR数字滤波器设计及应用一：实验目的加深理解IIR数字滤波器的特性，掌握IIR数字滤波器的设计原理与设计方法，以及IIR数字滤波器的应用。
二：实验原理N阶IIR数字滤波器的系统函数为:IIR数字滤波器的设计主要通过成熟的模拟滤波器设计方法来实现：将数字滤波器设计指标转换为模拟滤波器设计指标，设计出相应的模拟滤波器H(s)，再经过脉冲响应不变法或双线性变换法得到所需的IIR数字滤波器H(z)。
IIR数字滤波器设计的重要环节是模拟原型低通滤波器的设计，主要包括Butterworth、Chebyshev和椭圆等滤波器。
MATLAB 信号处理工具箱中提供了IIR滤波器设计的函数。
IIR 滤波器阶数选择buttord -巴特沃斯(Butterworth)滤波器阶数选择。
cheb1ord -切比雪夫(Chebyshev)I 型滤波器阶数选择。
cheb2ord -切比雪夫(Chebyshev)II 型滤波器阶数选择。
ellipord -椭圆(Elliptic)滤波器阶数选择。
IIR 滤波器设计butter -巴特沃斯(Butterworth)滤波器设计cheby1 -切比雪夫(Chebyshev)I 型滤波器设计cheby2 -切比雪夫(Chebyshev)II 型滤波器设计ellip -椭圆(Elliptic)滤波器设计maxflat -通用的巴特沃斯(Butterworth)低通滤波器设计yulewalk -Yule-Walker 滤波器设计(直接数字滤波器设计法)1. Butterworth滤波器设计Butterworth滤波器是通带、阻带都单调衰减的滤波器。
(1)调用buttord函数确定巴特沃斯滤波器的阶数，格式为 [N,Wc] = buttord(Wp,Ws,Ap,As)输入参数：Ap,As为通带最大衰减和阻带最小衰减，以dB为单位。
Wp,Ws为归一化通带截频和阻带截频，0<Wp,Ws<1 。
输出参数：N为滤波器的阶数；Wc为截频，0 < Wc < 1。
(2)调用butter函数设计出巴特沃斯滤波器，格式为[b,a] = butter(N,Wc,options)输入参数：N和Wc是buttord函数返回的参数，含义见上。
Options=’low’, ’high’, ’bandpass’, ’stop’, 分别对应低通、高通、带通、带阻，默认情况下为低通或带通。
输出参数：b和a为设计出的IIR数字滤波器H(s)的分子多项式和分母多项式的系数矩阵。
2. Chebyshev I型滤波器设计Chebyshev I型滤波器为通带纹波控制器：在通带呈现纹波特性，在阻带单调衰减。
[N,Wc] = cheb1ord(Wp, Ws, Ap, As)[b,a] = cheby1(N,Ap,Wc,options)参数含义与butter中参数一致。
2. Chebyshev II 型滤波器设计Chebyshev II 型滤波器为阻带纹波控制器：在阻带呈现纹波特性。
[N,Wc] = cheb2ord(Wp, Ws, Ap, As)[b,a] = cheby2(N,As,Wc,options)3. 椭圆滤波器设计椭圆滤波器在通阻带都呈现纹波特性。
[N,Wc] = ellipord(Wp,Ws,Ap,As)[b,a] = ellip(N,Ap,As,Wc,options)三：实验内容1（1）[N,Wc]=buttord(0.250,0.677,3,60)[b,a]=butter(N,Wc)freqz(b,a);axis([0,1,-120,0]);grid ontitle('巴特沃斯低通数字滤波器')（2）[N,Wc]=buttord(0.250,0.677,3,60)[b,a]=butter(N,Wc,'high')freqz(b,a);axis([0,1,-120,0]);grid ontitle('巴特沃斯高通数字滤波器')（3）Wp =[0.25 0.67]; Ws =[0.25-0.03 0.67+0.03];Rp = 3;Rs = 60;[N,Wc]=buttord(Wp,Ws,Rp,Rs)[b,a]=butter(N,Wc,'bandpass')freqz(b,a);axis([0,1,-120,0]);grid ontitle('巴特沃斯带通数字滤波器')N =40Wc =0.2499 0.6701b =Columns 1 through 90.0000 0 -0.0000 0 0.0000 0 -0.0000 0 0.0000Columns 10 through 180 -0.0000 0 0.0000 0 -0.0000 0 0.0000 0Columns 19 through 27-0.0001 0 0.0003 0 -0.0007 00.0017 0 -0.0037Columns 28 through 360 0.0072 0 -0.0125 0 0.0195 0 -0.0276 0Columns 37 through 450.0353 0 -0.0408 0 0.0429 0-0.0408 0 0.0353Columns 46 through 540 -0.0276 0 0.0195 0 -0.0125 0 0.0072 0Columns 55 through 63-0.0037 0 0.0017 0 -0.0007 00.0003 0 -0.0001Columns 64 through 720 0.0000 0 -0.0000 0 0.0000 0 -0.0000 0Columns 73 through 810.0000 0 -0.0000 0 0.0000 0-0.0000 0 0.0000a =1.0e+005 *Columns 1 through 90.0000 -0.0001 0.0003 -0.0011 0.0030 -0.0074 0.0160 -0.0318 0.0585Columns 10 through 18-0.1008 0.1637 -0.2519 0.3692 -0.5174 0.6952 -0.8980 1.1176 -1.3427Columns 19 through 271.5597 -1.7542 1.9125 -2.0234 2.0794 -2.07732.0188 -1.9098 1.7596Columns 28 through 36-1.5798 1.3828 -1.1803 0.9828 -0.7985 0.6332 -0.4902 0.3705 -0.2734Columns 37 through 450.1970 -0.1386 0.0953 -0.0639 0.0419 -0.0268 0.0167 -0.0102 0.0061Columns 46 through 54-0.0035 0.0020 -0.0011 0.0006 -0.0003 0.0002 -0.0001 0.0000 -0.0000Columns 55 through 630.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000Columns 64 through 72-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000Columns 73 through 810.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000（4）Wp =[0.25 0.67];Ws =[0.25-0.03 0.67+0.03];Rp = 3;Rs = 60;[N,Wc]=buttord(Wp,Ws,Rp,Rs)[b,a]=butter(N,Wc,'stop')freqz(b,a);axis([0,1,-120,0]);grid ontitle('巴特沃斯带阻数字滤波器')N =40Wc =0.2499 0.6701b =1.0e+005 *Columns 1 through 70.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000Columns 8 through 14-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000Columns 15 through 210.0001 -0.0001 0.0003 -0.0007 0.0015 -0.0029 0.0056Columns 22 through 28-0.0102 0.0179 -0.0305 0.0502 -0.0798 0.1227 -0.1828Columns 29 through 350.2636 -0.3686 0.4998 -0.6576 0.8396 -1.04091.2534Columns 36 through 42-1.4660 1.6661 -1.8401 1.9753 -2.0610 2.0904 -2.0610Columns 43 through 491.9753 -1.8401 1.6661 -1.4660 1.2534 -1.0409 0.8396Columns 50 through 56-0.6576 0.4998 -0.3686 0.2636 -0.1828 0.1227 -0.0798Columns 57 through 630.0502 -0.0305 0.0179 -0.0102 0.0056 -0.0029 0.0015Columns 64 through 70-0.0007 0.0003 -0.0001 0.0001 -0.0000 0.0000 -0.0000Columns 71 through 770.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000Columns 78 through 81-0.0000 0.0000 -0.0000 0.0000a =1.0e+005 *Columns 1 through 70.0000 -0.0001 0.0003 -0.0011 0.0030 -0.0074 0.0160Columns 8 through 14-0.0318 0.0585 -0.1008 0.1637 -0.2519 0.3692 -0.5174Columns 15 through 210.6952 -0.8980 1.1176 -1.3427 1.5597 -1.75421.9125Columns 22 through 28-2.0234 2.0794 -2.0773 2.0188 -1.9098 1.7596 -1.5798Columns 29 through 351.3828 -1.1803 0.9828 -0.7985 0.6332 -0.4902 0.3705Columns 36 through 42-0.2734 0.1970 -0.1386 0.0953 -0.0639 0.0419 -0.0268Columns 43 through 490.0167 -0.0102 0.0061 -0.0035 0.0020 -0.0011 0.0006Columns 50 through 56-0.0003 0.0002 -0.0001 0.0000 -0.0000 0.0000 -0.0000Columns 57 through 630.0000 -0.0000 0.0000 -0.0000 0.0000 -0.00000.0000Columns 64 through 70-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000Columns 71 through 770.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000Columns 78 through 81-0.0000 0.0000 -0.0000 0.00003（1）T0=204;N=205;T=1;k=0:T0;x=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k);subplot(2,1,1);stem(k,x);title('时域波形 ');Xm=fft(x,N)/N;f=(-(N-1)/2:(N-1)/2)/N/T;subplot(2,1,2);stem(f,abs(fftshift(Xm)));title('频谱图');（2）[N,Wc]=buttord(0.1925,0.30225,3,60)[b,a]=butter(N,Wc)freqz(b,a);axis([0,1,-120,0]);grid ontitle('巴特沃斯低通数字滤波器')T0=204;N=205;T=1;k=0:T0;x=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k); subplot(4,1,1);stem(k,x);title('时域波形');Xm=fft(x,N)/N;f=(-(N-1)/2:(N-1)/2)/N/T;subplot(4,1,2);stem(f,abs(fftshift(Xm)));title('频谱图');y=filter(b,a,x);stem(k,y);title('低通滤波后时域波形')ym=fft(y,N)/N;subplot(4,1,4);stem(f,abs(fftshift(ym)));title('低通滤波后频谱图')[N,Wc]=buttord(0.1925,0.30225,3,60)[b,a]=butter(N,Wc,'high')freqz(b,a);axis([0,1,-120,0]);grid onT0=204;N=205;T=1;k=0:T0;x=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k); subplot(4,1,1);stem(k,x);title('时域波形');Xm=fft(x,N)/N;f=(-(N-1)/2:(N-1)/2)/N/T;stem(f,abs(fftshift(Xm)));title('频谱图');y=filter(b,a,x);subplot(4,1,3);stem(k,y);title('高通滤波后时域波形')ym=fft(y,N)/N;subplot(4,1,4);stem(f,abs(fftshift(ym)));title('高通滤波后频谱图')（3）Wp1 =[680 720]/4000;Ws1=[650-20 720+20]/4000;Rp1 = 3;Rs1 = 40;[N1,Wn1]=cheb1ord(Wp1,Ws1,Rp1,Rs1); [b1,a1] =cheby1(N1,Rp1,Wn1); freqz(b1,a1,512,8000);title('Ⅰ型切比雪夫滤波器1');grid onWp2 =[750 790]/4000;Ws2 =[750-20 790+20]/4000;Rp2 = 3;Rs2 = 40;[N2,Wn2]=cheb1ord(Wp2,Ws2,Rp2,Rs2) [b2,a2] =cheby1(N2,Rp2,Wn2); figure;freqz(b2,a2,512,8000);title('Ⅰ型切比雪夫滤波器2'); grid onWp3 =[830 870]/4000;Ws3 =[830-20 870+20]/4000;Rp3 = 3;Rs3 = 40;[N3,Wn3]=cheb1ord(Wp3,Ws3,Rp3,Rs3) [b3,a3] =cheby1(N3,Rp3,Wn3); figure;freqz(b3,a3,512,8000);title('Ⅰ型切比雪夫滤波器3'); grid onWp4 =[920 960]/4000;Ws4 =[920-20 960+20]/4000;Rp4 = 3;Rs4 = 40;[N4,Wn4]=cheb1ord(Wp4,Ws4,Rp4,Rs4);[b4,a4] =cheby1(N4,Rp4,Wn4);figure;freqz(b4,a4,512,8000);title('Ⅰ型切比雪夫滤波器4');grid on;k=0:1:500;x=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k); y1=filter(b1,a1,x);y2=filter(b2,a2,x);y3=filter(b3,a3,x);y4=filter(b4,a4,x);figure;plot(k,y1,k,y2,'g--',k,y3,'r--',k,y4,'y--');title('滤波后4条输出曲线') ;legend('697HZ', '770HZ','852HZ','941HZ');（4）Wp1 =[1180 1220]/4000;Ws1 =[1180-30 1220+30]/4000;Rp1 = 3;Rs1 = 40;[N1,Wn1]=cheb1ord(Wp1,Ws1,Rp1,Rs1) [b1,a1] =cheby1(N1,Rp1,Wn1);freqz(b1,a1,512,8000);title('Ⅰ型切比雪夫滤波器1'); grid on;Wp2 =[1310 1350]/4000;Ws2 =[1310-30 1350+30]/4000;Rp2 = 3;Rs2 = 40;[N2,Wn2]=cheb1ord(Wp2,Ws2,Rp2,Rs2); [b2,a2] =cheby1(N2,Rp2,Wn2); figure;freqz(b2,a2,512,8000);title('Ⅰ型切比雪夫滤波器2'); grid on;Wp3 =[1460 1500]/4000;Ws3=[1460-30 1500+30]/4000;Rp3 = 3;Rs3 = 40;[N3,Wn3]=cheb1ord(Wp3,Ws3,Rp3,Rs3) [b3,a3] =cheby1(N3,Rp3,Wn3); figure;freqz(b3,a3,512,8000);title('Ⅰ型切比雪夫滤波器3');grid on;k=0:1:500;x=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k);y1=filter(b1,a1,x);y2=filter(b2,a2,x);y3=filter(b3,a3,x);figure;plot(k,y1,k,y2,'g--',k,y3,'y--'); title('输出曲线') ; legend('1209HZ','1336HZ','1477HZ');(5)k=0:1:500;x0=sin((2/8000)*941*pi*k)+sin((2/8000)*1336*pi*k);x1=sin((2/8000)*697*pi*k)+sin((2/8000)*1209*pi*k);x2=sin((2/8000)*697*pi*k)+sin((2/8000)*1336*pi*k);x3=sin((2/8000)*697*pi*k)+sin((2/8000)*1477*pi*k);x4=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k);x5=sin((2/8000)*770*pi*k)+sin((2/8000)*1336*pi*k);x6=sin((2/8000)*770*pi*k)+sin((2/8000)*1477*pi*k);x7=sin((2/8000)*852*pi*k)+sin((2/8000)*1209*pi*k);x8=sin((2/8000)*852*pi*k)+sin((2/8000)*1336*pi*k);x9=sin((2/8000)*852*pi*k)+sin((2/8000)*1477*pi*k); Wp1 =[680 720]/4000;Ws1=[650-20 720+20]/4000;Rp1 = 3;Rs1 = 40;[N1,Wn1]=cheb1ord(Wp1,Ws1,Rp1,Rs1);[B1,A1] =cheby1(N1,Rp1,Wn1);Wp2 =[750 790]/4000;Ws2 =[750-20 790+20]/4000;Rp2 = 3;Rs2 = 40;[N2,Wn2]=cheb1ord(Wp2,Ws2,Rp2,Rs2);[B2,A2] =cheby1(N2,Rp2,Wn2);Wp3 =[830 870]/4000;Ws3 =[830-20 870+20]/4000;Rp3 = 3;Rs3 = 40;[N3,Wn3]=cheb1ord(Wp3,Ws3,Rp3,Rs3);[B3,A3] =cheby1(N3,Rp3,Wn3);Wp4 =[920 960]/4000;Ws4 =[920-20 960+20]/4000;Rp4 = 3;Rs4 = 40;[N4,Wn4]=cheb1ord(Wp4,Ws4,Rp4,Rs4); [B4,A4] =cheby1(N4,Rp4,Wn4);wp1 =[1180 1220]/4000;ws1 =[1180-30 1220+30]/4000;rp1 = 3;rs1 = 40;[n1,wn1]=cheb1ord(wp1,ws1,rp1,rs1); [b1,a1] =cheby1(n1,rp1,wn1);wp2 =[1310 1350]/4000;ws2 =[1310-30 1350+30]/4000;rp2 = 3;rs2 = 40;[n2,wn2]=cheb1ord(wp2,ws2,rp2,rs2); [b2,a2] =cheby1(n2,rp2,wn2);wp3 =[1460 1500]/4000;ws3=[1460-30 1500+30]/4000;rp3 = 3;rs3 = 40;[n3,wn3]=cheb1ord(wp3,ws3,rp3,rs3); [b3,a3]=cheby1(n3,rp3,wn3);Y01=filter(B1,A1,x0);Y03=filter(B3,A3,x0);Y04=filter(B4,A4,x0);figure ; subplot(2,1,1);plot(k,Y01,k,Y02,'y--',k,Y03,'r--',k,Y04,'g--'); title('输出曲线1');legend('697HZ', '770HZ','852HZ','941HZ'); subplot(2,1,2);y01=filter(b1,a1,x0);y02=filter(b2,a2,x0);y03=filter(b3,a3,x0);plot(k,y01,k,y02,'g--',k,y03,'r--');legend('1209HZ', '1336HZ','1477HZ');Y11=filter(B1,A1,x1);Y12=filter(B2,A2,x1);Y13=filter(B3,A3,x1);Y14=filter(B4,A4,x1);figure;subplot(2,1,1);plot(k,Y11,k,Y12,'y--',k,Y13,'g--',k,Y14,'r--'); title('输出曲线2') ;legend('697HZ', '770HZ','852HZ','941HZ');subplot(2,1,2);y12=filter(b2,a2,x1);y13=filter(b3,a3,x1);plot(k,y11,k,y12,'y--',k,y13,'g--');legend('1209HZ', '1336HZ','1477HZ');(5)T0=500;N=501;T=1;k=0:T:T0;x0=sin((2/8000)*941*pi*k)+sin((2/8000)*1336*pi*k); x1=sin((2/8000)*697*pi*k)+sin((2/8000)*1209*pi*k);x2=sin((2/8000)*697*pi*k)+sin((2/8000)*1336*pi*k);x3=sin((2/8000)*697*pi*k)+sin((2/8000)*1477*pi*k);x4=sin((2/8000)*770*pi*k)+sin((2/8000)*1209*pi*k);x5=sin((2/8000)*770*pi*k)+sin((2/8000)*1336*pi*k);x6=sin((2/8000)*770*pi*k)+sin((2/8000)*1477*pi*k);x7=sin((2/8000)*852*pi*k)+sin((2/8000)*1209*pi*k);x8=sin((2/8000)*852*pi*k)+sin((2/8000)*1336*pi*k);x9=sin((2/8000)*852*pi*k)+sin((2/8000)*1477*pi*k);Wp1 =[680 720]/4000;Ws1=[650-20 720+20]/4000;Rp1 = 3;Rs1 = 40;[N1,Wn1]=cheb1ord(Wp1,Ws1,Rp1,Rs1); [B1,A1] =cheby1(N1,Rp1,Wn1);Wp2 =[750 790]/4000;Ws2 =[750-20 790+20]/4000;Rp2 = 3;Rs2 = 40;[N2,Wn2]=cheb1ord(Wp2,Ws2,Rp2,Rs2); [B2,A2] =cheby1(N2,Rp2,Wn2);Wp3 =[830 870]/4000;Ws3 =[830-20 870+20]/4000;Rp3 = 3;Rs3 = 40;[N3,Wn3]=cheb1ord(Wp3,Ws3,Rp3,Rs3); [B3,A3] =cheby1(N3,Rp3,Wn3);Wp4 =[920 960]/4000;Ws4 =[920-20 960+20]/4000;Rp4 = 3;Rs4 = 40;[N4,Wn4]=cheb1ord(Wp4,Ws4,Rp4,Rs4); [B4,A4] =cheby1(N4,Rp4,Wn4);wp1 =[1180 1220]/4000;ws1 =[1180-30 1220+30]/4000;rp1 = 3;rs1 = 40;[n1,wn1]=cheb1ord(wp1,ws1,rp1,rs1); [b1,a1] =cheby1(n1,rp1,wn1);wp2 =[1310 1350]/4000;ws2 =[1310-30 1350+30]/4000;rp2 = 3;rs2 = 40;[n2,wn2]=cheb1ord(wp2,ws2,rp2,rs2); [b2,a2] =cheby1(n2,rp2,wn2);wp3 =[1460 1500]/4000;ws3=[1460-30 1500+30]/4000;rp3 = 3;rs3 = 40;[n3,wn3]=cheb1ord(wp3,ws3,rp3,rs3); [b3,a3]=cheby1(n3,rp3,wn3);Y01=filter(B1,A1,x0);Y02=filter(B2,A2,x0);Y03=filter(B3,A3,x0);Y04=filter(B4,A4,x0);Ym01=fft(Y01,N)/N;Ym02=fft(Y02,N)/N;Ym04=fft(Y04,N)/N;figure;plot(k,abs(fftshift(Ym01)),k,abs(fftshift(Ym02)),'y--',k,abs(ff tshift(Ym03)),'g--',k,abs(fftshift(Ym04)),'r--');legend('697HZ', '770HZ','852HZ','941HZ'); title('第一组中"0"滤波后的频谱');y01=filter(b1,a1,x0);y02=filter(b2,a2,x0);y03=filter(b3,a3,x0);ym01=fft(y01,N)/N;ym02=fft(y02,N)/N;ym03=fft(y03,N)/N;figure;plot(k,abs(fftshift(ym01)),k,abs(fftshift(ym02)),'g--',k,abs(ff tshift(ym03)),'g--');legend('1209HZ', '1336HZ','1477HZ');title('第二组中"0"滤波后的频谱');Y11=filter(B1,A1,x1);Y12=filter(B2,A2,x1);Y13=filter(B3,A3,x1);Y14=filter(B4,A4,x1);Ym12=fft(Y12,N)/N;Ym13=fft(Y13,N)/N;Ym14=fft(Y14,N)/N;figure;plot(k,abs(fftshift(Ym11)),k,abs(fftshift(Ym12)),'y--',k,abs(ff tshift(Ym13)),'g--',k,abs(fftshift(Ym14)),'r--');legend('697HZ', '770HZ','852HZ','941HZ');title('第一组中"1"滤波后的频谱');y11=filter(b1,a1,x1);y12=filter(b2,a2,x1);y13=filter(b3,a3,x1);ym11=fft(y11,N)/N;ym12=fft(y12,N)/N;ym13=fft(y13,N)/N;figure;plot(k,abs(fftshift(ym11)),k,abs(fftshift(ym12)),'g--',k,abs(ff tshift(ym13)),'g--');legend('1209HZ', '1336HZ','1477HZ');title('第二组中"1"滤波后的频谱');四：思考题1. 哪些主要因素直接影响IIR数字滤波器的阶数？从工程概念进行定性解释。