用matlab设计数字滤波器

用matlab设计一个滤波器~

你自己整合吧，我没时间帮你整合，我给你提供一些程序：
绝对正确的代码：程序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.1 实验目的

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/sample;wp2=2*fp2/sample;

ws1=2*fs1/sample;ws2=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:

产品名称：K&L可调带阻滤波器3TNF-200/400-N/N
频率范围（MHz）：200-400
3dB带宽范围（MHz）：3-7
40dB带宽（Min KHz）：300
切槽深度：50 dB
长度（Inch/mm）：6.56/167
宽度（Inch/mm）：5.38/137
高度（Inch/mm）：2.75/70
驻波比：1.5:1
插入损耗：＜ 0.5 dB
阻抗：50 Ohms
功率：50 Watts CW
刻度精确度：+/- 1%
连接器形式：Type N请访问：上海欧桥微波网

摘要

《数字信号处理》课程是一门理论性和实践性都很强， 它具备高等代数、数值分析、概率统计、随机过程等计算学科的知识; 要求我们学生掌握扎实的基础知识和理论基础。 又是跟其他学科密切相关，即与通信理论、计算机、微电子技术不可分，又是人工智能、模式识别、神经网络等新兴学科的理论基础之一。 本次数字滤波器设计方法是基于MATLAB的数字滤波器的设计。此次设计的主要内容为：IIR数字滤波器和FIR数字滤波器的设计
关键词：IIR、FIR、低通、高通、带阻、带通

Abstract
"Digital Signal Processing" is a theoretical and practical nature are strong, and it has advanced algebra and numerical analysis, probability and statistics, random process such as calculation of discipline knowledge; requires students to acquire basic knowledge and a solid theoretical basis. Is closely related with other subjects, namely, and communication theory, computers, microelectronics can not be separated, but also in artificial intelligence, pattern recognition, neural network theory one of the emerging discipline. The digital filter design method is based on MATLAB for digital filter design. The main elements of design: IIR and FIR digital filter design of digital filter
Key Words: IIR, FIR, low pass, high pass, band stop, band pass

目录
一、 前言 3
二、 课程设计的目的 3
三、 数字信号处理课程设计说明及要求 3
四、 滤波器的设计原理 4
4.1 数字滤波器简介 4
4.2 IIR滤波器的设计原理 4
4.3 FIR滤波器的设计原理 5
4.4 FIR滤波器的窗函数设计法 6
五、 设计内容 6
5.1 设计题目： 6
5.2设计程序代码及结果： 7
六、 结束语 15
七、 参考文献 16

一、 前言

数字信号处理（Digital Signal Processing，简称DSP）是一门涉及许多学科而又广泛应用于许多领域的新兴学科。
随着信息时代和数字世界的到来，数字信号处理已成为今一门极其重要的学科和技术领域。数字信号处理在通信语音、图像、自动控制、雷达、军事、航空航天、医疗和家用电器等众多领域得到了广泛的应用。在数字信号处理应用中，数字滤波器十分重要并已获得广泛应用。

二、 课程设计的目的
1)

三、 数字信号处理课程设计说明及要求
所需硬件：PC机
四、 滤波器的设计原理
4.1 数字滤波器简介
数字滤波器是一种用来过滤时间离散信号的数字系统，通过对抽样数据进行数学处理来达到频域滤波的目的。可以设计系统的频率响应，让它满足一定的要求，从而对通过该系统的信号的某些特定的频率成分进行过滤，这就是滤波器的基本原理。如果系统是一个连续系统，则滤波器称为模拟滤波器。如果系统是一个离散系统，则滤波器称为数字滤波器。
信号 通过线性系统后，其输出 就是输入信号 和系统冲激响应 的卷积。除了 外， 的波形将不同于输入波形 。从频域分析来看，信号通过线性系统后，输出信号的频谱将是输入信号的频谱与系统传递函数的乘积。除非 为常数，否则输出信号的频谱将不同于输入信号的频谱，某些频率成分 较大的模，因此， 中这些频率成分将得到加强，而另外一些频率成分 的模很小甚至为零， 中这部分频率分量将被削弱或消失。因此，系统的作用相当于对输入信号的频谱进行加权。
4.2 IIR滤波器的设计原理
IIR数字滤波器的设计一般是利用目前已经很成熟的模拟滤波器的设计方法来进行设计，通常采用模拟滤波器原型有butterworth函数、chebyshev函数、bessel函数、椭圆滤波器函数等。
IIR数字滤波器的设计步骤：
（1） 按照一定规则把给定的滤波器技术指标转换为模拟低通滤波器的技术指标；
（2） 根据模拟滤波器技术指标设计为响应的模拟低通滤波器；
（3） 很据脉冲响应不变法和双线性不变法把模拟滤波器转换为数字滤波器；
（4） 如果要设计的滤波器是高通、带通或带阻滤波器，则首先把它们的技术指标转化为模拟低通滤波器的技术指标，设计为数字低通滤波器，最后通过频率转换的方法来得到所要的滤波器。
4.3 FIR滤波器的设计原理
FIR滤波器通常采用窗函数方法来设计。窗设计的基本思想是，首先选择一个适当的理想选频滤波器（它总是具有一个非因果，无限持续时间脉冲响应），然后街区（加窗）它的脉冲响应得到线性相位和因果FIR滤波器。我们用Hd（e^jw）表示理想的选频滤波器，它在通带上具有单位增益和线性相位，在阻带上具有零响应。一个带宽wc<pi的低通滤波器由下式给定：

为了从hd(n)得到一个FIR滤波器，必须同时在两边截取hd（n）。而要得到一个因果的线性相位滤波器，它的h(n)长度为N，必须有：

这种操作叫做加窗，h(n)可以看做是hd(n)与窗函数w（n）的乘积：
h(n)=hd(n)w(n)
其中

根据w（n）的不同定义，可以得到不同的窗结构。
在频域中，因果FIR滤波器响应H(e^jw)由Hd(e^jw)和窗响应W(e^jw)的周期卷积得到，即

常用的窗函数有矩形窗、巴特利特（BARTLETT）窗、汉宁（HANNING）窗、海明（HAMMING）窗、布莱克曼（BLACKMAN）窗、凯泽（KAISER）窗等。
4.4 FIR滤波器的窗函数设计法
FIR滤波器的设计方法有许多种，如窗函数设计法、频率采样设计法和最优化设计法等。窗函数设计法的基本原理是用一定宽度窗函数截取无限脉冲响应序列获得有限长的脉冲响应序列，主要设计步骤为：
(1) 通过傅里叶逆变换获得理想滤波器的单位脉冲响应hd(n)。

(2) 由性能指标确定窗函数W(n)和窗口长度N。
(3) 求得实际滤波器的单位脉冲响应h(n)， h(n)即为所设计FIR滤波器系数向量b(n)。

五、 设计内容
5.1 设计题目：
1-1.试用MATLAB设计一巴特沃斯低通数字滤波器，要求通带截至频率Wp=30HZ，主带截至频率为Ws=35HZ，通带衰减不大于0.5DB，主带衰减不小于40DB，抽样频Fs=100HZ。

1-2．基于Butterworth模拟滤波器原型，使用双线性状换设计数字滤波器：各参数值为：通带截止频率Omega=0.2*pi,阻带截止频率Omega=0.3*pi,通带波动值Rp=1dB,阻带波动值Rs=15dB,设Fs=20000Hz。

1-3设计一巴特沃斯高通数字滤波器，要求通带截止频率0.6*pi，通带衰减不大于1dB，阻带衰减15DB，抽样T=1。

1-4.设计一巴特沃斯带阻数字滤波器，要求通带上下截至频率为0.8*PI、0.2*PI，通带衰减不大于1DB，阻带上下截至频率0.7*PI、0.4*PI 阻带衰减不小于30DB，

2-1.用窗函数法设计一个线性相位FIR低通滤波器，并满足性能指标:通带边界频率
Wp=0.5*pi,阻带边界频率Ws=0.66*pi,阻带衰减不小于40dB,通带波纹不大于3dB。选择汉宁窗。

2-4．用海明窗设计一个FIR滤波器，其中Wp=0.2*pi，Ws=0.3*pi,通带衰减不大于0.25dB,阻带衰减不小于50dB。

5.2设计程序代码及结果：
1-1一.试用MATLAB设计一巴特沃斯低通数字滤波器，要求通带截至频率Wp=30HZ，阻带截至频率为Ws=35HZ，通带衰减不大于0.5DB，阻带衰减不小于40DB，抽样频Fs=100HZ。
代码为：
fp = 30;
fs = 35;
Fs = 100;
wp = 2*pi*fp/Fs;
ws = 2*pi*fs/Fs;
wp = tan(wp/2);
ws = tan(ws/2); % 通带最大衰减为0.5dB，阻带最小衰减为40dB
[N, wn] = buttord(wp, ws, 0.5, 40, 's'); % 模拟低通滤波器极零点
[z, p, k] = buttap(N); % 由极零点获得转移函数参数
[b, a] = zp2tf(z, p, k); % 由原型滤波器获得实际低通滤波器
[B, A] = lp2lp(b, a, wp);
[bz, az] = bilinear(B, A, .5);
[h, w] = freqz(bz, az, 256, Fs);
figure
plot(w, abs(h))
grid on

图1 巴特沃斯数字低通滤波器
1-2基于Butterworth模拟滤波器原型，使用双线性状换设计数字滤波器：各参数值为：通带截止频率Omega=0.2*pi,阻带截止频率Omega=0.3*pi,通带波动值Rp=1dB,阻带波动值Rs=15dB,设Fs=4000Hz。
代码：
wp=0.2*pi;ws=0.3*pi;
Fs=4000;T=1/Fs;
OmegaP=(2/T)*tan(wp/2);
OmegaS=(2/T)*tan(ws/2);
rp=1;rs=15;as=15;
ripple=10^(-rp/20);attn=10^(-rs/20);
[n,wn]=buttord(OmegaP,OmegaS,rp,rs,'s');
[z,p,k]=Buttap(n);
[b,a]=zp2tf(z,p,k);
[bt,at]=lp2lp(b,a,wn);
[b,a]=bilinear(bt,at,Fs);
[db,mag,pha,grd,w]=freqz_m(b,a);
%
%下面绘出各条曲线
subplot(2,2,1);plot(w/pi,mag);title('Magnitude Frequency幅频特性');
xlabel('w(/pi)');ylabel('|H(jw)|');
axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 1]);
set(gca,'YTickMode','manual','YTick',[0 attn ripple 1]);grid

subplot(2,2,2);plot(w/pi,db);title('Magnitude Frequency幅频特性(db)');
xlabel('w(/pi)');ylabel('dB');
axis([0,1,-30,5]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 1]);
set(gca,'YTickMode','manual','YTick',[-60 -as -rp 0]);grid

subplot(2,2,3);plot(w/pi,pha/pi);title('Phase Frequency相频特性');
xlabel('w(/pi)');ylabel('pha(/pi)');
axis([0,1,-1,1]);

subplot(2,2,4);plot(w/pi,grd);title('Group Delay群延时');
xlabel('w(/pi)');ylabel('Sample');
axis([0,1,0,15]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 1]);grid
运行结果：

图2巴特沃思数字低通滤波器幅频-相频特性

1-3设计一巴特沃斯高通数字滤波器，要求通带截止频率0.6*pi，通带衰减不大于1dB，阻带衰减15DB，抽样T=1。
Wp=0.6*pi;
Ws=0.4*pi;
Ap=1;
As=15;
[N,wn]=buttord(Wp/pi,Ws/pi,Ap,As); %计算巴特沃斯滤波器阶次和截止频率
%频率变换法设计巴特沃斯高通滤波器
[db,mag,pha,grd,w]=freqz_m(b,a); %数字滤波器响应
plot(w,mag);
title('数字滤波器幅频响应|H(ej\Omega)|')

图3巴特沃斯数字高通滤波器

2-1用窗函数法设计一个线性相位FIR低通滤波器，并满足性能指标:通带边界频率
Wp=0.5*pi,阻带边界频率Ws=0.66*pi,阻带衰减不小于40dB,通带波纹不大于3dB。选择汉宁窗。
代码：
wp =0.5*pi;
ws=0.66*pi;
wdelta =ws-wp;
N= ceil(8*pi/wdelta)
if rem(N,2)==0
N=N+1;
end

);
运行结果：

给分就给你个全的!

图6低通FIR滤波器

六、 结束语
本次数字滤波器设计方法是基于MATLAB的数字滤波器的设计，是用学过的数字信号理论为依据，用MATLAB代码来实现。课程设计过程中，通过IIR数字滤波器和FIR数字滤波器的设计实例，说明如何利用MATLAB来完成数字滤波器的设计。窗函数法中相位响应有严格的线性,不存在稳定性问题, 设计简单。双线性变换不会出现由于高频部分超过折叠频率而混淆到低频部分去的现象，但会产生频率混碟现象，使数字滤波器的频响偏移模拟滤波器的频响。由滤波器的频谱图和滤波前后的语音信号的频谱图对比可知本设计选用双线性变换法设计的IIR滤波器比较好。在同样的技术指标的要求下，IIR滤波器所要求的阶数N也比较小，实现起来比较容易。
通过综合运用数字信号处理的理论知识进行滤波器设计，通过理论推导得出相应结论，再利用 MATLAB 作为编程工具进行计算机实现，从而加深了对所学知识的理解，建立概念。对以前在课本上所学的东西有了更深入的理解和掌握。最后，无论做什么课程设计，都需要要有一定的理论知识作为基础，同时通过这次课程设计，我对于以前所学的数字信号处理知识有了更深的理解。

七、 参考文献
1. 程佩青《数字信号处理教程》北京清华大学出版社2007年2月.
2. 赵知劲、刘顺兰《数字信号处理实验》.浙江大学出版社.
3. S.K.Mitra．Digital Signal Processing:A Computer-Based Approach.
NewYork,NewYork:McGraw-Hill,thirded,2006．
4. 肖伟、刘忠等《 MATLAB程序设计与应用》清华大学出版社、北京交通大学出版社.
5. 胡良剑、孙晓君 《 MATLAB数学实验》.高等教育出版社.

我有做好的。要不

---

城北区13845724884： matlab 数字带通滤波器 - ？
井宰牛黄： 以前我做的一个完整的信号滤波前后的程序,可能数据有些不一样,但用MATLAB做滤波器的大体思路都差不多,希望对你有用.(完全自创,要给分哦) f=8*10^3; %信号频率 f1=30*10^3; %噪声频率 l=100; %信号的长度 fs=100*10^3; ...

城北区13845724884： 如何在MATLAB中设计CIC滤波器 - ？
井宰牛黄： 刚好最近在做这个项目.有好几种方法.滤波器Matlab有专门设计滤波器的工具,叫fdatool.你在命令栏输入 fdatool 就会打开设计的界面.点击左下角的第一个按钮,进入CIC滤波器的设计界面.填入参数指标以后,就可以生成滤波器了.

城北区13845724884： 如何用自己设计好的滤波器对频谱滤波 matlab - ？
井宰牛黄：(1)求出信号的频谱X(W)和滤波器的频率特性H(W),即把它们变换到频域去(2) 频域相乘,即Y(W)=X(W)*H(W),相当于时域做卷积.(可用分段卷积法)(3)把Y(W)反变换到时域为y(n)即可滤波器是一种用来消...

城北区13845724884： 利用matlab设计出数字滤波器后,怎么在现实中实现? - ？
井宰牛黄： 将随机信号输入到matlab中,与数字滤波器相卷积,得到输出信号

城北区13845724884： 数字切比雪夫滤波器MATLAB代码急求 - ？
井宰牛黄： 用Matlab设计切比雪夫滤波器是很简单的,比如设计I型的 N = 8; %阶数 wp=0.2; %归一化截止频率 Rp= 1; %通带波纹 [b a]=cheby1(N,Rp,wp,'low'); %设计低通滤波器 freqz(b,a) %频率响应曲线

城北区13845724884： 在MATLAB环境下设计带阻滤波器的程序 - ？
井宰牛黄： clear all;fp=[200 500];通带上下限 fs=[300 400];阻带 rp=3;rs=18;衰减db Fs=2000; wp=fp*2*pi/Fs; ws=fs*2*pi/Fs; % % Firstly to finish frequency prewarping; wap=2*Fs*tan(wp./2) was=2*Fs*tan(ws./2); [n,wn]=buttord(wap,was,rp,rs,'s'); % Note: ...

城北区13845724884： Matlab编程实现数字滤波器怎么设计 - ？
井宰牛黄： 你好好查点资料,这个有太多东西可以写,我给你个简单的例子:Wc=0.5*pi; N=21; window=hanning(N+1); %选择要用的窗函数 b=fir1(N,Wc/pi,window); %加窗的线性相位设计函数 [H,W]=freqz(b,1,512); %求幅频响应 subplot(3,1,1) plot(W/pi,20*log...

城北区13845724884： 关于matlab设计数字滤波器 - ？
井宰牛黄： T是采样周期,单位是秒,Fs是采样频率,是周期T的倒数

城北区13845724884： 用MATLAB设计数字滤波器,分析脑电波信号!急求!只要能滤波就行! - ？
井宰牛黄： 用simulink的digital filter design吧,自动设计很方便.

城北区13845724884： 联系matlab用双线性变换法设计Butterworth低通滤波器m - ？
井宰牛黄： %下面程序可能对你有用,需根据你的指标修改一下前5行相关e79fa5e98193e58685e5aeb931333264656133参数就行,至于最后显示部分,相应变化一下即可.%1.用Butterworth低通模拟滤波器原型%%%%%%%%%%%%%%%%%%%%%...