数字信号处理实验报告

数字信号处理实验报告

实验一 时域采样与频域采样

一、实验目的

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

的概念，以及频率域采样定理及其对频域采样点数选择的指导作用。

二、实验内容及步骤

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

给定模拟信号，xa(t)?Ae??tsin(?0t)u(t)

式中A=444.128，?=50特性曲线如图10.2.1

2π，?0=50

2πrad/s，它的幅频

图1 xa(t)的幅频特性曲线

现用DFT(FFT)求该模拟信号的幅频特性，以验证时域采样理论。

安照xa(t)的幅频特性曲线，选取三种采样频率，即Fs=1kHz，300Hz，200Hz。观测时间选T

p?50ms。

- 1 -

数字信号处理实验报告

为使用DFT，首先用下面公式产生时域离散信号，对三种采样频率，采样序列按顺序用x1(n)，x2(n)，x3(n)表示。 x(n)?xa(nT)?Ae??nTsin(?0nT)u(nT)

因为采样频率不同，得到的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?Mk。

要求： 编写实验程序，计算x1(n)、x2(n)和x3(n)的幅度特性，并绘图显示。观察分析频谱混叠失真。

（2）频域采样理论的验证。

给定信号如下：

?n?1? x(n)??27?n?0?0?n?1314?n?26其它

编写程序分别对频谱函数X(ej?)?FT[x(n)]在区间[0,2?]上等间隔采样32

和16点，得到X32(k)和X16(k)： X32(k)? X16(k)?X(eX(ej?)??2?32k , k?0,1,2,?31

j?)??2?16 , k?0,1,2,?15k

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

x32(n)?IFFT[X32(k)]32 , n?0,1,2,?,31

- 2 -

数字信号处理实验报告

x16(n)?IFFT[X16(k)]16 , n?0,1,2,?,15

分别画出

X(ej?)、

X32(k)和X16(k)的幅度谱，并绘图显示x(n)、

x32(n)和x16(n)的波形，进行对比和分析，验证总结频域采样理论。

三、实验程序清单

1 时域采样理论的验证程序清单

TSTEM:

function tstem(xnt,yn) %自制序列绘图函数

%yn代表xnt=y轴标示，x轴标示默认为‘n’; stem(xnt,'fill'); %只画一周期 xlabel('n'); ylabel(yn); ylim([-10 160]);

% 时域采样理论验证程序exp2a.m Tp=64/1000; %观察时间Tp=64微秒 %产生M长采样序列x(n) % Fs=1000;T=1/Fs; Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1;

A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5; xnt=A*exp(-alph*n*T).*sin(omega*n*T); Xk=T*fft(xnt,M); %M点FFT[xnt)]

- 3 -

数字信号处理实验报告

yn='xa(nT)';subplot(3,2,1);

tstem(xnt,yn); %调用自编绘图函数tstem绘制序列图 box on;title('(a) Fs=1000Hz'); k=0:M-1;fk=k/Tp;

subplot(3,2,2);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=1000Hz'); xlabel('f(Hz)');ylabel('

');axis([0,Fs,0,1.2*max(abs(Xk))])

%================================================= % Fs=300Hz和 Fs=200Hz的程序与上面Fs=1000Hz完全相同。

幅

度

2 频域采样理论的验证程序清单

%频域采样理论验证程序exp2b.m M=27;N=32;n=0:M; %产生M长三角波序列x(n)

xa=0:floor(M/2); xb= ceil(M/2)-1:-1:0; xn=[xa,xb]; Xk=fft(xn,1024); 24点FFT[x(n)], 用于近似序列x(n)的TF

X32k=fft(xn,32) ;2点FFT[x(n)]

x32n=ifft(X32k); 2点IFFT[X32(k)]得到x32(n) X16k=X32k(1:2:N); %隔点抽取X32k得到X16(K) x16n=ifft(X16k,N/2); 点IFFT[X16(k)]得到x16(n) subplot(3,2,2);stem(n,xn,'.');box on

- 4 -

数字信号处理实验报告

title('(b) 三角波序列

x(n)');xlabel('n');ylabel('x(n)');axis([0,32,0,20]) k=0:1023;wk=2*k/1024; %

subplot(3,2,1);plot(wk,abs(Xk));title('(a)FT[x(n)]'); xlabel('\\omega/\\pi');ylabel('|X(e^j^\\omega)|');axis([0,1,0,200]) k=0:N/2-1;

subplot(3,2,3);stem(k,abs(X16k),'.');box on title('(c)

16

点

频

域

采

样

');xlabel('k');ylabel('|X_1_6(k)|');axis([0,8,0,200]) n1=0:N/2-1;

subplot(3,2,4);stem(n1,x16n,'.');box on title('(d)

16

点

IDFT[X_1_6(k)]');xlabel('n');ylabel('x_1_6(n)');axis([0,32,0,20]) k=0:N-1;

subplot(3,2,5);stem(k,abs(X32k),'.');box on title('(e)

32

点

频

域

采

样

');xlabel('k');ylabel('|X_3_2(k)|');axis([0,16,0,200]) n1=0:N-1;

subplot(3,2,6);stem(n1,x32n,'.');box on title('(f)

- 5 -

32点

数字信号处理实验报告

IDFT[X_3_2(k)]');xlabel('n');ylabel('x_3_2(n)');axis([0,32,0,20])

四、实验程序运行结果

1、 时域采样理论的验证程序运行结果

exp2a.m如图10.3.2所示。由图可见，采样序列的频谱的确是以采样频率为周期对模拟信号频谱的周期延拓。当采样频率为1000Hz时频谱混叠很小；当采样频率为300Hz时，在折叠频率150Hz附近频谱混叠很严重；当采样频率为200Hz时，在折叠频率110Hz附近频谱混叠更很严重。

图2

2 时域采样理论的验证程序exp2b.m运行结果如图10.3.3所示。

- 6 -

数字信号处理实验报告

图10.3.3

该图验证了频域采样理论和频域采样定理。对信号x(n)的频谱函数X(ejω)在[0，2π]上等间隔采样N=16时， N点IDFT[XN(k)]得到的序列正是原序列x(n)以16为周期进行周期延拓后的主值区序列：

?xN(n)?IDFT[XN(k)]N?[?x(n?iN)]RN(n)i???

由于NM，频域采样定理，所以不存在时域混叠失真，因此。

xN(n)与

x(n)相同。

五、思考题

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

- 7 -

数字信号处理实验报告

采样？

答： 先对原序列x(n)以N为周期进行周期延拓后取主值区序列，

xN(n)?[?x(n?iN)]RN(n)

i????再计算N点DFT则得到N点频域采样

XN(k)?DFT[xN(n)]N =X(ej?)??2?Nk , k?0,1,2,?,N?1

六、实验心得

在实验前先要弄清实验原理，理解试验程序，熟悉matlab软件的使用方法，这样更有利于实验的进行以及对实验结果的了解。在实验时要认真操作，对出现的错误要认真思考做出及时修改。

本次试验让我对时域采样和频域采样加深了了解，让我对数字信号处理这门课加深了认识，同时通过理论与实验的结合让我对学习这门课更加有兴趣了。

- 8 -

数字信号处理实验报告

实验二、用FFT对信号作频谱分析

一．实验目的

学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析

误差及其原因，以便正确应用FFT。

二、 实验原理

用FFT对信号作频谱分析是学习数字信号处理的重要内容。经常需要进行谱分析的信号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是2?求2?/N?D/N，因此要

。可以根据此式选择FFT的变换区间N。误差主要来自于

用FFT作频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时离散谱的包络才能逼近于连续谱，因此N要适当选择大一些。

周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。

对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。如果是模拟周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分析进行。

三、实验步骤及内容

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

- 9 -

数字信号处理实验报告

x1(n)?R4(n)

?n?1,0?n?3?x2(n)??8?n,4?n?7??0,其它n?4?n,?x3(n)??n?3,??0,0?n?34?n?7

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

x4(n)?cos?4n

x5(n)?cos(?n/4)?cos(?n/8)

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

x6(t)?cos8?t?cos16?t?cos20?t 选择 采样频率Fs?64Hz，变换区间N=16,32,64 三种情况进行谱分

析。分别打印其幅频特性，并进行分析和讨论。

四、 实验程序

程序一

xn1=[1 1 1 1]; XK8=fft(xn1,8); XK16=fft(xn1,16); n1=0:15;

- 10 -

数字信号处理实验报告

subplot(3,2,1); stem(n1,abs(XK16),'.')

title('x1(n)的16点变换');xlabel('n');ylabel('|XK16|') grid; hold on; n2=0:7; subplot(3,2,2); stem(n2,abs(XK8),'.')

title('x1(n)的8点变换');xlabel('n');ylabel('|XK8|') grid; hold on;

xn2=[4 3 2 1 1 2 3 4]; XK8=fft(xn2,8); XK16=fft(xn2,16); n1=0:15; subplot(3,2,3); stem(n1,abs(XK16),'.')

title('x2(n)的16点变换');xlabel('n');ylabel('|XK16|') grid; hold on; n2=0:7; subplot(3,2,4);

- 11 -

数字信号处理实验报告

stem(n2,abs(XK8),'.')

title('x2(n)的8点变换');xlabel('n');ylabel('|XK8|') hold on; grid;

xn=[1 2 3 4 4 3 2 1]; XK8=fft(xn,8); XK16=fft(xn,16); n1=0:15; subplot(3,2,5); stem(n1,abs(XK16),'.')

title('x3(n)的16点变换');xlabel('n');ylabel('|XK16|') grid; hold on; n2=0:7; subplot(3,2,6); stem(n2,abs(XK8),'.')

title('x3(n)的8点变换');xlabel('n');ylabel('|XK8|') 实验结果：

- 12 -

数字信号处理实验报告

实验结果分析：对序列x（n)进行N点的DFT即对x（ejw）在{0 2pi）进行N点等间隔采样则第一行的两个图分别是对x（n）的16点和8点的采样

2由第二 三行的四个图可知由序列x2（n)=x（(n+3））8R8(n）则当序列进行8点的DFT时满足循环移位关系则，则8点的DFT的结果的X(K)的模相等，而进行16点的DFT时不满足循环位移所以16点的DFT的模不相等、 实验二程序 clc; clear all; n=0:15;

x4n=cos(pi*n/4); X4K16=fft(x4n,16); n1=0:15; subplot(3,2,1);

- 13 -

数字信号处理实验报告

stem(n1,abs(X4K16))

title('x4(n)的16点变换');xlabel('w/pi');ylabel('X4K16') n=0:7;

x4n=cos(pi*n/4); X4K8=fft(x4n,8); n2=0:7

subplot(3,2,2); stem(n2,abs(X4K8),'.')

title('x4(n)的8点变换');xlabel('w/pi');ylabel('X4K8') n=0:15;

x5n=cos(pi*n/4)+cos(pi*n/8); X5K16=fft(x5n,16); n1=0:15 subplot(3,2,3); stem(n1,X5K16,'.')

title('x5(n)的16点变换');xlabel('w/pi');ylabel('X5K16') n=0:7;

x5n=cos(pi*n/4)+cos(pi*n/8); X5K8=fft(x5n,8); n2=0:7

subplot(3,2,4); stem(n2,X5K8,'.')

- 14 -

数字信号处理实验报告

title('x5(n)的8点变换');xlabel('w/pi');ylabel('X5K8')

实验程序三 n=0:15; Fs=64; t=n/Fs;

x8t=cos(pi*8*t)+cos(pi*16*t)+cos(pi*20*t); XK16=fft(x8t,16); n1=0:15; subplot(3,1,1) stem(n1,abs(XK16),'.')

title('16点的变换');xlabel('n');ylabel('幅度') grid;

- 15 -

数字信号处理实验报告

n=0:31; Fs=64; t=n/Fs;

x8t=cos(pi*8*t)+cos(pi*16*t)+cos(pi*20*t); XK32=fft(x8t,32); n2=0:31; subplot(3,1,2) stem(n2,abs(XK32),'.')

title('32点的变换');xlabel('n');ylabel('幅度') n=0:63; Fs=64; t=n/Fs;

x8t=cos(pi*8*t)+cos(pi*16*t)+cos(pi*20*t); XK64=fft(x8t,64); n3=0:63; subplot(3,1,3); stem(n3,abs(XK64),'.')

title('64点的变换');xlabel('n');ylabel('幅度')

- 16 -

数字信号处理实验报告

实验三结果分析：

这是用DFT对模拟信号进行谱分析，

x6(t)?cos8?t?cos16?t?cos20?t

五．思考题

（1）对于周期序列，如果周期不知道，如何用FFT进行谱分析？ 答：尽量观察信号的时间长一些，以便得到周期。

（2）如何选择FFT的变换区间？（包括非周期信号和周期信号） 答：对非周期信号：N要选择大一些，离散谱的包络才能逼近于连续谱；对周期信号：只有用整数倍周期的长度做FFT，得到的离散谱才能代表周期信号的频谱。

（3）当N=8时，x2(n)和x3(n)的幅频特性会相同吗？为什么？N=16 呢？

答：由序列x2（n)=x（(n+3））8R8(n）则当序列进行8点的DFT时满

- 17 -

数字信号处理实验报告

足循环移位关系则，则8点的DFT的结果的X(K)的模相等，而进行16点的DFT时不满足循环位移所以16点的DFT的模不相等。

六、实验心得

通过本次试验我对FFT更加了解，让我对所学的章节有了更加深入的体会与了解。

- 18 -

数字信号处理实验报告

实验三、IIR数字滤波器设计及软件实现

一、实验目的

（1）熟悉用双线性变换法设计IIR数字滤波器的原理与方法； （2）学会调用MATLAB信号处理工具箱中滤波器设计函数（或滤波器设计分析工具fdatool）设计各种IIR数字滤波器，学会根据滤波需求确定滤波器指标参数。

（3）掌握IIR数字滤波器的MATLAB实现方法。

（3）通过观察滤波器输入输出信号的时域波形及其频谱，建立数字滤波的概念。

二、实验原理

设计IIR数字滤波器一般采用间接法（脉冲响应不变法和双线性变换法），应用最广泛的是双线性变换法。基本设计过程是：①先将给定的数字滤波器的指标转换成过渡模拟滤波器的指标； ②设计过渡模拟滤波器；③将过渡模拟滤波器系统函数转换成数字滤波器的系统函数。MATLAB信号处理工具箱中的各种IIR数字滤波器设计函数都是采用双线性变换法。第六章介绍的滤波器设计函数butter、cheby1 、cheby2 和ellip可以分别被调用来直接设计巴特沃斯、切比雪夫1、切比雪夫2和椭圆模拟和数字滤波器。本实验要求读者调用如上函数直接设计IIR数字滤波器。

本实验的数字滤波器的MATLAB实现是指调用MATLAB信号处理工具箱函数filter对给定的输入信号x(n)进行滤波，得到滤波后的输

- 19 -

数字信号处理实验报告

出信号y(n）。

三、 实验内容及步骤

（1）调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st，该函数还会自动绘图显示st的时域波形和幅频特性曲线，如图10.4.1所示。由图可见，三路信号时域混叠无法在时域分离。但频域是分离的，所以可以通过滤波的方法在频域分离，这就是本实验的目的。

图10.4.1 三路调幅信号st的时域波形和幅频特性曲线

（2）要求将st中三路调幅信号分离，通过观察st的幅频特性曲线，分别确定可以分离st中三路抑制载波单频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的通带截止频率和阻带截止频率。要求滤波器的通带最大衰减为0.1dB,阻带最小衰减为60dB。

提示：抑制载波单频调幅信号的数学表示式为

- 20 -

数字信号处理实验报告

s(t)?cos(2?f0t)cos(2?fct)?12[cos(2?(fc?f0)t)?cos(2?(fc?f0)t)]

其中，cos(2?fct)称为载波，fc为载波频率，cos(2?f0t)称为单频调制信

?f0。由上式可见，所谓抑

号，f0为调制正弦波信号频率，且满足fc制载波单频调幅信号，就是2个正弦信号相乘，它有2个频率成分：和频fc?f0和差频fc?f0，这2个频率成分关于载波频率fc对称。所以，

1路抑制载波单频调幅信号的频谱图是关于载波频率fc对称的2根谱线，其中没有载频成分，故取名为抑制载波单频调幅信号。容易看出，图10.4.1中三路调幅信号的载波频率分别为250Hz、500Hz、1000Hz。如果调制信号m(t)具有带限连续频谱，无直流成分，则

s(t)?m(t)cos(2?fct)就是一般的抑制载波调幅信号。其频谱图是关于载

波频率fc对称的2个边带（上下边带），在专业课通信原理中称为双边带抑制载波 (DSB-SC) 调幅信号,简称双边带 (DSB) 信号。如果调制信号m(t)有直流成分，则s(t)?m(t)cos(2?fct)就是一般的双边带调幅

信号。其频谱图是关于载波频率fc对称的2个边带（上下边带），并包含载频成分。

（3）编程序调用MATLAB滤波器设计函数ellipord和ellip分别设计这三个椭圆滤波器，并绘图显示其幅频响应特性曲线。 （4）调用滤波器实现函数filter，用三个滤波器分别对信号产生函数mstg产生的信号st进行滤波，分离出st中的三路不同载波频率的调幅信号y1(n)、y2(n)和y3(n)， 并绘图显示y1(n)、y2(n)和y3(n)的时域波形，观察分离效果。

四、信号产生函数mstg清单

- 21 -

数字信号处理实验报告

function st=mstg

%产生信号序列向量st,并显示st的时域波形和频谱

%st=mstg 返回三路调幅信号相加形成的混合信号，长度N=1600 N=1600 %N为信号st的长度。

Fs=10000;T=1/Fs;Tp=N*T; %采样频率Fs=10kHz，Tp为采样时间 t=0:T:(N-1)*T;k=0:N-1;f=k/Tp;

fc1=Fs/10; %第1路调幅信号的载波频率fc1=1000Hz, fm1=fc1/10; %第1路调幅信号的调制信号频率fm1=100Hz fc2=Fs/20; %第2路调幅信号的载波频率fc2=500Hz fm2=fc2/10; %第2路调幅信号的调制信号频率fm2=50Hz fc3=Fs/40; %第3路调幅信号的载波频率fc3=250Hz, fm3=fc3/10; %第3路调幅信号的调制信号频率fm3=25Hz xt1=cos(2*pi*fm1*t).*cos(2*pi*fc1*t); %产生第1路调幅信号 xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t); %产生第2路调幅信号 xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t); %产生第3路调幅信号 st=xt1+xt2+xt3; %三路调幅信号相加 fxt=fft(st,N); %计算信号st的频谱

%====以下为绘图部分，绘制st的时域波形和幅频特性曲线==================== subplot(3,1,1)

plot(t,st);grid;xlabel('t/s');ylabel('s(t)');

axis([0,Tp/8,min(st),max(st)]);title('(a) s(t)的波形')

- 22 -

数字信号处理实验报告

subplot(3,1,2)

stem(f,abs(fxt)/max(abs(fxt)),'.');grid;title('(b) s(t)的频谱')

axis([0,Fs/5,0,1.2]); xlabel('f/Hz');ylabel('幅度')

五、实验程序框图

调用函数mstg产生st，自动绘图 显示st的时域波形和幅频特性曲线 调用ellipord和ellip分别设计三个椭圆滤波器，并绘图显示其幅频响应特性曲线。 调用filter，用三个滤波器分别对信号st进行滤波，分离出三路不同载波频率的调幅信号y1(n)、y2(n)和y3(n) 绘图显示y1(n)、y2(n)和y3(n)的时域波形和幅频特性曲线 End

图10.4.2 实验4程序框图

六、滤波器参数及实验程序清单

1、滤波器参数选取

观察图10.4.1可知，三路调幅信号的载波频率分别为250Hz、500Hz、1000Hz。带宽（也可以由信号产生函数mstg清单看出）分别为50Hz、100Hz、200Hz。所以，分离混合信号st中三路抑制载波单

- 23 -

数字信号处理实验报告

频调幅信号的三个滤波器（低通滤波器、带通滤波器、高通滤波器）的指标参数选取如下：

对载波频率为250Hz的条幅信号，可以用低通滤波器分离，其指标为

带截止频率fp?280Hz，通带最大衰减?p?0.1dBdB； 阻带截止频率fs?450Hz，阻带最小衰减?s?60dBdB， 对载波频率为500Hz的条幅信号，可以用带通滤波器分离，其指标为

带截止频率fpl?440Hz，fpu?560Hz，通带最大衰减

?p?0.1dBdB；

阻带截止频率fsl?275Hz，fsu?900Hz，Hz，阻带最小衰减?s?60dBdB，

对载波频率为1000Hz的条幅信号，可以用高通滤波器分离，其指标为

带截止频率fp?890Hz，通带最大衰减?p?0.1dBdB； 阻带截止频率fs?550Hz，阻带最小衰减?s?60dBdB， 说明：（1）为了使滤波器阶数尽可能低，每个滤波器的边界频率选择原则是尽量使滤波器过渡带宽尽可能宽。

（2）与信号产生函数mstg相同，采样频率Fs=10kHz。 （3）为了滤波器阶数最低，选用椭圆滤波器。

按照图10.4.2 所示的程序框图编写的实验程序为exp4.m。

2、实验程序清单

- 24 -

数字信号处理实验报告

%实验4程序exp4.m

% IIR数字滤波器设计及软件实现 clear all;close all

Fs=10000;T=1/Fs; %采样频率

%调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st st=mstg; %

低

通

滤

波

器

设

计

与

实

现

========================================= fp=280;fs=450;

wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60; ?指标（低通滤波器的通、阻带边界频）

[N,wp]=ellipord(wp,ws,rp,rs); %调用ellipord计算椭圆DF阶数N和通带截止频率wp

[B,A]=ellip(N,rp,rs,wp); %调用ellip计算椭圆带通DF系统函数系数向量B和A

y1t=filter(B,A,st); %滤波器软件实现 % 低通滤波器设计与实现绘图部分 figure(2);subplot(3,1,1);

myplot(B,A); %调用绘图函数myplot绘制损耗函数曲线 yt='y_1(t)';

subplot(3,1,2);tplot(y1t,T,yt); %调用绘图函数tplot绘制滤

- 25 -

数字信号处理实验报告

波器输出波形 %

带

通

滤

波

器

设

计

与

实

现

==================================================== fpl=440;fpu=560;fsl=275;fsu=900;

wp=[2*fpl/Fs,2*fpu/Fs];ws=[2*fsl/Fs,2*fsu/Fs];rp=0.1;rs=60;

[N,wp]=ellipord(wp,ws,rp,rs); %调用ellipord计算椭圆DF阶数N和通带截止频率wp

[B,A]=ellip(N,rp,rs,wp); %调用ellip计算椭圆带通DF系统函数系数向量B和A

y2t=filter(B,A,st); %滤波器软件实现 % 带通滤波器设计与实现绘图部分（省略） %

高

通

滤

波

器

设

计

与

实

现

================================================ fp=890;fs=600;

wp=2*fp/Fs;ws=2*fs/Fs;rp=0.1;rs=60; ?指标（低通滤波器的通、阻带边界频）

[N,wp]=ellipord(wp,ws,rp,rs); %调用ellipord计算椭圆DF阶数N和通带截止频率wp

[B,A]=ellip(N,rp,rs,wp,'high'); %调用ellip计算椭圆带通DF系统函数系数向量B和A

y3t=filter(B,A,st); %滤波器软件实现

- 26 -

数字信号处理实验报告

% 高低通滤波器设计与实现绘图部分（省略）

七、 实验程序运行结果

实验4程序exp4.m运行结果如图104.2所示。由图可见，三个分离滤波器指标参数选取正确，算耗函数曲线达到所给指标。分离出的三路信号y1(n)，y2(n)和y3(n)的波形是抑制载波的单频调幅波。

(a) 低通滤波器损耗函数及其分离出的调幅信号y1(t)

(b) 带通滤波器损耗函数及其分离出的调幅信号y2(t)

- 27 -

数字信号处理实验报告

(c)高通滤波器损耗函数及其分离出的调幅信号y3(t)

图104. 实验4程序exp4.m运行结果

八、思考题

（1）请阅读信号产生函数mstg，确定三路调幅信号的载波频率和调制信号频率。

（2）信号产生函数mstg中采样点数N=800，对st进行N点FFT可以得到6根理想谱线。如果取N=1000，可否得到6根理想谱线？为什么？N=2000呢？请改变函数mstg中采样点数N的值，观察频谱图验证您的判断是否正确。

（3）修改信号产生函数mstg，给每路调幅信号加入载波成分，产生调幅（AM）信号，重复本实验，观察AM信号与抑制载波调幅信号的时域波形及其频谱的差别。 提示：AM信号表示式：s(t)?[1?cos(2?f0t)]cos(2?fct)。

答：（1）三路调幅信号的载波频率分别为250Hz、500Hz、1000Hz。调制信号分别为50Hz、100Hz、200Hz。

（2） 因为信号st是周期序列，谱分析时要求观察时间为整

- 28 -

数字信号处理实验报告

数倍周期。所以，本题的一般解答方法是，先确定信号st的周期，在判断所给采样点数N对应的观察时间Tp=NT是否为st的整数个周期。但信号产生函数mstg产生的信号st共有6个频率成分，求其周期比较麻烦，故采用下面的方法解答。

分析发现，st的每个频率成分都是25Hz的整数倍。采样频率Fs=10kHz=25×400Hz，即在25Hz的正弦波的1个周期中采样400点。所以，当N为400的整数倍时一定为st的整数个周期。因此，采样点数N=800和N=2000时，对st进行N点FFT可以得到6根理想谱线。如果取N=1000，不是400的整数倍，不能得到6根理想谱线。

九、实验心得

通过本次试验，我对数字滤波器有了更深入的了解，也对matlab的应用更加熟练，这次试验提高了我的工程素质。

- 29 -

数字信号处理实验报告

实验四、FIR数字滤波器设计与软件实现

一、实验目的

（1）掌握用窗函数法设计FIR数字滤波器的原理和方法。 （2）掌握用等波纹最佳逼近法设计FIR数字滤波器的原理和方法。

（3）掌握FIR滤波器的快速卷积实现原理。 （4）学会调用MATLAB函数设计与实现FIR滤波器。

二、实验内容及步骤

（1）、认真复习第七章中用窗函数法和等波纹最佳逼近法设计FIR数字滤波器的原理；

（2）、调用信号产生函数xtg产生具有加性噪声的信号xt，并自动显示xt及其频谱，如图10.5.1所示；

图10.5.1 具有加性噪声的信号x(t)及其频谱如图

（3）、请设计低通滤波器，从高频噪声中提取xt中的单频调幅信号，要求信号幅频失真小于0.1dB，将噪声频谱衰减60dB。先观察xt的频谱，确定滤波器指标参数。

- 30 -

数字信号处理实验报告

（4）、根据滤波器指标选择合适的窗函数，计算窗函数的长度N，调用MATLAB函数fir1设计一个FIR低通滤波器。并编写程序，调用MATLAB快速卷积函数fftfilt实现对xt的滤波。绘图显示滤波器的频响特性曲线、滤波器输出信号的幅频特性图和时域波形图。 （5）、重复（3），滤波器指标不变，但改用等波纹最佳逼近法，调用MATLAB函数remezord和remez设计FIR数字滤波器。并比较两种设计方法设计的滤波器阶数。

提示：① MATLAB函数fir1和fftfilt的功能及其调用格式请查阅本书第7章和第9章；

②采样频率Fs=1000Hz，采样周期T=1/Fs；

③根据图10.6.1(b)和实验要求，可选择滤波器指标参数：通带截止频率fp=120Hz，阻带截至频率fs=150Hz，换算成数字频率，通带截止频率?p?2?fp??0.24?，通带最大衰为0.1dB，阻带截至频率60dB。]

?s?2?fs??0.3?，阻带最小衰为

④实验程序框图如图10.5.2所示。

- 31 -

数字信号处理实验报告

Fs=1000，T=1/Fs xt=xtg 产生信号xt, 并显示xt及其频谱 用窗函数法或等波纹最佳逼近法 设计FIR滤波器hn 对信号xt滤波：yt=fftfilt(hn,xt) 1、计算并绘图显示滤波器损耗函数 2、绘图显示滤波器输出信号yt End

图10.5.2 实验程序框图

（6）实验报告要求

①对两种设计FIR滤波器的方法（窗函数法和等波纹最佳逼近法）进行分析比较，简述其优缺点。

②附程序清单、打印实验内容要求绘图显示的曲线图。 ③分析总结实验结果。 ④简要回答思考题。

三、信号产生函数xtg程序清单

function xt=xtg(N)

%实验五信号x(t)产生,并显示信号的幅频特性曲线

%xt=xtg(N) 产生一个长度为N,有加性高频噪声的单频调幅信号xt,采样频率Fs=1000Hz

%载波频率fc=Fs/10=100Hz,调制正弦波频率f0=fc/10=10Hz. Fs=1000;T=1/Fs;Tp=N*T;

- 32 -

数字信号处理实验报告

t=0:T:(N-1)*T;

fc=Fs/10;f0=fc/10; %载波频率fc=Fs/10，单频调制信号频率为f0=Fc/10;

mt=cos(2*pi*f0*t); %产生单频正弦波调制信号mt，频率为f0 ct=cos(2*pi*fc*t); %产生载波正弦波信号ct，频率为fc xt=mt.*ct; %相乘产生单频调制信号xt nt=2*rand(1,N)-1; %产生随机噪声nt

%=======设计高通滤波器hn,用于滤除噪声nt中的低频成分,生成高通噪声=======

fp=150; fs=200;Rp=0.1;As=70; % 滤波器指标

fb=[fp,fs];m=[0,1]; % 计算remezord函数所需参数f,m,dev

dev=[10^(-As/20),(10^(Rp/20)-1)/(10^(Rp/20)+1)]; [n,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez函数所需参数

hn=remez(n,fo,mo,W); % 调用remez函数进行设计,用于滤除噪声nt中的低频成分

yt=filter(hn,1,10*nt); %滤除随机噪声中低频成分，生成高通噪声yt

%================================================================

xt=xt+yt; %噪声加信号

- 33 -

数字信号处理实验报告

fst=fft(xt,N);k=0:N-1;f=k/Tp;

subplot(3,1,1);plot(t,xt);grid;xlabel('t/s');ylabel('x(t)');

axis([0,Tp/5,min(xt),max(xt)]);title('(a) 信号加噪声波形')

subplot(3,1,2);plot(f,abs(fst)/max(abs(fst)));grid;title('(b) 信号加噪声的频谱')

axis([0,Fs/2,0,1.2]);xlabel('f/Hz');ylabel('幅度')

四、 滤波器参数及实验程序清单

1、滤波器参数选取

根据10.5.1 节实验指导的提示③选择滤波器指标参数：通带截止频率fp=120Hz，阻带截至频率fs=150Hz。代入采样频率Fs=1000Hz，换算成数字频率，通带截止频率?阻带截至频率?sp?2?fp??0.24?，通带最大衰为0.1dB，

?2?fs??0.3?，阻带最小衰为60dB。所以选取blackman

窗函数。与信号产生函数xtg相同，采样频率Fs=1000Hz。 按照图10.5.2 所示的程序框图编写的实验程序为exp5.m。

2、实验程序清单

%《数字信号处理(第三版）学习指导》第10章实验5程序exp5.m % FIR数字滤波器设计及软件实现 clear all;close all;

%==调用xtg产生信号xt, xt长度N=1000,并显示xt及其频谱,=========

- 34 -

数字信号处理实验报告

N=1000;xt=xtg(N);

fp=120; fs=150;Rp=0.2;As=60;Fs=1000; % 输入给定指标 % (1) 用窗函数法设计滤波器

wc=(fp+fs)/Fs; %理想低通滤波器截止频率(关于pi归一化）

B=2*pi*(fs-fp)/Fs; %过渡带宽度指标 Nb=ceil(11*pi/B); %blackman窗的长度N hn=fir1(Nb-1,wc,blackman(Nb));

Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性 ywt=fftfilt(hn,xt,N); %调用函数fftfilt对xt滤波 %以下为用窗函数法设计法的绘图部分（滤波器损耗函数，滤波器输出信号波形）

%省略

% (2) 用等波纹最佳逼近法设计滤波器

fb=[fp,fs];m=[1,0]; % 确定remezord函数所需参数f,m,dev

dev=[(10^(Rp/20)-1)/(10^(Rp/20)+1),10^(-As/20)]; [Ne,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez函数所需参数

hn=remez(Ne,fo,mo,W); % 调用remez函数进行设计 Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性 yet=fftfilt(hn,xt,N); % 调用函数fftfilt对xt滤波

- 35 -

数字信号处理实验报告

%以下为用等波纹设计法的绘图部分（滤波器损耗函数，滤波器输出信号yw(nT)波形）

%省略

五、实验程序运行结果

用窗函数法设计滤波器，滤波器长度 Nb=184。滤波器损耗函数和滤波器输出yw(nT)分别如图10.5.3(a)和（b）所示。

用等波纹最佳逼近法设计滤波器，滤波器长度 Ne=83。滤波器损耗函数和滤波器输出ye(nT)分别如图10.5.3(c)和（d）所示。

两种方法设计的滤波器都能有效地从噪声中提取信号，但等波纹最佳逼近法设计的滤波器阶数低得多，当然滤波实现的运算量以及时延也小得多，从图10.5.3(b)和（d）可以直观地看出时延差别。

- 36 -

数字信号处理实验报告

图10.5.3

六、简答思考题

思考题

(1)如果给定通带截止频率和阻带截止频率以及阻带最小衰减,如何用窗函数法设计线性相位低通滤波器?请写出设计步骤.

答：用窗函数法设计FIR滤波器的步骤如下：

1、根据对阻带衰减及过渡带的指标要求，选择窗函数的类型，并估计窗口长度N。

2、构造希望逼近的频率响应函数 3、计算相应的单位脉冲响应 4、加窗得到设计结果。

(2)如果要求用窗函数法设计带通滤波器,且给定通带上、下截止频率为?pl和?pu，阻带上、下截止频率为?sl和?su，试求理想带通滤波器的截止频率?cl和?cu。

答： 希望逼近的理想带通滤波器的截止频率?cl和?cu分别为：

- 37 -

数字信号处理实验报告

?cl?(?sl??pl)/2, ?cu?(?su??pu)/2

（3）解释为什么对同样的技术指标，用等波纹最佳逼近法设计的滤波器阶数低？

答：①用窗函数法设计的滤波器，如果在阻带截止频率附近刚好满足，则离开阻带截止频率越远，阻带衰减富裕量越大，即存在资源浪费；

② 几种常用的典型窗函数的通带最大衰减和阻带最小衰减固定，且差别较大，又不能分别控制。所以设计的滤波器的通带最大衰减和阻带最小衰减通常都存在较大富裕。如本实验所选的blackman窗函数，其阻带最小衰减为74dB,而指标仅为60dB。

③ 用等波纹最佳逼近法设计的滤波器，其通带和阻带均为等波纹特性，且通带最大衰减和阻带最小衰减可以分别控制，所以其指标均匀分布，没有资源浪费，所以期阶数低得多。

七、实验心得

通过本次试验让我对FIR滤波器有了更深的了解，我熟悉了FIR滤波器的设计以及设计时要注意的事项与难点，增加了我对实验的兴趣。

- 38 -

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