7.1给定一理想低通FIR滤波器的频率特性
1,4jHde
0,4现希望用窗函数法设计该滤波器,要求具有线性相位。假定滤波器系数的长度为29
点,即M/2=14。
试计算并打印滤波器的系数、幅频响应及相频响应。滤波器系数的计算先用手算,然后再调用有关的MATLAB文件来计算。
解:线性相位理想低通滤波器的频率特性:
jwe,4jHde
0,4
1hdn2ejejn1d2ejnde2jj1ejn2jn441enj4n4n
sinn4,N114n2(1) 矩形窗:
1,0n28 窗函数:wrn0,其它sinn144,0n28 h1nhdnwrnn140,其它n
系数算法一:
clc clear
for n=0:1:28
x=sin(0.25*pi*(n-14))/(pi*(n-14)) n=n+1 end
系数算法二: clc clear n=0:1:28
b=sin(0.25*pi*(n-14))./(pi*(n-14))
求得滤波器系数为:h(0)——h(28)=
-0.0227 -0.0173 0.0000 0.0205 0.0318 0.0250 -0.0000 -0.0322 -0.0531 -0.0450 0.0000 0.0750 0.1592 0.2251 0.2500 0.2251 0.1592 0.0750 0.0000 -0.0450 -0.0531 -0.0322 -0.0000 0.0250 0.0318 0.0205 0.0000 -0.0173 -0.0227 (2) 汉明(Hamming)窗:
2n0.0.46cos,n0,1,,28
窗函数:whnN0,其它nsinn1440.0.46cos2n,0n28 h2nhdnwhnn14N0,其它n
Matlab程序:
clc clear n=0:1:28
b=sin(0.25*pi*(n-14))./(pi*(n-14)).*(0.-0.46*cos(2*pi*n/29))
求得滤波器系数为:h(0)——h(28)=
-0.0018 -0.0016 0.0000 0.0036 0.0077 0.0081 -0.0000 -0.0166 -0.0326 -0.0320 0.0000 0.0656 0.1487 0.2197 0.2493 0.2245 0.1553 0.0701 0.0000 -0.0359 -0.0377 -0.0198 -0.0000 0.0104 0.0103 0.0050 0.0000 -0.0021 -0.0021 (3) 布莱克曼(Blackman)窗:
2n4n0.420.5cos0.08cos,n0,1,,28
窗函数:wbnNN0,其它nh3nhdnwbnsinn1440.420.5cos2n0.08cos4n,0n28 n14NN0,其它n
Matlab程序:
clc clear n=0:1:28
b=sin(0.25*pi*(n-14))./(pi*(n-14)).*(0.42-0.5*cos(2*pi*n/29)+0.08*cos(4*pi*n/29))
求得滤波器系数为:h(0)——h(28)=
0.0000 -0.0001 0.0000 0.0009 0.0027 0.0035 -0.0000 -0.0101 -0.0226 -0.0246 0.0000 0.0591 0.1410 0.2155 0.2488 0.2240 0.1524 0.0665 0.0000 -0.0302 -0.0290 -0.0137 -0.0000 0.00 0.0045 0.0017 0.0000 -0.0003 -0.0001 (4) Matlab代码及运行结果:
clc; clear; N=28;
br=0.25*sinc((-14:1:14)/4);%计算滤波器系数
bb1=fir1(N,0.25,boxcar(N+1));%用fir1计算滤波器系数 [H1,wr]=freqz(br,1,100,2);%频域变换 [Hb1,wr]=freqz(bb1,1,100,2); subplot(331);
plot(wr,abs(Hb1),'-r',wr,abs(H1),'-b');grid;%画出幅频特性 title('29点矩形窗') subplot(334);
plot(wr,unwrap(angle(Hb1)),'-r',wr,unwrap(angle(H1)),'-b');grid;%画出相频特性 subplot(337);
plot(wr,20*log10(abs(Hb1)/sum(Hb1)),'-r',wr,20*log10(abs(H1)/sum(H1)),'-b');grid;
title('对数幅频特性')
bh=br.*(0.-0.46*cos(2*pi*(0:1:28)/29)); bb2=fir1(N,0.25,hamming(N+1)); [H2,wh]=freqz(bh,1,100,2); [Hb2,wh]=freqz(bb2,1,100,2); subplot(332);
plot(wh,abs(Hb2),'-r',wh,abs(H2),'-b');grid; title('29点汉明窗')
subplot(335);
plot(wh,unwrap(angle(Hb2)),'-r',wh,unwrap(angle(H2)),'-b');grid; subplot(338);
plot(wh,20*log10(abs(Hb2)/sum(Hb2)),'-r',wh,20*log10(abs(H2)/sum(H2)),'-b');grid;
title('对数幅频特性')
bb=br.*(0.42-0.5*cos(2*pi*(0:1:28)/29)+0.08*cos(4*pi*(0:1:28)/29)); [H3,wb]=freqz(bb,1,100,2); subplot(333);
plot(wb,abs(Hb2),'-r',wb,abs(H3),'-b');grid; title('29点布莱克曼窗') subplot(336);
plot(wb,unwrap(angle(Hb2)),'-r',wb,unwrap(angle(H3)),'-b');grid; subplot(339);
plot(wh,20*log10(abs(Hb2)/sum(Hb2)),'-r',wb,20*log10(abs(H3)/sum(H3)),'-b');grid;
title('对数幅频特性')
29点矩形窗1.511.5129点汉明窗1.5129点布莱克曼窗0.50.50.5000.20.40.60.81000.20.40.60.81000.20.40.60.810-50-5-100-10-20-3000.20.40.60.81-4000.20.40.60.81-10-1500.20.40.60.81-20-15对数幅频特性0-20-40-60-8000.20.40.60.81-1000-50对数幅频特性0-50对数幅频特性-100-15000.20.40.60.81-15000.20.40.60.81
—Matlab自动设计 —窗函数设计
图1 长度为29点时滤波器的频率特性
分析:
由图可知,矩形窗设计的滤波器,幅频接近Matlab自动设计的滤波器的幅频特性,相频特性曲线和对数幅频特性曲线几乎与之重合。汉明窗设计的滤波器,幅频特性曲线与Matlab自动设计的滤
波器的幅频特性重合,相频特性曲线和对数幅频特性曲线接近Matlab自动设计的滤波器的特性曲线。对于Hamming窗和Blackman窗来说,当长度为29点时,滤波器的相频特性不理想。
(5) 将滤波器长度增大为41点,重新用这三种窗函数设计滤波器:
1hdn2ejejn1ejn2jn441jnded2 sinn4,N120n21,0n40矩形窗:wrn
0,其它sinn204,0n40 h1nhdnwrnn200,其它n2n0.0.46cos,n0,1,,40
汉明(Hamming)窗:whnN0,其它nsinn2040.0.46cos2n,0n40 h2nhdnwhnn20N0,其它n2n4n0.420.5cos0.08cos,n0,1,,40
布莱克曼(Blackman)窗:wbnNN0,其它nh3nhdnwbnsinn2040.420.5cos2n0.08cos4n,0n40 n20NN0,其它nMatlab代码及运算结果: clc; clear; N=40;
br=0.25*sinc((-20:1:20)/4);%计算滤波器系数
bb1=fir1(N,0.25,boxcar(N+1));%用fir1计算滤波器系数
[H1,wr]=freqz(br,1,100,2);%频域变换 [Hb1,wr]=freqz(bb1,1,100,2); subplot(331);
plot(wr,abs(Hb1),'-r',wr,abs(H1),'-b');grid;%画出幅频特性 title('41点矩形窗') subplot(334);
plot(wr,unwrap(angle(Hb1)),'-r',wr,unwrap(angle(H1)),'-b');grid;%画出相频特性 subplot(337);
plot(wr,20*log10(abs(Hb1)/sum(Hb1)),'-r',wr,20*log10(abs(H1)/sum(H1)),'-b');grid; title('对数幅频特性')
bh=br.*(0.-0.46*cos(2*pi*(0:1:40)/41)); bb2=fir1(N,0.25,hamming(N+1)); [H2,wh]=freqz(bh,1,100,2); [Hb2,wh]=freqz(bb2,1,100,2); subplot(332);
plot(wh,abs(Hb2),'-r',wh,abs(H2),'-b');grid; title('41点汉明窗') subplot(335);
plot(wh,unwrap(angle(Hb2)),'-r',wh,unwrap(angle(H2)),'-b');grid; subplot(338);
plot(wh,20*log10(abs(Hb2)/sum(Hb2)),'-r',wh,20*log10(abs(H2)/sum(H2)),'-b');grid; title('对数幅频特性')
bb=br.*(0.42-0.5*cos(2*pi*(0:1:40)/41)+0.08*cos(4*pi*(0:1:40)/41)); [H3,wb]=freqz(bb,1,100,2); subplot(333);
plot(wb,abs(Hb2),'-r',wb,abs(H3),'-b');grid; title('41点布莱克曼窗') subplot(336);
plot(wb,unwrap(angle(Hb2)),'-r',wb,unwrap(angle(H3)),'-b');grid; subplot(339);
plot(wh,20*log10(abs(Hb2)/sum(Hb2)),'-r',wb,20*log10(abs(H3)/sum(H3)),'-b');grid; title('对数幅频特性')
41点矩形窗1.511.5141点汉明窗141点布莱克曼窗0.50.50.5000.20.40.60.81000.20.40.60.81000.20.40.60.810-5-10-15-2000.20.40.60.810-100-10-20-20-3000.20.40.60.81-3000.20.40.60.81对数幅频特性0-20-40-60-8000.20.40.60.81-10000.2-500对数幅频特性0-50对数幅频特性-1000.40.60.81-15000.20.40.60.81
—Matlab自动设计 —窗函数设计
图表 1长度为41点时滤波器的频率特性
与图1对比,当N=41时,滤波器的幅频特性的过度带变得更窄,相频特性有所改善。还可看出,对于矩形窗来说,N的大小对于其幅频谱的过冲大小没有影响。
7.2 给定一理想带阻滤波器的频率特性
1,6Hdej0,
361,3重复7.1题各项要求。
解:
为了使设计出的滤波器hd(n)的对称中心在M/2,并保证滤波器具有线性相位,我们按照下
面的方法规定滤波器的相位。
Mexp(-j),,Hdej263
0,其他对Hdej做IDTFT有
1hd(n)2Hd(ej)ejnd
MM12nn1/6j1j22ededed/6/322MMM sinnsinnsinn62232Mn2/3Mjn2对该序列加一个长度M=29的矩形窗,就得到了按矩形窗设计的FIR滤波器
sinhRnMMMnsinnsinn62232Mn2,n0,1,,28
clear all; close all; n=0:1:28
br=sin(pi/6*(n-14))./(pi*(n-14))+sin(pi*(n-14))./(pi*(n-14))-sin(pi/3*(n-14))./(pi*(n-14)) br(15)=0.833 stem(n,br)
求得滤波器系数为:h(0)——h(28)=
-0.0000 -0.0090 -0.0000 0.0106 -0.0000 -0.03 -0.06 -0.0621 -0.0000 0.0870 0.1378 0.1061 -0.0000 -0.1165 0.8330 -0.1165 -0.0000 0.1061 0.1378 0.0870 -0.0000 -0.0621 -0.06 -0.03 -0.0000 0.0106 -0.0000 -0.0090 -0.0000
若施加一个长度M=29的Hamming窗,就得到了按Hamming窗设计的FIR滤波器
sinhHnMn62MMsinnsinn232Mn20.0.46cos2n,n0,1,,28Mclear all; close all; n=0:1:28
br=sin(pi/6*(n-14))./(pi*(n-14))+sin(pi*(n-14))./(pi*(n-14))-sin(pi/3*(n-14))./(pi*(n-14)) br(15)=0.833
bh=br.*(0.-0.46*cos(2*pi*(0:1:28)/29)); stem(n,bh)
求得滤波器系数为:h(0)——h(28)=
-0.0000 -0.0090 -0.0000 0.0106 -0.0000 -0.03 -0.06 -0.0621 -0.0000 0.0870 0.1378 0.1061 -0.0000 -0.1165 0.8330 -0.1165 -0.0000 0.1061 0.1378 0.0870 -0.0000 -0.0621 -0.06 -0.03 -0.0000 0.0106 -0.0000 -0.0090 -0.0000
clear all; close all; n=0:1:28;
br=sin(pi/6*(n-14))./(pi*(n-14))+sin(pi*(n-14))./(pi*(n-14))-sin(pi/3*(n-14))./(pi*(n-14)) br(15)=0.833;
[H1,wr]=freqz(br,1,100,2);
f=[0 0.16 0.17 0.33 0.34 1]; % 给定频率轴分点; m=[1 1 0 0 1 1];% 给定在这些频率分点上理想的幅频响应 N1=28;
b1=fir2(N1,f,m);
[H2,wr]=freqz(b1,1,100,2)
subplot(321); stem(n,br,'-b') title('29点矩形窗系数')
subplot(323);
plot(wr,abs(H2),'-r',wr,abs(H1),'-b');grid; title('幅频特性') subplot(325);
plot(wr,unwrap(angle(H2)),'-r',wr,unwrap(angle(H1)),'-b');grid; title('相频特性')
bh=br.*(0.-0.46*cos(2*pi*(0:1:28)/29)); subplot(322); stem(n,bh)
title('29点汉明窗系数')
[H3,wr]=freqz(bh,1,100,2); subplot(324);
plot(wr,abs(H2),'-r',wr,abs(H3),'-b');grid; title('幅频特性') subplot(326);
plot(wr,unwrap(angle(H2)),'-r',wr,unwrap(angle(H3)),'-b');grid; title('相频特性')
29点矩形窗系数10.510.529点汉明窗系数00-0.5051015幅频特性202530-0.5051015幅频特性2025301.511.510.50.5000.20.40.60.81000.20.40.60.81相频特性0-200-20相频特性-40-40-6000.20.40.60.81-6000.20.40.60.81
—Matlab自动设计 —窗函数设计
图表 2 长度为29点的滤波器的频率特性
7.3 试用频域抽样法完成7.1及7.2题两个滤波器的设计。
解:
j利用频域抽样法设计数字滤波器时,需要将理想的连续频域相应Hde转换为理想的离
散频域响应Hdk,例如,对低通滤波器,
当N为偶数时,有
ejN1k/N,k0,1,,N/21 Hdk0,kN/2jN1k/N,kN/21,,N1ejN1k/N或令Hdke,k0,1,,N/21 jN1k/N则HdNke,k0,1,,N/21
Hdk0,kN/2
jN1k/N当N为奇数时,有Hdke,k0,1,,N1 jN1k/N,k0,1,,N1/2Hdke或令
HdNkHdk,k1,2,,N1/22jnk1N1根据hnnHdkeN,n0,1,,N1
Nk0即可求出所要设计的滤波器。上述赋值方法保证了所设计的滤波器具有线性相位。
j根据线性相位的条件,本题中N=29,为奇数,并且在02区间内,He关于偶
对称,选择第1种线性相位进行设计。
根据指标要求,在02内有29个采样点,所以第K点对应的频率为几个特殊的频率:位于
2k。 2942273和4之间,同理,位于第25个和第26个采样点之间,29294求得Hk:
ad =
1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1
kN122kN28k,0k28 292jnk1N1hnHkeN,n0,1,,N1Nk01HkeNk0
N1jN12ej2nkN 将H(k)代入,得到h(n):
,n0,1,,N1clc;
clear;
k1=[0:3];k2=[4:25];k3=[26:28]; a1=0*k1+1; a2=0*k2; a3=0*k3+1;
ad(1:29)=[a1,a2,a3] N=29; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) stem(k,h);
title('29点抽样h(n)图像');
h =
-0.0321 -0.0147 0.0114 0.0328 0.0376 0.0215 -0.0097 -0.0415 -0.0562 -0.0405 0.00 0.0823
0.1602 0.2193 0.2414 0.2193 0.1602 0.0823 0.00 -0.0405 -0.0562 -0.0415 -0.0097 0.0215 0.0376 0.0328 0.0114 -0.0147 -0.0321
clc; clear;
k1=[0:3];k2=[4:25];k3=[26:28]; a1=0*k1+1; a2=0*k2; a3=0*k3+1;
ad(1:29)=[a1,a2,a3] N=29; k=0:(N-1);
wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221); stem(k,h);
title('29点抽样h(n)的图像'); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222); plot(p,abs(Hd));
title('29点频率抽样幅频特性'); subplot(223); plot(p,angle(Hd));
title('29点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)));grid on; title('29点频率抽样对数幅频特性')
29点矩形窗2100-10-200-50-1002100-10-200-50-10029点汉明窗10.500-20-400-100-20029点布莱克曼窗00.5100.5100.5100.5对数幅频特性100.5对数幅频特性100.5对数幅频特性100.5100.5100.51
比较两图可以看出,用频域抽样法设计的滤波器的幅频响应没有用窗函数法设计的幅频响
应的曲线平滑,因此,对频域抽样法需要改进,改进的方法是在Hdk由1变为0的频率出加过度点。
插一个点:
插两个点:
29点抽样h(n)的图像0.30.20.10-0.101020300.51.5129点频率抽样幅频特性000.5129点频率抽样相频特性420-2-400.51-1000-5029点频率抽样对数幅频特性-15000.51
结论:差值要注意,不一定越多越好,规律有待总结。 clc; clear;
k1=[0:2];k2=[4:24];k3=[26:28]; a1=0*k1+1; a2=0.4 a3=0*k2; a4=0.4
a5=0*k3+1;
ad(1:29)=[a1,a2,a3,a4,a5] N=29; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221); stem(k,h,'g-');hold
title('29点抽样h(n)的图像'); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222);
plot(p,abs(Hd),'g-');hold
title('29点频率抽样幅频特性'); subplot(223);
plot(p,angle(Hd),'g-');hold
title('29点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)),'g-');grid on;hold title('29点频率抽样对数幅频特性') clc; clear;
k1=[0:2];k2=[5:24];k3=[27:28]; a1=0*k1+1; a2=0.4 a3=0.3 a4=0*k2; a5=0.3 a6=0.4
a7=0*k3+1;
ad(1:29)=[a1,a2,a3,a4,a5,a6,a7] N=29; k=0:(N-1);
wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221);
stem(k,h,'r-');grid on;hold title('29点抽样h(n)的图像'); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222);
plot(p,abs(Hd),'r-');grid on;hold title('29点频率抽样幅频特性'); subplot(223);
plot(p,angle(Hd),'r-');grid on;hold title('29点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)),'r-');grid on;hold title('29点频率抽样对数幅频特性')
clc; clear;
k1=[0:2];k2=[3:4];k3=[5:9];k4=[10:25];k5=[26:28]; a1=0*k1+1; a2=0*k2; a3=0*k3+1; a4=0*k4; a5=0*k5+1;
ad(1:29)=[a1,a2,a3,a4,a5]; N=29; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) stem(k,h);
title('29点抽样h(n)图像');
clc; clear;
k1=[0:2];k2=[3:4];k3=[5:9];k4=[10:25];k5=[26:28]; a1=0*k1+1; a2=0*k2; a3=0*k3+1; a4=0*k4; a5=0*k5+1;
ad(1:29)=[a1,a2,a3,a4,a5]; N=29; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221); stem(k,h);
title('29点抽样h(n)的图像'); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222); plot(p,abs(Hd));
title('29点频率抽样幅频特性'); subplot(223); plot(p,angle(Hd));
title('29点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)));grid on; title('29点频率抽样对数幅频特性')
比较两图可以看出,用频域抽样法设计的滤波器的幅频响应没有用窗函数法设计的幅频响应的曲线平滑,因此,对频域抽样法需要改进,改进的方法是在Hdk由1变为0的频率出加过度点。
7.4一滤波器的理想频率响应如图5所示:
图5
(1) 试用窗函数法设计该滤波器,要求具有线性相位,滤波器长度为33,用汉宁窗。 (2) 用频率抽样法设计,仍要求具有线性相位,滤波器长度为33,过渡点自行设置。 要求先用手算得出hn,然后上机求Hej。
解:线性相位理想滤波器的频率特性:
2j1e,22j1,0e2j Hde
21ej,0221ej,21hdn212212jjn1eed22jjn1ede202012jjn1eed22jjn1ede2
1144cosn8cosn,n222n21,n2N1其中,,当N33时,162(1) 窗函数法
(a) 采用33点汉宁窗设计滤波器:
2n0.50.5cos,n0,1,,32
窗函数:wh1nN0,其它nh1nhdnwh1n144cosn168cosn16211cos2n,n0,1,,15,17,18,,32222332n162 12n0.50.5cos,n162330,其它n
clc clear n=0:1:32
x=(4+4*cos((n-16)*pi)-8*cos(0.5*pi*(n-16)))./(2*pi*pi*(n-16).^2) b=x.*(0.5-0.5*cos(2*pi*n/33))
b(17)=0.5 stem(n,b);grid;
滤波器系数h1(0)——h1(32):
0 0 0.0003 0 0 0 0.0029 0 0 0 0.0163 0 0 0 0.1958 0 0.5000 0 0.1958 0 0 0 0.0163 0 0 0 0.0029 0 0 0 0.0003 0 0
滤波器的频率特性:
clc clear n=0:1:32
x=(4+4*cos((n-16)*pi)-8*cos(0.5*pi*(n-16)))./(2*pi*pi*(n-16).^2) b=x.*(0.5-0.5*cos(2*pi*n/33)) b(17)=0.5 stem(n,b);grid;
[H1,wr]=freqz(b,1,100,2); subplot(121);
plot(wr,abs(H1));grid; title('幅频特性') subplot(122);
plot(wr,unwrap(angle(H1)));grid; title('相频特性')
图6 33点汉宁窗设计的滤波器频率特性
(b) 采用33点矩形窗设计滤波器:
1,0n32 窗函数:wrn0,其它h2nhdnwrn144cosn168cosn162,n0,1,,15,17,18,,3222n162 12,n160,其它n
clc clear n=0:1:32
x=(2+2*cos((n-16)*pi)-4*cos(0.5*pi*(n-16)))./(pi*pi*(n-16).^2)
滤波器系数:h2(0)——h2(32) =
0 0 0.0041 0 0 0 0.0081 0 0 0 0.0225 0 0 0 0.2026 0 0.5000 0 0.2026 0 0 0 0.0225 0 0 0 0.0081 0 0 0 0.0041 0 0
(c)采用33点三角窗
2n,n0,1,,16 窗函数:wtn32w(32n),n17,18,,32h3nhdnwtn144cosn168cosn1622n,n0,1,,152322n1621 ,n1621n168cosn1644cos2232n,n17,18,,322322n1620,其它n滤波器系数:h3(0)——h3(32)=
0 0 0.0007 0 0 0 0.0033 0 0 0 0.0146 0 0 0 0.1788 0 0.5000 0 0.1788 0 0 0 0.0146 0 0 0 0.0033 0 0 0 0.0007 0 0
33点矩形窗设计滤波器10.80.60.90.80.70.60.50.40.200.40.30.200.20.40.60.810.100.20.40.60.8133点三角窗设计滤波器相频特性(wraped)43210-1-2-3-400.20.40.60.8143210-1-2-3-400.2相频特性(wraped)0.40.60.81
图7 33点矩形窗和三角窗设计的滤波器的频率特性
(d)采用65点矩形窗和三角窗设计滤波器 方法与33点相同,得到滤波器的频率特性:
65点矩形窗设计滤波器10.80.610.80.665点三角窗设计滤波器0.40.200.40.2000.20.40.60.8100.20.40.60.81相频特性(wraped)43210-1-2-3-400.20.40.60.8143210-1-2-3-400.2相频特性(wraped)0.40.60.81
图8 65点矩形窗和三角窗设计的滤波器的频率特性
结论:对于本题来说,采用多点矩形窗设计的FIR数字滤波器,其频率特性接近于理想滤波器。 (2) 用频率抽样法设计滤波器 (a)33点抽样
j根据线性相位的条件,本题中N=33,为奇数,并且在02区间内,He关于偶
对称,选择第1种线性相位进行设计。
根据指标要求,在02内有33个采样点,所以第K点对应的频率为几个特殊的频率:
2k。 33228和9之间,同理,位于第16个和第17个采样点之间,位于333323位于第24个和第26个采样点之间。根据对称性,求得Hk: 2 H0= 1.0000 H1=H32=0.8788 H2=H31=0.7576 H3=H30=0.63
H4=H29=0.5152 H5=H28=0.3939 H6=H27=0.2727 H7=H26=0.1515
H8=H25=0.0303 H9=H24=0.0909 H10=H23=0.2121 H11=H22=0.3333
H12=H21=0.45 H13=H20=0.5758 H14=H19=0.6970 H15=H18=0.8182 H16=H17=0.9394
kN122k32N33k,0k32 hn1N1j2nkNHkeN,n0,1,,N1k0N1 将H(k)代入,得到h(n):
1jN12j2NHkeeNnk,n0,1,,N1k0
k1=[0:8];k2=[9:16];k3=[17:24];k4=[25:32]; a1=(-4/33)*k1+1; a2=(4/33)*k2-1; a3=(-4/33)*k3+3; a4=(4/33)*k4-3;
ad(1:33)=[a1,a2,a3,a4]; N=33; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) stem(k,h);
title('33点抽样h(n)图像');
h(n)=
0.0005 0.0027 0.0043 0.0006 0.0005 0.0018 0.0083 0.0005 0.0014 0.0227 0.0007 0.0005 0.0011 0.2028 0.5005 0.0009 0.2028 0.0011 0.0005 0.0007 0.0227 0.0005 0.0006 0.0083 0.0018 0.0005 0.0006 0.0043 0.0005 h(n)的图像:
0.0006 0.0009 0.0014 0.0027
图9 33点频率抽样时h(n)的图像
MATLAB代码:
k1=[0:8];k2=[9:16];k3=[17:24];k4=[25:32]; a1=(-4/33)*k1+1; a2=(4/33)*k2-1; a3=(-4/33)*k3+3; a4=(4/33)*k4-3;
ad(1:33)=[a1,a2,a3,a4]; N=33; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221); stem(k,h);
title('33点抽样h(n)的图像 '); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222); plot(p,abs(Hd));
title('33点频率抽样幅频特性'); subplot(223); plot(p,angle(Hd));
title('33点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)));grid on; title('33点频率抽样对数幅频特性')
滤波器频率特性:
图10 33点频率抽样时滤波器的频率特性
(b)65点抽样
j 根据线性相位的条件, N=65,为奇数,并且在02区间内,He关于偶对
称,选择第1种线性相位进行设计。
根据指标要求,在02内有65个采样点,所以第K点对应的频率为几个特殊的频率:间,
2k。 652216和17之间,同理,位于第32个和第33个采样点之位于653323位于第48个和第49个采样点之间。根据对称性,求得Hk: 2H0——H=
1.0000 0.9385 0.8769 0.81 0.7538 0.6923 0.6308 0.5692 0.5077 0.4462 0.3846 0.3231 0.2615 0.2000 0.1385 0.0769 0.01 0.0462 0.1077 0.1692 0.2308 0.2923 0.3538 0.41 0.4769 0.5385 0.6000 0.6615 0.7231 0.7846 0.8462 0.9077 0.9692 0.9692 0.9077 0.8462 0.7846 0.7231 0.6615 0.6000 0.5385 0.4769 0.41 0.3538 0.2923 0.2308 0.1692 0.1077 0.0462 0.01 0.0769 0.1385 0.2000 0.2615 0.3231 0.3846 0.4462 0.5077 0.5692 0.6308 0.6923 0.7538 0.81 0.8769 0.9385
kN122kNk,0k 652jnk1N1hnHkeN,n0,1,,N1Nk01HkeNk0N1N1j2e2jnkN 将H(k)代入,得到h(n):
,n0,1,,N1k1=[0:16];k2=[17:32];k3=[33:48];k4=[49:]; a1=(-4/65)*k1+1; a2=(4/65)*k2-1; a3=(-4/65)*k3+3; a4=(4/65)*k4-3;
ad(1:65)=[a1,a2,a3,a4]; N=65; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) stem(k,h);
title('65点抽样h(n)图像');
h(0)—— h()=
0.0001 0.0007 0.0009 0.0001 0.0005 0.0017 0.0001 0.0004 0.0042 0.0001 0.0003 0.0226 0.5001 0.0002 0.2027 0.0001 0.0002 0.0081 0.0001 0.0002 0.0025 0.0001 0.0002 0.0012 0.0001
0.0001 0.0001 0.0002 0.0001 0.0002 0.0001 0.0002 0.0001 0.0003 0.0001 0.0003 0.0001 0.0004 0.0001 0.0006 0.0001 0.0006 0.0012 0.0004 0.0025 0.0003 0.0081 0.0003 0.2027 0.0002 0.0226 0.0002 0.0042 0.0002 0.0017 0.0001 0.0009 0.0002 0.0002 0.0002 0.0002 0.0003 0.0004 0.0005 0.0007
k1=[0:16];k2=[17:32];k3=[33:48];k4=[49:]; a1=(-4/65)*k1+1; a2=(4/65)*k2-1; a3=(-4/65)*k3+3; a4=(4/65)*k4-3;
ad(1:65)=[a1,a2,a3,a4]; N=65; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd)) subplot(221); stem(k,h);
title('65点抽样h(n)的图像 '); [Hd,p]=freqz(h,1,100,2); figure(1); subplot(222); plot(p,abs(Hd));
title('65点频率抽样幅频特性'); subplot(223); plot(p,angle(Hd));
title('65点频率抽样相频特性'); subplot(224);
plot(p,20*log10(abs(Hd)/sum(Hd)));grid on; title('65点频率抽样对数幅频特性')
滤波器的频率特性:
对比
结论:用频率抽样法设计数字滤波器时,增加采样点的数量,可以减小逼近误差。
clc; clear;
ka1=[0:8];ka2=[9:16];ka3=[17:24];ka4=[25:32];
aa1=(-4/33)*ka1+1; aa2=(4/33)*ka2-1; aa3=(-4/33)*ka3+3; aa4=(4/33)*ka4-3;
ada(1:33)=[aa1,aa2,aa3,aa4]; Na=33; ka=0:(Na-1); wma=2*pi*ka/Na;
hda=ada.*exp(-j*(Na-1)/2*wma); ha=real(ifft(hda))
k1=[0:16];k2=[17:32];k3=[33:48];k4=[49:]; a1=(-4/65)*k1+1; a2=(4/65)*k2-1; a3=(-4/65)*k3+3; a4=(4/65)*k4-3;
ad(1:65)=[a1,a2,a3,a4]; N=65; k=0:(N-1); wm=2*pi*k/N;
hd=ad.*exp(-j*(N-1)/2*wm); h=real(ifft(hd))
subplot(221);
stem(ka,ha,'b-');hold; stem(k,h,'g-');
title('抽样h(n)的图像')
[Hda,pa]=freqz(ha,1,100,2); [Hd,p]=freqz(h,1,100,2); subplot(222);
plot(pa,abs(Hda),'b-',p,abs(Hd),'g-'); title('频率抽样幅频特性'); subplot(223);
plot(pa,angle(Hda),'b-',p,angle(Hd),'g-'); title('频率抽样相频特性'); subplot(224);
plot(pa,20*log10(abs(Hda)/sum(Hda)),'b-',p,20*log10(abs(Hd)/sum(Hd)),'g-');grid on;
title('频率抽样对数幅频特性')
因篇幅问题不能全部显示,请点此查看更多更全内容
Copyright © 2019- awee.cn 版权所有 湘ICP备2023022495号-5
违法及侵权请联系:TEL:199 1889 7713 E-MAIL:2724546146@qq.com
本站由北京市万商天勤律师事务所王兴未律师提供法律服务