试验四IIR数字滤波器设计试验目的掌握用双线性变换法设计.docx

试验四IIR数字滤波器设计试验目的掌握用双线性变换法设计.docx

《试验四IIR数字滤波器设计试验目的掌握用双线性变换法设计.docx》由会员分享，可在线阅读，更多相关《试验四IIR数字滤波器设计试验目的掌握用双线性变换法设计.docx（8页珍藏版）》请在冰豆网上搜索。

试验四IIR数字滤波器设计试验目的掌握用双线性变换法设计

实验四IIR数字滤波器设计

一、实验目的

（1）掌握用双线性变换法设计IIR数字低通和高通滤波器。

（2）设计低通滤波器对实际心电图信号进行滤波。

（3）设计低通滤波器对含有啸叫噪声的音乐信号进行消噪。

*（4）设计IIR数字低通和高通滤波器对某个DTMF（双音多频）信号进行频带分离。

二、实验环境

1.Windows98以上操作系统

2.安装MATLAB6.0以上版本

三、实验原理

1．选频型数字滤波器的种类有低通、高通、带通和带阻滤波器。

2.从实现方法上，数字滤波器通常分为IIR和FIR滤波器。

3．IIR滤波器的设计目的是根据技术指标，找到H（z）分子/分母系数b,a；IIR滤波的MATLAB语句为y=filter（b,a,x）;

四、实验内容

1.人体心电图信号在测量过程中往往受到工业高频干扰，必须经过低通滤波处理后才能作为判断心脏功能的有用信息。

给出一实际心电图信号采样序列样本x（n），其中存在高频干扰。

试以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}。

低通滤波器

设计指标：

ωp=0.2πrad,ωs=0.3πrad,ap=1dB,as=15dB。

已设计出H（z）（p300）

其中

A=0.09036；B1=1.2686,C1=-0.7051；B2=1.0106,C2=-0.3583；B3=0.9044,C3=-0.2155

*IIR滤波的Matlab语句：

y=filter（b,a,x）

b,a----Hk（z）分子/分母系数；x---输入信号x（n）；y---滤波结果y（n）。

**求频响特性的Matlab语句：

[H,w]=freqz（b,a,N）

b,a----Hk（z）分子/分母系数；N---频率点数；H—滤波器的频响特性H（w）;

w---数字频率轴w。

（1）用IIR低通滤波器（原理设计）对实际的心电信号进行滤波（IIRECG.m）

%程序：

（IIRECG.m）

clear

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];%心电信号x（n）[含高频噪声]

A=0.09036;%IIR低通滤波器系数

B1=1.2686;C1=-0.7051;

B2=1.0106;C2=-0.3583;

B3=0.9044;C3=-0.2155;

b=[A2*AA];

a1=[1-B1-C1];

N=128;

[H1,w]=freqz（b,a1,N）;%频响特性H1（w）

y1=filter（b,a1,x）;%y1（n）是x（n）经H1（z）滤波的结果

a2=[1-B2-C2];

[H2,w]=freqz（b,a2,N）;%频响特性H2（w）

y2=filter（b,a2,y1）;%y2（n）是y1（n）经H2（z）滤波的结果

a3=[1-B3-C3];

[H3,w]=freqz（b,a3,N）;%频响特性H3（w）

H=H1.*H2.*H3;%总的频响特性H（w）=H1（w）H2（w）H3（w）

mag=abs（H）;%滤波器的幅频特性

db=20*log10（（mag+eps）/max（mag））;%幅频特性（dB）

y3=filter（b,a3,y2）;%y3（n）是y2（n）经H3（z）滤波的结果

X=fft（x,N）;%对x（n）做N点fft,得其频谱X

wx=2*pi*（0:

N/2-1）/N;%将坐标轴从频率点k转换为数字频率wx（wx=2*pi*k/N）

%%%绘图%%%

subplot（221）;plot（x）;gridon;title（'x（n）'）;

subplot（222）;plot（wx/pi,abs（X（1:

N/2）））;gridon;title（'|X（wx）|'）;

subplot（223）;plot（w/pi,db）;gridon;title（'|H（w）|（db）'）;

subplot（224）;plot（y3）;gridon;title（'y3（n）'）;

要求：

绘出并分析结果图形。

（2）用现有的Matlab函数设计上述IIR低通滤波器,对实际的心电信号进行滤波，同样绘出题1的结果图形。

Matlab函数:

*根据技术指标ωp,ωs,αp,αs计算巴特沃思滤波器的阶数N和3dB截止频率ωc

语句：

[N,wc]=buttord（wp/pi,ws/pi,ap,as）

注：

数字频率以π为单位。

**根据N,ωc确定数字滤波器H（z）分子/分母多项式的系数b,a

低通[b,a]=butter（N,wc）

高通[b,a]=butter（N,wc,’high’）

%程序IIRECGbt.m

clear

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];%心电信号x（n）[含高频噪声]

wp=0.2*pi;ws=0.3*pi;ap=1;as=15;%低通滤波器技术指标

N=128;

[N1,wc]=buttord（wp/pi,ws/pi,ap,as）%根据技术指标计算巴特沃思滤波器的阶数N和3dB截止频率ωc

[b,a]=butter（N1,wc）%根据N,ωc确定数字滤波器H（z）分子/分母多项式的系数b,a

[H,w]=freqz（b,a,N）;%频响特性H（w）

y3=filter（b,a,x）;%y3（n）是a经H（z）滤波的结果

mag=abs（H）;%滤波器的幅频特性

db=20*log10（（mag+eps）/max（mag））;%幅频特性（dB）

X=fft（x,N）;%对x做N点的傅里叶变换得到X

wx=2*pi*（0:

N/2-1）/N;%将坐标轴从频率点k转换为数字频率wx（wx=2*pi*k/N）

%%%绘图%%%

subplot（221）;plot（x）;gridon;title（'x（n）'）;

subplot（222）;plot（wx/pi,abs（X（1:

N/2）））;gridon;title（'|X（wx）|'）;

subplot（223）;plot（w/pi,db）;gridon;title（'|H（w）|（db）'）;

subplot（224）;plot（y3）;gridon;title（'y3（n）'）;

1）在空格中填写程序注释。

2）运行上述程序，得：

滤波器阶数N=6；3dB截止频率wc=0.2329；

b=0.00070.00440.01110.01480.01110.00440.0007；a=1.0000-3.18364.6222-3.77951.8136-0.48000.0544。

3）修改as=50,得：

N=15；wc=0.2127；

b=1.0e-004*

0.00000.00070.00510.02200.06590.14510.24180.31080.31080.24180.14510.06590.02200.00510.00070.0000；a=1.0000-8.610335.4878-92.5611170.4274-234.1787247.6705-205.0244133.7846-68.746527.5683-8.46621.9261-0.30630.0304-0.0014。

2.文件”yinyue.wav”是含有啸叫噪声的音乐信号，首先对其进行谱分析，观察信号和噪声的频带范围，再设计适当的IIR数字低通滤波器对其进行消噪处理恢复原信号，将结果存入音乐文件”yinyuexiaozao.wav”并用耳机监听消噪前后的音乐信号。

%程序:

IIRxiaozaoy.m

clear

clc

N1=81920;%做fft的点数

y1=wavread（'yinyue.wav'）;%读入含噪音乐信号

sound（y1）%播放y1

pause（10）%暂停10s

y1=y1*6;

Y1=fft（y1,N1）;%对y1（n）做N1点fft,得其频谱Y1

PSD1=Y1.*conj（Y1）/N1;%计算功率谱PSD1

fs=8192;

f=8192/N1*（0:

N1/2-1）;%将频率点转化为实际频率

%%%绘图%%%

figure

（1）

subplot（211）;plot（real（y1））;gridon;

subplot（212）;plot（f,PSD1（1:

N1/2））;axis（[050000500]）;gridon;

实验要求：

（1）运行上述程序，绘出结果波形并监听含噪音乐信号。

（2）从第二张子图可见音乐信号的频带范围在340-680Hz,啸叫噪声频率为4000Hz。

（3）根据%后的要求编写程序，将这些程序连在上述程序之后。

%%IIR低通滤波器%%

fp=1000;%通带截止频率fp=1000Hz

fs1=1200;%阻带截止频率fs1=1200Hz

wp=2*pi*fp/fs;%化为数字频率wp和ws

ws=2*pi*fs1/fs;

Rp=1;As=50;%Rp=1dB,As=50dB

[N,wc]=buttord（wp/pi,ws/pi,Rp,As）;

%根据技术指标计算巴特沃思滤波器的阶数N和3dB截止频率ωc

[b,a]=butter（N,wc）;

%根据N,ωc确定数字滤波器H（z）分子/分母多项式的系数b,a

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

[H,w]=freqz（b,a,fs,'whole'）;%频响特性H（w）

mag=abs（H）;%幅频特性

db=20*log10（（mag+eps）/max（mag））;%幅频特性（以dB为单位）

y=filter（b,a,y1）;%用所设计的滤波器对y1进行低通滤波消噪,得y

%%%%绘图%%%%

figure

（2）

subplot（211）;plot（w*fs/（2*pi）,db）;axis（[0fs/2