数字信号处理实验5

第一篇：数字信号处理实验5

实验五 FFT 算法的应用

1.进一步加深对离散信号 DFT 的理解。

2.运用其 FFT 算法解决实际问题。

1.微机。

2.Matlab 编程环境。

1.熟悉 Matlab 的编程环境和编程语言。

2.学习教材 P213-227，P249-263，掌握快速傅里叶变换（FFT）的原理。

1.实验重点、难点、特点

快速傅里叶变换（FFT）的原理及应用。难点在 FFT 的应用以及 Matlab 编程中矩阵乘和数乘的应用。

2.实验原理

一、实验目的

二、实验仪器设备

三、实验学时 学时

四、预习要求

五、实验特点及实验原理简介

1.已知 2N 点实数序列

N=64。理论计算 X(k)=DFT[x(n)]2N。

六、实验内容及步骤

N=64;n=0:1:N-1;n2=0:1:2*N-1;k=0:1:2*N-1;xn1=cos(2*pi/N*7*2*n)+1/2*cos(2*pi/N*19*2*n);xn2=cos(2*pi/N*7*(2*n+1))+1/2*cos(2*pi/N*19*(2*n+1));wn=xn1+i*xn2;[wk,kk]=dft_my(wn,N);xk1=1/2*(wk+conj(wk(mod(-kk,N)+1)));xk2=1/2*(wk-conj(wk(mod(-kk,N)+1)))/i;XK1=fft(xn1,N);XK2=fft(xn2,N);%subplot(2,2,1);%stem(kk,abs(xk1));grid on;%subplot(2,2,3);%stem(kk,abs(XK1));grid on;%subplot(2,2,2);%stem(kk,abs(xk2));grid on;%subplot(2,2,4);%stem(kk,abs(XK2));grid on;

xk=[xk1+(exp(-1i*pi/N).^kk).*xk2,xk1-(exp(-1i*pi/N).^kk).*xk2];subplot(2,1,1);stem(k,abs(xk),'r');grid on;xlabel('k');ylabel('2*NFFT');title('n点DFT完成2N点FFT');xn=cos(2*pi/N*7*n2)+1/2*cos(2*pi/N*19*n2);xkk=fft(xn,2*N);subplot(2,1,2);stem(k,abs(xkk),'g');grid on;xlabel('k');ylabel('2*NFFT');title('Matlab2N点FFT');

七、问题思考

1.两次编程计算与理论计算相比较，结果一致吗？说明原因。

2.两次编程计算在编程方面、运算量上的比较。

八、心得体会

第二篇：数字信号处理实验讲稿

邯 郸 学 院

讲 稿

2010 ～2011 学年 第 一 学期

分院(系、部)： 信息工程学院 教 研 室： 电子信息工程 课 程 名 称： 数字信号处理

授 课 班 级： 07级电子信息工程

主 讲 教 师： 王苗苗 职

称：

助教（研究生）

使 用 教 材： 《数字信号处理》

制 作 系 统：

Word2003

邯郸学院制

实验一..Matlab仿真软件介绍

一、实验目的

熟悉Matlab仿真软件

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、学习Matlab仿真软件的安装

2、熟悉Matlab仿真软件的操作环境

3、直接在Matlab仿真软件的命令窗口实现数值计算

4、编写M文件

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

1、熟悉Matlab仿真软件

2、参阅Matlab及在电子信息类课程中的应用（第2版）唐向宏 电子工业出版社

实验二 离散信号和系统分析的Matlab实现

一、实验目的

1、Matlab实现离散信号和系统分析

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab产生离散信号

2、利用Matlab计算离散卷积

3、利用Matlab求解离散LTI系统响应

4、利用Matlab计算DTFT

5、利用Matlab实现部分分式法

6、利用Matlab计算系统的零极点

7、利用Matlab进行简单数字滤波器设计

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验三 利用Matlab实现信号DFT的计算

一、实验目的

1、Matlab实现信号DFT的计算

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab计算信号的DFT

2、利用Matlab实现由DFT计算线性卷积

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验四 利用Matlab实现滤波器设计

一、实验目的

1、Matlab实现实现滤波器设计

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab实现模拟低通滤波器的设计

2、利用Matlab实现模拟域频率变换

3、利用Matlab实现脉冲响应不变法

4、利用Matlab实现双线性变换法

5、利用Matlab实现数字滤波器设计

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验五 利用Matlab实现FIR滤波器设计

一、实验目的

1、Matlab实现实现滤波器设计

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab实现窗函数法

2、利用Matlab实现频率取样法

3、利用Matlab实现优化设计

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验六..随机信号功率谱估计的Matlab实现

一、实验目的

1、Matlab实现实现滤波器设计

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab实现随机序列

2、利用Matlab计算相关函数的估计

3、利用Matlab进行非参数功率谱估计

4、利用Matlab进行AR模型功率谱估计

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验七..数字滤波器结构的Matlab实现

一、实验目的

1、Matlab实现实现滤波器设计

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、利用Matlab实现数字滤波器直接型设计

2、利用Matlab实现数字滤波器级联设计

3、利用Matlab实现数字滤波器并联型设计

4、利用Matlab实现数字滤波器格型设计

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

实验八....利用Matlab实现信号小波分析

一、实验目的

1、Matlab实现实现滤波器设计

2、进一步熟悉Matlab软件操作

二、实验设备和元器件

含Matlab仿真软件的计算机

三、实验内容和步骤

1、小波测试信号

2、分解与重构滤波器组

3、离散小波变换

4、离散小波反变换

5、基于小波的信号去噪

6、基于小波的信号压缩

四、实验报告要求

按照《Matlab程序设计》模板提交实验报告

五、预习要求

预习课本上的相关内容

第三篇：数字信号处理实验二

北京信息科技大学

数字信号处理实验报告

题 目：

数字信号处理课程设计实验

学 院： 信息与通信工程学院 专 业： 通信工程专业 姓 名： 班 级：

学 号： 指导老师：

实验目的

1、熟悉IIR数字滤波器的设计原理与方法。

2、掌握数字滤波器的计算机软件实现方法。

3、通过观察对实际心电图信号的滤波作用，学习数字滤波器在实际中的应用。

实验仪器及材料

计算机，MATLAB软件

实验内容及要求

1．设计巴特沃斯低通数字滤波器对人体心电信号进行滤波

(1)人体心电图信号在测量过程中会受到工业高频干扰，所以必须经过低通滤波处理，才能作为判断心脏功能的有用信息。以下为一个实际心电图信号采样序列x(n)，其中存在高频干扰，抽样周期Ts=1秒。在实验中，以x(n)作为输入序列，滤除其中干扰成分。

x(n)=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0] 对序列x(n)用FFT做频谱分析，生成x(n)的频谱图。

(2)设计一个巴特沃斯低通IIR数字滤波器H(z)。

设计指标参数为：在通带内频率低于0.2π时，最大衰减小于1dB； 在阻带内 [0.3π, π]频率区间上，最小衰减大于15dB。

j|H(e)|。写出数字滤波器H(z)的表达式，画出滤波器的幅频响应曲线

(3)用所设计的滤波器对实际心电图信号采样序列x(n)进行滤波处理，编写程序，求滤波后的序列y(n)，并分别画出滤波前后的心电图信号波形图和频谱图。

y(n)= [0,0,0,0, 0,0,0,0,-0.14025,0.40279,-0.56085 ,0.33328,0.023981,-0.18809,0.11843,-0.1038,0.11576,-0.1225,0.099815 ,-0.13769 ,0.095249,-0.0070273,0.018867,0.090543,-0.11257,-0.070884 ,0.17676,-0.55407,0.24813,-0.34732,-0.30428,0.59426,-0.29574,-0.063869,0.34018,-0.73334,1.0293,-0.57107,-0.2461,0.83605,-0.83026,0.45459,0.011551,-0.25667,0.23896,-0.17361,0.20829,-0.28417,0.28765 ,-0.2035,0.02865,0.066164,0.077916,-0.36052, 0.53517,-0.5571]

源程序

clear all，clc

x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];%未经滤波的心电图信号 L=length(x);l=0:L-1;

y=fft(x,L);Wp=0.2*pi;Ws=0.3*pi;Rp=1;Rs=15;

[N,Wn] = buttord(Wp,Ws,Rp,Rs,'s');[b,a] = butter(N,Wn,'s');[numa,dena]=impinvar(b,a,1);w=linspace(0,pi,1024);h=freqz(numa,dena,w);norm=max(abs(h));numa=numa/norm;[z,p]=tf2zp(b,a);figure(1)

plot(w,20*log10(abs(h)/norm));grid;

xlabel('数字频率');ylabel('幅度响应dB');figure(2)plot(w,abs(h));grid;

xlabel('数字频率');

ylabel('幅度响应|H(e^(jw))|');figure(3)zplane(z,p);xx=filter(b,a,x);yy=fft(xx,L);figure(4)subplot(2,1,1)stem(l,x);

title('未经滤波的心电图信号');xlabel('n');subplot(2,1,2)stem(l,xx);

title('经滤波之后的心电图信号');xlabel('n');figure(5)subplot(2,1,1)plot(l,abs(y));

title('未经滤波的心电图信号的频谱');subplot(2,1,2)plot(l,abs(yy));

title('经滤波处理的心电图信号的频谱');

2．用help查看内部函数cheb1ord.m及cheby1.m，了解调用格式。

编程设计教材习题6-2，求模拟滤波器Ha(s)的表达式。

源程序

close all clear all clc

Wp=2*pi*3000;Rp=2;Ws=2*pi*12000;Rs=50;Fs=24000;

w=linspace(0,pi,1024);

[N,Wn]=cheb1ord(Wp,Ws,Rp,Rs,'s');e=sqrt(10^(Rp/10)-1);[b,a]=cheby1(N,e,Wn,'s')[numa,dena]=impinvar(b,a,Fs);h=freqz(numa,dena,w);norm=max(abs(h));

plot(w*Fs/pi,20*log10(abs(h)/norm))title('幅度响应')xlabel('频率(Hz)')ylabel('幅度(dB)')grid

3．模拟滤波器的数字化

用内部函数impinvar及bilinear实现教材习题6-5，求数字滤波器H(z)的表达式。

源程序

close all clear all clc b1=[0 0 1];a1=[1 1 1];b2=[0 0 1];a2=[2 3 1];

[numa1,dena1]=impinvar(b1,a1,0.5)[numa2,dena2]=bilinear(b1,a1,0.5)[numa3,dena3]=impinvar(b2,a2,0.5)

[numa4,dena4]=bilinear(b2,a2,0.5)

本实验所用的部分MATLAB函数

 L=length(x)：求序列x长度。

 y=fft(x,L)：将序列x(n)做L点快速傅立叶变换，结果赋给序列y(n)。 [n,Wn] = buttord(Wp,Ws,Rp,Rs,'s')：计算模拟Butterworth滤波器的最小阶次n和截止频率为Wn。

 [b,a] = butter(n,Wn,'s')：设计模拟截止频率为Wn（rad/s）的n阶 Butterworth低通滤波器，返回值为模拟滤波器的系数。

 y=filter(b,a,x): 将序列x(n)通过滤波器滤波后生成序列y(n)，滤波器的分母多项式系数构成a向量，分子多项式系数构成b向量。

 [BZ,AZ] = impinvar(B,A,Fs)：冲激响应不变法，返回值为数字滤波器的系数。 [BZ,AZ] = bilinear(B,A,fs)：双线性变换，返回值为数字滤波器的系数。 [H w]=freqz(b,a)：由滤波器分母多项式系数构成的a向量和分子多项式系数构成的b向量求系统频响。

截图

实验体会

这次的实验让这让我看到了理论与实践相结合的优势与用处，让我受益匪浅。我认识到了自己理论知识的不足，也认识到了我们学习的基础知识究竟能运用于什么领域，如何运用。我们在老师的耐心指导下调试电路，直到得到要求的效果。让我们在学习电路、信号等理论知识的同时，明白如何把这些应用于实际。

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

实验一 自适应滤波器

一、实验目的

1、掌握功率谱估计方法

2、会用matlab对功率谱进行仿真

二、实验原理

功率谱估计方法有很多种，一般分成两大类，一类是经典谱估计；另一类是现代谱估计。经典谱估计可以分成两种，一种是BT法，另一种是周期法；BT法是先估计自相关函数，然后将相关函数进行傅里叶变换得到功率谱函数。相应公式如下所示：

1ˆxx(m)rNˆPBTmN|m|1n0x*(n)x(nm)(11)

(12)ˆxx(m)ejwnr周期图法是采用功率谱的另一种定义，但与BT法是等价的，相应的功率谱估计如下所示：

1jwˆ

Pxx(e)N其计算框图如下所示：

观测数据x(n)FFTjwnx(n)en0N120nN1(13)

取模的平方1/N Pxx(ejw)

图1.1周期图法计算用功率谱框图 由于观测数据有限，所以周期图法估计分辨率低，估计误差大。针对经典谱估计的缺点，一般有三种改进方法：平均周期图法、窗函数法和修正的周期图平均法。

三、实验要求

信号是正弦波加正态零均值白噪声，信噪比为10dB，信号频率为2kHZ，取样频率为100kHZ。

四、实验程序与实验结果（1）用周期图法进行谱估计 A、实验程序： %用周期法进行谱估计 clear all;N1=128;%数据长度 N2=256;N3=512;N4=1024;f=2;%正弦波频率，单位为kHZ fs=100;%抽样频率,单位为kHZ n1=0:N1-1;n2=0:N2-1;n3=0:N3-1;n4=0:N4-1;a=sqrt(20);%由信噪比为10dB计算正弦信号的幅度 wn1=randn(1,N1);xn1=a*sin(2*pi*f*n1./fs)+wn1;Pxx1=10*log10(abs(fft(xn1).^2)/N1);%周期法求功率谱 f1=((0:length(Pxx1)-1))/length(Pxx1);wn2=randn(1,N2);xn2=a*sin(2*pi*f*n2./fs)+wn2;Pxx2=10*log10(abs(fft(xn2).^2)/N2);f2=((0:length(Pxx2)-1))/length(Pxx2);wn3=randn(1,N3);xn3=a*sin(2*pi*f*n3./fs)+wn3;Pxx3=10*log10(abs(fft(xn3).^2)/N3);f3=((0:length(Pxx3)-1))/length(Pxx3);wn4=randn(1,N4);xn4=a*sin(2*pi*f*n4./fs)+wn4;Pxx4=10*log10(abs(fft(xn4).^2)/N4);f4=((0:length(Pxx4)-1))/length(Pxx4);subplot(2,2,1);plot(f1,Pxx1);xlabel('频率');ylabel('功率(dB)');title('功率谱Pxx，N=128');subplot(2,2,2);plot(f2,Pxx2);xlabel('频率');ylabel('功率(dB)');title('功率谱Pxx，N=256');subplot(2,2,3);plot(f3,Pxx3);xlabel('频率');ylabel('功率(dB)');title('功率谱Pxx，N=512');subplot(2,2,4);plot(f4,Pxx4);xlabel('频率');ylabel('功率(dB)');title('功率谱Pxx，N=1024');B、实验仿真结果：

（2）采用汉明窗，分段长度L=32，用修正的周期图求平均法进行谱估计 A：实验程序: clear all;N=512;%数据长度 Ns=32;%分段长度

f1=2;%正弦波频率，单位为kHZ fs=100;%抽样频率,单位为kHZ n=0:N-1;a=sqrt(20);%由信噪比为10dB计算正弦信号的幅度 wn=randn(1,N);xn=a*sin(2*pi*f1*n./fs)+wn;w=hamming(32)';%汉明窗

Pxx1=abs(fft(w.*xn(1:32),Ns).^2)/norm(w)^2;Pxx2=abs(fft(w.*xn(33:64),Ns).^2)/norm(w)^2;Pxx3=abs(fft(w.*xn(65:96),Ns).^2)/norm(w)^2;Pxx4=abs(fft(w.*xn(97:128),Ns).^2)/norm(w)^2;Pxx5=abs(fft(w.*xn(129:160),Ns).^2)/norm(w)^2;Pxx6=abs(fft(w.*xn(161:192),Ns).^2)/norm(w)^2;Pxx7=abs(fft(w.*xn(193:224),Ns).^2)/norm(w)^2;Pxx8=abs(fft(w.*xn(225:256),Ns).^2)/norm(w)^2;Pxx9=abs(fft(w.*xn(257:288),Ns).^2)/norm(w)^2;Pxx10=abs(fft(w.*xn(289:320),Ns).^2)/norm(w)^2;Pxx11=abs(fft(w.*xn(321:352),Ns).^2)/norm(w)^2;Pxx12=abs(fft(w.*xn(353:384),Ns).^2)/norm(w)^2;Pxx13=abs(fft(w.*xn(385:416),Ns).^2)/norm(w)^2;Pxx14=abs(fft(w.*xn(417:448),Ns).^2)/norm(w)^2;Pxx15=abs(fft(w.*xn(449:480),Ns).^2)/norm(w)^2;Pxx16=abs(fft(w.*xn(481:512),Ns).^2)/norm(w)^2;Pxx=10*log10((Pxx1+Pxx2+Pxx3+Pxx4+Pxx5+Pxx6+Pxx7+Pxx8+Pxx9+Pxx10+Pxx11+Pxx12+Pxx13+Pxx14+Pxx15+Pxx16)/16);f=((0:length(Pxx)-1))/length(Pxx);plot(f,Pxx);xlabel('频率');ylabel('功率(dB)');title('加窗平均周期图法功率谱Pxx，N=512');grid on;B：实验仿真结果：

五．参考文献：

[1] 丁玉美，阔永红，高新波.数字信号处理-时域离散随机信号处理[M].西安：西安电子科技大学出版社，2002.[2] 万建伟，王玲.信号处理仿真技术[M].长沙：国防科技大学出版社，2008.实验二

卡尔曼滤波器的设计

一．实验目的

1.熟悉并掌握卡尔曼滤波、自适应滤波和谱估计的原理。

2.可以仿真符合要求的卡尔曼滤波器、自适应滤波器和各种谱估计方法。

3.掌握卡尔曼滤波器的递推公式和仿真方法。4.熟悉matlab的用法。二．实验原理

卡尔曼滤波是用状态空间法描述系统的，由状态方程和测量方程所组成。卡尔曼滤波用前一个状态的估计值和最近一个观测数据来估计状态变量的当前值，并以状态变量的估计值的形式给出。其状态方程和量测方程如下所示：

xk1AkxkwkykCkxkvk(11)

(12)其中，k表示时间，输入信号wk是一白噪声，输出信号的观测噪声vk也是一个白噪声，输入信号到状态变量的支路增益等于1，即B=1；A表示状态变量之间的增益矩阵，可随时间变化，Ak表示第k次迭代的取值，C表示状态变量与输出信号之间的增益矩阵，可随时间变化，其信号模型如图1.1所示(k用k1代替)。

vkCk+ykWk-1+XkZ-1Xk-1Ak-1

图1.1 卡尔曼滤波器的信号模型

卡尔曼滤波是采用递推的算法实现的，其基本思想是先不考虑输入信号wk和观测噪声vk的影响，得到状态变量和输出信号的估计值，再用输出信号的估计误差加权矫正状态变量的估计值，使状态变量估计误差的均方值最小。其递推公式如下所示：

xˆke0.02xˆk1Hk(yke0.02xk1)1HP(P1)kkk

Pke0.04Pk11e0.04Pk(IHk)Pk(112a)(112b)(112c)(112d)假设初始条件Ak,Ck,Qk,Rk,yk,xk1,Pk1已知，其中x0E[x0],P0var[x0]，那么递推流程见图1.2所示。

Pk1式(1-5)Pk'式(1-4)Hk式(1-3)式(1-6)

图1.2 卡尔曼滤波递推流程图

三．实验要求

一连续平稳的随机信号x(t),自相关rx()e

xkPk，信号x(t)为加性噪声所干扰，噪声是白噪声，测量值的离散值y(k)为已知。Matlab仿真程序如下：

%编卡尔曼滤波递推程序，估计信号x(t)的波形 clear all;clc;Ak=exp(-0.02);%各系数由前面确定； Ck=1;Rk=0.1;p(1)=20;%各初值； Qk=1-exp(-0.04);

p1(1)=Ak*p(1)*Ak'+Qk;%由p1代表p'； x(1)=0;%设信号初值为0；

H(1)=p1(1)*Ck'*inv(Ck*p1(1)*Ck'+Rk);

zk=[-3.2,-0.8,-14,-16,-17,-18,-3.3,-2.4,-18,-0.3,-0.4,-0.8,-19,-2.0,-1.2,-11,-14,-0.9,0.8,10,0.2,0.5,-0.5,2.4,-0.5,0.5,-13,0.5,10,-12,0.5,-0.6,-15,-0.7,15,0.5,-0.7,-2.0,-19,-17,-11,-14] %zk为测量出来的离散值； N=length(zk);%要测量的点数； for k=2:N

p1(k)=Ak*p(k-1)*Ak'+Qk;%未考虑噪声时的均方误差阵； H(k)=p1(k)*Ck'*inv(Ck*p1(k)*Ck'+Rk);%增益方程； I=eye(size(H(k)));%产生和H（k）维数相同的单位矩阵； p(k)=(I-H(k)*Ck)*p1(k);%滤波的均方误差阵；

x(k)=Ak*x(k-1)+H(k)*(zk(k)-Ck*Ak*x(k-1));%递推公式； end, x %显示信号x(k)的数据； m=1:N;n=m*0.02;

plot(n,zk,'-r*',n,x,'-bo');%便于比较zk和x(k)在同一窗口输出； xlabel('t/s','Fontsize',16);ylabel('z(t),x(t)','fontsize',16);

title('卡尔曼滤波递推——x(t)的估计波形与z(t)波形','fontsize',16)legend('观测数据z(t)','信号估计值x(t)',2);grid;四．实验结果

五．实验小结

通过卡尔曼滤波估计信号与观测信号比较知，卡尔曼滤波输出的估计信号x(t)与实际观测到的离散值z(t)还是存在一定的误差，卡尔曼滤波是从初始状态 就采用递推方法进行滤波，那么在初值迭代后的一段时间内可能会出现较大的误差，随着迭代进行，各参数逐渐趋于稳定，后面的估计值与观察值的误差就减少了。

六．参考文献

[1] 丁玉美，阔永红，高新波.数字信号处理-时域离散随机信号处理[M].西安：西安电子科技大学出版社，2002.[2] 万建伟，王玲.信号处理仿真技术[M].长沙：国防科技大学出版社，2008.

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

数字信号处理实验教案

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

数字信号处理是一门理论和实际密切结合的课程，为深入掌握课程内容，最好在学习理论的同时，做习题和上机实验。上机实验不仅可以帮助读者深入的理解和消化基本理论，而且能锻炼同学们的独立解决问题的能力。本讲义在第三版的基础上编写了五个实验，前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的变换区间？（包括非周期信号和周期信号）