数字信号处理matlab实验第三章报告

1、第三章 离散时间信号的频域分析一、实验目的通过matlab仿真进一步研究离散时间系统的时域特性。二、实验要求1 学习并调试本章所给的例子。2 回答书后给出的问题。三、实验过程Q.2 运行程序P3.1,求离散时间傅里叶变换的实部、虚部以及幅度和相位谱。离散时间傅里叶变换是w的周期函数吗?若是周期是多少？描述这四个图形表示的对称性。clf;w = -4*pi:8*pi/511:4*pi;num = 2 1;den = 1 -0.6;h = freqz(num, den, w);% Plot the DTFTsubplot(2,2,1) plot(w/pi,real(h);gridtitle(H(e

2、jomega)的实部)xlabel(omega /pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,imag(h);gridtitle(H(ejomega)的虚部)xlabel(omega /pi);ylabel(振幅);subplot(2,2,3)plot(w/pi,abs(h);gridtitle(|H(ejomega)|幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,4)plot(w/pi,angle(h);gridtitle(H(ejomega)相位谱)xlabel(omega /pi);ylabel(以弧度为单位的

3、相位);由图得出是周期函数，周期为2，实部是和幅度谱是实偶对称，虚部和相位谱是实奇对称Q.4修改程序P3.1，计算如下有限长序列的离散时间傅里叶变换：gn=1 3 5 7 9 11 13 15 17并重做习题Q.2。讨论你的结果。你能解释相位谱中的跳变吗？clf;w = -4*pi:8*pi/511:4*pi;num = 1 3 5 7 9 11 13 15 17;den = 1;h = freqz(num, den, w);% Plot the DTFTsubplot(2,2,1)plot(w/pi,real(h);gridtitle(H(ejomega)的实部) xlabel(omega

4、/pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,imag(h);gridtitle(H(ejomega)的虚部)xlabel(omega /pi);ylabel(振幅);subplot(2,2,3)plot(w/pi,abs(h);gridtitle(|H(ejomega)|幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,4)plot(w/pi,angle(h);gridtitle(H(ejomega)相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);1.由图得出是周期函数，周期为2，实部是

5、和幅度谱是实偶对称，虚部和相位谱是实奇对称。2.要计算一个系统相频特性，就要用到反正切函数，计算机中反正切函数规定，在一、二象限中的角度为0pi，三四象限的角度为0-pi。若一个角度从0变到2pi，但实际得到的结果是0pi，再由-pi0，在w=pi处发生跳变。Q.6通过加入合适的注释语句和程序语句，修改程序P3.2，对程序生成的图形的两个轴加标记。哪个参数控制时移量？clf;w = -pi:2*pi/255:pi; wo = 0.4*pi; D = 10;num = 1 2 3 4 5 6 7 8 9;h1 = freqz(num, 1, w);h2 = freqz(zeros(1,D) nu

6、m, 1, w);subplot(2,2,1)plot(w/pi,abs(h1);gridtitle(原序列的幅度谱) xlabel(omega /pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,abs(h2);gridtitle(时移后序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,3)plot(w/pi,angle(h1);gridtitle(原序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);subplot(2,2,4)plot(w/pi,angle(h2);gridtitl

7、e(时移后序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);参数D控制位移量Q.8选取不同的时移值重做习题Q.6令Q.6中的位移量D为20，程序结果如图： Q.10通过加入合适的注释语句和程序语句，修改程序P3.3，对程序生成的图形的两个轴加标记。哪个参数控制频移量？clf;w = -pi:2*pi/255:pi; wo = 0.4*pi;num1 = 1 3 5 7 9 11 13 15 17;L = length(num1);h1 = freqz(num1, 1, w);n = 0:L-1; num2 = exp(wo*i*n).*num1;h2 = f

8、reqz(num2, 1, w);subplot(2,2,1)plot(w/pi,abs(h1);gridtitle(原序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,abs(h2);gridtitle(频移后序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,3)plot(w/pi,angle(h1);gridtitle(原序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);subplot(2,2,4)plot(w/pi,angle(h2);g

9、ridtitle(频移后序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);参数wo控制频移量Q.12选取不同的频移值重做习题Q.10令Q.10中wo值为，程序结果如图Q.14通过加入合适的注释语句和程序语句，修改程序P3.4，对程序生成的图形的两个轴加标记。clf;w = -pi:2*pi/255:pi;x1 = 1 3 5 7 9 11 13 15 17;x2 = 1 -2 3 -2 1;y = conv(x1,x2);h1 = freqz(x1, 1, w);h2 = freqz(x2, 1, w);hp = h1.*h2;h3 = freqz(y,1,

10、w); subplot(2,2,1)plot(w/pi,abs(hp);gridtitle(幅度谱的乘积)xlabel(omega /pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,abs(h3);gridtitle(卷积后序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,3)plot(w/pi,angle(hp);gridtitle(相位谱的和)xlabel(omega /pi);ylabel(以弧度为单位的相位);subplot(2,2,4)plot(w/pi,angle(h3);gridtitle(卷积后序列的相

11、位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);Q.16选取两个改变了长度的序列，重做习题Q.14更改Q.14中的两个序列为X1=1 1 2 4,X2=4 2 3 1,程序结果如图：Q.18通过加入合适的注释语句和程序语句，修改程序P3.5，对程序生成的图形的两个轴加标记。clf;w = -pi:2*pi/255:pi;x1 = 1 3 5 7 9 11 13 15 17;x2 = 1 -1 1 -1 1 -1 1 -1 1;y = x1.*x2;h1 = freqz(x1, 1, w);h2 = freqz(x2, 1, w);h3 = freqz(y,1,w)

12、; subplot(3,1,1)plot(w/pi,abs(h1);gridtitle(第一个序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(3,1,2)plot(w/pi,abs(h2);gridtitle(第二个序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(3,1,3)plot(w/pi,abs(h3);gridtitle(乘积序列的幅度谱)xlabel(omega /pi);ylabel(振幅);Q.20通过加入合适的注释语句和程序语句，修改程序P3.6，对程序生成的图形的两个轴加标记。试解释程序怎样进行时

13、间反转运算。clf;w = -pi:2*pi/255:pi;num = 1 2 3 4;L = length(num)-1;h1 = freqz(num, 1, w);h2 = freqz(fliplr(num), 1, w); h3 = exp(w*L*i).*h2;subplot(2,2,1)plot(w/pi,abs(h1);gridtitle(原序列的幅度谱)xlabel(omega /pi);ylabel(振幅);subplot(2,2,2)plot(w/pi,abs(h3);gridtitle(时间反转后序列的幅度谱)xlabel(omega /pi);ylabel(振幅);sub

14、plot(2,2,3)plot(w/pi,angle(h1);gridtitle(原序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);subplot(2,2,4)plot(w/pi,angle(h3);gridtitle(时间反转后序列的相位谱)xlabel(omega /pi);ylabel(以弧度为单位的相位);程序先将num沿时间原点翻转，再进行离散时间傅里叶变换，然后变换后的函数乘以时移分量exp(w*L*i)。Q.22选取两个不同长度的序列，重做习题Q.201.将Q.20中的序列改为num=1 1 2 4 5,程序结果如图：2. 将Q.20中的序列

15、改为num=1 6 8 3 1 2,程序结果如图： Q.24编写一个MATLAB程序，用一个N点复数离散傅里叶变换计算两个长度为N的实数序列的N点离散傅里叶变换，并将结果同直接使用两个N点离散傅里叶变换得到的结果进行比较。clf;g=1 1 2 4;h=2 3 2 1;x=g+i*h;N = length(x)-1; n = 0:N;gk=fft(g);hk=fft(h);xk=fft(x);xk1=fft(conj(x);gk1=(xk+xk1)/2;hk1=(xk-xk1)/2i;subplot(4,2,1)stem(n,abs(gk);gridtitle(实部序列gk的离散傅里叶变换的幅

16、度)xlabel(时间序号n);ylabel(振幅);subplot(4,2,2)stem(n,abs(hk);gridtitle(虚部序列gk的离散傅里叶变换的幅度)xlabel(时间序号n);ylabel(振幅);subplot(4,2,3)stem(n,abs(gk1);gridtitle(通过xk得到的gk1的离散傅里叶变换的幅度)xlabel(时间序号n);ylabel(振幅);subplot(4,2,4)stem(n,abs(hk1);gridtitle(通过xk得到的hk1的离散傅里叶变换的幅度)xlabel(时间序号n);ylabel(振幅);subplot(4,2,5)ste

17、m(n,angle(gk);gridtitle(实部序列gk的离散傅里叶变换的相位)xlabel(时间序号n);ylabel(以弧度为单位的相位);subplot(4,2,6)stem(n,angle(hk);gridtitle(虚部序列hk的离散傅里叶变换的相位)xlabel(时间序号n);ylabel(以弧度为单位的相位);subplot(4,2,7)stem(n,angle(gk1);gridtitle(通过xk得到的gk1的离散傅里叶变换的相位)xlabel(时间序号n);ylabel(以弧度为单位的相位);subplot(4,2,8)stem(n,angle(hk1);gridtit

18、le(通过xk得到的hk1的离散傅里叶变换的相位)xlabel(时间序号n);ylabel(以弧度为单位的相位);程序运行结果如下，说明两种方法得到的结果一样：Q.26在函数circshift中，命令rem的作用是什么？函数cirshift中rem是用来确定要补的零的个数。Q.28在函数circonv中，运算符=的作用是什么？函数circonv中=是用来确定两个序列的长度是否相等。Q.30通过加入合适的注释语句和程序语句，修改程序P3.7，对程序生成的图形的两个轴加标记。哪一个参数决定时移量？若时移量大于序列长度，将会发生什么？clf;M = 6;a = 0 1 2 3 4 5 6 7 8 9

19、;b = circshift(a,M);L = length(a)-1;n = 0:L;subplot(2,1,1);stem(n,a);axis(0,L,min(a),max(a);title(原序列)xlabel(时间序号n);ylabel(a);subplot(2,1,2);stem(n,b);axis(0,L,min(a),max(a);title(圆周移位,num2str(M),个样本得到的序列)xlabel(时间序号n);ylabel(a);M控制时移量，若时移量大于序列长度会得到时移M1=M-序列长度 的序列 下图是M=6时的结果，此时时移量M小于序列长度： 下图是M=13时的结

20、果，此时时移量M大于序列长度： Q.32通过加入合适的注释语句和程序语句，修改程序P3.8，对程序生成的图形的两个轴加标记。clf;x = 0 2 4 6 8 10 12 14 16;N = length(x)-1; n = 0:N;y = circshift(x,5);XF = fft(x);YF = fft(y);subplot(2,2,1)stem(n,abs(XF);gridtitle(原序列的离散傅里叶变换的幅度);xlabel(时间序号n);ylabel(振幅);subplot(2,2,2)stem(n,abs(YF);gridtitle(圆周移位后的序列的离散傅里叶变换的幅度);

21、xlabel(时间序号n);ylabel(振幅);subplot(2,2,3)stem(n,angle(XF);gridtitle(原序列的离散傅里叶变换的相位);xlabel(时间序号n);ylabel(以弧度为单位的相位);subplot(2,2,4)stem(n,angle(YF);gridtitle(圆周移位后的序列的离散傅里叶变换的相位);xlabel(时间序号n);ylabel(以弧度为单位的相位);时移量是5.结果如下图：Q.34选取两个不同的时移量，重做习题Q.321.时移量为10，结果如下图：2.时移量为15，结果如下图：Q.36运行程序P3.9并验证离散傅里叶变换的圆周卷积

22、性质clf;g1 = 1 2 3 4 5 6; g2 = 1 -2 3 3 -2 1;ycir = circonv(g1,g2);disp(圆周卷积的结果 = );disp(ycir)G1 = fft(g1); G2 = fft(g2);yc = real(ifft(G1.*G2);disp(离散傅里叶变换乘积的离散傅里叶逆变换的结果 = );disp(yc)结果如下图：两个结果一致，说明两个序列循环卷积等于这两个序列分别的DFT相乘后的IDFTQ.38运行程序P3.10并验证线性卷积可通过圆周卷积得到clf;g1 = 1 2 3 4 5;g2 = 2 2 0 1 1;g1e = g1 zer

23、os(1,length(g2)-1);g2e = g2 zeros(1,length(g1)-1);ylin = circonv(g1e,g2e);disp(通过圆周卷积的线性卷积 = );disp(ylin);y = conv(g1, g2);disp(直接线性卷积 = );disp(y)结果如下图：两个结果一致，说明线性卷积可通过圆周卷积得到Q.40编写一个MATLAB程序，对两个序列分别做离散傅里叶变换，以生成他们的线性卷积。用此程序验证习题Q.38的结果clf;g1 = 1 2 3 4 5;g2 = 2 2 0 1 1;g1e = g1 zeros(1,length(g2)-1);g2

24、e = g2 zeros(1,length(g1)-1);G1e = fft(g1e); G2e = fft(g2e);ylin = real(ifft(G1e.*G2e);disp(通过离散傅里叶变换得到的线性卷积 = );disp(ylin);y = conv(g1, g2);disp(直接线性卷积 = );disp(y)结果如下图：两者结果一致，验证了习题Q.38的结果Q.42运行程序P3.11。由于周期序列的偶数部分的离散傅里叶变换是原序列的XEF的实数部分，XEF的虚部应该为零。你能验证它们吗？你怎样解释仿真结果？x = 1 2 4 2 6 32 6 4 2 zeros(1,247)

25、;x1 = x(1) x(256:-1:2);xe = 0.5 *(x + x1);XF = fft(x); XEF = fft(xe);clf;k = 0:255;subplot(2,2,1);plot(k/128,real(XF); grid;ylabel(振幅);title(Re(DFTxn);subplot(2,2,2);plot(k/128,imag(XF); grid;ylabel(振幅);title(Im(DFTxn);subplot(2,2,3)plot(k/128,real(XEF); grid;xlabel(时间序号n);ylabel(振幅);title(Re(DFTx_e

26、n);subplot(2,2,4);plot(k/128,imag(XEF); grid;xlabel(时间序号n);ylabel(振幅);title(Im(DFTx_en);Q.44运行程序P3.12。你得到的a和b 的值相等吗？clf;x = (1:128) (128:-1:1);XF = fft(x);a = sum(x.*x)b = round(sum(abs(XF).2)/256)运行结果如图：a和b的值相等Q.46使用程序P3.1在单位圆上求下面的z变换：Gz=2+5z-1+9z-2+5z-3+3z-45+45z-1+2z-2+z-3+z-4clf;w = -pi:2*pi/511:pi;num = 2 5 9 5 3;den = 5 45 2 1 1;h = f