数字信号处理实验指导书汇总

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

一. 实验目的

（1）掌握求系统响应的方法。 （2）掌握时域离散系统的时域特性。

二.实验原理与方法

在时域中，描写系统特性的方法是差分方程和单位脉冲响应，在频域可以用系统函数描述系统特性。已知输入信号可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应，本实验仅在时域求解。在计算机上适合用递推法求差分方程的解，最简单的方法是采用MATLAB语言的工具箱函数filter函数。也可以用MATLAB语言的工具箱函数conv函数计算输入信号和系统的单位脉冲响应的线性卷积，求出系统的响应。

三．实验内容及结果

（1）编制程序，包括产生输入信号、单位脉冲响应序列的子程序，用filter函数或conv函数求解系统输出响应的主程序。程序中要有绘制信号波形的功能。 给定一个低通滤波器的差分方程为：

5n(?) y(n)?0.0x0.x0?n5(?1)y0.?n 9输入信号：

82 1

a) 分别求出系统对两个输入信号的响应序列，并画出其波形。 b) 求出系统的单位冲响应，画出其波形。

x(n)?R(n), x(n)?u(n)% A=[1,-0.9]; % B=[0.05,0.05];

% x1n=[1 1 1 1 1 1 1 1 zeros(1,25)]; % x2n=ones(1,32); % y1n=impz(B,A,33); % n=0:length(y1n)-1; % figure(1); % stem(n,y1n);

% title('1102303005 y1(n)'); % box on; % xlabel('n'); % ylabel('y1(n)'); % x2n=ones(1,33); % y2n=filter(B,A,x2n); % n=0:length(y2n)-1; % figure(2); % stem(n,y2n);

% title('1102303005 y2(n)'); % box on; % xlabel('n'); % ylabel('y2(n)'); % hn=filter(B,A,x1n); % n=0:length(hn)-1; % figure(3); % stem(n,hn);

% title('1102303005 h(n)'); % box on % xlabel('n'); %ylabel('h(n)');

（2）给定系统的单位脉冲响应为：

用线性卷积法分别求两个系统对输入信号x(n)?R8(n)的输出响应，并画出波形。

% x1n=[1 1 1 1 1 1 1 1 zeros(1,24)]; % h1n=[1 1 1 1 1 1 1 1 1 1 zeros(1,22)]; % h2n=[1 2.5 2.5 1 zeros(1,28)]; % y21n=conv(h1n,x1n); % y22n=conv(h2n,x1n); % n=0:length(y21n)-1; % figure(4); % stem(n,y22n);

% title('1102303005 x1(n)'); % box on; % xlabel('n'); % ylabel('x1(n)');

% n=0:length(y21n)-1; % figure(5); % stem(n,y21n);

% title('1102303005 y21(n)'); % box on; % xlabel('n');

% ylabel('y21(n)'); % n=0:length(y22n)-1; % figure(6); % stem(n,y22n);

% title('1102303005 y22(n)'); % box on; % xlabel('n');

% ylabel('y22(n)');

h1(n)?R10(n), h2(n)??(n)?2.5?(n?1)?2.5?(n?2)??(n?3)

（3）测量人耳辨别回声的最小时间——设计一个混合声音数字系统h(n)，利用

MATLAB语言实现声音的混合。

A=[1,-0.9]; B=[0.05,0.05];

x1n=[1 1 1 1 1 1 1 1 zeros(1,25)]; x2n=ones(1,32); xn=[1 zeros(1,32)]; y1n=filter(B,A,33); h=filter(B,A,xn); load mtlb.mat yn=conv(h,mtlb) n=0:length(yn)-1; figure(1); stem(n,yn);

title('1102303005 yn)'); box on; xlabel('n'); ylabel('yn');

四．思考题

（1）简述在时域求系统响应的方法。

有两种，一种是通过解差分方程求得系统的输出，注意要合理的选择初始条件，一种是已知系统的单位脉冲响应，通过球输出信号和系统的单位脉冲的线性卷积求得性输出。

（2）简述判断系统稳定性的方法。

如果信号经过低通滤波器，把信号的高频部分滤去，时域信号的的剧烈变化将被平滑。

实验二 离散系统的频域分析

一、实验目的

（1）熟悉对离散系统的频率响应分析方法。

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

二、实验原理

离散系统的时域方程为：

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

系统的频率响应为:

Z域:

系统函数为:

分解因式:

其中和称为零、极点。

，

三、实验内容及结果

（1）求差分方程

所对应的系统的频率响应。 程序：

k=256;%采样频率为256Hz 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(‘1102303005 实部’)；

xlabel('\\omega/\\pi');ylabel('幅度'); subplot(2,2,2);

plot(w/pi,imag(h));grid; title('1102303005 虚部‘)；

xlabel('\\omega/\\pi');ylabel('Amplitude'); subplot(2,2,3);

plot(w/pi,abs(h));grid;

title(‘1102303005 幅度谱’)；

xlabel('\\omega/\\pi');ylabel('幅值'); subplot(2,2,4);

plot(w/pi,angle(h));grid;

title(‘1102303005 相位谱’)；

xlabel('\\omega/\\pi');ylabel('弧度'); 结果：

（2）求系统

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

程序：

clc;

clear ;

num=[0.0528 0.0797 0.1295 0.1295 0.0797 0.0528]; den=[1 -1.8107 2.4947 -1.8801 0.9537 -0.2336]; figure(1)

zplane(num,den); grid;

[z,p,k]=tf2zp(num,den); disp('零点=');disp(z); disp('极点=');disp(p); disp('增益系数=');disp(k); k=256; %采样频率为K=256Hz w=0:pi/k:pi;

h=freqz(num,den,w); figure(2)

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('弧度'); 结果：

10.80.60.4Imaginary Part0.20-0.2-0.4-0.6-0.8-1-1-0.50Real Part0.51

虚部1实部2Amplitude1幅度0.50-0.5-100.5?/?相位谱10-100.5?/?幅度谱11.51幅值弧度420-200.5?/?1-400.5?/?10.50

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

实验三 FFT算法的应用

一、实验目的

（1）加深对离散信号的DFT的理解。

（2）熟悉FFT算法的应用。

二、实验原理

N点序列的DFT和IDFT变换定义式如下：

,

利用旋转因子具有周期性，可以得到快速算法（FFT）。 在MATLAB中，可以用函数和 计算N点序列的DFT正、反变换。

三、实验内容及结果

（1）观察一句话的频谱特点：将自己录制的一段话均匀分段，并用DFT计算每段信号的频谱。 程序： Fs=8000; T=4;

x=wavrecord(T*Fs,Fs); L=length(x); t=(1:L)/Fs; subplot(411); plot(t,x); xlabel('t/s'); ylabel('x(t)');

title('1102303005‘)； axis([0 4 -0.7 0.7]); subplot(413);

specgram(x,256,Fs,256,250); xlabel('t/s'); ylabel('f/Hz');

title('1102303005‘)； 结果：

（2）2N点实数序列

N=64。一次算出，

程序： N=64;

n=0:2*N-1;

x1=cos(2*pi*7*n/N);

x2=0.5*cos(2*pi*19*n/N); x=x1+x2;

subplot(2,1,1); stem(n,x);

xlabel('t/T');ylabel('x(n)'); title('1102303005‘)； xk=fft(x,N); disp(xk); xk=abs(xk); k=0:N-1;

subplot(2,1,2);

并绘出

的图形。

stem(k,xk);

xlabel('k');ylabel('X(k)'); title('1102303005‘)； 结果：

（3）已知某序列

在单位圆上的N=64等分样点的Z变换为:

。

用N点IFFT程序计算出

程序： N=64; k=0:N-1;

xk=1./(1-0.8*exp((2*pi*k/N)*-1*j)); subplot(2,1,1) stem(k,xk)

xlabel('k');ylabel('X(k)'); title('1102303005‘)； axis([0 63 0 6]) x=ifft(xk,N); subplot(2,1,2)

。

stem(k,x)

xlabel('t/T');ylabel('x(n)'); axis([0 63 0 1.5]) title('1102303005‘)； 结果：

四、思考题

（1）如果周期信号的周期预先不知道，如何用FFT进行谱分析？

答：周期信号的周期预先不知道时,可先截取M点进行DFT,再将截取长度扩大1倍截取,比较结果,如果二者的差别满足分析误差要求,则可以近似表示该信号的频谱,如果不满足误差要求就继续将截取长度加倍,重复比较,直到结果满足要求。

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

答：一、对于非周期信号：有频谱分辨率F，而频谱分辨率直接和FFT的变换区间有关，因为FFT能够实现的频率分辨率是2π/N...因此有最小的N>2π/F。就可以根据此式选择FFT的变换区间。

二、对于周期信号，周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。

实验四 IIR数字滤波器的设计

1、实验目的

（1）掌握双线性变换法及脉冲相应不变法设计IIR数字滤波器的具体设计方法; （2）熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。 2、实验原理

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

（1）利用buttord和cheb1ord确定低通原型巴特沃斯和切比雪夫滤波器的阶数和截止频率；

如：求阶数[N, Wn] = cheb1ord (Wp, Ws, Rp, Rs,’s’)

选择项说明：high－类别。缺省为low； s－模/数，缺省为数Rp即 ?p, Rs即?s ； Wn－Chebyshev自然频率（3dB频率），数字设计： Wp＝?p/? Ws＝?s/? 。 （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可以完成脉冲响应不变法的模拟滤波器到数字滤波器的转换。 3、实验内容

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

?p?0.3?，通带峰值起伏：；

阻带边缘频率：?s?0.8?，最小阻带衰减：的损耗函数和幅频特性。

%脉冲响应不变法 Fs=400; T=1/Fs;

。画出两种设计方法下

wp=0.3*pi; ws=0.8*pi; Wp=wp*Fs; Ws=ws*Fs; Ap=1; As=40;

[N,Wc]=buttord(Wp,Ws,Ap,As,'s'); [B,A]=butter(N,Wc,'s'); [D,C]=impinvar(B,A,Fs); Hz=freqz(D,C,W); f=(W*Fs)/(2*pi); figure(1)

plot(f,abs(Hz)); grid on;

title('1102303005脉冲响应不变法、幅频特性'); xlabel('频率/Hz'); ylabel('归一化幅值'); figure(2)

a=20*log10(abs(Hz)); plot(f,a); grid on;

title('1102303005脉冲响应不变法、衰减函数') xlabel('频率/Hz'); ylabel('幅度/dB');

%双线性变化 Fs=400; T=1/Fs; wp=0.3*pi; ws=0.8*pi;

Wp=(2/T)*tan(wp/2); Ws=(2/T)*tan(ws/2); Ap=1; As=40;

[N,Wc]=buttord(Wp,Ws,Ap,As,'s'); [B,A]=butter(N,Wc,'s'); [D,C]=bilinear(B,A,Fs); [Hz,W]=freqz(D,C); f=(W*Fs)/(2*pi); figure(1)

plot(f,abs(Hz)); grid on;

title('1102303005双线性变换法、幅频特性'); xlabel('频率/Hz'); ylabel('归一化幅值'); figure(2)

a=20*log10(abs(Hz));

plot(W/pi,a); grid on;

title('1102303005双线性变换法、幅频特性'); xlabel('频率/Hz'); ylabel('幅度/dB');

实验五 FIR数字滤波器的设计

1、实验目的：

（1）加深对数字滤波器的常用指标理解。 （2）学习数字滤波器的设计方法。 2、实验原理：

（1）低通滤波器的常用指标：

（a）通带边缘频率（d）通带峰值起伏（e）阻带起伏

；（b）阻带边缘频率；（c）通带起伏，

。

；

，最小阻带衰减

（2）在MATLAB中，熟悉函数fir1、kaiserord 、remezord、remez的使用； B = fir1(n,Wn,'high','noscale')设计滤波器；[n,Wn,beta,ftype] = kaiserord(f,a,dev)估计滤波器阶数；（默认情况下，过滤器的第一通带中心为精确的一个窗口。使用一个尾随'noscale'参数

来防止这种缩放remezord、remez一起使用设计等波纹线性相位FIR滤波器）

3、实验内容：

利用MATLAB编程,用窗函数法设计FIR数字滤波器，指标要求如下：

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

阻带边缘频率：

clear all

clc

wp1=0.45*pi; wp2=0.65*pi; ws1=0.3*pi; ws2=0.75*pi;

，最小阻带衰减：。

width1 =wp1-ws1; width2 =ws2-wp2;

width=min(width1,width2) N=ceil(6.62*pi/width) n =0: (N-1);

a =(N-1)/2+0.00001; wc1=(ws1+wp1)/2; wc2=(ws2+wp2)/2;

hd=(sin(wc2*(n-a))-sin(wc1*(n-a)))./(pi*(n-a)); win =(hanning(N))'; h=hd.*win; freqz(h);

subplot(2,1,1)

title('1102303005 幅频特性') xlabel('频率') ylabel('db') subplot(2,1,2)

title('1102303005 相频特性') xlabel('频率') ylabel('db')

实验六 数字信号处理在双音多频拨号系统中的应用

1、 DTMF的组成

在电话中，数字0~9的中每一个都用两个不同的单音频传输，所用的8个频率分成高频带和低频带两组，低频带有四个频率：679Hz,770Hz,852Hz和941Hz；高频带也有四个频率：1209Hz,1336Hz,1477Hz和1633Hz.。每一个数字均由高、低频带中各一个频率构成，例如1用697Hz和1209Hz两个频率，信号用sin(2?f1t)?sin(2?f2t)表示，其中f1?697Hz,f2?1209Hz。这样8个频率形成16种不同的双频信号。

2、 电话中的双音多频（DTMF）信号的产生与检测 （1）双音多频信号的产生

假设时间连续的DTMF信号用sin(2?f1t)?sin(2?f2t)表示，式中f1和f2是按照上表选择的两个频率，f1代表低频带中的一个频率，f2代表高频带中的一个频率。采用数字方法产生DTMF信号，规定用8KHz对DTMF信号进行采样，得到时域离散信号：

x(n)?sin(2?f1n/8000)?sin(2?f2n/8000)

因为采样频率是8000Hz，因此要求每125ms输出一个样本，得到的序列再送到D/A变换器和平滑滤波器，输出便是连续时间的DTMF信号。DTMF信号通过电话线路送到交换机。

（2）双音多频信号的检测

在接收端，要对收到的双音多频信号进行检测，检测两个正弦波的频率是

多少，以判断所对应的十进制数字或者符号。显然这里仍然要用数字方法进行检测，因此要将收到的时间连续 DTMF信号经过A/D变换，变成数字信号进行检测。检测的方法有两种，一种是用一组滤波器提取所关心的频率，根据有输出信号的2个滤波器判断相应的数字或符号。另一种是用DFT（FFT）对双音多频信号进行频谱分析，由信号的幅度谱，判断信号的两个频率，最后确定相应的数字或符号。

3、 根据实验原理利用数字信号处理方法进行双音多频信号的产生与解调。 clear all; clc;

tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68]; N=205;

K=[18,20,22,24,31,34,38,42]; f1=[697,770,852,941]; f2=[1209,1336,1477,1633]; TN=input('键入10学号='); TNr=0; for l=1:10;

d=fix(TN/10^(10-l)); TN=TN-d*10^(10-l); for p=1:4; for q=1:4;

if tm(p,q)==abs(d); break,end end

if tm(p,q)==abs(d); break,end end

n=0:1023;

x = sin(2*pi*n*f1(p)/8000) + sin(2*pi*n*f2(q)/8000); sound(x,8000); pause(0.1)

X=goertzel(x(1:205),K+1); val = abs(X); subplot(5,2,l); stem(K,val,'.'); grid;

title('1102303005'); xlabel('k'); ylabel('|H(k)|'); axis([10 50 0 120]) limit = 80; for s=5:8;

if val(s)limit, break, end end for r=1:4;

if val(r)limit, break, end end

TNr=TNr+tm(r,s-4)*10^(10-l); end

disp('接收端检测到的号码为：') disp(TNr) 运行结果：

发送学号：1102303005 检测学号：1102303005

对时域离散DTMF信号进行频率检测，幅度谱图如下：

实验七 滤除心电图信号中的电网干扰

心脏的健康状况可以通过采集心电图信号来观察，不过信号必须是没有环境污染的，才能提供给医生诊断。最常见的污染源来自50Hz的交流电源线，它辐射的电磁场会经过电容耦合还有电磁感应进入心电图信号。假设h(n)是心

电图的一个周期的信号，它的采样频率为50Hz；现将h(n)变为三个周期的心跳信号，并将采样频率提升到200Hz，加入50Hz的交流干扰。

h(n)?[78,?39.5,?70.7,?17,?17,?16.5,?16.5,?15,?15,?14.5,?10,?3.5,1,?3.5,?12,?16,?17,?18,?18,?18,?18,?18,?18,?18,?18,?16.5,?4,9,19,6.5,?11,?17,?17,?18, ?23,?18.3,49]请运用零极点设计法设计一个二阶的点阻滤波器，并用MATLAB的基本数学函数在计

算机上消除信号的交流干扰。

实验八 软件无线电的通信

软件无线电是利用计算程序来控制DSP工作的通信方式。任何信号通过电磁波传输时都是以模拟的方式进行，数字化的信号也不例外。下图描述了一种数字调幅通信，它将英文字符变为高频的电信号，进行所谓的无线电通信。

软件无线电发射机字符编码4元码内插h(m)基带信号脉冲形成调幅波幅度调制u(n)v(m)w(m)载波x(t)y(t)噪声软件无线电接收机字符解码抽取低通滤波幅度解调e(n)d(m)基带信号c(m)本机振荡b(t)a(t)

图8-1 软件无线电的通信原理

英文字符通常用ASCII编码，MATLAB能将字符表示成ASCII码。例如在命令窗口敲上字符串s =‘I am Yangyiming’，再敲上real(s)就能看到s的内在形式。字符在上图8-1中编码为4元码u(n)，4元码的基本符号是{-3，-1，1，3}，每个字符用四个基本符号表示。

4元码也是数字信号，它要想通过电波传送，就必须变成为某种脉冲。有些脉冲的频带比较宽，有些脉冲的频带比较窄，合理地选择脉冲形状可以提高信号传输效率。还有，电信号想要通过天线在空气中传播，天线的长度就必须大于被传输的电信号的波长的1/10。图8-1采用幅度调制的方法，降低电信号的波长。软件无线电发射机通过数字信号处理，将字符送上无线电传播的高频载波。

软件无线电接收机得到是携带字符含噪声高频信号，接下来的工作就是恢复或者是猜测发射机播出的内容。处理过程与发射的过程正好相反。

请用MATLAB模仿图8-1编写软件无线电程序。

本文来源：https://www.bwwdw.com/article/babp.html