《DSP using MATLAB》Problem 5.10

代码：

第1小题：

%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

%%            Output Info about this m-file

fprintf('\n***********************************************************\n');

fprintf('        <DSP using MATLAB> Problem 5.10 \n\n');

banner();

%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

% --------------------------------------------------------------

%              1  x(n) = 0.6^abs(n)[u(n+10)-u(n-11)]

% --------------------------------------------------------------

nn1 = [-15:16];

xx1 = 0.6 .^ abs(nn1) .* (stepseq(-10, min(nn1), max(nn1)) - stepseq(11, min(nn1), max(nn1)));

NN1 = length(xx1);                                % length is 32

m = mod_1(nn1, NN1);

x = [xx1 zeros(1, 8)];                        % padding zeros

n = [m max(m)+1:max(m)+8];

figure('NumberTitle', 'off', 'Name', 'P5.10.1 x(n)')

set(gcf,'Color','white');

subplot(2,1,1); stem(nn1, xx1);

xlabel('n'); ylabel('x(n)');

title('x(n) ori sequence');  grid on;

subplot(2,1,2); stem(n, x);

xlabel('n'); ylabel('x(n)');

title('x(n) padding zeros');  grid on;

%% =============================================================================

%%                 DTFT X(w) of xn sequence, w=[0:2pi],

%% =============================================================================

MM = 500;

[Xw_DTFT, w] = dtft1(x, n, MM);

 magXw_DTFT = abs(Xw_DTFT);   angXw_DTFT = angle(Xw_DTFT)/pi;

realXw_DTFT = real(Xw_DTFT); imagXw_DTFT = imag(Xw_DTFT);

%% --------------------------------------------------------------

%%        START X_DTFT's  mag ang real imag

%% --------------------------------------------------------------

figure('NumberTitle', 'off', 'Name', 'P5.10.1 X(w) DTFT of x(n)');

set(gcf,'Color','white');

subplot(2,2,1); plot(w/pi,magXw_DTFT); grid on;  % axis([-2,2,0,15]);

title('Magnitude Part');

xlabel('frequency in \pi units'); ylabel('Magnitude  |X\_DTFT|');

subplot(2,2,3); plot(w/pi, angXw_DTFT); grid on;  % axis([-2,2,-1,1]);

title('Angle Part');

xlabel('frequency in \pi units'); ylabel('Rad \pi'); %axis([-200,200,0,2]);

subplot('2,2,2'); plot(w/pi, realXw_DTFT); grid on;

title('Real Part');

xlabel('frequency in \pi units'); ylabel('Real');

subplot('2,2,4'); plot(w/pi, imagXw_DTFT); grid on;

title('Imaginary Part');

xlabel('frequency in \pi units'); ylabel('Imaginary');

%% --------------------------------------------------------------

%%        END X_DTFT's  mag ang real imag

%% --------------------------------------------------------------

%% ------------------------------------------------------------------

%%                 DFT(k) of xn sequence, k=[0:N-1]

%%                      w=2pi*k/N       k=Nw/(2pi)

%% ------------------------------------------------------------------

N1 = length(x);

k1 = [0 : N1-1];

%k2 = [-N : N-1];

%k3 = [-N/2 : N/2];

Xk_DFT = dft(x, N1);                                         % DFT

    magXk_DFT = abs( [ Xk_DFT ] );                          % DFT magnitude

    angXk_DFT = angle( [Xk_DFT] )/pi;                       % DFT angle

   realXk_DFT = real(Xk_DFT); imagXk_DFT = imag(Xk_DFT);

figure('NumberTitle', 'off', 'Name', 'P5.10.1 DFT(k) of x(n)')

set(gcf,'Color','white');

subplot(2,1,1); stem(k1, magXk_DFT); hold on; plot(N1*w/(2*pi), magXw_DTFT,'r--'); hold off;

%axis([-N/2, N/2, -0.5, 50.5]);

xlabel('k'); ylabel('magnitude(k)');

title('DFT magnitude of x(n), N=40');  grid on;

subplot(2,1,2); stem(k1, angXk_DFT);  hold on; plot(N1*w/(2*pi), angXw_DTFT,'r--'); hold off;

%axis([-N/2, N/2, -0.5, 50.5]);

xlabel('k'); ylabel('angle(k)');

title('DFT angle of x(n), N=40');  grid on;

　　运行结果：

第4小题：

%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

%%            Output Info about this m-file

fprintf('\n***********************************************************\n');

fprintf('        <DSP using MATLAB> Problem 5.10 \n\n');

banner();

%% ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

% ------------------------------------------------------------------

%            4  x(n) = [1 2 3 4 3 2 1]

%                             |

% ------------------------------------------------------------------

nn1 = [-3:3];

xx1 = [1 2 3 4 3 2 1];

NN1 = length(xx1);                                % length is 7

m = mod_1(nn1, NN1);

x = [xx1 zeros(1, 9)];                            % padding 9 zeros

n = [m max(m)+1:max(m)+9];                        % length is 16

figure('NumberTitle', 'off', 'Name', 'P5.10.4 x(n)')

set(gcf,'Color','white');

subplot(2,1,1); stem(nn1, xx1);

xlabel('n'); ylabel('x(n)');

title('x(n) ori sequence');  grid on;

subplot(2,1,2); stem(n, x);

xlabel('n'); ylabel('x(n)');

title('x(n) padding zeros');  grid on;

%% =============================================================================

%%                 DTFT X(w) of xn sequence, w=[0:2pi],

%% =============================================================================

MM = 500;

[Xw_DTFT, w] = dtft1(xx1, m, MM);

 magXw_DTFT = abs(Xw_DTFT);   angXw_DTFT = angle(Xw_DTFT)/pi;

realXw_DTFT = real(Xw_DTFT); imagXw_DTFT = imag(Xw_DTFT);

%% --------------------------------------------------------------

%%        START X_DTFT's  mag ang real imag

%% --------------------------------------------------------------

figure('NumberTitle', 'off', 'Name', 'P5.10.4 X(w) DTFT of x(n)');

set(gcf,'Color','white');

subplot(2,2,1); plot(w/pi,magXw_DTFT); grid on;  % axis([-2,2,0,15]);

title('Magnitude Part');

xlabel('frequency in \pi units'); ylabel('Magnitude  |X\_DTFT|');

subplot(2,2,3); plot(w/pi, angXw_DTFT); grid on;  % axis([-2,2,-1,1]);

title('Angle Part');

xlabel('frequency in \pi units'); ylabel('Rad \pi'); %axis([-200,200,0,2]);

subplot('2,2,2'); plot(w/pi, realXw_DTFT); grid on;

title('Real Part');

xlabel('frequency in \pi units'); ylabel('Real');

subplot('2,2,4'); plot(w/pi, imagXw_DTFT); grid on;

title('Imaginary Part');

xlabel('frequency in \pi units'); ylabel('Imaginary');

%% --------------------------------------------------------------

%%        END X_DTFT's  mag ang real imag

%% --------------------------------------------------------------

%% ------------------------------------------------------------------

%%                 DFT(k) of xn sequence, k=[0:N-1]

%%                      w=2pi*k/N       k=Nw/(2pi)

%% ------------------------------------------------------------------

N1 = length(x);

k1 = [0 : N1-1];

%k2 = [-N : N-1];

%k3 = [-N/2 : N/2];

Xk_DFT = dft(x, N1);                                         % DFT

    magXk_DFT = abs( [ Xk_DFT ] );                          % DFT magnitude

    angXk_DFT = angle( [Xk_DFT] )/pi;                       % DFT angle

   realXk_DFT = real(Xk_DFT); imagXk_DFT = imag(Xk_DFT);

figure('NumberTitle', 'off', 'Name', 'P5.10.4 DFT(k) of x(n)')

set(gcf,'Color','white');

subplot(2,1,1); stem(k1, magXk_DFT); hold on; plot(N1*w/(2*pi), magXw_DTFT,'r--'); hold off;

%axis([-N/2, N/2, -0.5, 50.5]);

xlabel('k'); ylabel('magnitude(k)');

title('DFT magnitude of x(n), N=16');  grid on;

subplot(2,1,2); stem(k1, angXk_DFT);  hold on; plot(N1*w/(2*pi), angXw_DTFT,'r--'); hold off;

%axis([-N/2, N/2, -0.5, 50.5]);

xlabel('k'); ylabel('angle(k)');

title('DFT angle of x(n), N=16');  grid on;

　　运行结果：

这里只放两个序列的，至于DFT的长度N选择的方法，我也没搞懂，欢迎博友赐教！