数字信号处理(西电上机实验)

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

数字信号处理实验报告

实验一： 信号、 系统及系统响应

一、实验目的：

(1) 熟悉连续信号经理想采样前后的频谱变化关系， 加深对时域采样定理的理解。

(2) 熟悉时域离散系统的时域特性。

(3) 利用卷积方法观察分析系统的时域特性。 (4) 掌握序列傅里叶变换的计算机实现方法， 利用序列的傅里叶变换对连续信号、 离散信号及系统响应进行频域分析。 二、实验原理与方法： (1) 时域采样。

(2) LTI系统的输入输出关系。 三、 实验内容、步骤

(1) 认真复习采样理论、 离散信号与系统、 线性卷积、 序列的傅里叶变换及性质等有关内容， 阅读本实验原理与方法。 (2) 编制实验用主程序及相应子程序。

① 信号产生子程序， 用于产生实验中要用到的下列信号序列： a. xa(t)=A*e^-at *sin(Ω0t)u(t)

A=444.128;a=50*sqrt(2)*pi; b. 单位脉冲序列：xb(n)=δ(n) c. 矩形序列： xc(n)=RN(n), N=10

② 系统单位脉冲响应序列产生子程序。 本实验要用到两种FIR系统。 a. ha(n)=R10(n);

b. hb(n)=δ(n)+2.5δ(n-1)+2.5δ(n-2)+δ(n-3) ③ 有限长序列线性卷积子程序

用于完成两个给定长度的序列的卷积。 可以直接调用MATLAB语言中的卷积函数conv。 conv用于两个有限长度序列的卷积， 它假定两个序列都从n=0 开始。 调用格式如下：y=conv (x, h) 四、实验内容

调通并运行实验程序， 完成下述实验内容： ① 分析采样序列的特性。

a. 取采样频率fs=1 kHz, 即T=1 ms。

b. 改变采样频率， fs=300 Hz， 观察|X(ejω)|的变化， 并做记录(打印曲线)； 进一步降低采样频率， fs=200 Hz， 观察频谱混叠是否明显存在， 说明原因， 并记录(打印)这时的|X(ejω)|曲线。

② 时域离散信号、 系统和系统响应分析。

a. 观察信号xb(n)和系统hb(n)的时域和频域特性； 利用线性卷积求信号xb(n)通过系统hb(n)的响应y(n)， 比较所求响应y(n)和hb(n)的时域及频域特性， 注意它们之间有无差别， 绘图说明， 并用所学理论解释所得结果。 b. 观察系统ha(n)对信号xc(n)的响应特性。

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

③ 卷积定理的验证。 五、思考题

(1) 在分析理想采样序列特性的实验中， 采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同? 它们所对应的模拟频率是否相同? 为什么?

(2) 在卷积定理验证的实验中， 如果选用不同的频域采样点数M值， 例如， 选M=10和M=20， 分别做序列的傅里叶变换， 求得的结果有无差异？为什么？ 六、参考资料

数字信号处理教程——Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年

七、实验程序及结果

本实验使用自定义函数的方法产生信号： 门函数：

function[y,n]=gate(np,ns,nf) if ns>np|ns>nf|np>nf

error('ÊäÈëλÖòÎÊý²»Âú×ãns<=np<=nf'); else

n=[ns:nf]; y=(n>=0&n<np); end

冲击函数：

function[y,n]=impseq(np,ns,nf) if ns>np|ns>nf|np>nf

error('ÊäÈëλÖòÎÊý²»Âú×ãns<=np<=nf'); else

n=[ns:nf]; y=[(n-np)==0]; end

时域函数：

function[y,ts]=sig(t0,tp,t1);

A=444.128;a=50.*sqrt(2).*pi;w0=50.*sqrt(2).*pi; ts=t0:tp:t1;

y=A.*exp(-a.*ts).*sin(w0.*ts).*(ts>=0);

;hb=impseq(0,0,10)+2.5.*impseq(1,0,10)+2.5.*impseq(2,0,10)+impseq(3,0,10);

;[ha,N]=gate(10,0,10);

第一部分主程序：

[y1,t1]=sig(0,0.001,1); %t1=0:0.001:1;

%y1=0.5.*sin(2.*pi.*15.*t1); figure; plot(t1,y1);

axis([0 0.15 -10 150]);

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

title('f=1000ʱÓòÐźÅ'); figure; f=-500:500;

yy1=fftshift(fft(y1)); plot(f,abs(yy1)); title('f=1000ƵÆ×'); [y2,t2]=sig(0,1/300,1); figure; plot(t2,y2);

axis([0 0.15 -10 150]); title('f=300ʱÓòÐźÅ'); figure;

yy2=fftshift(fft(y2)); plot(-150:150,abs(yy2)); title('f=300ƵÆ×'); [y3,t3]=sig(0,0.02,1); figure; plot(t3,y3);

title('f=200ʱÓòÐźÅ'); axis([0 0.15 -8 5]); figure;

yy3=fftshift(fft(y3)); plot(-25:25,abs(yy3)); title('f=200ƵÆ×');

第二部分主程序：

xb=impseq(0,0,100); figure;

stem(0:length(xb)-1,xb); axis([-1,3,0,2]); figure; y1=fft(xb);

plot(0:length(y1)-1,y1);

hb=impseq(0,0,100)+2.5.*impseq(1,0,100)+2.5.*impseq(2,0,100)+impseq(3,0,100); figure;

stem(0:length(hb)-1,hb); axis([-1,5,0,3]) figure

y2=fftshift(fft(hb));

plot(-(length(y2)-1)/2:(length(y2)-1)/2,abs(y2)); figure;

y3=conv(double(xb),double(hb)); plot(0:length(y3)-1,y3); y4=fftshift(fft(y3));

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

figure;

plot(-(length(y4)-1)/2:(length(y4)-1)/2,abs(y4)); xc=gate(10,0,100); ha=xc;

y5=fftshift(fft(conv(double(xc),double(ha)))); figure;

plot(-(length(y5)-1)/2:(length(y5)-1)/2,abs(y5));

实验结果： 第一部分：

采样频率为200Hz时时域恢复 200Hz 采样频率为300Hz时时域恢复 300Hz 采样频率为1000Hz是时域恢复 1000Hz

第二部分：

时频谱

时频谱

时频谱

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

Xb频域分析

hb频域分析 信号通过xb频域

卷积定理的验证：

%%%%%卷积定理的验证

y1=[0,1,2,3,4,5,6,7,8]; y2=[8,7,6,5,4,3,2,1,0]; y3=conv(y1,y2); f1=fft(y1,20); f2=fft(y2,20); f3=fft(y3,20); f=f1.*f2;

xb时域分析

hb时域分析

信号通过xc频域

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

subplot(2,1,1);

stem(0:length(f)-1,f,'.'); subplot(2,1,2);

stem(0:length(f3)-1,f3,'.');

时域卷积后求频谱和频域相乘

可见，f= conv(y1,y2)的频谱和y1,y2的频谱相乘后结果相同。即满足卷积定理

八、回答问题

(1) 在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？它们所对应的模拟频率是否相同？为什么？

答：数字频率度量不相同，但他们所对应的模拟频率相同。由w=Ω*Ts公式得，采样间隔变化时模拟频率对应的数字频率会有相应的变化，故其度量会有所变化。而且采样频率的大小直接关系到能否将能否将原始信号恢复出来。 (2) 在卷积定理验证的实验中，如果选用不同的频域采样点数M值，例如，选M=10和M=20，分别做序列的傅里叶变换，求得的结果有无差异？

答：有差异，所到的结果点数不同。 九、总结

本实验主要是后续实验的基础。涉及内容也比较浅显。 单位冲击序列与hb卷积后得到的频谱与hb原频谱相同。原因很简单，是因为单位冲击序列卷积任何函数仍然是原函数。

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

实验二： 用FFT作谱分析

一、实验目的

(1) 进一步加深DFT算法原理和基本性质的理解(因为FFT只是DFT的一种快速算法， 所以FFT的运算结果必然满足DFT的基本性质)。 (2) 熟悉FFT算法原理和FFT子程序的应用。

(3) 学习用FFT对连续信号和时域离散信号进行谱分析的方法， 了解可能出现的分析误差及其原因， 以便在实际中正确应用FFT。 二、实验步骤

(1) 复习DFT的定义、 性质和用DFT作谱分析的有关内容。 (2) 复习FFT算法原理与编程思想， 并对照DIT-FFT运算流图和程序框图， 读懂本实验提供的FFT子程序。

(3) 编制信号产生子程序， 产生以下典型信号供谱分析用： x1(n) R4(n)

n 1,0 n 3

x2(n) 8 n4 n 7

0

4 n0 n 3

x(n) n 34 n 7 3

0

x4(n) cosn 4

x(n) sinn5 8 x6(n) cos8 t cos16 t cos20 t(4) 编写主程序。

下图给出了主程序框图， 供参考。 本实验提供FFT子程序和通用绘图子程序。 (5) 按实验内容要求， 上机实验， 并写出实验报告。

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

三、实验内容

(1) 对 2 中所给出的信号逐个进行谱分析。

(2) 令x(n)=x4(n)+x5(n)， 用FFT计算 8 点和 16 点离散傅里叶变换， X(k)=DFT［x(n)］

(3) 令x(n)=x4(n)+jx5(n)， 重复(2)。 四、思考题

(1) 在N=8时， x2(n)和x3(n)的幅频特性会相同吗? 为什么? N=16呢? (2) 如果周期信号的周期预先不知道， 如何用FFT进行谱分析? 五、参考资料

数字信号处理教程——Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年

数字信号处理 第三版，高西泉等著，西安电子科技大学出版社，第三版 六、实验程序及结果

x1 = [1,1,1,1,0,0,0,0,0,0]; x2 = [1,2,3,4,4,3,2,1,0,0]; x3 = [4,3,2,1,1,2,3,4,0,0]; n=0:1000; x4=cos(n); x5=sin(n); t=0:0.001:1;

x6=cos(8*pi*t)+cos(16*pi*t)+cos(20*pi*t); % figure;

% title('Öð¸öÆ×·ÖÎö')

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

stem(0:length(x1)-1,x1); subplot(3,2,2);

stem(0:length(x2)-1,x2); subplot(3,2,3);

stem(0:length(x3)-1,x3); subplot(3,2,4); stem(0:20,x4(1:21)); subplot(3,2,5); stem(0:20,x5(1:21)); subplot(3,2,6); stem(0:20,x6(1:21)); figure;

stem(0:200,x6(1:201)); y1=abs(fft(x1,8)); y2=abs(fft(x2,8)); y3=abs(fft(x3,8)); y4=abs(fft(x4,8)); y5=abs(fft(x5,8)); y6=abs(fft(x6,8)); figure;

subplot(3,2,1);

stem(0:length(y1)-1,y1); subplot(3,2,2);

stem(0:length(y2)-1,y2); subplot(3,2,3);

stem(0:length(y3)-1,y3); subplot(3,2,4);

stem(0:length(y1)-1,y1); subplot(3,2,5);

stem(0:length(y5)-1,y5); subplot(3,2,6);

stem(0:length(y6)-1,y6); %%%figure;

%%%stem(0:15,y6(1:16));

%%%%%%%%%%%%%%% y1=abs(fft(x1,16)); y2=abs(fft(x2,16)); y3=abs(fft(x3,16)); y4=abs(fft(x4,16)); y5=abs(fft(x5,16)); y6=abs(fft(x6,16)); figure;

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

stem(0:length(y1)-1,y1); subplot(3,2,2);

stem(0:length(y2)-1,y2); subplot(3,2,3);

stem(0:length(y3)-1,y3); subplot(3,2,4);

stem(0:length(y6)-1,y4); subplot(3,2,5);

stem(0:length(y6)-1,y5); subplot(3,2,6);

stem(0:length(y6)-1,y6); %figure;

%stem(0:15,y6(1:16)); x7=x4+j*x5; y71=abs(fft(x7,8)); y72=abs(fft(x7,16)); figure;

subplot(1,2,1);

stem(0:length(y71)-1,y71); subplot(1,2,2);

stem(0:length(y72)-1,y72);

程序结果如下所示：

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

每个信号时域分析

因为第六个信号不易在上图中展示（包含信息不到一个周期）所以单独画出图形如下：

x6的时域波形

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

逐个进行8

点离散傅里叶变换

逐个进行16点离散傅里叶变换

x(n)8点和16点离散傅里叶变换比较

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

七、回答问题

(1) 在N=8时，x2(n)和x3(n)的幅频特性会相同吗? 为什么? N=16呢? 答：N=8时幅频特性一样，N=16时幅频特性不一样。Fn所能分辨到频率为为Fs/N，如果采样频率Fs为1024Hz，采样点数为1024点，则可以分辨到1Hz。1024Hz的采样率采样1024点，刚好是1秒，也就是说，采样1秒时间的信号并做FFT，则结果可以分析到1Hz，如果采样2秒时间的信号并做FFT，则结果可以分析到0.5Hz。如果要提高频率分辨力，则必须增加采样点数，也即采样时间。频率分辨率和采样时间是倒数关系。

总结：采样点N的不同，分辨率不同，在N=8，俩函数不同的部分没有分辨出来，因此特性相同，而N=16时已经足以分辨出不同的频率，因此频谱不同。 (2) 如果周期信号的周期预先不知道，如何用FFT进行谱分析?

答：设一个定长的m值,先取2m,看2m/m的误差是否大，如大的话再取4m,看4m/2m的误差是否大，如不大，4m(4倍的m值)则可近似原来点的谱分析。 八、总结 本实验主要掌握fft函数的用法，包括采样点的意义。这个意义已经在七、回答问题中描述，不再赘述。

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

实验三： 用窗函数法设计FIR数字滤波器

一、实验目的

(1) 掌握用窗函数法设计FIR数字滤波器的原理和方法。 (2) 熟悉线性相位FIR数字滤波器特性。 (3) 了解各种窗函数对滤波特性的影响。 二、实验内容及步骤

(1) 复习用窗函数法设计FIR数字滤波器一节内容， 阅读本实验原理， 掌握设计步骤。

(2) 编写程序。

① 编写能产生矩型窗、 升余弦窗、 改进升余弦窗和二阶升余弦窗的窗函数子程序。

② 编写主程序。 其中幅度特性要求用dB表示。

窗函数法设计滤波器主程序框图

三、实验内容 j a

e, c j

Hd(e)

c 0

N 1

a

2 1 1 c j a j nj j n

hd(n) H(e)e eed d c2 2 sin c(n a) (n a)

四、思考题

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

(2) 如果要求用窗函数法设计带通滤波器， 且给定上、 下边带截止频率为ω1

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

和ω2，试求理想带通的单位脉冲响应hd(n)。 五、参考资料

数字信号处理教程——Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年

数字信号处理 第三版，高西泉等著，西安电子科技大学出版社，第三版 百度知道

六、实验程序及结果

myfreqz(b,a)函数，是自编的改进freqz()的函数。

function [db,mag,pha,grd,w] =myfreqz(b,a); %advanced freqz

%[db,mag,pha,grd,w] =myfreqz(b,a); [H,w]=freqz(b,a,1000,'whole'); H=(H(1:1:501))';w=(w(1:1:501))'; mag=abs(H);

db=20*log10((mag+eps)/max(mag)); pha=angle(H);

grd=grpdelay(b,a,w); 以下是主程序部分

wp=0.2*pi;ws=0.3*pi;deltaw=ws-wp; N0=ceil(6.6*pi/deltaw); N=N0+mod(N0+1,2); wdham=(hamming(N))'; wc=(ws+wp)/2;

% hd=ideallp(wc,N); tao=(N-1)/2; n=[0:(N-1)]; m=n-tao+eps;

hd=sin(wc*m)./(pi*m); h=hd.*wdham;

[db,mag,pha,grd,w]=myfreqz(h,[1]); dw=2*pi/1000;

Rp=-(min(db(1:wp/dw+1)));

As=-round(max(db(ws/dw+1:501))); figure;

subplot(2,2,1);

stem(0:length(hd)-1,hd,'.'); subplot(2,2,2);

stem(0:length(wdham)-1,wdham,'.'); subplot(2,2,3);

stem(0:length(h)-1,h,'.'); subplot(2,2,4);

plot(1/length(db):1/length(db):1,db);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

wn1=boxcar(N); figure;

subplot(2,2,1); stem(0:N-1,wn1,'.'); title('boxcar'); wn2=bartlett(N); subplot(2,2,2);

stem(0:N-1,wn2,'.');title('bartlett'); wn3=hanning(N); subplot(2,2,3);

stem(0:N-1,wn3,'.');title('hanning'); % wn4=hamming(N); wn4=blackman(N); subplot(2,2,4);

stem(0:N-1,wn4,'.');title('blackman'); % wn=(N,beta);

得到的结果如下图：

其中，图一是理想脉冲响应，第一行图二是汉明窗；图三是实际脉冲响应，图四是幅度响应。一下是几个典型窗函数的图像。，汉明窗上例已画出，不再画。每个窗的名字已在每个子图上方标出。

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

典型窗函数

七、回答问题

(1) 给定通带截止频率和阻带截止频率以及阻带最小衰减，用窗函数法设计线性相位低通滤波器的设计步骤： 答：

技术指标 Wp=0.2*pi,Ws=0.4*pi,Ap=0.25dB,As=50dB 方法一 选择海明窗 clear all; Wp=0.2*pi; Ws=0.4*pi;

tr_wide=Ws-Wp; %过渡带宽度

N=ceil(6.6*pi/tr_wide)+1; %滤波器长度 n=0:1:N-1;

Wc=(Wp+Ws)/2; %理想低通滤波器的截止频率 hd=ideal_lp1(Wc,N); %理想滤波器的单位冲击响应 w_ham=(hamming(N))'; %海明窗

h=hd.*w_ham; %实际海明窗的响应

[db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;

Ap=-(min(db(1:1:Wp/delta_w+1))) %实际通带纹波

As=-round(max(db(Ws/delta_w+1:1:501))) %实际阻带纹波 subplot(221) stem(n,hd)

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_ham) title('海明窗') subplot(223) stem(n,h)

title('实际单位脉冲响应hd(n)') subplot(224) plot(wi/pi,db)

title('幅度响应(dB)') axis([0,1,-100,10])

方法二

Window=blackman(16);

b=fir1(15,0.3*pi ,'low',Window); freqz(b,128)

(2)用窗函数法设计带通滤波器， 且给定上、 下边带截止频率为ω1和ω2，理想带通的单位脉冲响应hd(n)的求解过程：

答：由hd(n)=1/2π∫w1woe-jwaejnwdw+1/2π∫w2w3e-jwaejnwdw （其中w0=-w0-wc,w1=-w0+wc,w2=w0-wc,w3=w0+wc） 计算整理后可得：

hd（n）=2/((n-a)*π)*sin[(n-a)wc]*cos[(n-a)w0] =2wc/π* sa[(n-a)wc]*cos[(n-a)w0]

包括1~4 4个程序,分别是 信号、 系统及系统响应,用FFT作谱分析,用窗函数法设计FIR数字滤波器,用双线性变换法设计IIR数字滤波器。本人参考了陈老师的书籍,最后一个程序参考了学长的代码。来自11级卓越

实验四：用双线性变换法设计IIR数字滤波器

一、实验目的

(1) 掌握用双线性变换法设计IIR数字滤波器的原理和方法。 (2) 熟悉IIR数字滤波器特性。 二、实验内容及步骤

(1) 复习用双线性变换法设计IIR数字滤波器一节内容， 阅读本实验原理， 掌握设计步骤。 (2) 编写程序。

① 设计产生巴特沃斯滤波器原型，采用双线性变换法设计巴特沃斯数字低通滤波器。

② 编写主程序。 其中幅度特性要求用dB表示。 三、上机实验内容

(1) 实验所需设计的IIR数字低通滤波器的指标为wp=0.2*pi，ws=0.3*pi，Rp=1dB，As=15dB，滤波器采样频率为100Hz。 (2) 采用巴特沃斯滤波器为原型进行设计。

(3) 将相同指标要求的IIR滤波器和FIR滤波器进行比较。 四、思考题

(1) 如果给定数字滤波器的通带截止频率和阻带截止频率以及阻带最小衰减， 如何用双线性变换法设计IIR低通滤波器? 写出设计步骤。 (2) 将IIR滤波器和FIR滤波器的比较，各自具有什么特点？ 五、参考资料

数字信号处理教程——Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年

数字信号处理 第三版，高西泉等著，西安电子科技大学出版社，第三版 六、实验程序及结果

求系数部分

wp=0.2*pi;ws=0.3*pi; Rp=1;As=15; T=0.01;Fs=1/T;

OmegaP=(2/T)*tan(wp/2); OmegaS=(2/T)*tan(ws/2); %

[N,OmegaC]=buttord(OmegaP,OmegaS,Rp,As,'s') [z0,p0,k0]=buttap(N);

p =p0*OmegaC;z=z0*OmegaC;%O0O0 k=k0*OmegaC^N; ba=k*real(poly(z)); aa=real(poly(p));

[bd,ad]=bilinear(ba,aa,Fs); [sos,G]=tf2sos(bd,ad)

得到如下截图所示结果

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