数字信号处理实验-2014

《数字信号处理》实验指导书

目录

实验一 FFT频谱分析及应用 (3)

实验二数字滤波器的设计 (5)

附录： MATLAB基本操作及常用命令 .................. 错误！未定义书签。

实验一 FFT频谱分析及应用

实验学时：2学时

实验类型：验证

实验要求：必修

一、实验目的：

1、学习MATLAB语言的编程和调试技巧；、

2、通过实验加深对FFT的理解；

3、熟悉应用FFT对典型信号进行频谱分析的方法。

二、实验内容

使用MATLAB程序实现信号频域特性的分析。涉及到离散傅立叶变换(DFT)、快速傅立叶变换(FFT)及信号频率分辨率等知识点。

三、实验原理与方法和手段

在各种信号序列中，有限长序列占重要地位。对有限长序列可以利用离散傅立叶变换(DFT)进行分析。DFT不但可以很好的反映序列的频谱特性，而且易于用快速算法(FFT)在计算机上进行分析。

有限长序列的DFT是其z变换在单位圆上的等距离采样，或者说是序列傅立叶的等距离采样，因此可以用于序列的谱分析。FFT是DFT的一种快速算法，它是对变换式进行一次次分解，使其成为若干小数据点的组合，从而减少运算量。

在MATLAB信号处理工具箱中的函数fft（x,n）,可以用来实现序列的N点快速傅立叶变换。

经函数fft求得的序列一般是复序列，通常要求出其幅值和相位。MATLAB中提供了求复数的幅值和相位的函数：abs、angle，这些函数一般和fft同时使用。

四、实验组织运行要求

1、学生在进行实验前必须进行充分的预习，熟悉实验内容；

2、学生根据实验要求，读懂并理解相应的程序；

3、学生严格遵守实验室的各项规章制度，注意人身和设备安全，配合和服从实验室人员管理；

4、教师在学生实验过程中予以必要的辅导，独立完成实验；

5、采用集中授课形式。

五、实验条件

1、具有WINDOWS 98/2000/NT/XP 操作系统的计算机一台； 2.、MATLAB 编程软件。 六、实验步骤

在“开始--程序”菜单中，找到MATLAB 程序，运行启动；

进入MATLAB 后 ，在Command Window 中输入实验程序，并执行； 记录运行结果图形，作分析。 具体步骤如下：

1、模拟信号)8cos(5)4sin(2)(t t t x ππ+=，以)1:0(01.0-==N n n t 进行采样，求： （1）N ＝40点FFT 的幅度频谱，从图中能否观察出信号的2个频谱分量? （2）提高采样点数，如N ＝128，再求该信号的幅度频谱，此时幅度频谱发生了什么变化？信号的2个模拟频率和数字频率各为多少？FFT 频谱分析结果与理论上是否一致？

程序：

N=40;n=0:N-1; t=0.01*n;

x=2*sin(4*pi*t)+5*cos(8*pi*t); k=0:N/2;w=2*pi/N*k; X=fft(x,N);

magX=abs(X(1:N/2+1));

subplot(2,1,1);stem(n,x,'.');title('signal x(n)'); subplot(2,1,2);plot(w/pi,magX);title('FFT N=40'); xlabel('f (unit :pi)');ylabel('|X|');grid

程序2：

N=128;n=0:N-1; t=0.01*n;

x=2*sin(4*pi*t)+5*cos(8*pi*t); k=0:N/2;w=2*pi/N*k; X=fft(x,N);

magX=abs(X(1:N/2+1));

subplot(2,1,1);stem(n,x,'.');title('signal x(n)'); subplot(2,1,2);plot(w/pi,magX);title('FFT N=128'); xlabel('f (unit :pi)');ylabel('|X|');grid

2、一个连续信号含三个频谱分量，经采样得以下序列：

))215.0(2cos())15.0(2cos()15.02sin()(n df n df n n x +?++?+?=πππ

（1）N ＝64，df 分别为161、1/64，观察其频谱；

（2）N ＝64、128，df 为1/64，做128点得FFT ，其结果有何不同？

3、被噪声污染得信号，比较难看出所包含得频率分量，如一个由50Hz 和120Hz 正弦信号构成的信号，受零均值随机噪声的干扰，数据采样率为1000Hz ，试用FFT 函数来分析其信号频率成分，要求：（1）画出时域波形；（2）分析信号功率谱密度。

注：在MATLAB 中，可用函数rand （1，N ）产生均值为0，方差为1，长度为N 的高斯随机序列。 程序： clear

t=0:0.001:1;

x=sin(2*pi*50*t)+sin(2*pi*120*t)+rand(size(t)); subplot(2,1,1) plot(x(1:50));

title('óé50Hzoí120Hz?y?òD?o?113é?òêü??éù??è?μ?D?o?í?:');grid on ; y=fft(x,512); subplot(2,1,2) plot(abs(y).^2);

title('êü??è?D?o?1|?ê?×?ü?èí?:');grid on ;

实验二 数字滤波器的设计

实验学时：2学时 实验类型：设计综合 实验要求：必修

知识点：IIR 、FIR 数字滤波器，脉冲响应不变法，双线性变换法，窗函数法 六、实验步骤

在“开始--程序”菜单中，找到MATLAB 程序，运行启动；

进入MATLAB 后 ，在Command Window 中输入自己编写的主程序，并执行；

记录运行结果图形，作分析对比。 具体步骤如下：

1、查看帮助文件，了解相关函数的调用格式。

2、KHz fp 1.0=，dB Rp 1=，KHz fs 3.0=，dB R 25s =，ms T 1=；分别用脉冲响应不变法和双线性变换法设计一个Butterworth 数字低通滤波器

wp=2*pi*100; ws=2*pi*300; rp=1; rs=25; fs=1000;

[n,wn]=buttord(wp,ws,rp,rs,'s');%??μ?μíí¨butterworth??2¨?÷?×êyoíí¨′????1?μ?ê

[z,p,k]=buttap(n);%éè???￡?aμíí¨butterworth??2¨?÷ [b0,a0]=zp2tf(z,p,k);%×a?ˉ?a′?í3oˉêy±íê???2¨?÷

[b,a]=lp2lp(b0,a0,wn);%????2¨?÷?-Dí±????aòa?óéè??μ?μíí¨??2¨?÷ [bz,az]=impinvar(b,a,fs);%?￡?a??2¨?÷×a?ˉ?aêy×???2¨?÷ [h,w]=freqz(bz,az);%?????μ?ê?ìó| subplot(2,1,1)

plot(w/(2*pi),20*log10(abs(h))); grid on

title('??3??ìó|2?±?·¨') xlabel('?μ?ê￡¨Hz￡?'); ylabel('·ù?μ?ìó|￡¨dB￡?'); Ap=-20*log10(abs(h(1))); As=-20*log10(abs(h(2))); disp(Ap) disp(As)

wp=2*pi*100; ws=2*pi*300; rp=1; rs=25; fs=1000;

[n,wc]=buttord(wp,ws,rp,rs,'s');%??μ?μíí¨butterworth??2¨?÷?×êyoíí¨′????1?μ?ê

[z,p,k]=buttap(n);%éè???￡?aμíí¨butterworth??2¨?÷ [b0,a0]=zp2tf(z,p,k);%×a?ˉ?a′?í3oˉêy±íê???2¨?÷

[b,a]=lp2lp(b0,a0,wc);%????2¨?÷?-Dí±????aòa?óéè??μ?μíí¨??2¨?÷ [bz,az]=bilinear(b,a,fs)%?￡?a??2¨?÷×a?ˉ?aêy×???2¨?÷

[h,w]=freqz(bz,az)%?????μ?ê?ìó| subplot(2,1,2)

plot(w/(2*pi),20*log10(abs(h))); grid on

title('????D?±???·¨') xlabel('?μ?ê￡¨Hz￡?'); ylabel('·ù?μ?ìó|￡¨dB￡?');

3、用窗函数法设计一线性相位FIR 低通滤波器，设计指标为：

π3.0=Wp ，π5.0=Ws ，dB Rp 25.0=，dB Rs 50=

程序：WP=0.3*pi; WS=0.5*pi; RP=0.5; RS=45; N=45; M=N-1;

w=boxcar(N)'; wc=(WP+WS)/2; k=0:M;

hd=(wc/pi)*sinc(wc*(k-0.5*M)/pi); h=hd.*w;

omega=linspace(0,pi,512);

mag=freqz(h,[1],omega);subplot(3,1,1) plot(omega/pi,20*log10(abs(mag)));

w=hanning(N)'; wc=(WP+WS)/2; k=0:M;

hd=(wc/pi)*sinc(wc*(k-0.5*M)/pi); h=hd.*w;

omega=linspace(0,pi,512);

mag=freqz(h,[1],omega);subplot(3,1,2) plot(omega/pi,20*log10(abs(mag))); w=hamming(N)'; wc=(WP+WS)/2; k=0:M;

hd=(wc/pi)*sinc(wc*(k-0.5*M)/pi); h=hd.*w;

omega=linspace(0,pi,512);

mag=freqz(h,[1],omega);subplot(3,1,3) plot(omega/pi,20*log10(abs(mag)));