DSP实验报告B13011025分析

clf;

x = [0 2 4 6 8 10 12]; 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));grid

title('Magnitude of DFT of Original Sequence'); xlabel('Frequency index k'); ylabel('|X[k]|'); subplot(2,2,2) stem(n,abs(YF));grid

title('Magnitude of DFT of Circularly Shifted Sequence'); xlabel('Frequency index k'); ylabel('|Y[k]|'); subplot(2,2,3) stem(n,angle(XF));grid

title('Phase of DFT of Original Sequence'); xlabel('Frequency index k'); ylabel('arg(X[k])'); subplot(2,2,4) stem(n,angle(YF));grid

title('Phase of DFT of Circularly Shifted Sequence'); xlabel('Frequency index k'); ylabel('arg(Y[k])');

Q3.36运行程序P3.9并验证离散傅里叶变换的圆周卷积性质

g1 = [1 2 3 4 5 6]; g2 = [1 -2 3 3 -2 1]; ycir = circonv(g1,g2);

disp('Result of circular convolution = ');disp(ycir) G1 = fft(g1); G2 = fft(g2); yc = real(ifft(G1.*G2));

disp('Result of IDFT of the DFT products = ');disp(yc)

结果：

x1 =

1 2 3 4 5 6 x2 =

1 -2 3 3 -2 1 Result of circular convolution = 12 28 14 0 16 14 Result of IDFT of the DFT products = 12 28 14 0 16 14

Q3.37选取另外两组等长序列重做习题Q3.36。 g1 = [2 3 7 1 6 -5]; g2 = [-4 6 7 0 2 9]; ycir = circonv(g1,g2);

disp('Result of circular convolution = ');disp(ycir) G1 = fft(g1); G2 = fft(g2); yc = real(ifft(G1.*G2));

disp('Result of IDFT of the DFT products = ');disp(yc) x1 =

2 3 7 1 6 -5 x2 =

-4 6 7 0 2 9 Result of circular convolution = 45 30 25 103 -10 87 Result of IDFT of the DFT products =

45.0000 30.0000 25.0000 103.0000 -10.0000 87.0000 Q3.38运行程序P3.10并验证线性卷积可通过圆周卷积得到。 g1 = [1 2 3 4 5];g2 = [2 2 0 1 1]; g1e = [g1 zeros(1,length(g2)-1)]; g2e = [g2 zeros(1,length(g1)-1)]; ylin = circonv(g1e,g2e);

disp('Linear convolution via circular convolution = ');disp(ylin); y = conv(g1, g2);

disp('Direct linear convolution = ');disp(y) x1 =

1 2 3 4 5 0 0 0 0 x2 =

2 2 0 1 1 0 0 0 0 Linear convolution via circular convolution =

2 6 10 15 21 15 7 9 5 Direct linear convolution =

2 6 10 15 21 15 7 9 5

Q3.39选取两组长的不等的序列重做习题Q3.38。 g1 = [1 2 3 4 5];g2 = [2 2 0 1 1 5]; g1e = [g1 zeros(1,length(g2)-1)]; g2e = [g2 zeros(1,length(g1)-1)]; ylin = circonv(g1e,g2e);

disp('Linear convolution via circular convolution = ');disp(ylin); y = conv(g1, g2);

disp('Direct linear convolution = ');disp(y) x1 =

1 2 3 4 5 0 0 0 0 0 x2 =

2 2 0 1 1 5 0 0 0 0 Linear convolution via circular convolution =

2 6 10 15 21 20 17 24 25 25 Direct linear convolution =

2 6 10 15 21 20 17 24 25 25

Q3.40编写一个MATLAB程序，对两个序列做离散傅里叶变换，以生成他们的线性卷积。用此程序验证习题Q3.38和习题Q3.39的结果。 验证Q3.38:

g1=[1 2 3 4 5]; g2=[2 2 0 1 1];

g1e = [g1 zeros(1,length(g2)-1)]; g2e = [g2 zeros(1,length(g1)-1)]; G1EF=fft(g1e); G2EF=fft(g2e);

ylin=real(ifft(G1EF.*G2EF));

disp('Linear convolution via DFT ='); disp(ylin);

Linear convolution via DFT =

2.0000 6.0000 10.0000 15.0000 21.0000 15.0000 7.0000 9.0000

5.0000 验证Q3.39:

g1=[1 2 3 4 5]; g2=[2 2 0 1 1 5];

g1e = [g1 zeros(1,length(g2)-1)]; g2e = [g2 zeros(1,length(g1)-1)]; G1EF=fft(g1e); G2EF=fft(g2e);

ylin=real(ifft(G1EF.*G2EF));

disp('Linear convolution via DFT ='); disp(ylin);

Linear convolution via DFT =

2.0000 6.0000 10.0000 15.0000 21.0000 20.0000 17.0000 24.0000 25.0000 25.0000

2?5z?1?9z?2?5z?3?3z?4Q3.46使用程序P3.1在单位圆生求下面的z变换：G(z)? ?1?2?3?45?45z?2z?z?zclf;

w = -4*pi:8*pi/511:4*pi;

num = [2 5 9 5 3];den = [5 45 2 1 1]; h = freqz(num, den, w); subplot(2,1,1)

plot(w/pi,real(h));grid

title('Real part of H(e^{j\\omega})') xlabel('\\omega /\\pi'); ylabel('Amplitude'); subplot(2,1,2)

plot(w/pi,imag(h));grid

title('Imaginary part of H(e^{j\\omega})') xlabel('\\omega /\\pi'); ylabel('Amplitude'); pause

subplot(2,1,1)

plot(w/pi,abs(h));grid

title('Magnitude Spectrum |H(e^{j\\omega})|') xlabel('\\omega /\\pi'); ylabel('Amplitude'); subplot(2,1,2)

plot(w/pi,angle(h));grid

title('Phase Spectrum arg[H(e^{j\\omega})]')