数字信号处理第七章
数字信号处理第七章习题解答
————第七章———— FIR 数字滤波器设计7.1 学 习 要 点7.1.1 线性相位FIR 数字滤波器特点归纳1. 线性相位概念设()()[]n h FT eH j =ω为FIR 滤波器的频响特性函数。
()ωj e H 可表示为()()()ωθωωj g j e H e H =()ωg H 称为幅度函数，为ω的实函数。
应注意()ωg H 与幅频特性函数()ωj e H 的区别，()ωj e H 为ω的正实函数，而()ωg H 可取负值。
()ωθ称为相位特性函数，当()ωτωθ-=时，称为第一类（A 类）线性相位特性；当()ωτθωθ-=0时，称为第二类（B 类）线性相位特性。
2. 具有线性相位的FIR 滤波器的特点（()n h长度为N ）1）时域特点A 类：()()()()⎪⎪⎩⎪⎪⎨⎧--=-=--=2121,1N N n n h n N h n h ωωθ偶对称关于 (7.1)B 类：()()()()⎪⎪⎩⎪⎪⎨⎧---=-=---=21221,1N N n n h n N h n h ωπωθ奇对称关于 (7.2)群延时：()21-==-N d d τωωθ为常数，所以将A 类和B 类线性相位特性统称为恒定群时延特性。
2）频域特点A 类：N 为奇数（情况1）：()ωg H 关于ππω2,,0=三点偶对称。
N 为偶数（情况2）：()ωg H 关于πω=奇对称（()0=πg H ）。
B 类：N 为奇数（情况3）：()ωg H 关于ππω2,,0=三点奇对称。
N 为偶数（情况4）：()ωg H 关于πω2,0=奇对称，关于πω=偶对称。
3. 要点（1）情况1：可以实现所有滤波特性（低通、高通、带通、带阻和点阻等）。
（2）情况2：()0=πg H ，不能实现高通、带通和点阻滤波器。
（3）情况3：只能实现带通滤波器。
（4）情况4：不能实现低通、带阻和点阻滤波器。
7.1.2 FIR 数字滤波器设计方法 FIR 滤波器设计方法： （1）窗函数法 （2）频率采样法 （3）切比雪夫逼近法1. 窗函数法的设计步骤与要点设()()[]n h FT eH d j d =ω为希望逼近的频响特性函数，()()[]n h FT e H j d =ω为用窗函数法设计的实际滤波器的频响函数。
数字信号处理第七章
数字信号处理第七章第七章数字滤波器设计7.1：无限脉冲响应滤波器的阶数估计q7.1用mattab确定一个数字无限冲激响应低通滤波器所有四种类型的最低阶数。
指标如下:40khz的抽样率，,4khz的通带边界频率，8khz的阻带边界频率，0.5db的通带波纹，40db的最小阻带衰减。
评论你的结果。
答：标准通带边缘角频率wp是：标准阻带边缘角频率WS为：理想通带波纹rp是0.5db理想阻带波纹rs是40db1.使用这些值，巴特沃斯低通滤波器的最低阶数为n=8，相应的标准通带边缘频率wn 为0.24692.使用这些值得到切比雪夫1型低通滤波器最低阶数n=5，相应的标准通带边缘频率wn是0.2000.3/使用这些值，切比雪夫2型低通滤波器n=5的最低阶数和相应的标准通带边缘频率wn为0.40004.使用这些值得到椭圆低通滤波器最低阶数n=8，相应的标准通带边缘频率wn是0.2000.从以上结果中观察到椭圆滤波器的阶数最低，并且符合要求。
问题7。
2.用MATLAB确定四种数字无限冲激响应高通滤波器的最低阶数。
指标如下：3500hz采样率、1050hz通带边界频率、600Hz阻带边界频率、1dB通带纹波和50dB最小阻带衰减。
对结果的评论a：标准通带边缘角频率WP为：标准阻带边缘角频率ws是：理想通带纹波RP为1dB，理想阻带纹波RS为50dB1.使用这些值得到巴特沃斯高通滤波器最低阶数n=8，相应的标准通带边缘频率wn是0.5646.2.使用这些值，切比雪夫1高通滤波器的最低阶数为n=5，相应的标准通带边缘频率wn为0.60003.使用这些值得到切比雪夫2型高通滤波器最低阶数n=5，相应的标准通带边缘频率wn是0.3429.4.使用这些值，椭圆低通滤波器的最低阶数n=4，相应的标准通带边缘频率wn为0.6000。
从上述结果可以看出，椭圆滤波器的阶数最低，满足要求。
q7.3用matlab确定一个数字无限冲激响应带通滤波器所有四种类型的最低阶数。
数字信号处理第七章有限单位冲激响应FIR数字滤波器的设计方法(共95张PPT)
线性相位分析
H (z)z (N 2 1 )N n 0 1h (n ) 1 2Z (n (N 2 1 )) 1 2Z (n (N 2 1 ))
H (ej)e e j j(( N )2 1) N n 0 1 h( n) c o s(n (N 2 1 ) ) (1) H ()
m 0
即 H (z) z (N 1 )H (z 1 )
H (z) z (N 1 )H (z 1 )
所以有： h (z) 1H (z) z (N 1 )H (z 1 ) 2
1N 1h (n )z nz (N 1 )zn 2n 0
z (N 2 1 )N n 0 1 h (n ) 1 2Z (n (N 2 1 )) 1 2Z (n (N 2 1 ))
m1
(N 1)/2a(n)con s)(
n0
其中: a ( 0 ) h (N 1 ),a ( n ) 2 h ( n N 1 ),( n 1 )
2
2
由于con s对 0,,2
是偶对称的。
因此，H()对0,,2
为偶对称。
线性相位滤波器的幅度特点
2、h(n)偶对称，N为偶数
对（1）式与如上合并项，注意到由于N为偶数， h(N 1) 项即为0，则
四种线性相位滤波器
偶对称单位冲激响应
h (n ) ＝h (N－ 1－n )
相位响应
( ) N 1 2
情
况
( )
1
o
－ N( － 1)
N为 奇 数 h (n )
0 a (n )
N－ 1 n
0
N 1
n
2
( N 1) / 2
H ( ) a (n) cos n
n0
精品课件-数字信号处理(第三版) 刘顺兰-第7章
第７章数字信号处理中的有限字长效应
7.1.2 定点制误差分析 1. 数的定点表示 定点制下，一旦确定了小数点在整个数码中的位置，在整个
运算过程中即保持不变。因此，根据系统设计要求、 数值范围来 确定小数点处于什么位置很重要，这就是数的定标。 数的定标有Q表示法和S表示法两种。Q表示法形如Qn，字母Q后的 数值n表示包含n位小数。如Q0表示小数点在第0位的后面，数为整 数；Q15 表示小数点在第15位的后面，0～14位都是小数位。S表 示法则形如Sm.n，m表示整数位，n表示小数位。以16位DSP为例， 通过设定小数点在16位数中的不同位置，可以表示不同大小和不 同精度的小数。表7.1列出了一个16位数的16种Q表示、 S表示及 它们所能表示的十进制数值范围。
小的正数： (01.000..0)2×2－127=1×2－127≈5.9×10－39
(4) 当S=1，E=－127，F的23位均为1时，表示的浮点数为绝 对值最小的负数：
(10.111..1)2×2－127=(－1－2－23)×2－127≈－5.9×10－39 双精度浮点数占用8个字节(64位)存储空间，包括1位符号位、 11位阶码、 52位尾数，数值范围为1.7E-308～1.7E+308。
第７章数字信号处理中的有限字长效应
乘除运算时，假设进行运算的两个数分别为x和y，它们的Q 值分别为Qx和Qy，则两者进行乘法运算的结果为xy，Q值为Qx+Qy， 除法运算的结果为x/y，Q值为Qx－Qy。
在程序或硬件实现中，上述定标值的调整可以直接通过寄存 器的左移或右移完成。若b>0，实现x×2b需将存储x的寄存器左 移b位；若b<0，实现x×2b则需将存储x的寄存器右移|b|位即可。
称为小数点位置。
精品课件-数字信号处理-第7章
xA (n) x(n) jxˆ(n)
(7-7)
式中 xˆ(n) 是时间离散信号x(n)
xˆ(n) x(n) h(n)
（7-8）
解析信号对实信号来说就是有一阶导数的连续信号。由此意 义来说，任何序列都不是解析信号，因为它是一个以整数为变量 的函数，但xA(n)是xA(t)的采样，如果xA(t)是解析的， 我们仍认 为xA(n)也是解析的， 这是对解析信号的修正。
第七章 离散希尔伯特变换
７
h(n) |H(j )|
－
n) π
2
－ π 2
(a)
－7 －6 －5 －4 －3 －2－1
1 2 3 4 5 6 7 n
(b)
图7.1 (a) 频域特性； (b) 时域特性
第七章 离散希尔伯特变换
７
7.4 因果序列傅里叶变换下的希尔伯特变换
当xA(n)是解析序列时，其实部和虚部成希尔伯特变换关系。 它对应的频谱则是单边的。如果把频谱看成解析的，即其实部与 虚部成希尔伯特变换关系，则对应的时域序列应是单边的， 即 因果的。本节主要讨论因果序列傅里叶变换的希尔伯特变换。
第七章 离散希尔伯特变换
７
7.2 时间连续信号的希尔伯特变换
给定一时间连续信号x(t)，其希尔伯特变换 xˆ(t)定义为
xˆ(t) 1 x( )d 1 x(t )d x(t) 1
π t
π
πt
（7-1）
xˆ(t) 可以看成是x(t)通过单位冲激响应 h(t) 1 滤波器
的输出。
第七章 离散希尔伯特变换
７
在时间连续信号处理中解析信号是一个重要的概念，本章我 们将其推广到时间离散信号。从形式上说不能把复时间离散信号 或复序列看成是解析函数，因为它是一个以整数为变量的函数， 但是也可以按照类似的处理方式，将复序列之实部和虚部联系起 来使复序列的频谱在单位圆上的-π≤ω<0范围内为零。用类似 的方法也可以将周期性（或有限时宽）序列的傅里叶变换之实部 和虚部联系起来，在这种情况下，“因果性”条件是，该周期序 列在各周期的后半部为零。根据对偶关系，对于时间序列呈单边 特性的因果序列，在频域(其实部与虚部)也应存在某种变换关系。 最小相位序列是一类很重要的信号， 其傅里叶变换幅度和相位 之间存在希尔伯特变换关系。
数字信号处理第三版第七章
对称，是满足式（7.1.9）的一组解，
因为cos［ω(n－τ)］关于n=τ偶对称，所以要求τ和h(n)满
足如下条件：
()
,
N1
2
2
h(n)h(N1n), 0≤ n≤ N1
（7.1.10）
2． 线性相位FIR滤波器幅度特性Hg(ω)的特点 实质上，幅度特性的特点就是线性相位FIR滤波
因为cos［ω(n-τ)］关于ω=0, π, 2π三点偶对称，所以由 式（7.1.11）可以看出，Hg(ω)关于ω=0, π, 2π三点偶对称。 因此情况1可以实现各种（低通、高通、带通、带阻）滤 波器。
情况2： h(n)=h(N－n－1), N为偶数。
仿照情况1的推导方法得到:
H ( e j ) H g () e j = N 1 h ( n ) e j n e j M 2 h ( n )c o s (( n ) )
第7章 有限脉冲响应数字滤波器的设计
7.1 线性相位FIR数字滤波器的条件和特点 7.2 利用窗函数法设计FIR滤波器 7.3 利用频率采样法设计FIR滤波器 7.4 利用等波纹最佳逼近法设计FIR滤波器 7.5 IIR和FIR数字滤波器的比较 7.6 几种特殊类型滤波器简介 7.7 滤波器分析设计工具FDATool
用情况3的推导过程可以得到:
M
Hg() 2h(n)sin[(n)] n0
（7.1.13）
N是偶数，τ=(N－1)/2=N/2－1/2。所以，当ω=0, 2π时，
sin［ω(n－τ)］=0；当ω=π时，sin［ω(n－τ)］=(－1)n－ N/2, 为峰值点。而且sin ［ω(n－τ)］关于过零点ω=0和
如何减少吉布斯效应的影响，设计一个满足要求的FIR滤波器呢？ 直观上，增加矩形窗口的宽度（即加大N）可以减少吉布斯效应 的影响。N 时， 在主瓣附近， WRg(ω)近似为：
第七章 数字信号处理中的有限字长效应
设系数采用b位量化长度和舍入方式进行量化，系数量化误
差为e(n)，其变化范围 ( / 2, / 2) ，均值为0，方差为 2 /12
则实际系数为：
ˆ h(n) h(n) e(n)
0 n ( N 1) / 2
ˆ 且量化后 h(n) 也一定满足偶对称，即
ˆ ˆ h(n) h( N 1 n)
2.有限字长效应对信号量化的影响;
3.有限字长效应对系统参数表示的影响
4.有限字长效应在运算过程中的影响
7.1
数字信号处理中的有限长效应
有限字长效应:
在实际的处理过程中，数字信号和系统都不是无限精度的，而是有 限精度，精度的大小则有字长的大小决定，正是由于有限精度，从而给 原有的数字信号处理系统带来了影响，这种影响称为数字信号处理中的 有限字长效应。
z1 0.85 j 0.15
求得a2对z1和z2的影响
z2 0.85 j 0.15
z1 1 j 900 3.3333e a2 z1 z2
z2 1 j 900 3.3333e a2 z2 z1
可见， a2对z1和z2的影响是相同的。因而
z2 z2 a2 a2
i 1 i 1
b
b1
i b 1
b1
ai 2 i
故截尾误差满足：
0 ET (2b 2b1 ), x 0
即
0 ET , x 0
②对于反码负数
b
x 1 2 b1 ai 2 i
i 1
b1
ET QT [ x] x 1 2 ai 2 (1 2
若采用截尾处理，试分别求出原码负数1.1001、反码负数1.1100
数字信号处理讲义第7章滤波器的设计方法
第７章滤波器的设计方法教学目的1．掌握由连续时间滤波器设计离散时间IIR滤波器的方法，包括冲激响应不变法，双线性变换法等；2．了解常用的窗函数，掌握低通IIR滤波器的频率变换法、用窗函数法设计FIR滤波器的方法；3．掌握FIR滤波器的逼近原理与设计方法。
教学重点与难点重点：本章是本课程的重中之重，滤波器的设计是核心内容之一。
1．连续时间滤波器设计离散时间IIR滤波器的方法，包括冲激响应不变法，双线性变换法等；2．常用的窗函数，掌握低通IIR滤波器的频率变换法、用窗函数法设计FIR滤波器的方法；3．掌握FIR滤波器的逼近原理与设计方法。
难点：1.冲激响应不变法，双线性变换法2.用窗函数法设计FIR滤波器FIR滤波器的逼近原理与设计方法基本概念7.0.1 选频滤波器的分类数字滤波器是数字信号处理的重要基础。
在对信号的过滤、检测与参数的估计等处理中, 数字滤波器是使用最广泛的线性系统。
数字滤波器是对数字信号实现滤波的线性时不变系统。
它将输入的数字序列通过特定运算转变为输出的数字序列。
因此，数字滤波器本质上是一台完成特定运算的数字计算机。
我们已经知道，一个输入序列x(n)，通过一个单位脉冲响应为h(n)的线性时不变系统后，其输出响应y(n)为∑∞-)(y))()()(n(nn=m*=xmhnhx将上式两边经过傅里叶变换，可得式中，Y (e j ω)、X (e j ω)分别为输出序列和输入序列的频谱函数， H (ejω)是系统的频率响应函数。
可以看出，输入序列的频谱X (e j ω)经过滤波后，变为X (e j ω)H (e j ω)。
如果|H (e j ω)|的值在某些频率上是比较小的，则输入信号中的这些频率分量在输出信号中将被抑制掉。
因此，只要按照输入信号频谱的特点和处理信号的目的，适当选择H (ej ω)，使得滤波后的X (e j ω)H (e j ω)符合人们的要求，这就是数字滤波器的滤波原理。
和模拟滤波器一样，线性数字滤波器按照频率响应的通带特性可划分为低通、高通、带通和带阻几种形式。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
成绩：《数字信号处理》作业与上机实验（第七章）班级：电信学号：姓名：任课老师：李宏民完成时间：信息与通信工程学院2015—2016学年第1 学期第7章 有限脉冲响应数字滤波器设计一、教材p238：19,20,21,25,26二、某信号()x t 为：123()0.5cos(2)0.7cos(20.1)0.4cos(2)x t f t f t f t ππππ=+++，其中121100,130,600.f Hz f Hz f Hz ===设计最低阶FIR 数字滤波器，按下图所示对()x t 进行数字滤波处理，实现：(x t ()y t 1）将3f 频率分量以高于50dB 的衰减抑制，同时以低于2dB 的衰减通过1f 和2f 频率分量；2）将1f 和2f 频率分量以高于50dB 的衰减抑制，同时以低于2dB 的衰减通过3f 频率分量；要求：按数字滤波器直接型结构图编写滤波程序，求得()y n ；1)中的FIR 滤波器采用窗函数法设计；2)中的FIR 滤波器采用频率采样法设计。
画出所设计的滤波器频率特性图、信号时域图；给出滤波器设计的MATLAB 代码与滤波器实现的代码；选择合适的信号采样周期T 。
3）与第6章作业2的IIR 滤波方法进行比较研究。
一、19、 Fs=80000; fp=15000;fs=20000;rs=40;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs; rp=-20*log10(1-0.02);rs=40; [N1,wpo]=ellipord(wp/pi,ws/pi,rp,rs); [B,A]=ellip(N1,rp,rs,wpo); [Hk,wk]=freqz(B,A,500);Bt=ws-wp;alph=0.5842*(rs-21)^0.4+0.07886*(rs-21); M=ceil((rs-8)/2.285/Bt) wc=(wp+ws)/2/pi;hn=fir1(M,wc,kaiser(M+1,alph)); [Hk1,wk1]=freqz(hn,1,500); figure(1);plot(wk1/pi,20*log10(abs(Hk1)),'k'); hold on plot(wk/pi,20*log10(abs(Hk)),'r--'); hold off legend('FIR 滤波器,'IIR 滤波器');axis([0,1,-80,5]);xlabel('w/\pi');ylabel('幅度/dB'); title('损耗函数'); figure(2)plot(wk1/pi,angle(Hk1)/pi,'k'); hold on plot(wk/pi,angle(Hk)/pi,'r--'); hold off legend('FIR 滤波器','IIR 滤波器');xlabel('w/\pi');ylabel('相位/\pi'); title('相频特性曲线');0.20.40.60.810w/幅度/d B损耗函数0.20.40.60.81-1-0.8-0.6-0.4-0.200.20.40.60.81w/π相位/π相频特性曲线20、N=21;n=0:20;wc=pi/4;hn1=fir1(N-1,wc,'s',hanning(N)); hn2=fir1(N-1,wc,'s',hamming(N)); hn3=fir1(N-1,wc,'s',boxcar(N)); hn4=fir1(N-1,wc,'s',blackman(N)); figure(1)plot(n,hn1,'*b');hold on ;plot(n,hn2,'--','linewidth',2); plot(n,hn3,'r:','linewidth',3); plot(n,hn4);hold off ;xlabel('n');ylabel('h(n)');legend('汉宁窗','哈明窗','矩形窗','布莱克曼窗'); title('单位冲击响应'); figure(2)[Hk1,wk1]=freqz(hn1,1,500);plot(wk1/pi,20*log10(abs(Hk1)),'*b');hold on [Hk2,wk2]=freqz(hn2,1,500);plot(wk2/pi,20*log10(abs(Hk2)),'--','linewidth',2); [Hk3,wk3]=freqz(hn3,1,500);plot(wk3/pi,20*log10(abs(Hk3)),'r:','linewidth',3);[Hk4,wk4]=freqz(hn4,1,500);plot(wk4/pi,20*log10(abs(Hk4)));hold offlegend('汉宁窗','哈明窗','矩形窗','布莱克曼窗'); axis([0,1,-80,30]);xlabel('w/\pi');ylabel('幅度'); title('四种低通滤波器的损耗函数');5101520nh (n )单位冲击响应0.20.40.60.81w/幅度四种低通滤波器的损耗函数21、N=21;n=0:20;wc=pi/4;hn1=fir1(N-1,wc,'high',hanning(N)); hn2=fir1(N-1,wc,'high',hamming(N)); hn3=fir1(N-1,wc,'high',boxcar(N)); hn4=fir1(N-1,wc,'high',blackman(N)); figure(1)plot(n,hn1,'*b');hold on ;plot(n,hn2,'--','linewidth',2); plot(n,hn3,'r:','linewidth',3); plot(n,hn4);hold off ;xlabel('n');ylabel('h(n)');legend('汉宁窗','哈明窗','矩形窗','布莱克曼窗'); title('单位冲击响应'); figure(2)[Hk1,wk1]=freqz(hn1,1,500);plot(wk1/pi,20*log10(abs(Hk1)),'*b');hold on [Hk2,wk2]=freqz(hn2,1,500);plot(wk2/pi,20*log10(abs(Hk2)),'--','linewidth',2); [Hk3,wk3]=freqz(hn3,1,500);plot(wk3/pi,20*log10(abs(Hk3)),'r:','linewidth',3); [Hk4,wk4]=freqz(hn4,1,500);plot(wk4/pi,20*log10(abs(Hk4)));hold offlegend('汉宁窗','哈明窗','矩形窗','布莱克曼窗'); axis([0,1,-80,30]);xlabel('w/\pi'); ylabel('幅度'); title('四种窗的损耗函数');5101520nh (n )单位冲击响应0102030w/幅度四种窗的损耗函数25、wp=0.6*pi;ws=0.45*pi;rp=0.2;rs=45;Bt=wp-ws;N1=ceil(6.2*pi/Bt); N11=N1+mod(N1+1,2); disp('N11='),disp(N11) N2=ceil(11*pi/Bt); N22=N2+mod(N2+1,2); disp('N22='),disp(N22),alph=0.5842*(rs-21)^0.4+0.07886*(rs-21); M=ceil((rs-8)/2.285/Bt); disp('M='),disp(M) wc=(wp+ws)/2/pi;hn1=fir1(N11-1,wc,'high',hamming(N11)); [Hk1,w]=freqz(hn1,1);hn2=fir1(N22-1,wc,'high',blackman(N22)); [Hk2,w]=freqz(hn2,1);hn3=fir1(M-1,wc,'high',kaiser(M,alph)); [Hk3,w]=freqz(hn3,1); figure(1)plot(w/pi,20*log10(Hk1),w/pi,20*log10(Hk2),'-.',w/pi,20*log10(Hk3),'--');grid onxlabel('w/\pi');ylabel('幅度');legend('哈明窗', '布莱克曼窗', '凯塞窗') title('三种窗FIR 高通损耗函数曲线') figure(2)stem(hn1)title('哈明窗单位脉冲响应h(n)') figure(3) stem(hn2)title('布莱克曼窗单位脉冲响应h(n)') figure(4) stem(hn3)title('凯塞窗单位脉冲响应h(n)')0.20.40.60.81-160-140-120-100020w/幅度三种窗FIR 高通损耗函数曲线51015202530354045-0.4-0.3-0.2-0.100.10.20.30.40.5哈明窗单位脉冲响应h(n)1020304050607080布莱克曼窗单位脉冲响应h(n)5101520253035-0.4-0.3-0.2-0.100.10.20.30.40.5凯塞窗单位脉冲响应h(n)26、wp1=0.55*pi;wp2=0.7*pi;ws1=0.45*pi;ws2=0.8*pi;Bt=wp2-wp1; N=ceil(6.2*pi/Bt);wc=[(wp1+ws1)/2/pi,(ws2+wp2)/2/pi]; hn=fir1(N-1,wc,hanning(N)); [Hk1,wk1]=freqz(hn,1,500); n=0:N-1; figure(1); stem(n,hn);xlabel('n');ylabel('h(n)'); title('单位冲击响应'); figure(2);plot(wk1/pi,20*log10(abs(Hk1))); axis([0,1,-100,5]);xlabel('w/\pi');ylabel('·幅度'); title('损耗函数');grid on051015202530354045nh (n )单位冲击响应-1000w/?幅度损耗函数二、(1) Fs=3500; fp=130; fs=600; rs=50;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs; Bt=ws-wp;alph=0.5842*(rs-21)^0.4+0.07886*(rs-21); N=ceil((rs-8)/2.285/Bt) wc=(wp+ws)/2/pi;hn=fir1(N,wc,kaiser(N+1,alph)) figure(1) freqz(hn,1); N=500;n=0:N-1; T=1/Fs;t=n*T;x=0.5*cos(2*pi*100*t)+0.7*cos(2*pi*130*t+0.1*pi)+0.4*cos(2*pi*600*t); figure(2),plot(n,x);title('信号x(n)');ylabel('·幅值');xlabel('n');m1=0;m2=0;m3=0;m4=0;m5=0;m6=0;m7=0;m8=0;m9=0;m10=0;m11=0;m12=0;m13=0;m14=0; m15=0;m16=0;m17=0;m18=0;m19=0;m20=0;m21=0;m22=0; for m=1:500y(m)=hn(1)*x(m)+m1*hn(2)+m2*hn(3)+m3*hn(4)+m4*hn(5)+m5*hn(6)+...m6*hn(7)+m7*hn(8)+m8*hn(9)+m9*hn(10)+m10*hn(11)+m11*hn(12)+m12*hn(13)+...m13*hn(14)+m14*hn(15)+m15*hn(16)+m16*hn(17)+m17*hn(18)+m18*hn(19)+... m19*hn(20)+m20*hn(21)+m21*hn(22)+m22*hn(23);m22=m21;m21=m20;m20=m19;m19=m18;m18=m17;m17=m16;m16=m15;m15=m14;m14=m13;m13=m12;m12=m11;m11=m10;m10=m9;m9=m8;m8=m7;m7=m6;m6=m5; m5=m4;m4=m3;m3=m2;m2=m1;m1=x(m);end figure(3) plot(n,y); xlabel('n'); ylabel('y(n)');title('直接网络型y(n)');-800-600-400-2000Normalized Frequency (⨯π rad/sample)P h a s e (d e g r e e s )-150-100050Normalized Frequency (⨯π rad/sample)M a g n i t u d e (d B )0100200300400500012信号x(n)?幅值n01002003004005001ny (n )直接网络型y(n)(2) T=0.3; Fs=3500;fp=600;fs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs; Bt=wp-ws;m=1;N=ceil((m+1)*2*pi/Bt+1)N=N+mod(N+1,2);Np=fix(wp/(2*pi/N));Ns=N-2*Np-1;Hw=[zeros(1,Np+1),ones(1,Ns),zeros(1,Np)]; Hw(Np+2)=T;Ak(N-Np)=T; thetak=-pi*(N-1)*(0:N-1)/N; Hdk=Hw.*exp(1j*thetak); hn=real(ifft(Hdk)) figure(1) freqz(hn,1); N=500;n=0:N-1; T=1/Fs;t=n*T;x=0.5*cos(2*pi*100*t)+0.7*cos(2*pi*130*t+0.1*pi)+0.4*cos(2*pi*600*t); figure(2),plot(n,x);title('x(n)');ylabel('·幅度');xlabel('n');m1=0;m2=0;m3=0;m4=0;m5=0;m6=0;m7=0;m8=0;m9=0;m10=0;m11=0;m12=0;m13=0;m14=0; for m=1:500y(m)=hn(1)*x(m)+m1*hn(2)+m2*hn(3)+m3*hn(4)+m4*hn(5)+m5*hn(6)+...m6*hn(7)+m7*hn(8)+m8*hn(9)+m9*hn(10)+m10*hn(11)+m11*hn(12)+m12*hn(13)+... m13*hn(14)+m14*hn(15);m14=m13;m13=m12;m12=m11;m11=m10;m10=m9;m9=m8;m8=m7;m7=m6;m6=m5; m5=m4;m4=m3;m3=m2;m2=m1;m1=x(m);end figure(3) plot(n,y); xlabel('n'); ylabel('y(n)');title('直接网络型y(n)');00.20.40.60.81-1000-5000500Normalized Frequency (⨯π rad/sample)P h a s e (d e g r e e s )0.20.40.60.81-400-300-200-1000100Normalized Frequency (⨯π rad/sample)M a g n i t u d e (d B )0100200300400500012x(n)?幅度n01002003004005000ny (n )直接网络型y(n)(3) T=1;Fs=3500;fp=130;fs=600;wp1=(2*pi*fp)/Fs/T;ws1=(2*pi*fs)/Fs/T; rp=2;rs=50;[N,wc]=ellipord(wp1,ws1,rp,rs,'s') [B,A]=ellip(N,rp,rs,wc,'s'); [Bz,Az]=impinvar(B,A,1/T) [Hk,wk]=freqz(Bz,Az,500);wp=(2*pi*fp)/Fs;ws=(2*pi*fs)/Fs; Bt=ws-wp;alph=0.5842*(rs-21)^0.4+0.07886*(rs-21); N=ceil((rs-8)/2.285/Bt) wc=(wp+ws)/2/pi;hn=fir1(N,wc,kaiser(N+1,alph)); [Hk1,wk1]=freqz(hn,1,500); figure(1);plot(wk1/pi,20*log10(abs(Hk1)),'k'); hold onplot(wk/pi,20*log10(abs(Hk)),'r--','linewidth',2); hold off legend('FIR 滤波器','IIR 滤波器');axis([0,1,-80,5]);xlabel('w/\pi');ylabel('·幅度'); title('损耗函数'); figure(2)plot(wk1/pi,angle(Hk1)/pi,'k'); hold onplot(wk/pi,angle(Hk)/pi,'--','linewidth',2); hold off legend('FIR 滤波器','IIR 滤波器');xlabel('w/\pi');ylabel('相位/\pi'); title('相频特性曲线');w/π?幅度损耗函数01w/π相位/π相频特性曲线。