什么滤波器设计软件好用
转来的内容:
Filter Solutions
Nuhertz公司出品,Nuhertz公司是滤波器设计软件的行业领军企业。
网上只有试用版下载,仅有20天的试用期,过了软件就基本瘫痪了,当然咯,细心点的话还是可以找到序列号的。
功能非常齐全,值得一提的是,Filter Solutions绘制的曲线可以与Protel相媲美,我基本就是用Filter Wiz Pro设计滤波器,但是使用Filter Solutions的曲线,通常是拿来就能用的。
Filter Wiz Pro
Schematica公司出品。
不多说了,本人使用次数最多的分立元件滤波器软件,功能非常齐全,基本能想到的问题它都替你考虑到了,唯一缺点是不注册的话对极点数和阻值作了一定的限制,注册费用为$199。
FilterCAD
Linear公司出品,在集成滤波器设计中应用非常广的软件,免费使用,无限制。
提供设计向导,方便快速设计集成滤波器,具有相频、幅频和群延迟曲线,频率轴可选择线性和对数两种模式。
在Linear Technology网站有免费下载。
FilterLab
Microchip公司出品,免费使用,无限制,非常实惠。
在Microchip主页有免费下载。
缺点是只有相频和幅频曲线。
FilterPro
Ti公司出品,大概我下载的安装包有文件损坏,在安装过程中出了点错误,使用过程中经常run time error,再者上述这些设计软件已经足够应付大多情况了,因此没有好好体会这款软件的功能,只是知道它的口碑一直很不错。
在滤波器设计窗口填写相关滤波器参数,即可出现窗口信息。
1、首先设计一个滤波器:打开Matlab,在命令行窗口键入fdatool等待软件打开滤波器设计工具窗口。
2、在滤波器设计窗口填写相关滤波器参数,即可出现窗口信息。
滤波器是由电容、电感和电阻组成的滤波电路。滤波器可以对电源线中特定频率的频段或该频段以外的频率进行有效滤除,得到一个特定频率的电源信号,或消除一个特定频率后的电源信号。
wp=0.25 %通带截止频率
ws=0.35 %阻带起始频率
deltaw=ws-wp %过渡带宽
N0=ceil(6.6/deltaw) %窗口长度
N=N0+mode(N0+1,2) %保证窗口长度为奇数
n=N-1 %滤波器阶数
wn=(ws+wp)/2 %滤波器截止频率
%根据通带波动和阻带最小衰减指标选择hamming窗
b=fir1(n,wn,'low',hamming(N))
freqz(b,1) %频率响应曲线
%计算验证通带波动和阻带最小衰减
[II,w]=freqz(b,1)
mag=abs(II)
db=20*log10(mag/max(mag)) %化为分贝值
dw=pi/512
Rp=-(min(db(1:wp*pi/dw+1))) %检验通带波动
As=-(max(db(ws*pi/dw+1:512))) %检验最小阻带衰减
>>Rp
Rp =
0.0279
>>As
As =
57.8095
符合要求
Order: 阶
filter circuits:滤波电路
frequency response:幅频响应
Passband :通频带、传输带宽
repeatedly cycle:重复周期
maximum signal to noise ratio:最大信噪比
gain constants:增益系数,放大常数
circuit topologies:电路拓扑结构
gain shortfall:增益不足
maximum output:最大输出功率
last stage:末级
preceding stage:前级
stage filter:分级过滤器
Gain Stage:增益级
voltage amplitude:电压振幅
Component values: 元件值
maximum valued: 最大值
minimum valued: 最小值
standard value:标准值
resistors: 电阻器
capacitors:电容器
operational amplifiers:运算放大器(OA)
circuit board:(实验用)电路板
active filters:有源滤波器
supply currents:源电流
power supplies:电源
bypassing capacitors:旁路电容
optimal:最佳的最理想的
Gain Bandwidth:带宽增益
passive component:无源元件
active component: 有源元件
overall spread:全局总范围
Component characteristics:组件特性
Modification:修改更改
data book:数据手册
typical values:标准值典型值
default values:省略补充
program execution:程序执行
Reset button:复原按钮
positive temperature coefficient:正温度系数
variable resistors:可变电阻器
cermet resistor:金属陶瓷电阻器
output resistance:输出电阻
distortion:失真
single amplifier:单级放大器
voltage follower:电压输出跟随器
troubleshooting:发现并修理故障
control panel,:控制面板
1、打开crack的软件后,根据滤波器的设计要求,在filter type中选择滤波器的类型(Gaussian:高斯滤波器、Bessel:贝塞尔滤波器、butterworth:巴特沃斯;Chebyshev1切比雪夫1;Chebyshev2切比雪夫2;Hourglass:对三角滤波器、Elliptic:椭圆滤波器、Custom:自定义滤波器、Raised Cos:升余弦滤波器、Matche:匹配滤波器 、Delay:延迟滤波器);
2、在filter class中选择滤波器的种类(低通、高通、带通、带阻);
3、在filter Attributes中设置滤波器的阶数(Order)、通频带频率(Passband frequency);
4、在Implementation中选择有源滤波器(active )、无源滤波器(passive)和数字滤波器(Digital);
5、在Freq Scale中选择Hertz和Log,如果选择了Rad/Sec(一种单位),则要注意Rad/Sec=6.28*Hertz;
6、在Graph Limits中设置好图像的最大频率和最小频率,最大频率要大于通频带的截止频率;在Passive Design/Ideal Filter Response中观察传输函数(Transfer Function)、时域响应(Time Response)、零极点图(Pole Zero Plots)、频域响应(Frequency Response)的图像;
7、在Circuit Parmaters中设置源电阻(Source Res)和负载电阻(Load Res)最后点击Circuits观察滤波器电路图;
8、在设计有缘滤波器的时候还要注意在Active Implementation 中选择滤波器的电路布局形式一般有源滤波器选择Pos SAB型的,在Circuit Parmaters中设置增益大小(gain)。
1低通滤波器的基本参数和设计特点为获得好的数据采集效果,获得稳定的数据,在ADC前应该加有源滤波电路,滤去不需要的干扰。
1.1低通滤波器的基本参数低通滤波幅频特性曲线。
中:Passband通带,可以通过频率的区域,Stopband阻带,不允许通过频率区域。
在通带与阻带之间是
绝对正确的代码:程序1:
fs=22050 %语音信号采样频率为22050
x1=wavread('Windows Critical Stop.wav')%读取语音信号的数据,赋给变量x1
sound(x1,22050) %播放语音信号
y1=fft(x1,1024) %对信号做1024点FFT变换
f=fs*(0:511)/1024
figure(1)
plot(x1) %做原始语音信号的时域图形
title('原始语音信号')
xlabel('time n')
ylabel('fuzhi n')
figure(2)
freqz(x1) %绘制原始语音信号的频率响应图
title('频率响应图')
figure(3)
subplot(2,1,1)
plot(abs(y1(1:512))) %做原始语音信号的FFT频谱图
title('原始语音信号FFT频谱')
subplot(2,1,2)
plot(f,abs(y1(1:512)))
title('原始语音信号频谱')
xlabel('Hz')
ylabel('fuzhi')
程序2:
fs=22050 %语音信号采样频率为22050
x1=wavread('Windows Critical Stop.wav')%读取语音信号的数据,赋给变量x1
t=0:1/22050:(size(x1)-1)/22050
y1=fft(x1,1024) %对信号做1024点FFT变换
f=fs*(0:511)/1024
x2=randn(1,length(x1)) %产生一与x长度一致的随机信号
sound(x2,22050)
figure(1)
plot(x2) %做原始语音信号的时域图形
title('高斯随机噪声')
xlabel('time n')
ylabel('fuzhi n')
randn('state',0)
m=randn(size(x1))
x2=0.1*m+x1
sound(x2,22050)%播放加噪声后的语音信号
y2=fft(x2,1024)
figure(2)
plot(t,x2)
title('加噪后的语音信号')
xlabel('time n')
ylabel('fuzhi n')
figure(3)
subplot(2,1,1)
plot(f,abs(y2(1:512)))
title('原始语音信号频谱')
xlabel('Hz')
ylabel('fuzhi')
subplot(2,1,2)
plot(f,abs(y2(1:512)))
title('加噪后的语音信号频谱')
xlabel('Hz')
ylabel('fuzhi')
根据以上代码,你可以修改下面有错误的代码
程序3:双线性变换法设计Butterworth滤波器
fs=22050
x1=wavread('h:\课程设计2\shuzi.wav')
t=0:1/22050:(size(x1)-1)/22050
Au=0.03
d=[Au*cos(2*pi*5000*t)]'
x2=x1+d
wp=0.25*pi
ws=0.3*pi
Rp=1
Rs=15
Fs=22050
Ts=1/Fs
wp1=2/Ts*tan(wp/2)%将模拟指标转换成数字指标
ws1=2/Ts*tan(ws/2)
[N,Wn]=buttord(wp1,ws1,Rp,Rs,'s') %选择滤波器的最小阶数
[Z,P,K]=buttap(N) %创建butterworth模拟滤波器
[Bap,Aap]=zp2tf(Z,P,K)
[b,a]=lp2lp(Bap,Aap,Wn)
[bz,az]=bilinear(b,a,Fs) %用双线性变换法实现模拟滤波器到数字滤波器的转换
[H,W]=freqz(bz,az)%绘制频率响应曲线
figure(1)
plot(W*Fs/(2*pi),abs(H))
grid
xlabel('频率/Hz')
ylabel('频率响应幅度')
title('Butterworth')
f1=filter(bz,az,x2)
figure(2)
subplot(2,1,1)
plot(t,x2) %画出滤波前的时域图
title('滤波前的时域波形')
subplot(2,1,2)
plot(t,f1)%画出滤波后的时域图
title('滤波后的时域波形')
sound(f1,22050) %播放滤波后的信号
F0=fft(f1,1024)
f=fs*(0:511)/1024
figure(3)
y2=fft(x2,1024)
subplot(2,1,1)
plot(f,abs(y2(1:512)))%画出滤波前的频谱图
title('滤波前的频谱')
xlabel('Hz')
ylabel('fuzhi')
subplot(2,1,2)
F1=plot(f,abs(F0(1:512))) %画出滤波后的频谱图
title('滤波后的频谱')
xlabel('Hz')
ylabel('fuzhi')
程序4:窗函数法设计滤波器:
fs=22050
x1=wavread('h:\课程设计2\shuzi.wav')
t=0:1/22050:(size(x1)-1)/22050
Au=0.03
d=[Au*cos(2*pi*5000*t)]'
x2=x1+d
wp=0.25*pi
ws=0.3*pi
wdelta=ws-wp
N=ceil(6.6*pi/wdelta) %取整
wn=(0.2+0.3)*pi/2
b=fir1(N,wn/pi,hamming(N+1)) %选择窗函数,并归一化截止频率
figure(1)
freqz(b,1,512)
f2=filter(bz,az,x2)
figure(2)
subplot(2,1,1)
plot(t,x2)
title('滤波前的时域波形')
subplot(2,1,2)
plot(t,f2)
title('滤波后的时域波形')
sound(f2,22050) %播放滤波后的语音信号
F0=fft(f2,1024)
f=fs*(0:511)/1024
figure(3)
y2=fft(x2,1024)
subplot(2,1,1)
plot(f,abs(y2(1:512)))
title('滤波前的频谱')
xlabel('Hz')
ylabel('fuzhi')
subplot(2,1,2)
F2=plot(f,abs(F0(1:512)))
title('滤波后的频谱')
xlabel('Hz')
ylabel('fuzhi')
1.了解数字信号处理系统的一般构成;
2.掌握奈奎斯特抽样定理。
1.2 实验仪器
1.YBLD智能综合信号源测试仪 1台
2.双踪示波器 1台
3.MCOM-TG305数字信号处理与现代通信技术实验箱 1台
4.PC机(装有MATLAB、MCOM-TG305配套实验软件) 1台
1.3 实验原理
一个典型的DSP系统除了数字信号处理部分外,还包括A/D和D/A两部分。这是因为自然界的信号,如声音、图像等大多是模拟信号,因此需要将其数字化后进行数字信号处理,模拟信号的数字化即称为A/D转换。数字信号处理后的数据可能需还原为模拟信号,这就需要进行D/A转换。一个仅包括A/D和D/A两部分的简化数字信号处理系统功能如图1所示。
A/D转换包括三个紧密相关的过程,即抽样、量化和编码。A/D转换中需解决的以下几个重要问题:抽样后输出信号中还有没有原始信号的信息?如果有能不能把它取出来?抽样频率应该如何选择?
奈奎斯特抽样定理(即低通信号的均匀抽样定理)告诉我们,一个频带限制在0至fx以内的低通信号x(t),如果以fs≥2fx的抽样速率进行均匀抽样,则x(t)可以由抽样后的信号xs(t)完全地确定,即xs(t)包含有x(t)的成分,可以通过适当的低通滤波器不失真地恢复出x(t)。最小抽样速率fs=2fx称为奈奎斯特速率。
低通
译码
编码
量化
抽样
输入信号 样点输出 滤波输出
A/D(模数转换) D/A(数模转换)
图1 低通采样定理演示
为方便实现,实验中更换了一种表现形式,即抽样频率固定(10KHz),通过改变输入模拟信号的频率来展示低通抽样定理。我们可以通过研究抽样频率和模拟信号最高频率分量的频率之间的关系,来验证低通抽样定理。
1.4 实验内容
1.软件仿真实验:编写并调试MATLAB程序,分析有关参数,记录有关波形。
2.硬件实验:输入不同频率的正弦信号,观察采样时钟波形、输入信号波形、样点输出波形和滤波输出波形。
1.5 MATLAB参考程序和仿真内容
%*******************************************************************%
%f—余弦信号的频率
% M—基2 FFT幂次数 N=2^M为采样点数,这样取值是为了便于作基2的FFT分析
%2. 采样频率Fs
%*******************************************************************%
function samples(f,Fs,M)
N=2^M% fft点数=取样总点数
Ts=1/Fs% 取样时间间隔
T=N*Ts% 取样总时间=取样总点数*取样时间间隔
n=0:N-1
t=n*Ts
Xn=cos(2*f*pi*t)
subplot(2,1,1)
stem(t,Xn)
axis([0 T 1.1*min(Xn) 1.1*max(Xn)])
xlabel('t -->')
ylabel('Xn')
Xk=abs(fft(Xn,N))
subplot(2,1,2)
stem(n,Xk)
axis([0 N 1.1*min(Xk) 1.1*max(Xk)])
xlabel('frequency -->')
ylabel('!Xk!')
%*******************************************************************%
假如有一个1Hz的余弦信号y=cos(2*π*t),对其用4Hz的采样频率进行采样,共采样32点,只需执行samples(1,4,5),即可得到仿真结果。
软件仿真实验内容如下表所示:
仿真参数
f
Fs
Wo(计算)
Xn(图形)
Xk(图形)
(1,4,5)
另外记录图形,并标图号
(1,8,5)
(2,8,6)
自 选
1.6 硬件实验步骤
本实验箱采样频率fs固定为10KHz,低通滤波器的截止频率约为4.5KHz。
1、用低频信号源产生正弦信号,正弦信号源频率f自定,并将其接至2TP2(模拟输入)端,将示波器通道一探头接至2TP6(采样时钟)端观察采样时钟波形,示波器通道二探头接至2TP2观察并记录输入信号波形。
2、将示波器通道二探头接至2TP3观察并记录样点输出波形。
3、将示波器通道二探头接至2TP4观察并记录滤波输出波形。
4、根据采样定理,分f=fs /8、f=fs/4、f=fs/2等3种情况更改正弦信号频率,重复步骤2至步骤3。
5、用低频信号源产生方波信号,重复步骤1至步骤4。
1.7 思考题
1、 讨论在仿真实验中所计算的数字域频率Wo和Xk的图形中非零谱线位置之间的对应关系。
2、 讨论在仿真实验中自选参数的意义。
3、将在2TP2端加方波信号后的恢复波形,与相同频率的正弦信号的恢复波形相比,能够得出哪些结论?
2 FFT频谱分析实验
2.1 实验目的
1.通过实验加深对快速傅立叶变换(FFT)基本原理的理解。
2.了解FFT点数与频谱分辨率的关系,以及两种加长序列FFT与原序列FFT的关系。
2.2 实验仪器
1.YBLD智能综合信号源测试仪 1台
2.双踪示波器 1台
3.MCOM-TG305数字信号处理与现代通信技术实验箱 1台
4.PC机(装有MATLAB、MCOM-TG305配套实验软件) 1台
2.3 实验原理
离散傅里叶变换(DFT)和卷积是信号处理中两个最基本也是最常用的运算,它们涉及到信号与系统的分析与综合这一广泛的信号处理领域。实际上卷积与DFT之间有着互通的联系:卷积可化为DFT来实现,其它的许多算法,如相关、滤波和谱估计等都可化为DFT来实现,DFT也可化为卷积来实现。
对N点序列x(n),其DFT变换对定义为:
在DFT运算中包含大量的重复运算。FFT算法利用了蝶形因子WN的周期性和对称性,从而加快了运算的速度。FFT算法将长序列的DFT分解为短序列的DFT。N点的DFT先分解为2个N/2点的DFT,每个N/2点的DFT又分解为2个N/4点的DFT。按照此规律,最小变换的点数即所谓的“基数(radix)。”因此,基数为2的FFT算法的最小变换(或称蝶形)是2点DFT。一般地,对N点FFT,对应于N个输入样值,有N个频域样值与之对应。一般而言,FFT算法可以分为时间抽取(DIT)FFT和频率抽取(DIF)两大类。
在实际计算中,可以采用在原来序列后面补0的加长方法来提高FFT的分辨率;可以采用在原来序列后面重复的加长方法来增加FFT的幅度。
2.4 实验内容
1.软件仿真实验:分别观察并记录正弦序列、方波序列及改变FFT的点数后的频谱;分别观察并记录正弦序列、方波序列及2种加长序列等信号的频谱。
2.硬件实验:分别观察并记录正弦信号、方波信号及改变FFT的点数后的频谱。
2.5 MATLAB参考程序和仿真内容
%*******************************************************************%
function[x]=ffts(mode,M)
Nfft=2^M
x=zeros(1,Nfft)%定义一个长度为Nfft的一维全0数组
if mode= =1 for n=0:Nfft-1 x(n+1)=sin(2*pi*n/Nfft)end
end %定义一个长度为Nfft的单周期正弦序列
if mode= =2 for n=0:Nfft-1 x(n+1)=sin(4*pi*n/Nfft)end
end %定义一个长度为Nfft的双周期正弦序列
if mode= =3 for n=0:Nfft/2-1 x(n+1)=sin(4*pi*n/Nfft)end
end %定义一个长度为Nfft/2的正弦序列,后面一半为0序列。
if mode= =4 for n=0:Nfft-1 x(n+1)=square(2*pi*n/Nfft)end
end
if mode= =5 for n=0:Nfft-1 x(n+1)=square(2*pi*n/Nfft)end
end
if mode= =6 for n=0:Nfft/2-1 x(n+1)=square(4*pi*n/Nfft)end
end
n=0:Nfft-1
subplot(2,1,1)
stem(n,x)
axis([0 Nfft-1 1.1*min(x) 1.1*max(x)])
xlabel('Points-->')
ylabel('x(n)')
X=abs(fft(x,Nfft))
subplot(2,1,2)
stem(n,X)
axis([0 Nfft-1 1.1*min(X) 1.1*max(X)])
xlabel('frequency-->')
ylabel('!X(k)!')
%*******************************************************************%
假设需观察方波信号的频谱,对一个周期的方波信号作32点的FFT,则只需在MATLAB的命令窗口下键入:[x]=ffts(21,5) ,程序进行模拟,并且输出FFT的结果。
关于软件仿真实验内容,建议在完成大量仿真例子的基础上,选择能够体现实验要求的4个以上的例子进行记录。例如要观察后面补0的加长方法来提高FFT的分辨率的现象,可以仿真ffts(4,5)和ffts(6,6)两个例子。
2.6 硬件实验步骤
1.将低频信号源输出加到实验箱模拟通道1输入端,将示波器探头接至模拟通道1输出端。
2.在保证实验箱正确加电且串口电缆连接正常的情况下,运行数字信号处理与DSP应用实验开发软件,在“数字信号处理实验”菜单下选择“FFT频谱分析”子菜单,出现显示FFT频谱分析功能提示信息的窗口。
3.用低频信号产生器产生一个1KHz的正弦信号。
4.选择FFT频谱分析与显示的点数为64点,开始进行FFT运算。此后,计算机将周期性地取回DSP运算后的FFT数据并绘图显示
5.改信号源频率,观察并记录频谱图的变化。
6.选择FFT的点数为128点,观察并记录频谱图的变化。
7.更改正弦信号的频率,重复步骤4 ~步骤6。
8.用低频信号产生器产生一个1KHz的方波信号,重复步骤4 ~步骤7。注意:应根据实验箱采样频率fs为10KHz和方波信号的频带宽度选择方波信号的频率。
本硬件实验要进行两种信号,每个信号两种频率,每个信号两种点数等共8次具体实验内容,性质能够体现实验要求的4个以上的例子进行记录。
2.7 思考题
1.对同一个信号,不同点数FFT观察到的频谱图有何区别?
2.序列加长后FFT与原序列FFT的关系是什么,试推导其中一种关系。
3.用傅立叶级数理论,试说明正弦信号频谱和方波信号频谱之间的关系。
3 IIR滤波器设计实验
3.1 实验目的
1.通过实验加深对IIR滤波器基本原理的理解。
2.学习编写IIR滤波器的MATLAB仿真程序。
3.2 实验仪器
1.YBLD智能综合信号源测试仪 1台
2.双踪示波器 1台
3.MCOM-TG305数字信号处理与现代通信技术实验箱 1台
4.PC机(装有MATLAB、MCOM-TG305配套实验软件) 1台
3.3 实验原理
IIR滤波器有以下几个特点:
1.IIR数字滤波器的系统函数可以写成封闭函数的形式。
2.IIR数字滤波器采用递归型结构,即结构上带有反馈环路。IIR滤波器运算结构通常由延时、乘以系数和相加等基本运算组成,可以组合成直接型、正准型、级联型、并联型四种结构形式,都具有反馈回路。由于运算中的舍入处理,使误差不断累积,有时会产生微弱的寄生振荡。
3.IIR数字滤波器在设计上可以借助成熟的模拟滤波器的成果,如巴特沃斯、契比雪夫和椭圆滤波器等,有现成的设计数据或图表可查,其设计工作量比较小,对计算工具的要求不高。在设计一个IIR数字滤波器时,我们根据指标先写出模拟滤波器的公式,然后通过一定的变换,将模拟滤波器的公式转换成数字滤波器的公式。
4.IIR数字滤波器的相位特性不好控制,对相位要求较高时,需加相位校准网络。
在MATLAB下设计IIR滤波器可使用Butterworth函数设计出巴特沃斯滤波器,使用Cheby1函数设计出契比雪夫I型滤波器,使用Cheby2设计出契比雪夫II型滤波器,使用ellipord函数设计出椭圆滤波器。下面主要介绍前两个函数的使用。
与FIR滤波器的设计不同,IIR滤波器设计时的阶数不是由设计者指定,而是根据设计者输入的各个滤波器参数(截止频率、通带滤纹、阻带衰减等),由软件设计出满足这些参数的最低滤波器阶数。在MATLAB下设计不同类型IIR滤波器均有与之对应的函数用于阶数的选择。
一、巴特沃斯IIR滤波器的设计
在MATLAB下,设计巴特沃斯IIR滤波器可使用butter函数。
Butter函数可设计低通、高通、带通和带阻的数字和模拟IIR滤波器,其特性为使通带内的幅度响应最大限度地平坦,但同时损失截止频率处的下降斜度。在期望通带平滑的情况下,可使用butter函数。
butter函数的用法为:
[b,a]=butter(n,Wn,/ftype/)
其中n代表滤波器阶数,Wn代表滤波器的截止频率,这两个参数可使用buttord函数来确定。buttord函数可在给定滤波器性能的情况下,求出巴特沃斯滤波器的最小阶数n,同时给出对应的截止频率Wn。buttord函数的用法为:
[n,Wn]= buttord(Wp,Ws,Rp,Rs)
其中Wp和Ws分别是通带和阻带的拐角频率(截止频率),其取值范围为0至1之间。当其值为1时代表采样频率的一半。Rp和Rs分别是通带和阻带区的波纹系数。
不同类型(高通、低通、带通和带阻)滤波器对应的Wp和Ws值遵循以下规则:
1.高通滤波器:Wp和Ws为一元矢量且Wp>Ws;
2.低通滤波器:Wp和Ws为一元矢量且Wp<Ws;
3.带通滤波器:Wp和Ws为二元矢量且Wp<Ws,如Wp=[0.2,0.7],Ws=[0.1,0.8]
4.带阻滤波器:Wp和Ws为二元矢量且Wp>Ws,如Wp=[0.1,0.8],Ws=[0.2,0.7]。
二、契比雪夫I型IIR滤波器的设计
在期望通带下降斜率大的场合,应使用椭圆滤波器或契比雪夫滤波器。在MATLAB下可使用cheby1函数设计出契比雪夫I型IIR滤波器。
cheby1函数可设计低通、高通、带通和带阻契比雪夫I型滤IIR波器,其通带内为等波纹,阻带内为单调。契比雪夫I型的下降斜度比II型大,但其代价是通带内波纹较大。
cheby1函数的用法为:
[b,a]=cheby1(n,Rp,Wn,/ftype/)
在使用cheby1函数设计IIR滤波器之前,可使用cheblord函数求出滤波器阶数n和截止频率Wn。cheblord函数可在给定滤波器性能的情况下,选择契比雪夫I型滤波器的最小阶和截止频率Wn。
cheblord函数的用法为:
[n,Wn]=cheblord(Wp,Ws,Rp,Rs)
其中Wp和Ws分别是通带和阻带的拐角频率(截止频率),其取值范围为0至1之间。当其值为1时代表采样频率的一半。Rp和Rs分别是通带和阻带区的波纹系数。
3.4 实验内容
1.软件仿真实验:编写并调试MATLAB程序,选择不同形式,不同类型的4种滤波器进行仿真,记录幅频和相频特性,对比巴特沃斯滤波器和契比雪夫滤波器。
2.硬件实验:设计IIR滤波器,在计算机上观察冲激响应、幅频特性和相频特性,然后下载到实验箱。用示波器观察输入输出波形,测试滤波器的幅频响应特性。
3.5 MATLAB参考程序和仿真内容
%*******************************************************************%
%mode: 1--巴特沃斯低通;2--巴特沃斯高通;3--巴特沃斯带通;4--巴特沃斯带阻
% 5--契比雪夫低通;6--契比雪夫高通;7--契比雪夫带通;8--契比雪夫带阻
%fp1,fp2: 通带截止频率,当高通或低通时只有fp1有效
%fs1, fs2: 阻带截止频率,当高通或低通时只有fs1有效
%rp: 通带波纹系数
%as: 阻带衰减系数
%sample: 采样率
%h: 返回设计好的滤波器系数
%*******************************************************************%
function[b,a]=iirfilt(mode,fp1,fp2,fs1,fs2,rp,as,sample)
wp1=2*fp1/samplewp2=2*fp2/sample
ws1=2*fs1/samplews2=2*fs2/sample
%得到巴特沃斯滤波器的最小阶数N和3bd频率wn
if mode<3[N,wn]=buttord(wp1,ws1,rp,as)
elseif mode<5[N,wn]=buttord([wp1 wp2],[ws1 ws2],rp,as)
%得到契比雪夫滤波器的最小阶数N和3bd频率wn
elseif mode<7[N,wn]=cheb1ord(wp1,ws1,rp,as)
else[N,wn]=cheblord([wp1 wp2],[ws1 ws2],rp,as)
end
%得到滤波器系数的分子b和分母a
if mode= =1[b,a]=butter(N,wn)end
if mode= =2[b,a]=butter(N,wn,/high/)end
if mode= =3[b,a]=butter(N,wn)end
if mode= =4[b,a]=butter(N,wn,/stop/)end
if mode= =5[b,a]=cheby1(N,rp,wn)end
if mode= =6[b,a]=cheby1(N,rp,wn,/high/)end
if mode= =7[b,a]=cheby1(N,rp,wn)end
if mode= =8[b,a]=cheby1(N,rp,wn,/stop/)end
set(gcf,/menubar/,menubar)
freq_response=freqz(b,a)
magnitude=20*log10(abs(freq_response))
m=0:511
f=m*sample/(2*511)
subplot(3,1,1)plot(f,magnitude)grid%幅频特性
axis([0 sample/2 1.1*min(magnitude) 1.1*max(magnitude)])
ylabel('Magnitude')xlabel('Frequency-->')
phase=angle(freq_response)
subplot(3,1,2)plot(f,phase)grid%相频特性
axis([0 sample/2 1.1*min(phase) 1.1*max(phase)])
ylabel('Phase')xlabel('Frequency-->')
h=impz(b,a,32)%32点的单位函数响应
t=1:32
subplot(3,1,3)stem(t,h)grid
axis([0 32 1.2*min(h) 1.1*max(h)])
ylabel('h(n)')xlabel('n-->')
%*******************************************************************%
假设需设计一个巴特沃斯低通IIR滤波器,通带截止频率为2KHz,阻带截止频率为3KHz,通带波纹系数为1,阻带衰减系数为20,采样频率为10KHz,则只需在MATLAB的命令窗口下键入:
[b,a]=iirfilt(1,2000,3000,2400,2600,1,20,10000)
程序进行模拟,并且按照如下顺序输出数字滤波器系统函数
的系数
b= b0 b1 ……bn
a= a0 a1 ……an
关于软件仿真实验内容,建议在完成大量仿真例子的基础上,选择能够体现实验要求的4个例子进行记录,系统函数只要记录系统的阶数。
3.6 硬件实验步骤
1.根据实验箱采样频率fs为10KHz的条件,用低频信号发生器产生一个频率合适的低频正弦信号,将其加到实验箱模拟通道1输入端,将示波器通道1探头接至模拟通道1输入端,通道2探头接至模拟通道2输出端。
2.在保证实验箱正确加电且串口电缆连接正常的情况下,运行数字信号处理与DSP应用实验开发软件,在“数字信号处理实验”菜单下选择“IIR滤波器”子菜单,出现提示信息。
3.输入滤波器类型、滤波器截止频率等参数后,分别点击“幅频特性”和“相频特性”按钮,在窗口右侧观察IIR滤波器的幅频特性和相频特性。此时提示信息将消失,如需查看提示信息,可点击“设计说明”按钮。
4.点击“下载实现”按钮,IIR滤波器开始工作,此时窗口右侧将显示IIR滤波器的幅频特性。
5.根据输入滤波器类型,更改低频信号源的频率,观察示波器上输入输出波形幅度的变化情况,测量IIR滤波器的幅频响应特性,看其是否与设计的幅频特性一致。
6.更改滤波器类型、滤波器截止频率等参数(共4种),重复步骤3至步骤5。所选择的例子参数最好和MATLAB仿真程序的例子一样。
7.用低频信号产生器产生一个500Hz的方波信号,分别设计3种滤波器,完成如下表要求的功能,并且记录参数和波形。
功 能
滤波器类型
参 数
输出波形
fp1
fp2
fs1
fs2
通过3次及以下次数的谐波
另外记录图形,并标图号
滤除5次及以下次数的谐波
通过3次到5次的谐波
3.7 思考题
1.在实验箱采样频率fs固定为10KHz的条件下,要观察方波信号频带宽度内的各个谐波分量,方波信号的频率最高不能超过多少,为什么?
2.硬件实验内容7中输出信号各个谐波分量,与原来方波信号同样谐波分量相比,有没有发生失真?主要发生了什么类型的失真?为什么?
4 窗函数法FIR滤波器设计实验
4.1 实验目的
1.通过实验加深对FIR滤波器基本原理的理解。
2.学习使用窗函数法设计FIR滤波器,了解窗函数的形式和长度对滤波器性能的影响。
4.2 实验仪器
1.YBLD智能综合信号源测试仪 1台
2.双踪示波器 1台
3.MCOM-TG305数字信号处理与现代通信技术实验箱 1台
4.PC机(装有MATLAB、MCOM-TG305配套实验软件) 1台
4.3 实验原理
数字滤波器的设计是数字信号处理中的一个重要内容。数字滤波器设计包括FIR(有限单位脉冲响应)滤波器与IIR(无限单位脉冲响应)滤波器两种。
与IIR滤波器相比,FIR滤波器在保证幅度特性满足技术要求的同时,很容易做到严格的线性相位特性。设FIR滤波器单位脉冲响应h(n)长度为N,其系统函数H(z)为:
H(z)是z-1的N-1次多项式,它在z平面上有N-1个零点,原点z=0是N-1阶重极点,因此H(z)是永远稳定的。稳定和线性相位特性是FIR滤波器突出的优点。
FIR滤波器的设计任务是选择有限长度的h(n)。使传输函数H( )满足技术要求。FIR滤波器的设计方法有多种,如窗函数法、频率采样法及其它各种优化设计方法,本实验介绍窗函数法的FIR滤波器设计。
窗函数法是使用矩形窗、三角窗、巴特利特窗、汉明窗、汉宁窗和布莱克曼窗等设计出标准响应的高通、低通、带通和带阻FIR滤波器。
一、firl函数的使用
在MATLAB下设计标准响应FIR滤波器可使用firl函数。firl函数以经典方法实现加窗线性相位FIR滤波器设计,它可以设计出标准的低通、带通、高通和带阻滤波器。firl函数的用法为:
b=firl(n,Wn,/ftype/,Window)
各个参数的含义如下:
b—滤波器系数。对于一个n阶的FIR滤波器,其n+1个滤波器系数可表示为:b(z)=b(1)+b(2)z-1+…+b(n+1)z-n。
n—滤波器阶数。
Wn—截止频率,0≤Wn≤1,Wn=1对应于采样频率的一半。当设计带通和带阻滤波器时,Wn=[W1 W2],W1≤ω≤W2。
ftype—当指定ftype时,可设计高通和带阻滤波器。Ftype=high时,设计高通FIR滤波器;ftype=stop时设计带阻FIR滤波器。低通和带通FIR滤波器无需输入ftype参数。
Window—窗函数。窗函数的长度应等于FIR滤波器系数个数,即阶数n+1。
二、窗函数的使用
在MATLAB下,这些窗函数分别为:
1.矩形窗:w=boxcar(n),产生一个n点的矩形窗函数。
2.三角窗:w=triang(n),产生一个n点的三角窗函数。
当n为奇数时,三角窗系数为w(k)=
当n为偶数时,三角窗系数为w(k)=
3.巴特利特窗:w=Bartlett(n),产生一个n点的巴特利特窗函数。
巴特利特窗系数为w(k)=
巴特利特窗与三角窗非常相似。巴特利特窗在取样点1和n上总以零结束,而三角窗在这些点上并不为零。实际上,当n为奇数时bartlett(n)的中心n-2个点等效于triang(n-2)。
4.汉明窗:w=hamming(n),产生一个n点的汉明窗函数。
汉明窗系数为w(k+1)=0.54-0.46cos( ) k=0,…,n-1
5.汉宁窗:w=hanning(n),产生一个n点的汉宁窗函数。
汉宁窗系数为w(k)=0.5[1-cos( )] k=1,…,n
6.布莱克曼窗:w=Blackman(n),产生一个n点的布莱克曼窗函数。
布莱克曼窗系数为w(k)=0.42-0.5cos(2π )+0.8cos(4π )] k=1,…,n
与等长度的汉明窗和汉宁窗相比,布莱克曼窗的主瓣稍宽,旁瓣稍低。
7.凯泽窗:w=Kaiser(n,beta),产生一个n点的凯泽窗数,其中beta为影响窗函数旁瓣的β参数,其最小的旁瓣抑制α与β的关系为:
0.1102(α-0.87) α>50
β= 0.5842(α-21)0.4+0.07886(α-21) 21≤α≤50
0 α<21
增加β可使主瓣变宽,旁瓣的幅度降低。
8.契比雪夫窗:w=chebwin(n,r)产生一个n点的契比雪夫窗函数。其傅里叶变换后的旁瓣波纹低于主瓣r个db数。
4.4 实验内容
1.软件仿真实验:编写并调试MATLAB程序,观察不同窗,不同类型滤波器不同点数等共4种FIR滤波器的h(n),并记录幅频特性和相频特性。
2.硬件实验:用窗函数法设计标准响应的FIR滤波器,在计算机上观察窗函数幅频特性、幅频特性和相频特性,然后下载到实验箱。用示波器观察输入输出波形,测试滤波器的幅频响应特性。
4.5 MATLAB参考程序和仿真内容
%*******************************************************************%
%mode: 模式(1--高通;2--低通;3--带通;4--带阻)
%n: 阶数,加窗的点数为阶数加1
%fp: 高通和低通时指示截止频率,带通和带阻时指示下限频率
%fs: 带通和带阻时指示上限频率
%window:加窗(1--矩形窗;2--三角窗;3--巴特利特窗;4--汉明窗;
% 5--汉宁窗;6--布莱克曼窗;7--凯泽窗;8--契比雪夫窗)
%r: 代表加chebyshev窗的r值和加kaiser窗时的beta值
%sample: 采样率
%h: 返回设计好的FIR滤波器系数
%*******************************************************************%
%mode: 模式(1--高通;2--低通;3--带通;4--带阻)
%n: 阶数,加窗的点数为阶数加1
%fp: 高通和低通时指示截止频率,带通和带阻时指示下限频率
%fs:
摘要:传统的数字滤波器的设计过程复杂,计算工作量大,滤波特性调整困难,影响了它的应用。本文介绍了一种利用MATLAB信号处理工具箱(Signal Processing Toolbox)快速有效的设计由软件组成的常规数字滤波器的设计方法。给出了使用MATLAB语言进行程序设计和利用信号处理工具箱的FDATool工具进行界面设计的详细步骤。利用MATLAB设计滤波器,可以随时对比设计要求和滤波器特性调整参数,直观简便,极大的减轻了工作量,有利于滤波器设计的最优化。本文还介绍了如何利用MATLAB环境下的仿真软件Simulink对所设计的滤波器进行模拟仿真。
关键词:数字滤波器 MATLAB FIR IIR
引言:
在电力系统微机保护和二次控制中,很多信号的处理与分析都是基于对正弦基波和某些整次谐波的分析,而系统电压电流信号(尤其是故障瞬变过程)中混有各种复杂成分,所以滤波器一直是电力系统二次装置的关键部件【1】。目前微机保护和二次信号处理软件主要采用数字滤波器。传统的数字滤波器设计使用繁琐的公式计算,改变参数后需要重新计算,在设计滤波器尤其是高阶滤波器时工作量很大。利用MATLAB信号处理工具箱(Signal Processing Toolbox)可以快速有效的实现数字滤波器的设计与仿真。
1 数字滤波器及传统设计方法
数字滤波器可以理解为是一个计算程序或算法,将代表输入信号的数字时间序列转化为代表输出信号的数字时间序列,并在转化过程中,使信号按预定的形式变化。数字滤波器有多种分类,根据数字滤波器冲激响应的时域特征,可将数字滤波器分为两种,即无限长冲激响应(IIR)滤波器和有限长冲激响应(FIR)滤波器。
IIR数字滤波器具有无限宽的冲激响应,与模拟滤波器相匹配。所以IIR滤波器的设计可以采取在模拟滤波器设计的基础上进一步变换的方法。FIR数字滤波器的单位脉冲响应是有限长序列。它的设计问题实质上是确定能满足所要求的转移序列或脉冲响应的常数问题,设计方法主要有窗函数法、频率采样法和等波纹最佳逼近法等。
在对滤波器实际设计时,整个过程的运算量是很大的。例如利用窗函数法【2】设计M阶FIR低通滤波器时,首先要根据(1)式计算出理想低通滤波器的单位冲激响应序列,然后根据(2)式计算出M个滤波器系数。当滤波器阶数比较高时,计算量比较大,设计过程中改变参数或滤波器类型时都要重新计算。
(1)
(2)
设计完成后对已设计的滤波器的频率响应要进行校核,要得到幅频相频响应特性,运算量也是很大的。我们平时所要设计的数字滤波器,阶数和类型并不一定是完全给定的,很多时候都是要根据设计要求和滤波效果不断的调整,以达到设计的最优化。在这种情况下,滤波器的设计就要进行大量复杂的运算,单纯的靠公式计算和编制简单的程序很难在短时间内完成设计。利用MATLAB强大的计算功能进行计算机辅助设计,可以快速有效的设计数字滤波器,大大的简化了计算量,直观简便。
2数字滤波器的MATLAB设计
2.1 FDATool界面设计
2.1.1 FDATool的介绍
FDATool(Filter Design Analysis Tool)是MATLAB信号处理工具箱里专用的滤波器设计分析工具,MATLAB6.0以上的版本还专门增加了滤波器设计工具箱(Filter Design Toolbox)。FDATool可以设计几乎所有的基本的常规滤波器,包括FIR和IIR的各种设计方法。它操作简单,方便灵活。
FDATool界面总共分两大部分,一部分是Design Filter,在界面的下半部,用来设置滤波器的设计参数,另一部分则是特性区,在界面的上半部分,用来显示滤波器的各种特性。Design Filter部分主要分为:
Filter Type(滤波器类型)选项,包括Lowpass(低通)、Highpass(高通)、Bandpass(带通)、Bandstop(带阻)和特殊的FIR滤波器。
Design Method(设计方法)选项,包括IIR滤波器的Butterworth(巴特沃思)法、Chebyshev Type I(切比雪夫I型)法、 Chebyshev Type II(切比雪夫II型) 法、Elliptic(椭圆滤波器)法和FIR滤波器的Equiripple法、Least-Squares(最小乘方)法、Window(窗函数)法。
Filter Order(滤波器阶数)选项,定义滤波器的阶数,包括Specify Order(指定阶数)和Minimum Order(最小阶数)。在Specify Order中填入所要设计的滤波器的阶数(N阶滤波器,Specify Order=N-1),如果选择Minimum Order则MATLAB根据所选择的滤波器类型自动使用最小阶数。
Frenquency Specifications选项,可以详细定义频带的各参数,包括采样频率Fs和频带的截止频率。它的具体选项由Filter Type选项和Design Method选项决定,例如Bandpass(带通)滤波器需要定义Fstop1(下阻带截止频率)、Fpass1(通带下限截止频率)、Fpass2(通带上限截止频率)、Fstop2(上阻带截止频率),而Lowpass(低通)滤波器只需要定义Fstop1、Fpass1。采用窗函数设计滤波器时,由于过渡带是由窗函数的类型和阶数所决定的,所以只需要定义通带截止频率,而不必定义阻带参数。
Magnitude Specifications选项,可以定义幅值衰减的情况。例如设计带通滤波器时,可以定义Wstop1(频率Fstop1处的幅值衰减)、Wpass(通带范围内的幅值衰减)、Wstop2(频率Fstop2处的幅值衰减)。当采用窗函数设计时,通带截止频率处的幅值衰减固定为6db,所以不必定义。
Window Specifications选项,当选取采用窗函数设计时,该选项可定义,它包含了各种窗函数。
2.1.2 带通滤波器设计实例
本文将以一个FIR 滤波器的设计为例来说明如何使用MATLAB设计数字滤波器:在小电流接地系统中注入83.3Hz的正弦信号,对其进行跟踪分析,要求设计一带通数字滤波器,滤除工频及整次谐波,以便在非常复杂的信号中分离出该注入信号。参数要求:96阶FIR数字滤波器,采样频率1000Hz,采用Hamming窗函数设计。
本例中,首先在Filter Type中选择Bandpass(带通滤波器);在Design Method选项中选择FIR Window(FIR滤波器窗函数法),接着在Window Specifications选项中选取Hamming;指定Filter Order项中的Specify Order=95;由于采用窗函数法设计,只要给出通带下限截止频率Fc1和通带上限截止频率Fc2,选取Fc1=70Hz,Fc2=84Hz。设置完以后点击Design Filter即可得到所设计的FIR滤波器。通过菜单选项Analysis可以在特性区看到所设计滤波器的幅频响应、相频响应、零极点配置和滤波器系数等各种特性。设计完成后将结果保存为1.fda文件。
在设计过程中,可以对比滤波器幅频相频特性和设计要求,随时调整参数和滤波器类型,
以便得到最佳效果。其它类型的FIR滤波器和IIR滤波器也都可以使用FDATool来设计。
图1 滤波器幅频和相频响应(特性区)
Fig.1 Magnitude Response and Phase Response of the filter
2.2 程序设计法
在MATLAB中,对各种滤波器的设计都有相应的计算振幅响应的函数【3】,可以用来做滤波器的程序设计。
上例的带通滤波器可以用程序设计:
c=95%定义滤波器阶数96阶
w1=2*pi*fc1/fs
w2=2*pi*fc2/fs%参数转换,将模拟滤波器的技术指标转换为数字滤波器的技术指标
window=hamming(c+1)%使用hamming窗函数
h=fir1(c,[w1/pi w2/pi],window)%使用标准响应的加窗设计函数fir1
freqz(h,1,512)%数字滤波器频率响应
在MATLAB环境下运行该程序即可得到滤波器幅频相频响应曲线和滤波器系数h。篇幅所限,这里不再将源程序详细列出。
3 Simulink仿真
本文通过调用Simulink中的功能模块构成数字滤波器的仿真框图,在仿真过程中,可以双击各功能模块,随时改变参数,获得不同状态下的仿真结果。例如构造以基波为主的原始信号,,通过Simulink环境下的Digital Filter Design(数字滤波器设计)模块导入2.1.2中FDATool所设计的滤波器文件1.fda。仿真图和滤波效果图如图2所示。
图2 Simulink仿真图及滤波效果图
Fig.2 Simulated connections and waveform
可以看到经过离散采样、数字滤波后分离出了83.3Hz的频率分量(scope1)。之所以选取上面的叠加信号作为原始信号,是由于在实际工作中是要对已经经过差分滤波的信号进一步做带通滤波,信号的各分量基本同一致,可以反映实际的情况。本例设计的滤波器已在实际工作中应用,取得了不错的效果。
4 结论
利用MATLAB的强大运算功能,基于MATLAB信号处理工具箱(Signal Processing Toolbox)的数字滤波器设计法可以快速有效的设计由软件组成的常规数字滤波器,设计方便、快捷,极大的减轻了工作量。在设计过程中可以对比滤波器特性,随时更改参数,以达到滤波器设计的最优化。利用MATLAB设计数字滤波器在电力系统二次信号处理软件和微机保护中,有着广泛的应用前景。
参考文献
1. 陈德树. 计算机继电保护原理与技术【M】北京:水利电力出版社,1992.
2. 蒋志凯. 数字滤波与卡尔曼滤波【M】北京:中国科学技术出版社,1993
3. 楼顺天、李博菡. 基于MATLAB的系统分析与设计-信号处理【M】西安:西安电子科技大学出版社,1998.
4. 胡广书. 数字信号处理:理论、算法与实现【M】.北京:清华大学出版社,1997.
5. 蒙以正. MATLAB5.X应用与技巧【M】北京:科学出版社,1999.
在工程设计中,高通、带通和带阻滤波器的设计通常是利用低通滤波器的原型,经过频率变换得到。这样就转变为一个低通滤波器的设计。
不过现在滤波器的设计一般都是用软件了,参数的计算很方便,不用懂太多的理论,稍微有点基础就行,仿真出结果,在做板。当然,最重要的还是最后的调试,也就是这里需要些理论基础。不过这都是建立在 滤波器设计软件 仿真软件 PCB用好的情况下。还是建议你多学些理论,容易找出不足,完善自己的设计。Filter Solutions 2009一个不错的滤波软件,你如果用,我有软件和中文资料