《频率归一化与转换》PPT课件.ppt

数字信号处理,2/95,8 频率归一化与滤波器转换,模拟滤波器频率转换 巴特沃斯滤波器设计 数字滤波器设计 实例分析,3/95,FIR滤波器 单位抽样响应是有限长的，滤波器一定是稳定的。只要经过一定的延时，任何非因果有限长序列都变成因果的有限序列。 FIR数字滤波器具有严格的线性相位、任意的幅度特性。 因此： FIR数字滤波可用FFT算法、时域卷积定理等来实现。,回顾,稳定系统的条件？,4/95,回顾,第I类线性相位： 严格线性,h(n)必须以(N-1)/2为中心偶对称,5/95,回顾,N又分两种情况： 偶数 奇数,此时，幅度为：,6/95,回顾,第II类线性相位,h(n)必须以(N-1)/2为中心奇对称,7/95,回顾,N分两种情况： 偶数 奇数,此时，幅度为：,回顾,设计任务 选择有限长 ，得到系统函数 ，使得： 幅频特性 满足技术指标要求； 相频特性 达到线性相位。,设计目标 给定：拟设计理想滤波器的频率响应： 设计：一个FIR滤波器频率响应： 去逼近理想滤波器的频率响应。,逼近,9/95,回顾,其中：,设计思路,无限长、非因果,有限长、因果,截断,逼近,10/95,设计步骤： 给定希望逼近的频率响应函数 ：,回顾,若没有90相移的特殊要求，一般选择滤波器具有第I类线性相位：,其中：,11/95,求单位脉冲响应 ：,如果 很复杂或不能直接计算积分，则必须用求和代替积分，以便在计算机上计算，也就是要计算离散傅里叶逆变换，一般都采用FFT来计算。,回顾,12/95,由过渡带宽及阻带最小衰减的要求，可选定窗形状，并估计窗口长度N。原则是在保证阻带衰减满足要求的情况下，尽量选择主瓣窄的窗函数： 计算所设计的FIR滤波器的单位脉冲响应（时域乘积）：,回顾,由h(n)求FIR滤波器的系统函数H(z) ：,13/95,设计一个线性相位的FIR低通滤波器，指标如下： 采样频率 fs=15kHz 通带截止频率 p=21.5103 rad/s 阻带起始频率 s=23103 rad/s 阻带最小衰减 -50 dB,例 1,14/95,【解】 （1） 首先由给出的模拟频率， 求出对应的数字频率： 通带截止频率,阻带截止频率,例 1,最后，得到理想低通滤波器的频响：,再求3dB截止频率：,其对应的数字频率:,15/95,例 1,（2）进一步得到时域表达式：,（3）由于2=50dB，查表7.2.2，可选哈明窗，其阻带最小衰减-53dB满足要求。所要求的过渡带宽：,由于哈明窗过渡带宽满足 ，因此，得到：,确定,16/95,（4）确定FIR滤波器的h(n)：,例 1,（5）由h(n)求H(ej)，检验各项指标是否满足要求。 如不满足要求，可再改变N，或改变窗函数形状（或两者都改变）来重新计算。,17/95,回顾 fir1 函数,MATLAB提供了fir1函数，以实现线性相位FIR滤波器。阅读P213 调用格式如下： hn=fir1(N, wc, ftype, window) 参数： N：阶数 wc：归一化的数字频率，0wc1。 ftype：滤波器类型，如高通、带阻等。当ftype=high时，设计高通FIR滤波器；当ftype=stop时，设计带阻FIR滤波器。注意，在设计高通和带阻滤波器时，阶数N只能取偶数(h(n)长度N+1为奇数)。不过，当用户将N设置为奇数时，fir1会自动对N加1。 window：应用的窗函数类型，默认为hamming窗。,18/95,回顾 fir1 函数,例如： hn=fir1(N,wc,bartlett(N+1) %使用Bartlett窗设计 hn=fir1(N,wc,chebwin(N+1,R) %使用Chebyshev窗设计 下面的实例代码： wc=29/66; N=32; %N=h(n)的长度减1 hn=fir1(N,wc,high),19/95,例 2 用函数设计滤波器,用窗函数法设计线性相位的FIR低通滤波器，实现对模拟信号采样后进行数字低通滤波，对模拟信号的滤波要求如下： 通带截止频率： =2 kHz 阻带截止频率： =3 kHz 阻带最小衰减： 40 dB 采样频率： =10 kHz 选择合适的窗函数及其长度，求出 ，并画出幅频特性衰减曲线和相频特性曲线。,模拟参数,20/95,例 2 用函数设计滤波器,【解】 （1）根据题意确定相应的数字滤波器指标： 通带截止频率： 阻带截止频率： 阻带最小衰减： 过渡带宽度：,转换为数字参数角频率,21/95,例 2 用函数设计滤波器,（2）选择窗函数w(n)，计算窗函数长度N 由阻带衰减为40 dB知道，hanning和hamming窗函数都满足要求，这两种窗函数过渡带宽度相同。这里选hamming窗： hamming窗函数的过渡带宽度为： 所以应当满足： 解之得到：,22/95,例 2 用函数设计滤波器,（3）确定希望逼近的理想低通滤波器截止频率 我们知道， 频率点 (该点为过渡带中心)处滤波器频响幅度为0.5(等价于衰减6dB)， 所以：,代码如下： fp=2000; %初始化 fs=3000; Fs=10000; wp=2*pi*fp/Fs; %参数转化 ws=2*pi*fs/Fs; B=ws-wp; %计算过渡带的宽度 N=ceil(8*pi/B); % 计算阶数N，ceil(x)取大于等于 x 的整数 wc=(wp+B/2)/pi; %设置理想带通截止频率(关于归一化),23/95,例 2 用函数设计滤波器,hn=fir1(N-1,wc); % 求h(n) M=1024; % 以下为计算频率相应函数和绘图部分 hk=fft(hn,M); % 求h(n)的FFT变换 n=0:N-1; figure(1);subplot(1,2,1); stem(n,hn,.); % 绘制h(n)的火柴杆图 xlabel(n); ylabel( h(n); grid on k=1:M/2+1; w=2*(0:M/2)/M; subplot(1,2,2); plot(w,20*log10(abs(hk(k); % 绘制h(n)的角频谱图(rad/s) axis(0,1,-100,5); xlabel(/); ylabel(20lg|Hg()|); grid on,24/95,例 2 用函数设计滤波器,所设计的滤波器如下：,25/95,例 2 用函数设计滤波器,绘制频谱（Hz）： figure(2); HK,F=freqz(hn，1，1024，Fs); % 没有分母，相当于其系数为1 plot(F，20*log10(abs(HK); xlabel(频率(Hz); ylabel(幅度);,26/95,例 3 信号提取,已知某信号由三个子信号叠加组成，表达式如下： 用窗函数法设计一个FIR滤波器，提取出第1个信号； 绘制滤波前、后的信号波形与频谱。,27/95,例 3 信号提取,第1步：初始化 clear all; A1=3; % 信号1的幅度 A2=1.5; % 信号2的幅度 A3=5; % 信号3的幅度 f1=100; % 信号1的频率(Hz) f2=250; % 信号2的频率(Hz) f3=270; % 信号3的频率(Hz) Fs=1024; % 采样频率(Hz) P1=-30; % 信号1的相位(度) P2=90; % 信号2的相位(度) P3=0; % 信号3的相位(度) N=1024; % 采样点数 t=0:1/Fs:2*pi; % 采样时刻 x=A1*sin(2*pi*f1*t+pi*P1/180)+A2*cos(2*pi*f2*t+pi*P2/180)+A3*sin(2*pi*f3*t+pi*P3/180);,28/95,例 3 信号提取,第2步：使用窗函数法设计FIR滤波器 fp=100; % 通带截止频率 实际实验时根据需要调整 fs=150; % 阻带截止频率 wp=2*pi*fp/Fs; % 转换为数字频率： ws=2*pi*fs/Fs; B=ws-wp; % 计算过渡带的宽度 M=ceil(8*pi/B); % 计算阶数M wc=(wp+B/2)/pi; % 设置理想带通截止频率(归一化) hn=fir1(M-1,wc); % 得到h(n),29/95,例 3 信号提取,第3步：绘制相关曲线 subplot(2,2,1); f=0:2*pi/100:2*pi; plot(f,x(1:101);title(滤波前的信号); Xk = fft(x,N); % FFT变换 Axx = (abs(Xk); % 取模 Axx=Axx/(N/2); % 换算成实际幅度:An=A/(N/2) Axx(1)=Axx(1)/2; % 换算第1个点模值：A0=A/N F=(1:N-1)*Fs/N; % 换算成实际频率值：F=(n-1)*Fs/N subplot(2,2,2); plot(F(1:N/2),Axx(1:N/2); % 显示换算后结果,仅显示左半部分 title(滤波前信号的频谱);,30/95,例 3 信号提取,Hk=fft(hn,N); Yk=Xk.*Hk; y=ifft(Yk); subplot(2,2,3); plot(f,y(1:101);title(滤波后的信号); Ayy = (abs(Yk); Ayy=Ayy/(N/2); Ayy(1)=Ayy(1)/2; F=(1:N-1)*Fs/N; subplot(2,2,4); plot(F(1:N/2),Ayy(1:N/2); title(滤波后信号的频谱);,进行滤波,31/95,课后思考题,在例3上述基础上，查阅有关知识，实现： 提取第3个信号 高通滤波器 提取第2个信号 带通滤波器 提取第1、3个信号 带阻滤波器,下次课请同学上台演示,32/95,回顾 噪声基础,一、白噪声 功率谱密度函数在整个频域内服从均匀分布。 之所以称为“白”噪声，是因为其类似于光学中包括全部可见光频率在内的白光。 凡是不符合上述条件的噪声就称为有色噪声。,白噪声是一种理想化模型，白噪声在数学处理上比较方便。 一般地，只要一个噪声过程所具有的频谱宽度远远大于它所作用系统的带宽，并且在该带宽中其频谱密度基本上可以作为常数来考虑，就可以把它作为白噪声来处理。 例如，热噪声和散弹噪声在很宽的频率范围内具有均匀的功率谱密度，通常均认为是白噪声.,33/95,回顾 噪声基础,通信系统中的某些噪声的统计特性或数字特征符合高斯过程（正态随机过程）的统计特性或数字特征，则称该噪声为高斯噪声： 经大量观察表明，高斯噪声始终存在于任何一种信道中，因而，对它的研究具有特别重要的实际意义。,二、高斯噪声,34/95,回顾 噪声基础,概率密度函数分布的正态性 功率谱密度函数分布的均匀性,高斯白噪声,当高斯噪声通过窄带系统时，即形成窄带高斯噪声 。,窄带高斯噪声,35/95,回顾 噪声基础,1. wgn()函数：产生高斯白噪声 y = wgn(m,n,p) 产生m行n列高斯白噪声矩阵，p指定输出噪声的强度。 y = wgn(m,n,p,imp) 以欧姆(Ohm)为单位指定负载阻抗。 y = wgn(m,n,p,imp,state) 重置RANDN的状态。 在数值变量后还可附加一些标志性参数： y = wgn(,POWERTYPE) 指定p的单位。POWERTYPE可以是dBW,dBm或linear。线性强度(linear power)以瓦特(Watt)为单位。 y = wgn(,OUTPUTTYPE) 指定输出类型。OUTPUTTYPE可以是real或complex。,三、生成噪声,36/95,回顾 噪声基础,例如： x = wgn(100,1,0); plot(x);,再如： x = wgn(100,1,0,complex); plot(x);,37/95,回顾 噪声基础,2. awgn()函数：在信号x中加入高斯白噪声 y = awgn(x,SNR) 信噪比SNR以dB为单位。x的强度假定为0dBW。 y = awgn(x,SNR,SIGPOWER) 如果SIGPOWER是数值，则其代表以dBW为单位的信号强度；如果SIGPOWER为measured，则函数将在加入噪声之前测定信号强度。 y = awgn(x,SNR,SIGPOWER,STATE) 重置RANDN的状态。 y = awgn(,POWERTYPE) 指定SNR和SIGPOWER单位。POWERTYPE可以是dB或linear。如果POWERTYPE是dB，那么SNR以dB为单 位，而SIGPOWER以dBW为单位。如果POWERTYPE是linear，那么SNR作为比值来度量，而SIGPOWER以瓦特为单位。,38/95,回顾 噪声基础,例如下面的代码： N=1000; t=linspace(0,50,N); x=3*sin(2*pi*0.2*t); 下面为x信号加噪声： y=awgn(x,15,measured); subplot(2,1,1); plot(t,x); subplot(2,1,2); plot(t,y);,39/95,回顾 噪声基础,MATLAB提供了高斯噪声的滤除函数： h = gaussfir(bt) h = gaussfir(bt,n) h = gaussfir(bt,n,o),四、滤除高斯噪声,40/95,回顾 噪声基础,bt = . 3; o = 8; n = 2； h = gaussfir(bt,n,o); % 产生高斯滤波器 t=linspace(0,50,1000); x=3*sin(2*pi*0.2*t); y=awgn(x,15,measured); subplot(3,1,1); plot(t,x); % 原信号 subplot(3,1,2); plot(t,y); % 加噪后的信号 Signal=conv(y,h); % 时域滤波（通过卷积运算） subplot(3,1,3); plot(Signal); % 滤波后信号,41/95,例 4 滤波器法消除噪声,（1）原始信号 Fs=256; N=256; t=0:1/N:1; % 产生原始信号 tt=0:2*pi/100:2*pi; ut=2*sin(2*pi*30*t); figure(1)；subplot(3,2,1); plot(tt,ut(1:101);grid on; title(原始信号波形); ylabel(u(t); Uk = fft(ut,N); Axx = (abs(Uk); Axx=Axx/(N/2); Axx(1)=Axx(1)/2; F=(1:N-1)*Fs/N; subplot(3,2,2); plot(F(1:N/2),Axx(1:N/2); ylabel(FFT(u); title(原始信号频谱);,42/95,例 4 滤波器法消除噪声,（2）加噪信号 xt=awgn(ut,10,0); %加入高斯白噪声 subplot(3,2,3); plot(tt,xt(1:101);grid on; title(加噪后的信号波形); ylabel(x(t); Xk=fft(xt,N); Axx = (abs(Xk); Axx=Axx/(N/2); Axx(1)=Axx(1)/2; F=(1:N-1)*Fs/N; subplot(3,2,4); plot(F(1:N/2),Axx(1:N/2); title(加噪后的信号频谱); ylabel(FFT(x);,43/95,例 4 滤波器法消除噪声,（3）使用FIR低通滤波器进行滤波 fp=30; %初始化 fs=40; wp=2*pi*fp/Fs; %参数转化 ws=2*pi*fs/Fs; B=ws-wp; %计算过渡带的宽度 M=ceil(8*pi/B); % 计算阶数N，ceil(x)取大于等于x的整数 wc=(wp+B/2)/pi; %设置理想截止频率(关于归一化) hn=fir1(M-1,wc); Hk=fft(hn,N); Tempk=Xk.*Hk; yt=ifft(Tempk); % 信号滤波,44/95,例 4 滤波器法消除噪声,（4）消噪后的信号 subplot(3,2,5); plot(tt,yt(1:101); grid on; title(消噪后的信号波形) ylabel(y(t); Yk=fft(yt,N); Axx = (abs(Yk); Axx=Axx/(N/2); Axx(1)=Axx(1)/2; F=(1:N-1)*Fs/N; subplot(3,2,6); plot(F(1:N/2),Axx(1:N/2); title(消噪后的信号频谱); ylabel(FFT(y);,45/95,例 4 滤波器法消除噪声,结果如下：,46/95,引言,已知某信号由三个子信号叠加组成，表达式如下： 设计一个高通滤波器，提取出第3个正弦信号； 并绘制滤波前、后的信号波形与频谱。,47/95,引言,1：原信号 A1=3; % 信号1的幅度 A2=1.5; % 信号2的幅度 A3=5; % 信号3的幅度 f1=100; % 信号1的频率(Hz) f2=180; % 信号2的频率(Hz) f3=270; % 信号3的频率(Hz) P1=-30; % 信号1的相位(度) P2=90; % 信号2的相位(度) P3=0; % 信号3的相位(度) Fs=1024; % 采样频率(Hz) N=Fs; % 采样点数 t=0:1/Fs:2*pi; % 采样时刻 x=A1*sin(2*pi*f1*t+pi*P1/180)+A2*cos(2*pi*f2*t+pi*P2/180)+A3*sin(2*pi*f3*t+pi*P3/180); figure(1); f=0:2*pi/100:2*pi; subplot(2,2,1);plot(f,x(1:101);title(原信号的波形); Xk = fft(x,N); F=(1:N-1)*Fs/N; subplot(2,2,2);plot(F(1:N/2),abs(Xk(1:N/2)/(N/2); title(原信号的频谱);,48/95,引言,2：双线性变换法直接设计高通滤波器 fp=250; % 通带截止频率 fs=190; % 阻带截止频率 Rp=1;Rs=25; wp=(fp/Fs)*2*pi; % 转换为数字的角频率 ws=(fs/Fs)*2*pi; OmegaP=2*Fs*tan(wp/2); %频率预畸 OmegaS=2*Fs*tan(ws/2); M,Wc=buttord(OmegaP,OmegaS,Rp,Rs,s); B1,A1=butter(M,Wc,high,s); Bz,Az=bilinear(B1,A1,Fs); y=filter(Bz,Az,x); subplot(2,2,3); plot(f,y(1:101); title(滤波后信号波形); Yk=fft(y,N); subplot(2,2,4); plot(F(1:N/2),abs(Yk(1:N/2)/(N/2); title(滤波后信号频谱);,49/95,引言,二种方法设计（模拟）高通滤波器 Fs=512; fp=250; fs=190； Rp=1;Rs=25; (1)：直接设计 wp=(fp/Fs)*2*pi; ws=(fs/Fs)*2*pi; N,Wc=buttord(wp,ws,Rp,Rs,s); B1,A1=butter(N,Wc,high,s);,(2)：间接设计 z,p,k=buttap(N); b,a = zp2tf(z, p, k); B2,A2 = lp2hp(b, a, Wc); % 低通高通 h1,w = freqs(B1,A1); subplot(2,1,1);plot(w*Fs/(2*pi),abs(h1);axis(0,Fs,0,1.2); title(直接设计高通滤波器); h2,w = freqs(B2, A2); subplot(2,1,2);plot(w*Fs/(2*pi), abs(h2); axis(0,Fs,0,1.2); title(频率变换法设计高通滤波器);,这就是 滤波器转换,50/95,因此，为了使设计规范化，通常将滤波器的频率参数进行归一化处理。,一、频率归一化,由于实际滤波器的频率范围直接取决于应用，因此必然千差万别。,1 频率归一化与转换,借助低通滤波器的系统函数，经合适的频率变换，可以得到高通、带通、带阻滤波器的系统函数。 因此不论那一种滤波器的设计，都可以： 先将该滤波器的技术指标转化为频率归一化的低通滤波器的技术指标， 按照该指标先设计低通滤波器的系统函数， 再通过频率变换，得到所需类型系统函数。,51/95,归一化复变量 q 与归一化频率的关系为：,1 频率归一化与转换,归一化原型模拟低通滤波器,52/95,若用通带截止频率 进行归一化处理，得归一化频率：,二、原型低通低通,通带、阻带截止频率分别为：,1 频率归一化与转换,若令归一化复变量为 p，则有：,前面用过,53/95,三、原型低通高通,用模拟高通滤波器的通带截止频率 进行归一化：,1 频率归一化与转换,由于滤波器的幅频响应是频率的偶函数，所以只要将低通的右（或左）半边曲线对应于高通的左（或右）半边曲线，就可以实现低通到高通的变换。,54/95,归一化低通频率与高通的对应关系如下：,显然，与之间的频率变换关系为：,1 频率归一化与转换,55/95,通带的上限、下限 阻带的上限、下限,通带带宽,通带的中心频率,用 对频率进行归一化处理：,1 频率归一化与转换,且：,四、原型低通带通,56/95,只要将 的整个曲线对应于 的右半边曲线，就可以实现低通到带通的变换 ：,1 频率归一化与转换,57/95,与之间的频率变换关系为：,的归一化复变量 p 与 的复变量 s 之间的映射关系为：,1 频率归一化与转换,58/95,边界频率参数、通带带宽 、阻带的中心频率的定义与同模拟带通滤波器相同。,其归一化频率为：,幅频特性为：,1 频率归一化与转换,五、原型低通带阻,59/95,与 的映射关系为：,与之间的频率变换关系为：,边界频率 ， 满足：,1 频率归一化与转换,60/95,按下式：,一、归一化模拟低通滤波器设计,即可求出 对应的阶次 N 和 3dB 截止频率 。,2 模拟滤波器设计,查表可确定巴特沃斯归一化原型系统函数 ，则有：,61/95,二、模拟高通、带通和带阻滤波器设计, 确定所要设计的模拟滤波器的技术指标,技术指标：,高通：,带通、带阻：,归一化：,带通、带阻：,高通：,2 模拟滤波器设计, 按照频率变换关系，将上面的技术指标转换为归一化原型模拟低通滤波器 的归一化技术指标要求：, 设计归一化原型模拟低通滤波器,62/95,【例】设计一个模拟高通滤波器，指标如下：, 由题意，得到模拟高通滤波器的技术指标（角频率）：,【解】,例 5,63/95, 归一化原型模拟低通滤波器 的技术指标要求为：, 设计归一化原型模拟低通滤波器,查P157表6.2.1得：,例 5,64/95, 确定模拟高通滤波器,例 5,主要代码如下（直接设计高通滤波器）： Wp=2*pi*100; % 转换为角频率 Ws=2*pi*50; Ap=3; As=30; N,Wc=buttord(Wp,Ws,Ap,As,s) ; B,A=butter(N,Wc,high, s),65/95,（1）buttord ( ) 函数,N,wc=buttord(wp，ws，Rp，Rs，s),输入参数：wp，ws分别为通带和阻带的实际截止频率，单位为 rad/s。 当wpws时，此时为低通滤波器； 当wswp时，为高通滤波器； 当wp和ws为二元时，为带通或带阻滤波器，此时wc也是二元。 输入参数：Rp，Rs为通带最大衰减和阻带最小衰减； 输出参数：N，wc将作为另一个函数butter的输入参数。,2 模拟滤波器设计,用于求解巴特沃斯模拟滤波器的阶数N和3dB截止频率wc。格式如下：,注意： buttord函数使用阻带指标计算3dB截止频率，这样阻带会刚好满足要求，而通带会有富余。,66/95,（2）buttap( ) 函数,z,p,k =buttap(N) % N：滤波器阶数,返回参数z、p 和 k 分别为零点、极点向量和增益常数。 由于没有零点，所以一般z为空。 p为N长的向量，表示极点。 k为增益因子。,2 模拟滤波器设计,用于设计巴特沃斯归一化（ ）模拟低通原型系统的函数:,67/95,用于将由z、p和k表示的系统函数转换为有理分式的形式：,2 模拟滤波器设计,（3）zp2tf( ) 函数,b,a = zp2tf(z,p,k),此时仍然是归一化频率的原型滤波器，即截止频率为1 rad/s。,68/95,（4）lp2lp( )、lp2hp( )、lp2bp( )和lp2bs( )函数,用于将截止频率为1 rad/s的归一化原型模拟低通滤波器的系统函数分别转换为指定频率Wo（及带宽Bw）的非归一化的模拟低通、高通、带通和带阻滤波器的系统函数。,低通：B,A=lp2lp(b,a,Wo) 高通：B,A=lp2hp(b,a,Wo) 带通：B,A=lp2bp(b,a,Wo,Bw) 带阻：B,A=lp2bs(b,a,Wo,Bw),2 模拟滤波器设计,Wo是目标滤波器的截止频率（或通带、阻带的中心频率），单位为rad/s Bw是目标滤波器的带宽,69/95,（5）freqs ( ) 函数,h,w = freqs(B,A,n) h=freqs(B,A,w) 实质：进行拉普拉斯变换，返回S域的频响特性。 B：分子多项式系数； A：分母多项式系数； w：一个实向量，一般表示为角频 rad/s。如作为输入参数向量，则表示需要计算频响的对应频率点。 n：为需计算的点数。此时自动产生一个长度为n的向量w。 如果省略n，则自动取n=200个点。 h：返回的复频响应。,2 模拟滤波器设计,70/95,例 6 低通滤波器,clear all; fp = 100;fs = 150;Rp=3;Rs=40; wp = 2*pi*fp;ws = 2*pi*fs; N, wc = buttord(wp, ws, Rp, Rs, s); z, p, k = buttap(N); % 模拟低通滤波器极零点 b, a = zp2tf(z, p, k); % 由极零点获得转移函数参数 B, A = lp2lp(b, a, wc); % 原型低通实际低通 h1,w1 = freqs(b, a); % 此时的w1是归一化频率 subplot(2,1,1); plot(w1*wc/(2*pi),abs(h1); %由于w1是归一化频率，需要转换为实 % 际角频率，再转换为线频率（Hz） title(原型滤波器);xlim(0 300); h2,w2 = freqs(B, A); subplot(2,1,2);plot(w2/(2*pi),abs(h2); title(低通滤波器);xlim(0 300);,有的书上要求归一化是对Wp进行的，此时可以改为w1*Wp。,71/95,例 6 低通滤波器,clear all; fp = 260;fs = 220; Rp=3;Rs=40; wp = 2*pi*fp; ws = 2*pi*fs; N, wc = buttord(wp, ws, Rp, Rs, s); z, p, k = buttap(N); b, a = zp2tf(z, p, k); B, A = lp2hp(b, a, wc); % 原型低通实际高通 h1,w1 = freqs(b, a); subplot(2,1,1);plot(w1*wc/(2*pi),abs(h1); title(原型滤波器);xlim(0 600); h2,w2 = freqs(B, A); subplot(2,1,2);plot(w2/(2*pi),abs(h2); title(高通滤波器);axis(0 600 0 1.5);,此时wpws，表示为高通滤波器,72/95,例 6 低通滤波器,clear all; fp1 = 160;fs1 = 140; fp2 = 200;fs2 = 220; Rp=3;Rs=40; wp1 = 2*pi*fp1;ws1 = 2*pi*fs1; wp2 = 2*pi*fp2;ws2 = 2*pi*fs2; Bw=wp2-wp1;w0=sqrt(wp1*wp1); N, wc = buttord(wp1,wp2,ws1,ws2, Rp, Rs, s); z, p, k = buttap(N); b, a = zp2tf(z, p, k); B, A = lp2bp(b, a, w0,Bw); % 原型低通实际带通 h1,w1 = freqs(b, a); subplot(2,1,1);plot(w1*w0/(2*pi),abs(h1); title(原型滤波器);xlim(0 600); h2,w2 = freqs(B, A); subplot(2,1,2);plot(w2/(2*pi),abs(h2); title(带通滤波器);axis(0 600 0 1.5);,此时的通带、阻带截止频率均为二元，表示为带通滤波器。,73/95,例 6 低通滤波器,clear all; fp1 = 140;fs1 = 160; fp2 = 220;fs2 = 200; Rp=3;Rs=40; wp1 = 2*pi*fp1;ws1 = 2*pi*fs1; wp2 = 2*pi*fp2;ws2 = 2*pi*fs2; Bw=ws2-ws1;w0=sqrt(ws1*ws1); N, wc = buttord(wp1,wp2,ws1,ws2, Rp, Rs, s); z, p, k = buttap(N); b, a = zp2tf(z, p, k); B, A = lp2bs(b, a, w0,Bw); % 原型低通实际带阻 h1,w1 = freqs(b, a); subplot(2,1,1);plot(w1*w0/(2*pi),abs(h1); title(原型滤波器);xlim(0 600); h2,w2 = freqs(B, A); subplot(2,1,2);plot(w2/(2*pi),abs(h2); title(带阻滤波器);axis(0 600 0 1.5);,74/95,3 数字滤波器设计,确定拟设计的数字滤波器（如数字高通）的技术指标； 将上述指标转换成对应的模拟滤波器（如模拟高通）的技术指标： 将转换后的模拟滤波器指标再转换成模拟低通滤波器指标； 按前面所学的方法设计模拟低通滤波器； 将得到的模拟低通通过频率变换，转换成所需类型的模拟滤波器。 采用双线性变换法等，将所设计的模拟滤波器转换成所对应类型的数字滤波器。,基本步骤,双线性,75/95,3 数字滤波器设计,设计的数字 滤波器指标,对应的模拟 滤波器指标,模拟低通 滤波器指标,设计 模拟 低通 滤波器,需要的 数字滤波器,对应的 模拟滤波器,双线性法等,76/95,综合实例（1） 高通滤波器,设计数字高通滤波器，要求如下： 通带截止频率 通带最大衰减 阻带截止频率 阻带最小衰减 【解】 （1）将上述指标转换为模拟高通的技术指标（令T=1）：,采用双线性设计法,77/95,（2）再转换为模拟低通滤波器的技术指标：,将 和 归一化（这里取 ）：,综合实例（1） 高通滤波器,78/95,（3）设计归一化模拟低通滤波器 ： 模拟低通滤波器的阶数N计算如下：,这里，取：,综合实例（1） 高通滤波器,79/95,通过查表，可轻易得到归一化模拟低通传输函数 为,下面进行去归一化。将 代入上式得到：,综合实例（1） 高通滤波器,80/95,综合实例（1） 高通滤波器,（4）将模拟低通转换成模拟高通： 将 中的变量进行替换，即： 即可得到模拟高通滤波器 ：,81/95,（5）用双线性变换法将模拟高通转换成数字高通,综合实例（1） 高通滤波器,82/95,综合实例（2） 带通滤波器,设计一个数字带通滤波器，技术指标如下： 通带上截止频率 通带下截止频率 阻带上截止频率 阻带下截止频率 通带内最大衰减 阻带内最小衰减,83/95,【解】（1）将上述指标转换为模拟带通滤波器的指标（设T=1）：,通带中心,通带带宽,综合实例（2） 带通滤波器,衰减指标不变,84/95,（2）得到模拟归一化低通滤波器技术指标 归一化阻带截止频率：,归一化通带截止频率：,综合实例（2） 带通滤波器,（3）设计模拟低通滤波器,查表，得到归一化低通传输函数：,85/95,综合实例（2） 带通滤波器,（4） 将归一化模拟低通转换成模拟带通： （5）应用双线性变换法转换成数字带通滤波器,86/95,综合实例（3） 带阻滤波器,设计一个数字带阻滤波器，指标如下： 通带上限频率 通带下限频率 阻带上截止频率 阻带下截止频率 通带内最大衰减 阻带内最小衰减,87/95,【解】（1）转换为模拟带阻滤波器的技术指标（设T=1）：,阻带中心频率的平方： 阻带的带宽为：,综合实例（3） 带阻滤波器,88/95,将以上边界频率对B归一化：,（2）得到归一化的模拟低通滤波器的技术指标：,综合实例（3） 带阻滤波器,89/95,（3）设计模拟低通滤波器：,（4）转换成模拟阻带滤波器,综合实例（3） 带阻滤波器,查表得到 。,90/95,（5）转换为数字阻带滤波器,综合实例（3） 带阻滤波器,91/95,综合实例（4） 代码设计,再看引言中的例题： 在例2的基础上，设计数字低通、高通、带通、带阻滤波器。 （一）原信号处理代码 如前面所述。 （二）绘制滤波器频响 figure(2); h, w = freqz(Bz,Az,Fs); plot(w*Fs/(2*pi),abs(h);,92/95,综合实例（4） 代码设计,（三）低通滤波：分离信号1 fp = 150; fs = 200; Rp = 3;Rs = 40; wp = 2*pi*fp/Fs;ws = 2*pi*fs/Fs; OmegaP=2*Fs*tan(wp/2); OmegaS=2*Fs*tan(ws/2); M,Wc=buttord(OmegaP,OmegaS,Rp,Rs,s); z, p, k = buttap(M); b, a = zp2tf(z, p, k); B, A = lp2lp(b, a, Wc); Bz,Az=bilinear(B,A,Fs); y=filter(Bz,Az,x); subplot(2,2,3);plot(f,y(1:101);title(滤波后信号波形); Yk=fft(y,N); subplot(2,2,4);plot(F(1:N/2),abs(Yk(1:N/2)/(N/2);title(滤波后信号频谱);,93/95,综合实例