DSP课程大作业英文版_图文(精)

Lab Report of DSP

姚志浩 5090309481 F0903018

Lab 1 Discrete-Time Fourier Transform

(1 Generate and plot the discrete-time signal

Solution:

It’s a simple sequence ,we can get the discrete-time signal using Matlab.

Chat 1-1 simple sequence sample

(2 Determine and plot the DTFT of the signal

S olution:

DTFT is the Z-Transform on the Unit Circle, we can easily find the Z-Transform ：

W ith ‘freqz’ function，we can get the amplitude-frequency and phase frequency response curve.

Chart 1-2 DTFT Of the sequence

(3 A linear and time-invariant system is described by the difference equation

where abcd is equal to the last four digits of your student ID number.

(i Is the system BIBO-stable?

(ii Determine and sketch the impulse response , of the system. Determine the stability of the system by

observing .

(iii Determine and plot the output of the system, if the input is

.

S olution:

(1 My student ID number is 5090309481，corresponding a b c d will be9 4 8 1,and the difference equation will be

；

Characteristic equation will be ，we can get three pole points : ，

All the pole points are in the unit circle, so the system is stable.

(2 C reate a function “impseq” to get a impulse signal.

H ere we use the function “filter”, we can get the impulse response:

Chart 1-3 Unit Step Response

W e get ans = 22.1234, so the system is stable.

(3 if the input is ,the input sequence and output sequence are showed :

A ddendum：Matlab Program List 1-1.

n=[0:30]

x=(0.93.^n.*cos(0.14*pi*n+pi/3;

stem(n,x;title('x(n';xlabel('n';ylabel('x(n'

1-2 .

n=[0:30];x=(0.89.^n.*(stepseq(0,30,0

w=[0:500]*2*pi/500;

X=exp(j*w./(exp(j*w-0.89*ones(1,501;

mag=abs(X;ang=angle(X

figure(1;plot(w/pi,mag;grid

ylabel('|X(w|';xlabel('w(in pi unit';title('DTFT of x(n-magnitude'

figure(2;plot(w/pi,ang;grid

ylabel('argX(w/pi';xlabel('w (in pi unit';title('DTFT of x(n-angle'

1-3.

(2 num=[9,-4,8,-1];den=[1,-0.65,0.35,0]

x=impseq(0,0,100;n=[0:100]

h=filter(num,den,x

stem(n,h;xlabel('n';ylabel('y(n';title('Impluse Response' ;sum(abs(h

(3num=[9,-4,8,-1];den=[1,-0.65,0.35,0]

n=[0:200];x=(2+5*cos(0.27*pi*n+3*sin(0.77*pi*n.*stepseq(0,200,0

figure(1;stem(n,x;xlabel('n';ylabel('x(n';title('Input Sequence';

h=filter(num,den,x;figure(2;stem(n,h;xlabel('n';ylabel('y(n';title('Output Response' Functions :

(1function [x,n]=impseq(n0,n1,n2

n=[n1:n2];x=[(n-n0==0];

(2function [x, n] = stepseq(n1,n2,n0

if nargin ~=3

disp('input error;';

return;

elseif ((n0 n2|(n1>n2

error(' input error?'

end

n = n1:n2;

x = ((n-n0>=0;

Lab 2Spectral Analysis

(1 A signal , where and

, is plagued by a Gaussian white noise . Suppose that the signal-to-noise ratio is -3 or 0 dB and the sampling rate is

. Plot the spectra of and

, 。

S olution:

Without loss of generality, suppose that the amplitudes of the two frequency components are different, they are 2 and 1.5. So the signal without noise is :.we can use the function “awgn” to add a Gaussian white Noise, and here we assume the radio is -2dB. So we can get the charts:

With the function “fft”, we can get the DFT and the frequency spectrum. Com paring these two pictures, we can find the primary signal only has responses with the frequencies of 100Hz and 50Hz. After adding noise, we can see there are lots of responses with small amplitude of different frequency.

(2An analog signal is sampled with a sampling rate of . Suppose that the data record length is or . Determine the 200-point

DFT of . Plot .

Solution:

First, we create two sample sequences with L=20 and L=100, then we can use the function …fft? to get the DFT. From the charts, we c an find that L=20 is not long enough, so the sample sequence cannot describe the true characters of the analog signal. However, L=100 is long enough, so the sample sequence is close to the analog signal. So we can say that if L is larger, the DFT will be closer to the true analog signal.

Addendum ：Matlab Program List

2-1

t=0:0.001:1.023 ; f = (0:1023/(1024*0.001;

x =2*sin(100*pi*t + 1.5*sin(220*pi*t;

figure(1; plot(t,x ; title('原始信号';xlabel('t';

y=awgn(x,-2;

figure(2;plot(t,y;title('加噪信号';xlabel('t';

X=fft(x;figure(3;plot(f,abs(X;title('原始信号频谱';xlabel('Frequency (Hz';ylabel('|X|'

Y=fft(y;figure(4;plot(f,abs(Y;title('加噪信号频谱';xlabel('Frequency (Hz';ylabel('|Y|'

2-2

fs=20;

n=0:19;x1=exp(-n/fs;

m=0:99;x2=exp(-m/fs;

X1=fft(x1,200;

X2=fft(x2,200;

figure(1;

subplot(2,1,1;stem(n,x1;title('Sample sequence of L=20 ';xlabel('n';

subplot(2,1,2;k=0:1:199;stem(k,abs(X1;title('L=20 N=200 DFT';xlabel('frequency index k';ylabel('Magnitude Spectrum';

figure(2;

subplot(2,1,1;stem(m,x2;title(' Sample sequence of L=100 ';xlabel('n';

subplot(2,1,2;stem(k,abs(X2;title('L=100 N=200 DFT';xlabel('frequency index

k';ylabel('Magnitude Spectrum';

Lab 3 Digital Filter Design

(1 The transfer function of a discrete-time system is

where abcd is equal to the last four digits of your stu dent?s ID number. Determine the impulse response and step response of the system.

S olution:

M y student ID number is 5090309481, so a=9, b=4,c=8,d=1.In Matlab, we can directly create two vectors b=[1,-9,4,-8,1,0] and a=[1,0.22,0.037,0.142,-0.107,-0.013]. With the function …Filter?, we can get the output responses of unit impulse and step signal. T he output responses are shown：

(2 Plot the magnitude squared frequency response of 5th-, 10th-, 20th-, and x th-order

lowpass Butterworth filter, where

S olution:

My class number is F0903018, my student ID number is 5090309481.

So X=（9+3+1+8+5+9+3+9+4+8+1）mod (20=0.

T he magnitude squared frequency response of a Nth-order lowpass Butterworth filter is:

B ecause we do n?t know , so here we use w=as horizontal ordinate.

(3 Design an FIR bandstop filter with the following specifications:

Determine and plot and .

S olution:

A ccording to the problem, Ac=41dB, and the stopband attenuation of Hanning Window is 44dB. We can see 41 is very close to 44, S o here we choose Hanning Window to solve the problem.According to the design requirements of Hanning Window, we get Ap =0.1017dB, As =44d

B and M=125. T he following four charts show the ideal unit impulse response, Hanning Window, actual unit impulse response, amplitude frequency characteristics curve

of Band elimination filter:

Addendum ：Matlab Program List

3-1

b=[1,-9,4,-8,1,0];a=[1,0.22,0.037,0.142,-0.107,-0.013];n=0:50;

x1=impseq(0,0,50;y1=filter(b,a,x1;

figure(1;stem(n,y1;title('Unit Impulse Response';xlabel('n';

x2=stepseq(0,50,0;y2=filter(b,a,x2;

figure(2;stem(n,y2;title('Unit Step Response';xlabel('n';

3-2

w=0:0.01:2;N1=5;N2=10;N3=20;N4=0;

h1=1./(1+(w.^(2*N1;h2=1./(1+(w.^(2*N2;

h3=1./(1+(w.^(3*N3;h4=1./(1+(w.^(4*N4;

subplot(2,2,1;plot(w,h1;title('5th-order ';ylabel('〖|H（jw）|〗^2';xlabel('w/wc'; subplot(2,2,2;plot(w,h2;title('10th-order ';ylabel('〖|H（jw）|〗^2';xlabel('w/wc'; subplot(2,2,3;plot(w,h3;title('20th-order ';ylabel('〖|H（jw）|〗^2';xlabel('w/wc'; subplot(2,2,4;plot(w,h4;title('0th-order ';ylabel('〖|H（jw）|〗^2';xlabel('w/wc'; 3-3

wp1=0.4*pi;ws1=0.45*pi;wp2=0.7*pi;ws2=0.65*pi;

tr_width=min((ws1-wp1,(wp2-ws2;M=ceil(6.2*pi/tr_width+1

n=0:1:M-1;wc1=(ws1+wp1/2;wc2=(ws2+wp2/2;

hd=ideal_bs(wc1,wc2,M;

w_hanning=(hanning(M';

h=hd.*w_hanning;

[db,mag,pha,w]=freqz_m(h,[1];delta_w=2*pi/1000;

Ap=-(min(db(1:1:wp1/delta_w+1

As=-round(max(db(ws1/delta_w+1:1:ws2/delta_w+1

subplot(2,2,1;stem(n,hd, '. ';title('Ideal Impulse Response';

subplot(2,2,2;stem(n,w_hanning, '. ';title('Hanning Window';

subplot(2,2,3;stem(n,h, '. ';title('Actual Impulse response';

subplot(2,2,4;plot(w/pi,db;title('Manitude Response in dB';axis([0 1 -100 10]; Functions:

function [db,mag,pha,w]=freqz_m(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;

end

function hd=ideal_bs(wc1,wc2,M;

alpha=(M-1/2;

n=[0:1:(M-1];

m=n-alpha+eps;

hd=[sin(wc1*m+sin(pi*m-sin(wc2*m]./(pi*m; end