数字信号处理西电上机实验.docx

1、数字信号处理西电上机实验数字信号处理实验报告实验一： 信号、 系统及系统响应一、实验目的：(1) 熟悉连续信号经理想采样前后的频谱变化关系， 加深对时域采样定理的理解。(2) 熟悉时域离散系统的时域特性。(3) 利用卷积方法观察分析系统的时域特性。(4) 掌握序列傅里叶变换的计算机实现方法， 利用序列的傅里叶变换对连续信号、 离散信号及系统响应进行频域分析。 二、实验原理与方法：(1) 时域采样。(2) LTI系统的输入输出关系。三、 实验内容、步骤 (1) 认真复习采样理论、 离散信号与系统、 线性卷积、 序列的傅里叶变换及性质等有关内容， 阅读本实验原理与方法。 (2) 编制实验用主程序及

2、相应子程序。 信号产生子程序， 用于产生实验中要用到的下列信号序列： a. xa(t)=A*e-at *sin(0t)u(t) A=444.128;a=50*sqrt(2)*pi; b. 单位脉冲序列：xb(n)=(n) c. 矩形序列： xc(n)=RN(n), N=10 系统单位脉冲响应序列产生子程序。 本实验要用到两种FIR系统。 a. ha(n)=R10(n); b. hb(n)=(n)+2.5(n-1)+2.5(n-2)+(n-3) 有限长序列线性卷积子程序 用于完成两个给定长度的序列的卷积。 可以直接调用MATLAB语言中的卷积函数conv。 conv用于两个有限长度序列的卷积，

3、它假定两个序列都从n=0 开始。 调用格式如下：y=conv (x, h)四、实验内容调通并运行实验程序， 完成下述实验内容： 分析采样序列的特性。 a. 取采样频率fs=1 kHz, 即T=1 ms。 b. 改变采样频率， fs=300 Hz， 观察|X(ej)|的变化， 并做记录(打印曲线)； 进一步降低采样频率， fs=200 Hz， 观察频谱混叠是否明显存在， 说明原因， 并记录(打印)这时的|X(ej)|曲线。 时域离散信号、 系统和系统响应分析。a. 观察信号xb(n)和系统hb(n)的时域和频域特性； 利用线性卷积求信号xb(n)通过系统hb(n)的响应y(n)， 比较所求响应y

4、(n)和hb(n)的时域及频域特性， 注意它们之间有无差别， 绘图说明， 并用所学理论解释所得结果。b. 观察系统ha(n)对信号xc(n)的响应特性。 卷积定理的验证。五、思考题 (1) 在分析理想采样序列特性的实验中， 采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同? 它们所对应的模拟频率是否相同? 为什么? (2) 在卷积定理验证的实验中， 如果选用不同的频域采样点数M值， 例如， 选M=10和M=20， 分别做序列的傅里叶变换， 求得的结果有无差异？为什么？六、参考资料数字信号处理教程Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年七、实验

5、程序及结果本实验使用自定义函数的方法产生信号：门函数：functiony,n=gate(np,ns,nf)if nsnp|nsnf|npnf error(ns=np=0&nnp|nsnf|npnf error(ns=np=0);hb=impseq(0,0,10)+2.5.*impseq(1,0,10)+2.5.*impseq(2,0,10)+impseq(3,0,10);ha,N=gate(10,0,10);第一部分主程序：y1,t1=sig(0,0.001,1);%t1=0:0.001:1;%y1=0.5.*sin(2.*pi.*15.*t1);figure;plot(t1,y1);axis

6、(0 0.15 -10 150);title(f=1000);figure;f=-500:500;yy1=fftshift(fft(y1);plot(f,abs(yy1);title(f=1000);y2,t2=sig(0,1/300,1);figure;plot(t2,y2);axis(0 0.15 -10 150);title(f=300);figure;yy2=fftshift(fft(y2);plot(-150:150,abs(yy2);title(f=300);y3,t3=sig(0,0.02,1);figure;plot(t3,y3);title(f=200);axis(0 0.1

7、5 -8 5);figure;yy3=fftshift(fft(y3);plot(-25:25,abs(yy3);title(f=200);第二部分主程序：xb=impseq(0,0,100);figure;stem(0:length(xb)-1,xb);axis(-1,3,0,2);figure;y1=fft(xb);plot(0:length(y1)-1,y1);hb=impseq(0,0,100)+2.5.*impseq(1,0,100)+2.5.*impseq(2,0,100)+impseq(3,0,100);figure;stem(0:length(hb)-1,hb);axis(-1

8、,5,0,3)figurey2=fftshift(fft(hb);plot(-(length(y2)-1)/2:(length(y2)-1)/2,abs(y2);figure;y3=conv(double(xb),double(hb);plot(0:length(y3)-1,y3);y4=fftshift(fft(y3);figure;plot(-(length(y4)-1)/2:(length(y4)-1)/2,abs(y4);xc=gate(10,0,100);ha=xc;y5=fftshift(fft(conv(double(xc),double(ha);figure;plot(-(le

9、ngth(y5)-1)/2:(length(y5)-1)/2,abs(y5);实验结果：第一部分： 采样频率为200Hz时时域恢复 200Hz时频谱 采样频率为300Hz时时域恢复 300Hz时频谱 采样频率为1000Hz是时域恢复 1000Hz时频谱第二部分： Xb频域分析 xb时域分析 hb频域分析 hb时域分析 信号通过xb频域 信号通过xc频域卷积定理的验证：%卷积定理的验证y1=0,1,2,3,4,5,6,7,8;y2=8,7,6,5,4,3,2,1,0;y3=conv(y1,y2);f1=fft(y1,20);f2=fft(y2,20);f3=fft(y3,20);f=f1.*f2

10、;subplot(2,1,1);stem(0:length(f)-1,f,.);subplot(2,1,2);stem(0:length(f3)-1,f3,.);时域卷积后求频谱和频域相乘可见，f= conv(y1,y2)的频谱和y1,y2的频谱相乘后结果相同。即满足卷积定理八、回答问题(1) 在分析理想采样序列特性的实验中，采样频率不同时，相应理想采样序列的傅里叶变换频谱的数字频率度量是否都相同？它们所对应的模拟频率是否相同？为什么？ 答：数字频率度量不相同，但他们所对应的模拟频率相同。由w=*Ts公式得，采样间隔变化时模拟频率对应的数字频率会有相应的变化，故其度量会有所变化。而且采样频率的

11、大小直接关系到能否将能否将原始信号恢复出来。(2) 在卷积定理验证的实验中，如果选用不同的频域采样点数M值，例如，选M=10和M=20，分别做序列的傅里叶变换，求得的结果有无差异？ 答：有差异，所到的结果点数不同。九、总结本实验主要是后续实验的基础。涉及内容也比较浅显。单位冲击序列与hb卷积后得到的频谱与hb原频谱相同。原因很简单，是因为单位冲击序列卷积任何函数仍然是原函数。实验二： 用FFT作谱分析一、实验目的(1) 进一步加深DFT算法原理和基本性质的理解(因为FFT只是DFT的一种快速算法， 所以FFT的运算结果必然满足DFT的基本性质)。(2) 熟悉FFT算法原理和FFT子程序的应用。

12、(3) 学习用FFT对连续信号和时域离散信号进行谱分析的方法， 了解可能出现的分析误差及其原因， 以便在实际中正确应用FFT。二、实验步骤(1) 复习DFT的定义、 性质和用DFT作谱分析的有关内容。(2) 复习FFT算法原理与编程思想， 并对照DIT-FFT运算流图和程序框图， 读懂本实验提供的FFT子程序。 (3) 编制信号产生子程序， 产生以下典型信号供谱分析用：(4) 编写主程序。下图给出了主程序框图， 供参考。 本实验提供FFT子程序和通用绘图子程序。(5) 按实验内容要求， 上机实验， 并写出实验报告。 三、实验内容(1) 对 2 中所给出的信号逐个进行谱分析。 (2) 令x(n)

13、=x4(n)+x5(n)， 用FFT计算 8 点和 16 点离散傅里叶变换， X(k)=DFTx(n)(3) 令x(n)=x4(n)+jx5(n)， 重复(2)。四、思考题 (1) 在N=8时， x2(n)和x3(n)的幅频特性会相同吗? 为什么? N=16呢? (2) 如果周期信号的周期预先不知道， 如何用FFT进行谱分析?五、参考资料数字信号处理教程Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年数字信号处理 第三版，高西泉等著，西安电子科技大学出版社，第三版六、实验程序及结果x1 = 1,1,1,1,0,0,0,0,0,0;x2 = 1,2,3,4,4,3,2,1,0

14、,0;x3 = 4,3,2,1,1,2,3,4,0,0;n=0:1000;x4=cos(n);x5=sin(n);t=0:0.001:1;x6=cos(8*pi*t)+cos(16*pi*t)+cos(20*pi*t);% figure;% title()subplot(3,2,1);stem(0:length(x1)-1,x1);subplot(3,2,2);stem(0:length(x2)-1,x2);subplot(3,2,3);stem(0:length(x3)-1,x3);subplot(3,2,4);stem(0:20,x4(1:21);subplot(3,2,5);stem(0

15、:20,x5(1:21);subplot(3,2,6);stem(0:20,x6(1:21);figure;stem(0:200,x6(1:201);y1=abs(fft(x1,8);y2=abs(fft(x2,8);y3=abs(fft(x3,8);y4=abs(fft(x4,8);y5=abs(fft(x5,8);y6=abs(fft(x6,8);figure;subplot(3,2,1);stem(0:length(y1)-1,y1);subplot(3,2,2);stem(0:length(y2)-1,y2);subplot(3,2,3);stem(0:length(y3)-1,y3)

16、;subplot(3,2,4);stem(0:length(y1)-1,y1);subplot(3,2,5);stem(0:length(y5)-1,y5);subplot(3,2,6);stem(0:length(y6)-1,y6);%figure;%stem(0:15,y6(1:16); %y1=abs(fft(x1,16);y2=abs(fft(x2,16);y3=abs(fft(x3,16);y4=abs(fft(x4,16);y5=abs(fft(x5,16);y6=abs(fft(x6,16);figure;subplot(3,2,1);stem(0:length(y1)-1,y1

17、);subplot(3,2,2);stem(0:length(y2)-1,y2);subplot(3,2,3);stem(0:length(y3)-1,y3);subplot(3,2,4);stem(0:length(y6)-1,y4);subplot(3,2,5);stem(0:length(y6)-1,y5);subplot(3,2,6);stem(0:length(y6)-1,y6);%figure;%stem(0:15,y6(1:16);x7=x4+j*x5;y71=abs(fft(x7,8);y72=abs(fft(x7,16);figure;subplot(1,2,1);stem(

18、0:length(y71)-1,y71);subplot(1,2,2);stem(0:length(y72)-1,y72); 程序结果如下所示：每个信号时域分析因为第六个信号不易在上图中展示（包含信息不到一个周期）所以单独画出图形如下：x6的时域波形逐个进行8点离散傅里叶变换逐个进行16点离散傅里叶变换x(n)8点和16点离散傅里叶变换比较 七、回答问题(1) 在N=8时，x2(n)和x3(n)的幅频特性会相同吗? 为什么? N=16呢?答：N=8时幅频特性一样，N=16时幅频特性不一样。Fn所能分辨到频率为为Fs/N，如果采样频率Fs为1024Hz，采样点数为1024点，则可以分辨到1Hz。

19、1024Hz的采样率采样1024点，刚好是1秒，也就是说，采样1秒时间的信号并做FFT，则结果可以分析到1Hz，如果采样2秒时间的信号并做FFT，则结果可以分析到0.5Hz。如果要提高频率分辨力，则必须增加采样点数，也即采样时间。频率分辨率和采样时间是倒数关系。总结：采样点N的不同，分辨率不同，在N=8，俩函数不同的部分没有分辨出来，因此特性相同，而N=16时已经足以分辨出不同的频率，因此频谱不同。(2) 如果周期信号的周期预先不知道，如何用FFT进行谱分析?答：设一个定长的m值,先取2m,看2m/m的误差是否大，如大的话再取4m,看4m/2m的误差是否大，如不大，4m(4倍的m值)则可近似原

20、来点的谱分析。八、总结 本实验主要掌握fft函数的用法，包括采样点的意义。这个意义已经在七、回答问题中描述，不再赘述。实验三： 用窗函数法设计FIR数字滤波器一、实验目的(1) 掌握用窗函数法设计FIR数字滤波器的原理和方法。(2) 熟悉线性相位FIR数字滤波器特性。(3) 了解各种窗函数对滤波特性的影响。二、实验内容及步骤 (1) 复习用窗函数法设计FIR数字滤波器一节内容， 阅读本实验原理， 掌握设计步骤。 (2) 编写程序。 编写能产生矩型窗、 升余弦窗、 改进升余弦窗和二阶升余弦窗的窗函数子程序。 编写主程序。 其中幅度特性要求用dB表示。窗函数法设计滤波器主程序框图三、实验内容四、思

21、考题 (1) 如果给定通带截止频率和阻带截止频率以及阻带最小衰减， 如何用窗函数法设计线性相位低通滤波器? 写出设计步骤。 (2) 如果要求用窗函数法设计带通滤波器， 且给定上、 下边带截止频率为1和2，试求理想带通的单位脉冲响应hd(n)。 五、参考资料 数字信号处理教程Matlab释义与实现 第二版，陈怀琛著，电子工业出版社，2008年数字信号处理 第三版，高西泉等著，西安电子科技大学出版社，第三版XX知道六、实验程序及结果myfreqz(b,a)函数，是自编的改进freqz()的函数。function db,mag,pha,grd,w =myfreqz(b,a);%advanced fr

22、eqz%db,mag,pha,grd,w =myfreqz(b,a);H,w=freqz(b,a,1000,whole);H=(H(1:1:501);w=(w(1:1:501);mag=abs(H);db=20*log10(mag+eps)/max(mag);pha=angle(H);grd=grpdelay(b,a,w);以下是主程序部分wp=0.2*pi;ws=0.3*pi;deltaw=ws-wp;N0=ceil(6.6*pi/deltaw);N=N0+mod(N0+1,2);wdham=(hamming(N);wc=(ws+wp)/2;% hd=ideallp(wc,N);tao=(N

23、-1)/2;n=0:(N-1);m=n-tao+eps;hd=sin(wc*m)./(pi*m);h=hd.*wdham;db,mag,pha,grd,w=myfreqz(h,1);dw=2*pi/1000;Rp=-(min(db(1:wp/dw+1);As=-round(max(db(ws/dw+1:501);figure;subplot(2,2,1);stem(0:length(hd)-1,hd,.);subplot(2,2,2);stem(0:length(wdham)-1,wdham,.);subplot(2,2,3);stem(0:length(h)-1,h,.);subplot(2

24、,2,4);plot(1/length(db):1/length(db):1,db); %wn1=boxcar(N);figure;subplot(2,2,1);stem(0:N-1,wn1,.);title(boxcar);wn2=bartlett(N);subplot(2,2,2);stem(0:N-1,wn2,.);title(bartlett);wn3=hanning(N);subplot(2,2,3);stem(0:N-1,wn3,.);title(hanning);% wn4=hamming(N);wn4=blackman(N);subplot(2,2,4);stem(0:N-1,

25、wn4,.);title(blackman);% wn=(N,beta);得到的结果如下图：其中，图一是理想脉冲响应，第一行图二是汉明窗；图三是实际脉冲响应，图四是幅度响应。一下是几个典型窗函数的图像。，汉明窗上例已画出，不再画。每个窗的名字已在每个子图上方标出。典型窗函数七、回答问题(1) 给定通带截止频率和阻带截止频率以及阻带最小衰减，用窗函数法设计线性相位低通滤波器的设计步骤： 答： 技术指标 Wp=0.2*pi,Ws=0.4*pi,Ap=0.25dB,As=50dB方法一选择海明窗clear all;Wp=0.2*pi;Ws=0.4*pi;tr_wide=Ws-Wp; %过渡带宽度N=

26、ceil(6.6*pi/tr_wide)+1; %滤波器长度n=0:1:N-1;Wc=(Wp+Ws)/2; %理想低通滤波器的截止频率hd=ideal_lp1(Wc,N); %理想滤波器的单位冲击响应w_ham=(hamming(N); %海明窗h=hd.*w_ham; %实际海明窗的响应db,mag,pha,w=freqz_m2(h,1); %计算实际滤波器的幅度响应delta_w=2*pi/1000;Ap=-(min(db(1:1:Wp/delta_w+1) %实际通带纹波As=-round(max(db(Ws/delta_w+1:1:501) %实际阻带纹波subplot(221)ste

27、m(n,hd)title(理想单位脉冲响应hd(n)subplot(222)stem(n,w_ham)title(海明窗)subplot(223)stem(n,h)title(实际单位脉冲响应hd(n)subplot(224)plot(wi/pi,db)title(幅度响应(dB)axis(0,1,-100,10)方法二Window=blackman(16);b=fir1(15,0.3*pi ,low,Window);freqz(b,128) (2)用窗函数法设计带通滤波器， 且给定上、 下边带截止频率为1和2，理想带通的单位脉冲响应hd(n)的求解过程： 答：由hd(n)=1/2w1woe-

28、jwaejnwdw+1/2w2w3e-jwaejnwdw （其中w0=-w0-wc,w1=-w0+wc,w2=w0-wc,w3=w0+wc） 计算整理后可得： hd（n）=2/(n-a)*)*sin(n-a)wc*cos(n-a)w0 =2wc/* sa(n-a)wc*cos(n-a)w0实验四：用双线性变换法设计IIR数字滤波器一、实验目的(1) 掌握用双线性变换法设计IIR数字滤波器的原理和方法。(2) 熟悉IIR数字滤波器特性。二、实验内容及步骤 (1) 复习用双线性变换法设计IIR数字滤波器一节内容， 阅读本实验原理， 掌握设计步骤。 (2) 编写程序。 设计产生巴特沃斯滤波器原型，采用双线性变换法设计巴特沃斯数字低通滤波器。 编写主程序。 其中幅度特性要求用dB表示。三、上机实验内容 (1) 实验所需设计的IIR数字低通滤波器的指标为wp=0.2*pi，ws=0.3*pi，Rp=1dB，As=15dB，滤波器