《数字信号处理》实验指导书(版)

数字信号处理实验指导书电子信息工程学院2015年4月目录TOC\o"1-1"\h\z\u实验一离散信号产生和基本运算 3实验二基于MATLAB的离散系统时域分析 7实验三基于ICETEK-F2812-A教学系统软件的离散系统时域分析 9实验四基于MATLAB的FFT算法的应用 16实验五基于ICETEK-F2812-A的FFT算法分析 18实验六基于ICETEK-F2812-A的数字滤波器设计 20实验七基于ICETEK-F2812-A的交通灯综合控制 24实验八基于BWDSP100的步进电机控制 26实验一离散信号产生和基本运算实验目的1．掌握MATLAB最基本的矩阵运算语句。2．掌握对常用离散信号的理解与运算实现。二、实验原理向量的生成（1）利用冒号“：”运算生成向量，其语句格式有两种：A=m:nB=m:p:n第一种格式用于生成不长为1的均匀等分向量，m和n分别代表向量的起始值和终止值，n>m。第二种格式用于生成步长为p的均匀等分的向量。（2）利用函数linspace()生成向量，linspace()的调用格式为：A=linspace(m,n)B=linspace(m,n,s)第一种格式生成从起始值m开始到终止值n之间的线性等分的100元素的行向量。第二种格式生成从起始值m开始到终止值n之间的s个线性等分点的行向量。矩阵的算术运算（1）加法和减法对于同维矩阵指令的A+BA-B对于矩阵和标量（一个数）的加减运算，指令为：A+3A-9（2）乘法和除法运算 A*B是数学中的矩阵乘法，遵循矩阵乘法规则 A．*B是同维矩阵对应位置元素做乘法 B=inv（A）是求矩阵的逆 A/B是数学中的矩阵除法，遵循矩阵除法规则 A．/B是同维矩阵对应位置元素相除 另A’表示矩阵的转置运算数组函数下面列举一些基本函数，他们的用法和格式都相同。sin(A),cos(A),exp(A),log(A)(相当于ln)sqrt(A)开平方abs(A)求模real(A)求实部imag(A)求虚部式中A可以是标量也可以是矩阵例：利用等差向量产生一个正弦值向量t=0:0.1:10A=sin(t)plot(A)这时候即可看到一个绘有正弦曲线的窗口弹出另：每条语句后面加“；”表示不要显示当前语句的执行结果不加“；”表示要显示当前语句的执行结果。二维曲线的绘制plot()函数plot()函数是将各个数据点通过连折线的方式来绘制二维图形的，其命令格式有以下几种：plot(y)当y为向量时，以y的序号作为X轴，按向量y的值绘制曲线。plot(x,y)x,y均为向量时，以x向量作为X轴，向量y作为Y轴绘制曲线。（3）plot(x,y1,’option1’,x,y2,’option2’,……)以公共的x向量作为X轴，分别以向量y1，y2……的数据绘制多条曲线，每条曲线的属性由相应的‘option’来确定。（4）plot(x1,y1,’option1’,x2,y2,’option2’,……)分别以向量x1，x2，……作为X轴，以y1，y2，……的数据绘制多条曲线，每条曲线的属性由相应的选项‘option表1plot函数中option选项的取值和含义符号颜色符号线型符号标记符号标记‘b’

蓝色‘-’实线‘.’点‘d’

‘g’绿色‘:’虚线‘o’圆圈‘ˇ’▽‘m’品红-.点画线‘x’叉号‘^’△‘c’青色‘--’双画线‘+’加号‘<’左三角‘k’黑色‘none’无线‘*’星号‘>’右三角‘r’红色‘s’□‘p’五角星‘y’黄色‘h’六角星5.常用离散信号的MATLAB实现函数（1）.单位抽样序列在MATLAB中利用zeros()函数实现如果在时间轴上延迟了k个单位，得到即：（2）．单位阶跃序列在MATLAB中可以利用ones()函数实现。（3）.正弦序列，在MATLAB（4）．复指数序列，在MATLAB中，（5）．指数序列，在MATLAB中，三、实验内容（1）熟悉MATLAB的使用环境和方法；（2）练习使用基本的向量生成、矩阵运算、绘图等语句；利用冒号（：）生成向量X1=[12345]X2=[1.0001.5002.0002.500]X3=[54321]分别生成3*3，3*4的全0矩阵，全1矩阵和随机矩阵分别输入矩阵分别计算A+B，A-B，A+3，A-4，A*B，A.*B，C=inv（A）,A/B,A./B分别计算sin(x1),cos(x1),exp(x1),log(x2),sqrt(x2)（3）生成以上五种基本离散信号函数；（4）绘出信号，当、时、、、时的信号实部和虚部图；（5）绘出信号的频率是多少？周期是多少？产生一个数字频率为0.9的正弦序列，并显示该信号，说明其周期?四、实验预习要求1.预习实验原理。2.熟悉实验程序。3.思考课程设计实验部分程序的编写。五、实验报告要求1.简述实验目的与步骤2.写出各步输入语句和输出结果六、思考题1.矩阵的乘法和除法有几种，它们有什么区别？2.离散正弦序列的性质。实验二基于MATLAB的离散系统时域分析一、实验目的1.学习MATLAB语言的编程和调试技巧；2.掌握离散卷积方法和MATLAB语言实现。二、实验原理时域中，离散时间系统对输入信号或延迟信号进行运算处理，生成具有所需特性的输出信号。本实验通过MATLAB仿真一些简单的离散时间信号和系统，并研究其时域特性。涉及到离散时间信号、离散时间系统、系统性质及线性卷积等知识点。一个离散时间系统，输入信号为x(n)，输出信号为y(n)，运算关系用T[﹒]表示，则输入与输出的关系可表示为y(n)＝T[x(n)]。（1）线性时不变系统的输入输出关系可通过单位脉冲响应h(n)表示：式中*表示卷积运算。（2）线性时不变系统的实现：可物理实现的线性时不变系统是稳定的、因果的。这种系统的单位脉冲响应是因果的（单边）且绝对可和的，即：:，；。在MATLAB语言中采用conv实现卷积运算即：y=conv(x,h)，它默认从n=0开始。三、实验内容1.编制程序求解下列两个系统的单位冲激响应和阶跃响应，并绘出其图形。要求分别用filter、conv、impz三种函数完成。给出理论计算结果和程序计算结果并讨论。四、实验预习要求1.预习实验原理。2.熟悉实验程序。3.思考程序设计实验部分程序的编写。五、实验报告要求1.报告中要给出实验的MATLAB程序，并对每个语句给出注释，说明语句作用；2.简述实验目的和原理；3.给出用笔算时卷积和conv计算线性卷积对照图；4.给出收获和体会。六、思考题如何用matlab实现卷积函数？实验三离散系统的变换域分析一、实验目的：加深对离散系统的频率响应分析和零、极点分布的概念理解。二、实验原理：离散系统的时域方程为其变换域分析方法如下：频域系统的频率响应为Z域系统的转移函数为分解因式，其中和称为零、极点。在MATLAB中，可以用函数[z，p，K]=tf2zp（num，den）求得有理分式形式的系统转移函数的零、极点，用函数zplane（z，p）绘出零、极点分布图；也可以用函数zplane（num，den）直接绘出有理分式形式的系统转移函数的零、极点分布图。使h=freqz(num,den,w)函数可求系统的频率响应，w是频率的计算点，如w=0:pi/255:pi,h是复数，abs(h)为幅度响应，angle(h)为相位响应。另外，在MATLAB中，可以用函数[r，p，k]=residuez（num，den）完成部分分式展开计算；可以用函数sos=zp2sos（z，p，K）完成将高阶系统分解为2阶系统的串联。例1求下列直接型系统函数的零、极点，并将它转换成二阶节形式 解用MATLAB计算程序如下： num=[1-0.1-0.3-0.3-0.2]; den=[10.10.20.20.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.6997i0.5276-0.6997i-0.5776+0.5635i-0.5776-0.5635i

增益系数1

二阶节1.0000-0.3885-0.55091.00001.15520.65111.00000.28850.36301.0000-1.05520.7679

系统函数的二阶节形式为： 极点图如右图。

例2差分方程 所对应的系统的频率响应。解：差分方程所对应的系统函数为用MATLAB计算的程序如下： k=256; num=[0.8-0.440.360.02]; den=[10.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('弧度')三、实验内容：1、求系统的零、极点和幅度频率响应和相位响应。2、对于传输函数，求脉冲响应和阶跃相应3、滤波器的传输函数为：，不计算，画出幅度响应曲线四、实验预习要求1.预习实验原理。2.熟悉实验程序。3.思考课程设计实验部分程序的编写。五、实验报告要求1.简述实验目的与原理。2.写出实验内容的程序和结果。实验四基于MATLAB的FFT算法的应用一、实验目的1.进一步加深DFT算法原理和基本性质的理解(因为FFT只是DFT的一种快速算法，所以FPT的运算结果必然满足DFT的基本性质。2.熟悉FFT算法原理和FFT子程序的应用。3.学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。二、实验原理N点序列的DFT和IDFT变换定义式如下： ，利用旋转因子具有周期性，可以得到快速算法（FFT）。在MATLAB中，可以用函数X=fft（x，N）和x=ifft（X，N）计算N点序列的DFT正、反变换。三、实验内容与方法1．对连续的单一频率周期信号按采样频率采样，截取长度N分别选N=20和N=16，观察其DFT结果的幅度谱。2．2N点实数序列N=64。用一个64点的复数FFT程序，一次算出，并绘出。3．已知某序列在单位圆上的N=64等分样点的Z变换为用N点IFFT程序计算，绘出和。四、实验预习要求1.预习实验原理。2.熟悉实验程序。3.思考课程设计实验部分程序的编写。五、实验报告要求1.简述实验原理及目的。2.结合实验中所得给定典型序列幅频特性曲线，与理论结果比较，并分析说明误差产生的原因以及用FFT作谱分析时有关参数的选择方法。3.总结实验所得主要结论。4.简要回答思考题。六、思考题如果周期信号的周期预先不知道，如何使用FFT进行谱分析？实验五线性卷积与循环卷积的计算一、实验目的1、进一步加深对线性卷积的理解和分析能力；2、通过编程，上机调试程序，进一步增强使用计算机解决问题的能力；3、掌握线性卷积与循环卷积软件实现的方法，并验证二者之间的关系。二、实验原理1、线性卷积线性时不变系统（LinearTime-InvariantSystem,orL.T.I系统）输入、输出间的关系为：当系统输入序列为，系统的单位脉冲响应为，输出序列为，则系统输出为：或上式称为离散卷积或线性卷积。图6.1示出线性时不变系统的输入、输出关系。L.T.Ih(n)

L.T.Ih(n)→L.T.I—→—→—→

图6.1线性时不变系统的输入、输出关系2、循环卷积DFT设两个有限长序列和，均为点长DFTDFTDFT如果则N上式称为循环卷积或圆周卷积注：为序列的周期化序列；为的主值序列。上机编程计算时，可表示如下：3、两个有限长序列的线性卷积序列为点长，序列为点长，为这两个序列的线性卷积，则为且线性卷积的最大长，也就是说当和时。4、循环卷积与线性卷积的关系序列为点长，序列为点长，若序列和进行N点的循环卷积，其结果是否等于该两序列的线性卷积，完全取决于循环卷积的长度：当时循环（圆周）卷积等于线性卷积，即N当时，循环卷积等于两个序列的线性卷积加上相当于下式的时间混叠，即三、实验方法对于无限长序列不能用MATLAB直接计算线性卷积，在MATLAB内部只提供了一个conv函数计算两个有限长序列的线性卷积。对于循环卷积MATLAB内部没有提供现成的函数，我们可以按照定义式直接编程计算。例6.1：已知两序列：求它们的线性卷积yl(n)=h(n)*x(n)和N点的循环卷积yc=[h(n)*x(n)]n，并研究两者之间的关系.MATLAB实现程序：循环卷积的函数functionyc=circonv(x1,x2,N)%realizecircularconvolutionusedirectmethod%y=circonv(x1,x2,N)%y:outputsequences%x1,x2:inputsequences%N:circulationlengthiflength(x1)>N error(‘Nmustnotbelessthanlengthofx1’);endiflength(x2)>N error(‘Nmustnotbelessthanlengthofx2’);end%以上语句判断两个序列的长度是否小于Nx1=[x1,zeros(1,N-length(x1))]; %填充序列x1(n)使其长度为N1+N2-1(序列%h(n)的长度为N1，序列x(n)的长度为N2)x2=[x2,zeros(1,N-length(x2))]; %填充序列x2(n)使其长度为N1+N2-1n=[0:1:N-1];x2=x2(mod(-n,N)+1); %生成序列x2((-n))NH=zeros(N,N);forn=1:1:N H(n,:)=cirshiftd(x2,n-1,N); %该矩阵的k行为x2((k-1-n))Nendyc=x1*H’; %计算循环卷functiony=cirshiftd(x,m,N)%directlyrealizecircularshiftforsequencex%y=cirshiftd(x,m,N);%x:inputsequencewhoselengthislessthanN%m:howmuchtoshift%N:circularlength%y:outputshiftedsequenceiflength(x)>Nerror('lengthofxmustbelessthanN');endx=[x,zeros(1,N-length(x))];n=[0:1:N-1];y=x(mod(n-m,N)+1);

研究两者之间的关系clearall;n=[0:1:11];m=[0:1:5];N1=length(n);N2=length(m);xn=0.8.^n; %生成x(n)hn=ones(1,N2); %生成h(n)yln=conv(xn,hn); %直接用函数conv计算线性卷积ycn=circonv(xn,hn,N1); %用函数circonv计算N1点循环卷积ny1=[0:1:length(yln)-1];ny2=[0:1:length(ycn)-1];subplot(2,1,1); %画图stem(ny1,yln);subplot(2,1,2);stem(ny2,ycn);axis([0,16,0,4]);

运行结果：图4.2线性卷积和循环卷积的比较四、实验内容及要求

已知两个有限长序列实验前，预先笔算好这两个序列的线性卷积及下列几种情况的循环卷积eq\o\ac(○,*)eq\o\ac(○,*)eq\o\ac(○,*)eq\o\ac(○,*)编制一个计算两个序列线性卷积的通用程序，计算。编制一个计算循环卷积的通用程序，计算上述4种情况下两个序列与的循环卷积。上机调试并打印或记录实验结果。注：可在一个程序中用菜单形式实现上述两种卷积的计算。将实验结果与预先笔算的结果比较，验证其正确性。五、实验总结

思考题若系统脉冲响应长为，系统输入长为，且，应使用什么方法进行编程上机计算？

实验报告要求列出计算两种卷积的公式，画出程序框图，并列出实验程序清单（可略）（包括必要的程序说明）。记录调试运行情况及所遇问题的解决方法。给出实验结果，并对结果作出分析。验证循环卷积两者之间的关系。

简要回答思考题。实验六基于matlab的模拟Butterworth低通滤波器设计实验目的熟悉模拟滤波器设计的基本过程；熟悉用matlab实现模拟Butterworth低通滤波器设计原理、方法。实验原理模拟低通滤波器设计的基本过程和原理如下。给定技术指标模拟低通滤波器的技术指标有四个，通带截止频率OmegaP(弧度/秒)，阻带截止频率OmegaS(弧度/秒)，通带衰减Rp（dB），阻带衰减As（dB）。2．选定滤波器的类型模拟低通滤波器的类型有巴特沃思（Butterworth）滤波器、切比雪夫（Chebyshev）滤波器、椭圆（Ellipse）滤波器、贝塞尔（Bessel）滤波器等。本实验中选用巴特沃思（Butterworth）滤波器。3．求滤波器阶次N及设计用的通带截止频率OmegaC.N代表滤波器阶数，ΩC为巴特沃思低通滤波器频率响应幅度衰减到3dB的通带截止频率。可调用matlab函数实现。[N，OmegaC]=buttord(OmegaP,OmegaS,Rp,As,’s’)求归一化原型模拟低通Butterworth滤波器[z0,p0,k0]=buttap(N)；b0=k0*real(poly(z0));a0=real(poly(p0));Butterworth滤波器是全极点型的，没有零点。即归一化原型系统函数Han(s)可表示为求非归一化（即去归一化）的低通Butterworth滤波器。Butterworth非归一化低通滤波器，其极点由位于单位圆（归一化的）转换成半径为ΩC的圆上，故应该将极点乘以ΩC，而将增益乘以ΩCN，利用Ha(s)=Han(s/ΩC)可得上述结论，Han(s)为归一化的，Ha(s)为非归一化的。[z0,p0,k0]=buttap(N)；p=p0*OmegaC;a=real(poly(p));k=k0*OmegaC^N;b0=real(poly(z0));b=k*b0;非归一化系统函数可写成上面步骤4和5，可以合并用下列步骤6的matlab来实现。模拟低通Butterworth滤波器的设计的matlab函数[b,a]=butter(N,Omegac,’s’)其中b,a分别为设计出的滤波器的分子多项式与分母多项式的系数向量。所设计出的模拟低通滤波器的系统函数的为N阶的形式如下。求模拟低通Butterworth滤波器的零极点及增益的matlab函数表示[z,p,k]=butter(N,OmegaC,’s’)【例题1】设计一个模拟巴特沃思低通滤波器，其通带截止频率fp=3000HZ,通带最大衰减RP=2dB;阻带截止频率为fs=6000Hz；阻带最小衰减As=30dB。解（1）程序如下clc;clearOmegaP=2*pi*3000;OmegaS=2*pi*6000;%Butterworth模拟低通的各项指标Rp=2;As=30;[N,OmegaC]=buttord(OmegaP,OmegaS,Rp,As,'s');%Butterworth模拟低通的阶次[b,a]=butter(N,OmegaC,'s');%设计Butterworth模拟低通w0=[OmegaP,OmegaS];%以下是检验衰减指标[H,w]=freqs(b,a);%求模拟滤波器的频率响应自动取200个频率点Hx=freqs(b,a,w0);%求OmegaP和OmegaS的两个频率点的频率响应值。dbHx=-20*log10(abs(Hx)/max(abs(H)));%用分贝表示两频率点的频率响应值plot(w/(2*pi)/1000,20*log10(abs(H)));xlabel('f(kHz)');ylabel('dB');axis([-1,12,-55,1]);set(gca,'xtickmode','manual','xtick',[0,1,2,3,4,5,6,7,8,9,10]);set(gca,'xtickmode','manual','ytick',[-50,-40,-30,-20,-10,0]);grid;（2）实验结果N=6;OmegaC=21202,dbHx=0.94784 30;b=0 0 0 0 0 0 9.0825e+025a=1 81916 3.3552e+009 8.7121e+013 1.5082e+018 1.6552e+022 9.0825e+025系统函数可以写成三、实验内容1.设计一个模拟巴特沃思低通滤波器，其通带截止频率fp=6000HZ,通带最大衰减RP=3.5dB;阻带截止频率为fs=12000Hz；阻带最小衰减As=35dB。（1）写出matlab设计程序并运行；（2）绘出巴特沃思滤波器的幅度特性图；四、实验报告要求1.写出matlab的设计程序；2.绘出巴特沃思滤波器的幅度特性图；3.写出设计出的巴特沃思滤波器的系统函数表达式。实验八双线性变换法设计Butterworth低通数字滤波器实验目的1.熟悉利用模拟滤波器设计数字滤波器的基本过程；2.熟悉双线性变换法设计数字滤波器的基本原理；3.熟悉基于matlab的实现双线性变换法设计Butterworth低通数字滤波器。实验原理1.利用模拟滤波器设计数字滤波器利用模拟滤波器设计数字滤波器是数字滤波器设计的间接方法。其基本原理是先设计出模拟滤波器，然后通过频带变换和数字化的方法把模拟低通滤波器变换成数字滤波器。主要方法有冲激响应不变法，阶跃响应不变法和双线性变换法。其中双线性变换法由于没有频率响应的混叠失真，可适用于低通、带通、高通和带阻各种滤波器设计，比其它两种方法适用范围更广。双线性变换法（1）基本思路模拟滤波器数字化成数字滤波器的方法，就是要把s平面映射到z平面，使模拟系统函数Ha(s)变换成所需数字滤波器的系统函数H(z)。双线性变换是使数字滤波器的频率响应与模拟滤波器的频率响应相似的一种变换，它使得Ω和ω之间是单值映射关系可以避免频率响应的混叠失真。（2）变换关系使用双线性变换法，模拟滤波器的s域系统函数与数字滤波器的z域系统函数存在如下变换关系（1）T为抽样间隔，fs=1/T。模拟滤波器的频率响应与数字滤波器的频率响应存在如下变换关系（2）模拟滤波器的频率Ω到数字滤波器的频率ω是非线性频率变换关系（3）双线性变换法主要用于分段常数的频率响应滤波器中，例如低通、高通、带通、带阻等，这是大多数滤波器都具有的幅度特性。分段常数的模拟滤波器AF，经变换后仍然为分段常数的数字滤波器DF，但是各分段的临界频率点由于非线性频率变换而产生变化，这种变化可以用频率“预畸”来加以克服。频率预畸即若给定数字滤波器的截止频率为ωi，则根据上面公式（3）将它预畸为，以此Ωi来设计“样本”AF，将设计好的“样本”AF经双线性变换后，就得到所需的DF，它的截止频率正是原先要求的ωi，预畸是双线性变换法设计必须要的一步。双线性变换法实现模拟到数字的滤波器的matlab实现[bz,az]=bilinear(b,a,Fs)式中，b、a分别为模拟滤波器系统函数的分子、分母多项式系数；Fs为抽样频率；bz、az分别为数字滤波器系统函数的分子分母多项式系数。设计时必须有模拟原型滤波器的频率预畸。模拟低通Butterworth滤波器的设计的matlab函数[b,a]=butter(N,Omegac,’s’)其中b,a分别为设计出的滤波器的分子多项式与分母多项式的系数向量。所设计出的模拟低通滤波器的系统函数的为N阶的形式如下。Butterworth滤波器是全极点型的，没有零点。求模拟低通Butterworth滤波器的零极点及增益的matlab函数表示[z,p,k]=butter(N,OmegaC,’s’)【例题1】试用双线性变换法设计一个巴特沃思低通数字滤波器，低通滤波器的技术指标用模拟频率给出，其通带截止频率fp=1000HZ；通带最大衰减