数字信号处理实验1序列的离散傅立叶变换 2 用双线性变换法设计IIR数字滤波器Word文档格式.docx
《数字信号处理实验1序列的离散傅立叶变换 2 用双线性变换法设计IIR数字滤波器Word文档格式.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验1序列的离散傅立叶变换 2 用双线性变换法设计IIR数字滤波器Word文档格式.docx（13页珍藏版）》请在冰点文库上搜索。
窗函数的选择原则
1.主瓣带宽要尽可能小，以获得较陡的过渡带
2.与主瓣的幅度相比，旁瓣应尽可能的小，以减少通带内、阻带内波动的最大振幅
常用的窗函数
窗口函数
主瓣宽度
第一旁瓣比主瓣低（dB）
矩形窗
4π/N
13
三角形窗
8π/N
27
汉宁窗
32
哈明窗
43
布莱克曼窗
12π/N
58
实验截图：
分析：
分辨率下降：
频谱泄露：
若信号x（n）具有能量较强的频率分量，那么由它所产生的旁瓣可能会掩盖能量较弱的频率分量
如上实验矩形窗显示的幅频响应显示能量集中在主瓣但是由于频谱泄漏使得旁瓣0.4~0.6π之间的频率分量被掩盖；
而汉宁窗弱化了频谱泄漏使得在0.4~0.6π之间的频率表现出来，所以汉宁窗势必矩形窗更优的窗函数；
从窗函数的采样点对比来说，采样点数越高，其说明幅频响应分辨率越高上图尖峰的尖细程度可以说明
源代码：
n1=0:
1:
511;
n2=0:
2047;
xn1=0.001*cos（0.45*pi*n1）+sin（0.3*pi*n1）-cos（0.302*pi*n1-pi/4）;
%512点DFT原始信号方程
xn2=0.001*cos（0.45*pi*n2）+sin（0.3*pi*n2）-cos（0.302*pi*n2-pi/4）;
%2048点DFT原始信号方程
subplot（4,2,1）
stem（n2,xn2,'
.'
）
title（'
（2048点原始DFT波形图）'
）;
xlabel（'
n'
ylabel（'
x（n）'
subplot（4,2,2）
stem（n1,xn1,'
（512点原始DFT波形图）'
%boxcarwindow
Xk1=fft（xn1,512）;
Xk1=20*log10（abs（Xk1））;
%分贝表示
wk=0:
wk=2*wk/512;
%关于π的归一化坐标
subplot（4,2,3）
plot（wk,Xk1）;
grid;
n=512boccarwindow幅频响应波形'
|H（ejw）|'
w/pi'
subplot（4,2,4）
（加512矩形窗后下x（n）DFT波形）'
%hanningwindow（0：
511）
N=512;
wn=hanning（N）;
xn1=xn1.*wn'
;
subplot（4,2,5）
n=512hanningwindow幅频响应波形'
subplot（4,2,6）
（加512汉林窗下x（n）DFT波形'
2047）
N=2048;
xn2=xn2.*wn'
Xk1=fft（xn2,2048）;
wk=2*wk/2048;
%关于π归一化坐标
subplot（4,2,7）
n=2048hanningwindow幅频响应波形'
subplot（4,2,8）
加2048汉林窗下x（n）DFT波形'
实验二用双线性变换法设计IIR数字滤波器
实验目的
1.熟悉用双线性变换法设计IIR数字滤波器的原理与方法
2.掌握IIR滤波器的计算机仿真方法
3.通过观察对三路抑制载波调幅信号的滤波作用，获得数字滤波器的感性知识
实验内容
1.产生三路抑制载波调幅复合信号，绘出时域波形和幅频特性曲线：
三路的指标如下：
第一路：
载波频率fc=250Hz，单频信号频率fo=25Hz;
第二路：
载波频率fc=500Hz，单频信号频率fo=100Hz;
第三路：
载波频率fc=1000Hz，单频信号频率fo=200Hz;
采样频率Fs=1/Ts=10kHz.信号长度N=1600.
2.用双线性变换法设计一个椭圆低通IIR滤波器，滤波器设计指标参数为：
通带截止频率fp=280Hz,通带最大波动Rp=0.1dB
阻带截止频率fs=450Hz,阻带最小衰减As=60dB
采样频率Fs=1/Ts=10kHz.信号长度N=1600.
3.以0.05π为采样间隔，绘出数字滤波器在频率区间[0,π]上的幅值响应曲线。
4.应用所设计的椭圆低通IIR滤波器滤波分离出第一路抑制载波调幅信号。
绘出滤波器分离出的第一路抑制载波调幅信号的时域波形和幅频特性曲线。
观察总结滤波作用与效果。
5.实验步骤
1．参见教材实验四
2．用MATLAB语言编写仿真滤波程序，完成实验内容1~4
3．写出完整的实验报告
有波形图可知本实验设计的双线性数字滤波器能有效地的将分离出第一路波形，并且有效地抑制其他两路波形；
从波形图上可以看出滤波的效果是滤波前后的波形形状未发生变化，但是时域波形上滤波后信号发生了延迟，表现在频域上就是信号的整个相位发生的变化。
信号发生函数：
functionst=mstg
N=1600;
Fs=10000;
T=1/Fs;
Tp=N*T;
t=0:
T:
（N-1）*T;
k=0:
N-1;
f=k/Tp;
%采样频率10Khz，Tp为采样时间
fc1=Fs/40;
fm1=fc1/10;
fc2=Fs/20;
fm2=fc2/5;
fc3=Fs/10;
fm3=fc3/5;
xt1=cos（2*pi*fm1*t）.*cos（2*pi*fc1*t）;
xt2=cos（2*pi*fm2*t）.*cos（2*pi*fc2*t）;
xt3=cos（2*pi*fm3*t）.*cos（2*pi*fc3*t）;
fxt1=fft（xt1）;
st=xt1+xt2+xt3;
fxt=fft（st,N）;
plot（t,st）;
s/t'
s（t）'
axis（[0,Tp/8,min（st）,max（st）]）;
复合信号s（t）的波形'
stem（f,abs（fxt）/max（abs（fxt））,'
（复合信号s（t）的幅频响应波形'
axis（[0,Fs/5,0,1.2]）;
f/Hz'
幅度'
plot（t,xt1）;
xt1/t'
xt1（t）'
axis（[0,Tp/8,min（xt1）,max（xt1）]）;
参考波xt1（t）的波形'
滤波函数：
function[Bz,Az]=myfft
fp=280;
Rp=0.1;
fs=450;
As=60;
%模拟频率最高输入频率280Hz
wpz=2*fp/Fs;
wsz=2*fs/Fs;
%关于π归一化转换成数字指标wp=2pi*fp*T
Wp=2*Fs*tan（wpz*pi/2）;
Ws=2*Fs*tan（wsz*pi/2）;
%预畸校正确定相应模拟滤波器指标
[N,Wc]=ellipord（Wp,Ws,Rp,As,'
s'
%确定椭圆滤波器阶数N和截止频率Wc
[B,A]=ellip（N,Rp,As,Wc,'
%根据阶数和截止频率确定滤波器分子A,分母B系数
[Bz,Az]=bilinear（B,A,Fs）;
%有模拟滤波器和信号采样频率Fs利用双线性变换为数字滤波器的系数
[h,w]=freqz（Bz,Az,20,Fs）;
%以0.05π为采样间隔，绘出数字滤波器在频率区间[0,π]上的幅值响应曲线20分单位圆。
w=2*w/Fs;
%关于π归一化
Hk=20*log10（abs（h））;
plot（w,Hk）;
n=5EllipticBandpassFilter'
归一化频率w/π'
幅度/dB'
axis（[0,1,-100,100]）;
主函数：
%采样频率10Khz，Tp为采样时间
st=mstg;
%载入数字信号
[Bz,Az]=myfft;
%载入椭圆低通滤波器
xn=filter（Bz,Az,st）;
%利用生成的双线性椭圆滤波器对st行进滤波得到输出滤波数字信号点
Y=fft（xn,1600）;
%对滤波后信号进行频率响应求解
plot（t,xn）;
xn/t'
axis（[0,Tp/8,min（xn）,max（xn）]）;
分离出的第一路抑制载波调幅信号x（n）的时域波形'
subplot（4,2,8）%滤波信号的幅频响应图
stem（f,abs（Y）/max（abs（Y））,'
分离出的第一路抑制载波调幅信号幅频响应波形'
Matlab知识点备注：
IIR滤波器设计的常用函数
1.ellipord
功能：
椭圆滤波器阶的选择
格式：
[N,wc]=ellipord（wp,ws,rp,rs）----数字域
[N,wc]=ellipord（wp,ws,rp,rs,’s’）----模拟域
说明：
ellipord可在给定滤波器性能的情况下，选择模拟或数字椭圆滤波器的最小阶，其中wp和ws分别是通带和阻带的截止频率，其值0≤wp（ws）≤1，当其值为1时，表示0.5fs，rp,rs分别是通带和阻带区的波纹系数。
2.ellip
椭圆低通滤波器设计
[B,A]=ellip（N,rp,rs,wc）;
设计出截止频率为wc的N阶低通椭圆滤波器，其滤波器为
截止频率是滤波器幅度下降至处的频率。
3.freqz
[h,w]=freqz（B,A,N,Fs,'
whole'
B：
系统函数分子系数
A：
系统函数分母系数
N：
整个单位圆上分为多少个点用以计算频率响应
Fs：
采样频率
‘whole’：
在整个单位圆上计算频率响应，对应的角频率范围是[0,2π]，不加whole只在上半圆计算，对应角频率范围是[0,π]
3.freqz（）--matlab函数
MATLAB提供了专门用于求离散系统频响特性的函数freqz（），调用freqz（）的格式有以下两种：
（1）[H,w]=freqz（B,A,N）
（2）[H,w]=freqz（B,A,N,’whole’）
（1）中B和A分别为离散系统的系统函数分子、分母多项式的系数向量，返回量H则包含了离散系统频响在0~pi范围内N个频率等分点的值（其中N为正整数），w则包含了范围内N个频率等分点。
调用默认的N时，其值是512。
（2）中调用格式将计算离散系统在0~pi范内的N个频率等分店的频率响应的值。
因此，可以先调用freqz（）函数计算系统的频率响应，然后利用abs（）和angle（）函数及plot（）函数，即可绘制出系统在或范围内的频响曲线。
绘制如下系统的频响曲线
H（z）=（z-0.5）/z
MATLAB命令如下：
B=[1-0.5];
A=[10];
[H,w]=freqz（B,A,400,'
Hf=abs（H）;
Hx=angle（H）;
clf
figure
（1）
plot（w,Hf）
离散系统幅频特性曲线'
figure
（2）
plot（w,Hx）
离散系统相频特性曲线'
4filter是一维数字滤波器
使用方法:
Y=filter（B,A,X），输入X为滤波前序列，Y为滤波结果序列，B/A提供滤波器系数，B为分子，A为分母
整个滤波过程是通过下面差分方程实现的：
a
（1）*y（n）=b
（1）*x（n）+b
（2）*x（n-1）+...+b（nb+1）*x（n-nb）-a
（2）*y（n-1）-...-a（na+1）*y（n-na）
[Y,Zf]=filter（B,A,X,Zi），输入X为滤波前序列，Y为滤波结果序列，B/A提供滤波器系数，B为分子，A为分母，
《Simulink与信号处理》
并输入Zi指定X的初始状态，Zf为最终状态矢量《Simulink与信号处理》
filter（B,A,X,[],DIM）或filter（B,A,X,Zi,DIM）指定X的维数DIM进行操作
举例
k=-100:
100;
>
uk=[zeros（1,100）,ones（1,101）];
a=[1-56];
b=[2-1];
x=uk;
y=filter（b,a,x）
5.%"
'
"
运算符:
在矩阵或向量之后表示复共轭转置
%"
."
运算符：
设x和y是n*n的矩阵。
x*y是按照矩阵乘法；
%x.*y.是对应元素相乘（x（i,j）*y（i,j））。
%x^2是x*x;
x.^2是每个元素的平方。
升级会员 