第一篇:数字信号处理实验报告完整版
实验 1
利用 T DFT 分析信号频谱
一、实验目的
1.加深对 DFT 原理的理解。
2.应用 DFT 分析信号的频谱。
3.深刻理解利用 DFT 分析信号频谱的原理,分析实现过程中出现的现象及解决方法。
二、实验设备与环境
计算机、MATLAB 软件环境 三、实验基础理论
T 1.DFT 与 与 T DTFT 的关系
有限长序列 的离散时间傅里叶变换 在频率区间的 N 个等间隔分布的点 上的 N 个取样值可以由下式表示:
212 /0()|()()0 1Nj knjNk NkX e x n e X k k N 由上式可知,序列 的 N 点 DFT ,实际上就是 序列的 DTFT 在 N 个等间隔频率点 上样本。
2.利用 T DFT 求 求 DTFT
方法 1 1:由恢复出的方法如下:
由图 2.1 所示流程可知:
101()()()Nj j n kn j nNn n kX e x n e X k W eN 由上式可以得到:
IDFT DTFT
第二篇:数字信号处理实验报告
南京邮电大学
实 验 报 告
实验名称_____熟悉MATLAB环境 ___ 快速傅里叶变换及其应用 ____IIR数字滤波器的设计_ FIR数字滤波器的设计
课程名称 数字信号处理A
班级学号_______09002111___________ 姓 名 王都超
开课时间 2011/2012学年,第 二 学期
实验一
熟悉MATLAB环境
一、实验目的
(1)熟悉MATLAB的主要操作命令。(2)学会简单的矩阵输入和数据读写。(3)掌握简单的绘图命令。
(4)用MATLAB编程并学会创建函数。(5)观察离散系统的频率响应。
二、实验内容
(1)数组的加、减、乘、除和乘方运算。输入A=[1 2 3 4],B=[3,4,5,6],求
C=A+B,D=A-B,E=A.*B,F=A./B,G=A.^B。并用stem语句画出A、B、C、D、E、F、G。
D =
-2 E =
F =
0.3333
0.5000
0.6000
0.6667 G =
243
4096(2)用MATLAB实现下列序列: a)x(n)0.8n 0n1
5n=0:1:15;x1=0.8.^n;a=(0.2+3*i)*n;stem(x1)b)x(n)e(0.23j)n 0n15
n=0:1:15;x2=exp(a);a=(0.2+3*i)*n;stem(x2)
c)x(n)3cos(0.125n0.2)2sin(0.25n0.1)
0n15
(4)绘出下列时间函数的图形,对x轴、y轴以及图形上方均须加上适当的标注: a)x(t)sin(2t)0t10s
b)x(t)cos(100t)sin(t)0t4s t=0:0.01:4;x=cos(100*pi*t).*sin(pi*t);plot(t,x, 'r-');xlabel('t'),ylabel('x(t)'),title('cos')
(6)给定一因果系统H(z)(1频响应和相频响应。
2z1z2)/(10.67z10.9z2),求出并绘制H(z)的幅
(7)计算序列{8-2-1 2 3}和序列{2 3-1-3}的离散卷积,并作图表示卷积结果。
(8)求以下差分方程所描述系统的单位脉冲响应h(n), 0n50
y(n)0.1y(n1)0.06y(n2)x(n)2x(n1)
实验过程与结果(含实验程序、运行的数据结果和图形); clear all;N=50;a=[1-2];b=[1 0.1-0.06];x1=[1 zeros(1,N-1)];n=0:1:N-1;h=filter(a,b,x1);stem(n,h)axis([-1 53-2.5 1.2])
实验二
快速傅里叶变换及其应用
一、实验目的
(1)在理论学习的基础上,通过本实验,加深对FFT的理解,熟悉MATLAB中的有关函数。(2)应用FFT对典型信号进行频谱分析。
(3)了解应用FFT进行信号频谱分析过程中可能出现的问题,以便在实际中正确应用FFT。(4)应用FFT实现序列的线性卷积和相关。
二、实验内容
实验中用到的信号序列 a)高斯序列
(np)qxa(n)e020n15 其他
b)衰减正弦序列
eansin(2fn)xb(n)00n15其他
c)三角波序列 nxc(n)8n00n34n7 其他
d)反三角波序列
4nxd(n)n400n34n7 其他
(1)观察高斯序列的时域和幅频特性,固定信号xa(n)中参数p=8,改变q的值,使q分别等于2,4,8,观察它们的时域和幅频特性,了解当q取不同值时,对信号序列的时域幅频特性的影响;固定q=8,改变p,使p分别等于8,13,14,观察参数p变化对信号序列的时域及幅频特性的影响,观察p等于多少时,会发生明显的泄漏现象,混叠是否也随之出现?记录实验中观察到的现象,绘出相应的时域序列和幅频特性曲线。
(3)观察三角波和反三角波序列的时域和幅频特性,用N=8点FFT分析信号序列xc(n)和观察两者的序列形状和频谱曲线有什么异同?绘出两序列及其幅频特性xd(n)的幅频特性,曲线。
在xc(n)和xd(n)末尾补零,用N=32点FFT分析这两个信号的幅频特性,观察幅频特性发生了什么变化?两种情况的FFT频谱还有相同之处吗?这些变化说明了什么?
(5)用FFT分别实现xa(n)(p=8,q=2)和xb(n)(a=0.1,f=0.0625)的16点循环卷积和线性卷积。
n=0:15;p=8;q=2;
xa=exp(-(n-p).^2/q);subplot(2,3,1);stem(n,xa,'.');title('xa波形');
Xa=fft(xa,16);subplot(2,3,4);stem(abs(Xa),'.');
title('Xa(k)=FFT[xa(n)]的波形 ');A=1;f=0.0625;a=0.1;
xb=exp(-a*n).*sin(2*pi*f*n);subplot(2,3,2);stem(n,xb,'.');title('xb波形');Xb=fft(xb,16);subplot(2,3,5);stem(abs(Xb),'.');
title('Xb(k)=FFT[xb(n)]的波形 ');
实验过程与结果(含实验程序、运行的数据结果和图形);
实验三 IIR数字滤波器的设计
一、实验目的
(1)掌握双线性变换法及脉冲响应不变法设计IIR数字滤波器的具体设计方法及其原理,熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。
(2)观察双线性变换及脉冲响应不变法设计的滤波器的频域特性,了解双线性变换法及脉冲响应不变法的特点。
(3)熟悉巴特沃思滤波器、切比雪夫滤波器和椭圆滤波器的频率特性。
二、实验内容(1)P162 例4.4 设采样周期T=250s(采样频率fs=4kHz),分别用脉冲响应不变法和双线性变换法设计一个三阶巴特沃思低通滤波器,其3dB边界频率为fc=1kHz。
脉冲响应不变法: fc=1000;fs=4000;OmegaC=2*pi*fc;[B,A]=butter(3, OmegaC,'s');[num1,den1]=impinvar(B,A,fs);[h1,w]=freqz(num1,den1);f = w/pi*fs/2;plot(f,abs(h1));
双线性变换法: fc=1000;fs=4000;
OmegaC=2*fs*tan(pi*fc/fs);[B,A]=butter(3, OmegaC,'s');[num2,den2]=bilinear(B,A,fs);[h2,w]=freqz(num2,den2);f = w/pi*fs/2;plot(f,abs(h2));
同一图中画两条曲线: fc=1000;fs=4000;OmegaC=2*pi*fc;[B,A]=butter(3, OmegaC,'s');[num1,den1]=impinvar(B,A,fs);[h1,w]=freqz(num1,den1);f = w/pi*fs/2;
OmegaC=2*fs*tan(pi*fc/fs);[B,A]=butter(3, OmegaC,'s');[num2,den2]=bilinear(B,A,fs);[h2,w]=freqz(num2,den2);f = w/pi*fs/2;plot(f,abs(h1),'r-.');hold on;plot(f,abs(h2),'g-');
(选做)(2)fc=0.2kHz,=1dB,fr=0.3kHz,At=25dB,T=1ms;分别用脉冲响应不变法及双线性变换法设计一巴特沃思数字低通滤波器,观察所设计数字滤波器的幅频特性曲线,记录带宽和衰减量,检查是否满足要求。比较这两种方法的优缺点。
实验过程与结果(含实验程序、运行的数据结果和图形);
实验四
FIR数字滤波器的设计
一、实验目的
(1)掌握用窗函数法,频率采样法及优化设计法设计FIR滤波器的原理及方法,熟悉相应的计算机编程;
(2)熟悉线性相位FIR滤波器的幅频特性和相频特性;
(3)了解各种不同窗函数对滤波器性能的影响。
二、实验内容
(1)生成一个长度为20的矩形窗,画出其时域和幅频特性曲线。n=0:1:19;N=20;win(1:20)=1;[H,w]=freqz(win,1);subplot(2,1,1);stem(n,win)subplot(2,1,2);plot(w,abs(H));
(2)用矩形窗设计一个21阶的线性相位低通FIR数字滤波器,截止频率Wc=0.25π,求出滤波器系数,并绘出滤波器的幅频特性。修改程序,分别得到阶次为N=41,61的滤波器,并显示其各自的幅频曲线。
a)在上面所得的几幅图中,在截止频率两边可以观察到幅频响应的摆动行为。请问波纹的数量与滤波器脉冲响应的长度之间有什么关系?
b)最大波纹的高度与滤波器脉冲响应的长度之间有什么关系?
实验过程与结果(含实验程序、运行的数据结果和图形); 21阶的线性相位低通FIR数字滤波器: Wc=0.25*pi;N=21;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函数
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
41阶的线性相位低通FIR数字滤波器: Wc=0.25*pi;N=41;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函数
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
61阶的线性相位低通FIR数字滤波器: Wc=0.25*pi;N=61;M=(N-1)/2;
%位移量
for n=0:(N-1)
if(n== fix(M))
%中间的点单独算
hd(n+1)=Wc/pi;
else
hd(n+1)=sin(Wc*(n-M))/(pi*(n-M));end;end;win=boxcar(N);%%%不同窗函数
h=hd.*win';[H,w]=freqz(h,1);n=0:1:N-1;subplot(3,1,1);stem(n,h)subplot(3,1,2);plot(w,abs(H));subplot(3,1,3);plot(w,angle(H));
数字信号处理实验小结及心得体会:
通过这次实验,我对MATLAB语言有了一定的认识,虽然还不能完全用MATLAB独立编写程序,但对这种语言环境有了新的了解。我知道了一般的加减乘除在MATLAB中不同的意义。知道输入、输出语句怎么形成。通过快速傅里叶变换及其应用的实验,加深了我对FFT的理解,还有对各典型信号的频谱分析,改变参数后时域和幅频特性的变化。IIR数字滤波器的设计让我知道了巴特沃思滤波器和切比雪夫滤波器的频率特性,还有双线性变换及脉冲响应不变法设计的滤波器的频率特性。做这个实验的时候程序有点困难,很多细节问题不能考虑清楚,导致图形出不来。FIR数字滤波器的设计出来的是三种窗的图形,通过三种窗的比较,我了解了他们各自的特点,幅频和相频特性。我在这次实验中的收获很大,接触了很多新的知识,但在实验写程序时,我发现自己还有很多不足。很多程序写不完全。这是自己今后要加强的地方。
第三篇:数字信号处理实验报告
JIANGSU
UNIVERSITY OF TECHNOLOGY
数字信号处理实验报告
学院名称: 电气信息工程学院
专 业:
班 级: 姓 名: 学 号: 指导老师: 张维玺(教授)
2013年12月20日
实验一 离散时间信号的产生
一、实验目的
数字信号处理系统中的信号都是以离散时间形态存在的,所以对离散时间信号的研究是数字信号的基本所在。而要研究离散时间信号,首先需要产生出各种离散时间信号。使用MATLAB软件可以很方便地产生各种常见的离散时间信号,而且它还具有强大绘图功能,便于用户直观地处理输出结果。
通过本实验,学生将学习如何用MATLAB产生一些常见的离散时间信号,实现信号的卷积运算,并通过MATLAB中的绘图工具对产生的信号进行观察,加深对常用离散信号和信号卷积和运算的理解。
二、实验原理
离散时间信号是指在离散时刻才有定义的信号,简称离散信号,或者序列。离散序列通常用x(n)来表示,自变量必须是整数。常见的离散信号如下:(1)单位冲激序列δ(n)
如果δ(n)在时间轴上延迟了k个单位,得到δ(n-k),即长度为N的单位冲激序列δ(n)可以通过下面的MATLAB命令获得。
n=-(N-1):N-1 x=[zeros(1,N-1)1 zeros(1,N-1)]; stem(n,x)延迟K个采样点的长度为N的单位冲激序列δ(n-k)(k n=0:N-1 y=[zeros(1,M)1 zeros(1,N-M-1)]; stem(n,y) (2)单位阶跃序列u(n) 如果u(n)在时间轴上延迟了k个单位,得到u(n-k),即长度为N的单位阶跃序列u(n)可以通过下面的MATLAB命令获得。 n=-(N-1):N-1 x=[zeros(1,N-1)ones(1,N)]; stem(n,x)延迟的单位阶跃序列可以使用类似于单位冲激序列的方法获得。(3)矩形序列 矩形序列有一个重要的参数,就是序列的宽度N。矩形序列与u(n)之间的关系为矩形序列等= u(n)— u(n-N)。 因此,用MATLAB表示矩形序列可利用上面的单位阶跃序列组合而成。(4)正弦序列x(n) 这里,正弦序列的参数都是实数。与连续的正弦信号不同,正弦序列的自变量n必须为整数。可以证明,只有当2π/w为有理数时,正弦序列具有周期性。 长度为N的正弦序列x(n)可以通过下面的MATLAB命令获得。n=0:N-1 x=A*cos(2*pi*f*n/Fs+phase)(5)单边实指数序列x(n) 长度为N的实指数序列x(n)可以通过下面的MATLAB命令实现。n=0:N-1 x=a.^n stem(n,x)单边指数序列n的取值范围为n>=0。当|a|>1时,单边指数序列发散;当|a|<1时,单边指数序列收敛。当a>0时,该序列均取正值;当a<0时,序列在正负摆动。 (6)负指数序列x(n) 当a=0时,得到虚指数序列x(n)。 与连续负指数信号一样,我们将负指数序列实部和虚部的波形分开讨论,得到如下结论: 1)当a>0时,负指数序列x(n)的实部和虚部分别是按指数规律增长的正弦振荡序列; 2)当a<0时,负指数序列x(n)的实部和虚部分别是按指数规律衰减的正弦振荡序列; 3)当a=0时,负指数序列x(n)即为虚指数序列,其实部和虚部分别是等幅的正弦振荡序列; 长度为N的实指数序列x(n)可以通过下面的MATLAB命令实现。n=0:N-1 x=exp((a.+j*w)*n)stem(n,real(x))或 stem(n,imag(x)) 三、实验内容及分析 1n01、编制程序产生单位冲激序列n“并绘出其图及n”学号后两位0n0形。程序:(1)N=4; n=-(N-1):N-1; x=[zeros(1,N-1)1 zeros(1,N-1)];stem(n,x); title('单位冲激序列'); grid on; (2)N=6; M=1;%学号01 n=-(N-1):N-1; y=[zeros(1,N-M+1)1 zeros(1,N-M-1)];stem(n,y); title('单位冲激序列');grid on; 分析:在上图的基础上向右平移了1个单位。 1n02、编制程序产生单位阶跃序列un、un“学号后两位”及 0n0unun“学号后两位”,并绘出其图形。程序: 4 (1)N=5; n=-(N-1):N-1; x=[zeros(1,N-1)ones(1,N)];stem(n,x); title('单位阶跃序列');grid on; (2)N=6; M=1;%学号01 n=-(N-1):N-1; x=[zeros(1,N-M+1)ones(1,N-M)];stem(n,x); title('单位阶跃序列');grid on; 分析:在上图的基础上平移了1个单位.(3)N=6; M=1;%学号01 n=-(N-1):N-1; x=[zeros(1,N-1)ones(1,N)];y=[zeros(1,N-M+1)ones(1,N-M)];z=x-y;stem(n,z); title('单位阶跃序列');grid on; 2 3、编制程序产生正弦序列xncos2n、xncosn及 学号后两位xnsin2n并绘出其图形。 程序:(1)N=5; A=1; w=2*pi;phi=0;n=0:0.05:N-1;x=A*cos(w*n+phi);stem(n,x);title('余弦信号');grid on; 分析:该序列具有周期性,且输出为余弦信号.(2)N=5; A=1; w=2*pi/1;%学号01 phi=0;n=0:0.05:N-1;x=A*cos(w*n+phi);stem(n,x);title('余弦信号');grid on; ; 分析:该序列具有周期性,且输出为余弦信号.(3)N=5; A=1; w=2*pi;phi=0; n=0:0.05:N-1;x=A*sin(w*n+phi);stem(n,x);title('正弦信号');grid on; 分析:该序列具有周期性,且输出为正弦信号.4、编制程序产生复正弦序列xne(2j学号后两位)n,并绘出其图形。N=3; n=0:0.2:N-1; w=1;%学号01 x=exp((2+j*w)*n);subplot(2,1,1) stem(n,real(x)),title('实部');grid on;subplot(2,1,2) stem(n,imag(x)),title('虚部');grid on; 5、编制程序产生指数序列xnan,并绘出其图形。其中a=学号后两位、a=1/“学号后两位”。 (1)N=10; n=0:N-1; a=1;%学号01 x=a.^n;stem(n,x);title('指数序列');grid on; (2)N=10; n=0:N-1; a=1;%学号01 x=a.^(-n);stem(n,x);title('指数序列');grid on; 实验三 离散时间信号的频域分析 一、实验目的 信号的频域分析是信号处理中一种有效的工具。在离散信号的频域分析中,通常将信号表示成单位采样序列的线性组合,而在频域中,将信号表示成复变量或的线性组合。通过这样的表示,可以将时域的离散序列映射到频域以便于进一步的处理。 在本实验中,将学习利用MATLAB计算离散时间信号的DTFT和DFT,并加深对其相互关系的理解。 二、实验原理 (1)DTFT和DFT的定义及其相互关系。 (2)使用到的MATLAB命令有基于DTFT离散时间信号分析函数以及求解序列的DFT函数。 三、实验内容及分析 (1)编程计算并画出下面DTFT的实部、虚部、幅度和相位谱。 X(e)jw0.05180.1553e11.2828ex(n)cosjwjw0.1553ej2w1.0388ej2w0.0518ej3w0.3418ej3w (2)计算32点序列 5n16,0≦n≦31的32点和64点DFT,分别绘出幅度谱图形,并绘出该序列的DTFT图形。 3-1 clear; x=[0.0518,-0.1553,0.1553,0.0518];y=[1,1.2828,1.0388,0.3418];w=[0:500]*pi/500 H=freqz(x,y,w); magX=abs(H);angX=angle(H);realX=real(H);imagX=imag(H);subplot(221);plot(w/pi,magX);grid; xlabel('frequency in pi unit');ylabel('magnitude');title('幅度 part');axis([0 0.9 0 1.1]); subplot(223);plot(w/pi,angX);grid; xlabel('frequency in pi unit');ylabel('radians');title('相位 part');axis([0 1-3.2 3.2]); subplot(222);plot(w/pi,realX);grid; xlabel('frequency in pi unit');ylabel('real part');title('实部 part');axis([0 1-1 1]); subplot(224);plot(w/pi,imagX);grid; xlabel('frequency in pi unit');ylabel('imaginary');title('虚部 part');axis([0 1-1 1.1]); 3-2 N=32;n=0:N-1; xn=cos(5*pi*n/16);k=0:1:N-1;Xk=fft(xn,N);subplot(2,1,1);stem(n,xn);subplot(2,1,2);stem(k,abs(Xk));title('32点');figure N=64;n=0:N-1; xn=cos(5*pi*n/16);k=0:1:N-1;Xk=fft(xn,N);subplot(2,1,1);stem(n,xn);subplot(2,1,2);stem(k,abs(Xk));title('64点'); (1) (2) 实验四 离散时间LTI系统的Z域分析 一、实验目的 本实验通过使用MATLAB函数对离散时间系统的一些特性进行仿真分析,以加深对离散时间系统的零极点、稳定性,频率响应等概念的理解。学会运用MATLAB分析离散时间系统的系统函数的零极点;学会运用MATLAB分析系统函数的零极点分布与其时域特性的关系;学会运用MATLAB进行离散时间系统的频率特性分析。 二、实验原理 离散时间系统的系统函数定义为系统零状态响应的Z变化与激励的Z变化之比。 在MATLAB中系统函数的零极点可通过函数roots得到,也可借助函数tf2zp得到,tf2zp的语句格式为 [Z,P,K]=tf2zp(B,A)其中,B与A分别表示H(z)的分子与分母多项式的系数向量。它的作用是将H(z)的有理分式表示式转换为零极点增益形式。 若要获得系统函数H(z)的零极点分布图,可直接应用zplane函数,其语句格式为 Zplane(B,A) 其中,B与A分别表示H(z)的分子和分母多项式的系数向量。它的作用是在z平面上画出单位圆、零点与极点。 离散系统中z变化建立了时域函数h(n)与z域函数H(z)之间的对应关系。因此,z变化的函数H(z)从形式可以反映h(n)的部分内在性质。可根据系统的传递函数H(z)求单位冲激响应h(n)的函数impz、filter等。 利用系统的频率响应,可以分析系统对各种频率成分的响应特性,并推出系统的特性(高通、低通、带通、带阻等)。 MATLAB提供了求离散时间系统频响特性的函数freqz,调用freqz的格式主要有两种。一种形式为 [H,w]= reqz(B,A,N)其中,B与A分别表示H(z)分子和分母多项式的系数向量;N为正整数,默认值为512;返回值w包含[0,π]范围内的N个频率等分点;返回值H则是离散时间系统频率响应在0~π范围内N个频率处的值。另一种形式为 [H,w]= freqz(B,A,N,‘whole’) 与第一种方式不同之处在于角频率的范围由[0,π]扩展到[0,2π]。 三、实验内容与结果分析 已知LTI离散时间系统,要求由键盘实现系统参数输入,并绘出幅频和相频响应曲线和零极点分布图,进而分析系统的滤波特性和稳定性。 (一)程序 b=[0.0528,0.797,0.1295,0.1295,0.797,0.0528]; a=[1,-1.8107,2.4947,-1.8801,0.9537,-0.2336];w=[0:20:500]*pi/500; x1=0.0528+0.797*exp(-1*j*w)+0.1295*exp(-2*j*w)+0.1295*exp(-3*j*w)+0.797*exp(-4*j*w)+0.0528*exp(-5*j*w); x2=1-1.8107*exp(-1*j*w)+2.4947*exp(-2*j*w)+1.8801*exp(-3*j*w)+0.9537*exp(-4*j*w)+0.2336*exp(-5*j*w);x22=x2+(x2==0)*eps;x=x1./x22;magx=abs(x); angx=angle(x).*180/pi; subplot(2,2,3);zplane(b,a);title('零极点图');subplot(2,2,2);stem(w/pi,magx);title('幅度部分');ylabel('振幅');subplot(2,2,4);stem(w/pi,angx); xlabel('以pi为单位的频率');title('相位部分');ylabel('相位'); (二)波形图 图4-1 幅频、相频响应曲线、零极点分布图 实验六 IIR数字滤波器的设计 一、实验目的 从理论上讲,任何的线性是不变(LTI)离散时间系统都可以看做一个数字滤波器,因此设计数字滤波器实际就是设计离散时间系统。数字滤波器你包括IIR(无限冲激响应)和FIR(有限冲激响应)型,在设计时通常采用不同的方法。 本实验通过使用MATLAB函数对数字滤波器进行设计和和实现,要求掌握IIR数字巴特沃斯滤波器、数字切比雪夫滤波器的设计原理、设计方法和设计步骤;能根据给定的滤波器指标进行滤波器设计;同时也加深学生对数字滤波器的常用指标和设计过程的理解。 二、实验原理 在IIR滤波器的设计中,常用的方法是:先根据设计要求寻找一个合适的模拟原型滤波器,然后根据一定的准则将此模拟原型滤波器转换为数字滤波器。 IIR滤波器的阶数就等于所选的模拟原型滤波器的阶数,所以其阶数确定主要是在模拟原型滤波器中进行的。 IIR数字滤波器的设计方法如下:(1)冲激响应不变法。(2)双线性变化法。 一般来说,在要求时域冲激响应能模仿模拟滤波器的场合,一般使用冲激响应不变法。冲激响应不变法一个重要特点是频率坐标的变化是线性的,因此如果模拟滤波器的频率响应带限于折叠频率的话,则通过变换后滤波器的频率响应可不失真地反映原响应与频率的关系。 与冲激响应不变法比较,双线性变化的主要优点是靠频率的非线性关系得到s平面与z平面的单值一一对应关系,整个值对应于单位圆一周。所以从模拟传递函数可直接通过代数置换得到数字滤波器的传递函数。 MATLAB提供了一组标准的数字滤波器设计函数,大大简化了滤波器的设计工程。 (1)butter。 (2)cheby1、cheby2。 三、实验内容及分析 利用MATLAB编程方法或利用MATLAB中fdatool工具设计不同功能的IIR数字滤波器。 1、基于chebyshev I型模拟滤波器原型使用冲激不变转换方法设计数字滤波器,要求参数为通带截止频率p0.4;通带最大衰减Ap1dB;阻带截止频率s0.4;阻带最小衰减As35dB。 程序: wp=0.2*pi; %通带边界频率 ws=0.4*pi; %阻带截止频率 rp=1; %通带最大衰减 rs=35; %阻带最小衰减 Fs=1000; %¼ÙÉè³éÑùÂö³å1000hz [N,Wn]=cheb1ord(wp,ws,rp,rs,'s'); [Z,P,K]=cheby1(N,rp,Wn,'s');[H,W]=zp2tf(Z,P,K); figure(1);freqs(H,W);[P,Q]=freqs(H,W);figure(2);plot(Q*Fs/(2*pi),abs(P));grid on; xlabel('频率/Hz');ylabel('幅度'); 2、基于Butterworth型模拟滤波器原型使用双线性变换方法设计数字滤波器的,要求参数为截止频率p0.4;通带最大衰减Ap1dB;阻带截止频率s0.25;阻带最小衰减AS40dB。程序: wp=0.4*pi;ws=0.25*pi;rp=1;rs=40;fs=500;ts=1/fs;wp1=wp*ts;ws1=ws*ts; wp2=2*fs*tan(wp1/2);ws2=2*fs*tan(ws1/2); [N,Wn]=buttord(wp2,ws2,rp,rs,'s');[Z,P,K]=buttap(N);[Bap,Aap]=zp2tf(Z,P,K);[b,a]=lp2lp(Bap,Aap,Wn);[bz,az]=bilinear(b,a,fs);[H,W]=freqz(bz,az);subplot(2,1,1);plot(W/pi,abs(H));grid on;xlabel('频率')ylabel('幅度')subplot(2,1,2); plot(W/pi,20*log10(abs(H)));grid on;xlabel('频率');ylabel('幅度(dB)'); 实验七 FIR数字滤波器的设计 一、实验目的 掌握用窗函数设计FIR数字滤波的原理及其设计步骤;熟悉线性相位数字滤波器的特性。学习编写数字滤波器的设计程序的方法,并能进行正确编程;根据给定的滤波器指标,给出设计步骤。 二、实验原理 如果系统的冲激响应h(n)为已知,则系统的输入输出关系为 y(n)=x(n)*h(n) 对于低通滤波器,只要设计出低通滤波器的冲激响应函数,就可以由式得到系统的输出了。 但是将h(n)作为滤波器的脉冲响应有两个问题:一是它是无限长的;二是它是非因果的。对此,采取两项措施:一是将h(n)截短;二是将其右移。 设计时,要根据阻带的最小衰减和过渡带宽度来选择恰当的窗函数类型和窗口长度N。常用的窗函数有矩形窗、海明窗和布莱克曼窗等。 窗函数设计FIR滤波器步骤如下: (1)给定理想频率响应的幅频特性和相频特性; (2)求理想单位脉冲响应,在实际计算中,可对理想频率响应采样。(3)根据过渡带宽度和阻带最小衰减,确定窗函数类型和窗口长度N;(4)求FIR滤波器单位脉冲响应; (5)分析幅频特性,若不满足要求,可适当改变窗函数形式或长度N,重复上述设计过程,以得到满意的结果。 三、实验内容及分析 1、分别用海明窗和布莱克曼窗设计一个48阶的FIR带通滤波器,通带为Wn0.450.55。程序1:海明窗设计 N=48; Window=hamming(N+1);w1=0.45;w2=0.55;ws=[w1,w2]; b=fir1(N,ws/pi,Window);freqz(b,1,512);title('海明窗');grid on; 程序2:莱克曼窗设计 N=48; Window=blackman(N+1);w1=0.45;w2=0.55;ws=[w1,w2]; b=fir1(N,ws/pi,Window);freqz(b,1,512);title('布莱克曼窗');grid on; 2、用矩形窗设计一个线性相位高通滤波器。其中Hejwej00.3 00.3程序: N=9; alpha=(N-1)/2;Wc=0.7*pi;n=(0:8);i=n-alpha;i=i+(i==0)*eps; h=(-1).^n.*sin((i).*Wc)./((i).*pi);%矩形窗函数设计的系统脉冲响应 w=(0:1:500)*2*pi/500; H=h*exp(-j*n'*w);%矩形窗函数设计的频响 magH=abs(H);% 矩形窗函数设计的振幅 subplot(211);stem(n,h); axis([0,8,-0.4,0.4]);title('矩形窗设计h(n)');line([0,10],[0,0]);xlabel('n');ylabel('h');subplot(212);plot(w/pi,magH); xlabel('以pi为单位的频率');ylabel('H振幅');axis([0,2,0,1.7]);title('矩形窗设计振幅谱'); 实验心得体会: 这次实验使我进一步加深了对MATLAB软件的使用。从上次的信号系统实验的初步使用到这一次的深入了解,有了更深刻的认识。对这种语言环境也有了新的了解。 在实验的过程中,我对数字滤波器的整个过程有了很好的理解和掌握。IIR数字滤波器的设计让我知道了巴特沃思滤波器和切比雪夫滤波器的频率特性,还有双线性变换及脉冲响应不变法设计的滤波器的频率特性。做这两个实验的时候程序有点困难,但经过细心的改写图形最终出来了。FIR数字滤波器的设计出来的是两种窗的图形,通过两种窗的比较,我了解了他们各自的特点,幅频和相频特性。 最后,感谢张老师对我的谆谆教导! 数字信号处理 实验报告 实验一 序列的傅立叶变换 一、实验目的 1.进一步加深理解DFS,DFT算法的原理;2.研究补零问题;3.快速傅立叶变换(FFT)的应用。 二、实验步骤 1.复习DFS和DFT的定义,性质和应用; 2熟悉MATLAB语言的命令窗口、编程窗口和图形窗口的使用;3利用提供的程序例子编写实验用程序;4.按实验内容上机实验,并进行实验结果分析;5.写出完整的实验报告,并将程序附在后面。 三、实验内容 1.周期方波序列的频谱 试画出下面四种情况下的的幅度频谱, 并分析补零后,对信号频谱的影响。 x(n)cos(0.48n)cos(0.52n)2.有限长序列x(n)的DFT(1)取x(n)(n=0:10)时,画出x(n)的频谱X(k)的幅度;(2)将(1)中的x(n)以补零的方式,使x(n)加长到(n:0~100)时,画出x(n)的频谱X(k)的幅度; (3)取x(n)(n:0~100)时,画出x(n)的频谱X(k)的幅度。利用FFT进行谱分析x(t)2sin(4t)5cos(8t)3.已知:模拟信号 以t=0.01n(n=0:N-1)进行采样,求N点DFT的幅值谱。请分别画出N=45;N=50;N=55;N=60时的幅值曲线。 四、实验数据分析 1.周期方波序列的频谱分析 首先定义一个功能函数dfs function[Xk]=dfs(xn,N)n=[0:1:N-1];k=[0:1:N-1];WN=exp(-j*2*pi/N);nk=n'*k;WNnk=WN.^nk;Xk=xn*WNnk;(1)L=5,N=20;%题1.(1)L=5;N=20;%对于(2),(3),(4)问,只要修改L,N的数值就好。n=1:N;xn=[ones(1,L),zeros(1,N-L)];Xk=dfs(xn,N);magXk=abs([Xk(N/2+1:N)Xk(1:N/2+1)]);k=[-N/2:N/2];figure(1)subplot(2,1,1);stem(n,xn);xlabel('n');ylabel('xtide(n)');title('DFS of SQ.wave:L=5,N=20');subplot(2,1,2);stem(k,magXk);axis([-N/2,N/2,0,16]);xlabel('k');ylabel('Xtide(k)'); (2)L=5,N=40; (3).L=5,N=60 (4)L=7,N=60; 结果分析:虽然周期序列不存在FT,但是一个周期序列可以利用其DFS系数X(k)表示它的频谱分布规律,从以上各频谱图可以看出,随着补零点数的增加,周期序列的谐波次数越来越多,其频谱的包络线越来越平滑连续,更能反映幅度值随时间的变化。 2.有限长序列的DFT(1) %题2-(1)n=0:10;xn=cos(0.48*pi*n)+cos(0.52*pi*n);N=11;Xk=fft(xn,N); %序列x(n)的N点DFT k=0:N-1;wk=2*k/N; subplot(1,1,1);stem(wk,abs(Xk),'.');title('频谱X(K)的幅度');xlabel('ω/π');ylabel('幅度'); (2)%题2-2 M=10;N=100;n=1:M;xn=cos(0.48*pi*n)+cos(0.52*pi*n);n1=[0:1:N-1];y1=[xn(1:1:M),zeros(1,N-M)];figure(1)subplot(2,1,1);stem(n1,y1);xlabel('n');ylabel('x(n)');title('序列x(n),0<=n<=100');axis([0,N,-2.5,2.5]);Y1=fft(y1);magY1=abs(Y1(1:1:N/2+1));k1=0:1:N/2;w1=2*pi/N*k1;subplot(2,1,2);title('x(n)的幅频特性曲线');stem(w1/pi,magY1);axis([0,1,0,60]);xlabel('omega/pi');ylabel('|X(K)|'); (3) %Example2-3 M=10;N=100;n=0:M;xn=cos(0.48*pi*n)+cos(0.52*pi*n);n1=[0:1:N-1];y1=[xn(1:1:M),zeros(1,N-M)];figure(1)subplot(2,2,1);stem(n1,y1);xlabel('n');ylabel('x(n)');title('序列x(n),0<=n<=100');axis([0,N,-2.5,2.5]);YK=fft(y1);Y=abs(Y1(1:1:N/2+1));k1=0:1:N/2;w1=2*pi/N*k1;subplot(2,2,3);stem(w1/pi,Y);title('x(n)的幅频特性曲线');axis([0,1,0,60]);xlabel('omega/pi');ylabel('|X(K)|');subplot(2,2,4);plot(angle(Y1));title('x(n)的相频特性曲线');xlabel('omega/pi');ylabel('phi(omega)'); 结果分析:由上述仿真图可得,随着n取值范围的增大,其频谱在[0,2π]上的采样间隔越来越小,采样点越来越多。采样点越多,其DFS频谱越接近FT的频谱。其相频特性曲线呈现周期性变化。 3.问题三 %题3 %N=45 figure(1)subplot(2,2,1)N=45;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=45')%N=50 subplot(2,2,2)N=50;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=50')%N=55 subplot(2,2,3)N=55;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=55')%N=60 subplot(2,2,4)N=60;n=0:N-1;t=0.01*n;q=n*2*pi/N;x=2*sin(4*pi*t)+5*cos(8*pi*t);y=fft(x,N);plot(q,abs(y))stem(q,abs(y))title('DFT N=60') 结果分析:由上述仿真图可得,随着N取值的增大,其频谱在[0,2π]上的采样间隔越来越小,采样点越来越多。采样点越多,其DFS频谱越接近FT的频谱,幅值曲线越来越清晰,更能准确反应幅值随时间的变化规律。 五、心得体会 对于周期序列的离散傅里叶变换,通过matlab的模拟,可以更好的了解掌握序列采样间隔对其采样频谱的影响。在实验过程中,学习如何使用matlab程序语言解决问题,是很有价值的。让我对matlab的使用更为熟练。 实验二 用双线性变换法设计IIR数字滤波器 一、实验目的 1.熟悉用双线性变换法设计IIR数字滤波器的原理与方法; 2. 掌握数字滤波器的计算机仿真方法; 3、通过观察对实际心电图的滤波作用,获得数字滤波器的感性知识。 二、实验内容 1.用双线性变换法设计一个巴特沃斯低通IIR滤波器,设计指标参数为:在通带内频率低于0.2π时,最大衰减小于1dB;在阻带内[0.3π,π]频率区间上,最小衰减大于15dB 2.以0.2π为采样间隔,打印出数字滤波器在频率区 间[0, 0.2π]上的幅值响应曲线。 3.用所设计的滤波器对实际的心电图信号采样序列 x(n)=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6, 6,6, 4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];)进行仿真滤波处理,并分别打印出滤波前后的心电图信号波形图,观察总结滤波作用与效果。三.实验步骤 1.复习有关巴特沃斯模拟滤波器设计和双线性变换法设计IIR滤波器的内容 2.参考例子程序用MATLAB语言编写仿真实验用程序; 3.在通用计算机上运行仿真程序 4.写出完整的实验报告并回答思考题。四.实验数据分析 1.巴特沃斯低通IIR滤波器的设计 由题可得,数字低通技术指标为 wp=0.2πrad, αp=1dB,ws=0.3πrad, αs=15dB 若T=1s,预畸变校正计算相应的模拟低通的技术指标为 Ωp=1dB,Ωs=15dB 通过计算可得阶数N=5.3056,向上取整的N=6.Ωc=0.7663rad/s,这样保证阻带技术指标满足要求,通带指标有富余。 以下是通过matlab实现的巴特沃斯低通IIR滤波器的设计 %1 T=1;Fs=1/T;wpz=0.2;wsz=0.3;wp=2*tan(wpz*pi/2);ws=2*tan(wsz*pi/2);rp=1;rs=15;%预畸变校正转换指标 [N,wc]=buttord(wp,ws,rp,rs,'s');%设计过渡模拟滤波器 [B,A]=butter(N,wc,'s');[Bz,Az]=bilinear(B,A,Fs);%用双线性变换法转换成数字滤波器 fk=0:1/512:1;wk=2*pi*fk;Hk=freqs(B,A,wk);figure(1);subplot(2,1,1);plot(fk,20*log10(abs(Hk)));grid on;title('模拟滤波器幅值响应曲线');xlabel('omega/pi');ylabel('幅度(dB)');axis([0,1,-100,5]);[Nd,wdc]=buttord(wpz,wsz,rp,rs);%调用buttord和butter直接设计数字滤波器 [Bz,Az]=butter(N,wdc);wk=0:pi/512:pi;Hz=freqz(Bz,Az,wk);subplot(2,1,2);plot(wk/pi,20*log10(abs(Hz)));grid on;title('数字滤波器幅值响应曲线');xlabel('omega/pi');ylabel('幅度(dB)');axis([0,1,-100,5]); 2.以0.2π为采样间隔,打印出数字滤波器在频率区 间[0, 0.2π]上的幅值响应曲线。 在第一问的基础上加上下面的程序 %2 figure(2);freqz(Bz,Az,[0:0.02*pi:0.2*pi]) 3.在第一问程序的基础上加上如下程序即可 figure(3);x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];subplot(2,2,1);n=0:55;stem(n,x,'.');title('x(n)的脉冲响应');xlabel('n');ylabel('x(n)');A=0.09036;b1=[A,2*A,A];a1=[1,-1.2686,0.7051];h1=filter(b1,a1,x);[H1,w]=freqz(b1,a1,100);b2=[A,2*A,A];a2=[1,-1.0106,0.3583];h2=filter(b2,a2,h1);[H2,w]=freqz(b2,a2,100);b3=[A,2*A,A];a3=[1,-0.9044,0.2155];h3=filter(b3,a3,h2);[H3,w]=freqz(b3,a3,100);subplot(2,2,2);stem(n,h3,'.');xlabel('n');ylabel('y(n)');title('通过滤波器H1(z),H2(z),H3(z)后的y3(n)函数');subplot(2,2,3);H4=H1.*(H2);H=H4.*(H3);mag=abs(H);db=20*log10((mag+eps)/max(mag));plot(w/pi,db);xlabel('ω/π');ylabel('20log[Ha3(ejw)]');title('通过滤波器H1(z),H2(z),H3(z)后的对数频率响应20log[Ha3(ejw)]函数');grid;figure(4);N=1024;n=0:N/2-1;Xk=fft(x,N);AXk=abs(Xk(1:N/2));f=(0:N/2-1)*Fs/N;f=f/Fs;subplot(211);plot(f,AXk);title('x(n)的频谱');xlabel('f');ylabel('| X(k)|');axis([0,0.5,0,400]); 五、思考题 用双线性变换设计数字滤波器的过程中,下面变换公式的T值的取值,对设计结果是否有影响?为什么? 21z1s T1z1 虽然采用双线性变换法设计数字滤波器不会产生混叠现象,T得取值可以任选。双线性变换法具有非线性,T小一些,非线性的影响也就少一些。 六、心得体会 通过这次实验,我学会了如何使用matlab语言来实现IIR低通滤波器的设计,同时对于其数字低通技术指标的就算有了进一步的认识体会。 在做实验的过程中,面对种种困难,但是却在解决困难的过程中收获了很多。是一个很有意义的经历。 实验三 用窗函数法设计FIR数字滤波器 一.实验目的 1.掌握用窗函数法设计FIR数字滤波器的原理和方法。2.熟悉线性相位FIR数字滤波器特性。3.了解各种窗函数对滤波特性的影响。二.实验原理 如果所希望的滤波器的理想频率响应函数为 Hd(e jω),则其对应的单位脉冲响应为 用窗函数w(n)将hd(n)截断,并进行加权处理,得到: h(n)就作为实际设计的FIR数字滤波器的单位脉冲响应序列,其频率响应函数H()为 如果要求线性相位特性,则h(n)还必须满足: 根据上式中的正、负号和长度N的奇偶性又将线性相位FIR滤波器分成四类。要根据所设计的滤波特性正确选择其中一类。例如,要设计线性相位低通特性,可选择h(n)=h(N-1-n)一类,而不能选h(n)=-h(N-1-n)一类。三.实验内容 1.用MATLAB产生各种窗函数 %窗函数 subplot(4,2,1)m=200;a=boxcar(m);%矩形窗 m=1:200;plot(m,a)title('矩形窗');subplot(4,2,2)m=200;b=bartlett(m)%三角窗 m=1:200;plot(m,b)title(' 三角窗');subplot(4,2,3)m=200;c=hanning(m);%汉宁窗 m=1:200;plot(m,c)title('汉宁窗');subplot(4,2,4)m=200;d=hamming(m);%哈明窗 m=1:200;title(' 哈明窗');plot(m,d)subplot(4,2,5)m=200;e=blackman(m);%布莱克曼窗 m=1:200;plot(m,e)title('布莱克曼窗');subplot(4,2,6)m=200;f=kaiser(m,7.865);% 凯塞窗 m=1:200;plot(m,f)title(' 凯塞窗');subplot(4,2,7)plot(m,a,'r*',m,b,'g+',m,c,'y*',m,d,'b.',m,e,'y.',m,f,'k.')title(' 各种窗函数'); 2.利用窗函数设计FIR滤波器 设计具有下列指标p=0.25,Rp=0.25dB,s=0.3,Rp=50dB的低通数字滤波器。 由已知条件Rp=50dB,通过计算可知应该选择哈明窗,哈明窗的Rps=53dB,选择窗函数时应该选择Rps>=50的,因而选择离50dB的哈明窗。 通过Bt=6.6π/N,Bt=ws-wp得到阶数N。 以下是利用matlab来实现FIR滤波器的设计。 %窗函数法设计FIR数字滤波器 wp=pi*0.25;ws=pi*0.3;%Rs=50dB,所以选择哈明窗 DB=ws-wp;%计算过渡带宽度 N=ceil(6.6*pi/DB);%计算哈明窗所需长度N wc=(wp+ws)/2/pi;%计算理想低通滤波器通带截止频率(关于π归一化)hn=fir1(N,wc);%调用fir1计算低通数字滤波器 %以下是绘图部分 figure(2)M=1024;hk=fft(hn,M);n=0:N;subplot(1,2,1);stem(n,hn,'.');%绘制序列h(n)xlabel('n');ylabel('h(n)');title('数字滤波器h(n)');k=1:M/2;w=2*(0:M/2-1)/M;subplot(1,2,2);plot(w,20*log10(abs(hk(k))));%绘制损耗函数曲线 axis([0,1,-80,5]);xlabel('ω/π');ylabel('20lg|Hg(ω)|');title('h(n)的损耗函数曲线');grid on; 结果分析:根据损耗函数可得,利用哈明窗实现了题目所要求的FIR低通滤波器。满足p=0.25,Rp=0.25dB,s=0.3,Rp=50dB。因而选择哈明窗可以很好的满足所需要求。 四、心得体会 本次实验是使用窗函数来设计FIR滤波器,在MATLAB中只需要将给定的参数输入到函数中即可马上得到处结果。通过本次实验,我对使用MATLAB快速设计滤波器的流程更为熟练,同时,也由衷地感叹MATLAB的功能强大之处,它让我们在设计时能节约大量的时间。 根据输入的参数以及结果,使我对课本上的设计滤波器的知识更为了解。更加巩固了理论知识。 怀化学院数学系实验报告 实验项目名称:IIR数字滤波器的设计(1) 指 导老 师: 欧卫华 学 姓 实验项目制定人:实验项目审批人: 年月日 一、实验目的掌脉冲相应不变法设计IIR-Butterworth数字滤波器的具体设计方法及原理。 二、实验原理与方法 1.确定数字滤波器的性能指标:通带临界频率fp、阻带临界频率fs;通带内的最大衰减Ap;阻带内的最小衰减As;采样周期T; 2.确定相应的数字角频率,ωp=2πfp;ωr=2πfr; 3.根据Ωp和Ωs计算模拟低通原型滤波器的阶数N,并求得低通原型的传递函 数Ha(s); 4.用上面的脉冲响应不变法公式代入Ha(s),求出所设计的传递函数H(z); 5.分析滤波器特性,检查其是否满足指标要求。 三、实验内容及步骤 冲激响应不变法设计数字Butterworth低通滤波器 (1)、模拟滤波器的最小阶数[N,wn]=buttord(wp,ws,rp,rs,'s'); (2)、设计模拟低通滤波器原型,[z,p,k]=buttap(N); (3)、将零极点形式转换为传递函数形式,[Bap,Aap]=zp2tf(z,p,k); (4)、进行频率变换,[b,a]=lp2lp(Bap,Aap,wn); (5)用脉冲相应不变法得到数字滤波器的系统函数[bz,az]=impinvar(b,a,fs); 四、实验范例 用脉冲相应不变法设计一个Butterworth低通数字滤波器,使其特征逼近一个低通Butterworth模拟滤波器的下列性能指标,通带截止频率Wp=2*pi*2000rad/s,通带波纹Rp小于3dB,阻带边界频率为Ws=2*pi*3000rad/s阻带衰减大于15dB,采样频率Fs=10000;z,假设一个信号x(t)=sin(2*pi*f1*t) +0.5*cos(2*pi*f2*t),其中f1=1000Hz,f2=4000Hz,试将原信号与通过该滤波器的输出信号进行比较。 wp=2000*2*pi;%滤波器截止频率 ws=3000*2*pi; rp=3;rs=15;%通带波纹和阻带衰减 fs=10000;%采样频率 Nn=128; [N,wn]=buttord(wp,ws,rp,rs,'s');%模拟滤波器的最小阶数 [z,p,k]=buttap(N);%设计模拟低通滤波器原型 [Bap,Aap]=zp2tf(z,p,k);%将零极点形式转换为传递函数形式 [b,a]=lp2lp(Bap,Aap,wn);%进行频率变换 [bz,az]=impinvar(b,a,fs);%应用脉冲相应不变法得到数字滤波器的系统函数 figure(1); [h,f]=freqz(bz,az,Nn,fs);%画出数字滤波器的幅频特性和相频特性 subplot(2,1,1),plot(f,20*log10(abs(h))); xlabel('频率/Hz');ylabel('振幅/dB');grid on; subplot(2,1,2),plot(f,180/pi*unwrap(angle(h))); xlabel('频率/Hz');ylabel('振幅/^o');grid on; figure(2); f1=1000;f2=4000;%输入信号的频率 N=100;%数据长度 dt=1/fs;n=0:N-1;t=n*dt;%采样间隔和时间序列 x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t);%滤波器输入信号 subplot(2,1,1),plot(t,x),title('输入信号')%画出输入信号 %y=filtfilt(bz,az,x); y1=filter(bz,az,x);%用上面设计的滤波器对输入信号滤波 subplot(2,1,2),plot(t,y1,'r-'),title('输出信号'),xlabel('时间/s');legend('filter') 五、实验习题 用脉冲相应不变法设计一个Butterworth低通数字滤波器,通带频率为0= 六,实验结果第四篇:数字信号处理实验报告
第五篇:六 数字信号处理实验报告--IIR数字滤波器设计