Dsp Lab -15ecl57 Part_i

DSP Lab Manual VISVESVARAYA TECHNOLOGICAL UNIVERSITY “JNANASANGAMA” BELAGAVI – 590018

KARNATAKA STATE, INDIA Page | 1

DIGITAL SIGNAL PROCESSING LAB MANUAL SUBJECT CODE: 15ECL57 V SEMESTER E&C ENGINEERING

Manual authored by

Manual Approved by

ANIL B.V. Associate Professor

Dr.Ravikumar M.S. Professor and Head

DEPT. OF ELECTRONICS AND COMMUNICATION ENGINEERING

K.V.G. COLLEGE OF ENGINEERING SULLIA, D.K. - 574327 KARNATAKA STATE, INDIA

K.V.G. College of Engineering, Sullia D.K.

DSP Lab Manual DSP Lab B.E., V Semester, EC/TC [As per Choice Based Credit System (CBCS) scheme] Page | 2

Subject Code:15ECL57 Code:15ECL57

IA Marks 20

Number of Lecture Hours/Week 01Hr 01Hr Tutorial (Instructions)+02 (Instructions)+ 02 Hours Hours Laboratory=03 Laboratory= 03 Exam Marks 80

Exam Hours 03

CREDITS – 02 Laboratory Experiments Following Experiments to be done using MATLAB MATLAB / / SCILAB / OCTAVE or Equivalent: Equivalent: 1. Verification of sampling theorem. 2. Linear and circular convolution of two given g iven sequences, Commutative, Commutative, Distributive and Associative property of convolution. 3. Auto and cross correlation of two sequences and verification of their properties 4. Solving a given difference equation. equation. 5. Computation of N point DFT of a given sequence and to plot magnitude and phase spectrum (using DFT equation and verify it by built-in routine). 6. (i) Verification of DFT properties (Linearity (L inearity and Parseval’s theorem, etc.)(ii) DFT computation of square pulse and sinc si nc function etc. 7. Design and implementation of FIR filter to meet given specifications (using different window techniques). 8. Design and implementation of IIR filter to meet given specifications. Following Experiments Experiments to be done using DSP kit 9. Linear convolution of two sequences 10. Circular convolution of two sequences 11. N-point DFT of a given sequence s equence 12. Impulse response of first order and second order system 13. Implementation of FIR filter Conduct of Practical Practical Examination: 1. All laboratory experiments are to be included for practical examination. 2. Strictly follow the instructions as printed on the cover page of answer script for breakup of marks. 3. Change of experiment is allowed only once and Marks allotted to the procedure part to be made zero. K.V.G. College of Engineering, Sullia D.K.

DSP Lab Manual DSP Lab B.E., V Semester, EC/TC [As per Choice Based Credit System (CBCS) scheme] Page | 2

Subject Code:15ECL57 Code:15ECL57

IA Marks 20

Number of Lecture Hours/Week 01Hr 01Hr Tutorial (Instructions)+02 (Instructions)+ 02 Hours Hours Laboratory=03 Laboratory= 03 Exam Marks 80

Exam Hours 03

CREDITS – 02 Laboratory Experiments Following Experiments to be done using MATLAB MATLAB / / SCILAB / OCTAVE or Equivalent: Equivalent: 1. Verification of sampling theorem. 2. Linear and circular convolution of two given g iven sequences, Commutative, Commutative, Distributive and Associative property of convolution. 3. Auto and cross correlation of two sequences and verification of their properties 4. Solving a given difference equation. equation. 5. Computation of N point DFT of a given sequence and to plot magnitude and phase spectrum (using DFT equation and verify it by built-in routine). 6. (i) Verification of DFT properties (Linearity (L inearity and Parseval’s theorem, etc.)(ii) DFT computation of square pulse and sinc si nc function etc. 7. Design and implementation of FIR filter to meet given specifications (using different window techniques). 8. Design and implementation of IIR filter to meet given specifications. Following Experiments Experiments to be done using DSP kit 9. Linear convolution of two sequences 10. Circular convolution of two sequences 11. N-point DFT of a given sequence s equence 12. Impulse response of first order and second order system 13. Implementation of FIR filter Conduct of Practical Practical Examination: 1. All laboratory experiments are to be included for practical examination. 2. Strictly follow the instructions as printed on the cover page of answer script for breakup of marks. 3. Change of experiment is allowed only once and Marks allotted to the procedure part to be made zero. K.V.G. College of Engineering, Sullia D.K.

DSP Lab Manual

Page | 3

I - CYCLE EXPERIMENTS

Expt. No. 1 2

3

4

Experiment Name

Verification of Sampling Theorem a) Linear convolutionof two sequences i) Causal ii) Noncausal b) Circular convolution of two sequences Commutative, Distributive and Associative properties of convolution. a) Auto a) Auto correlation of two sequences b) Cross correlation of two sequences Verification of the properties of Correlation. a) Solving a difference equation without initial conditions b) Solving a difference equation with initial conditions


DSP Lab Manual

Page | 4

PART - A MATLAB® PROGRAMS


DSP Lab Manual EXPT No: 1

DATE:

VERIFICATION OF SAMPLING THEOREM AIM: To execute MATLAB m-file program for the sampling of a composite Page | 5 continuous time signal and to reconstruct it from the samples based on Nyquist Sampling criteria.

THEORY:In signal processing, sampling is the mapping of a continuous-time continuous-valued signal to a discrete-time continuous-valued signal. A common example is the conversion of audio signal (a continuous signal) to a sequence of samples (a discrete-time signal).

Figure 1.1 The sampled continuous time signal

Nyquist Sampling Theorem: If the highest frequency contained in an analog signal xa(t) is Fmax = B Hz and the signal is sampled at a rate Fs> 2 Fmax = 2B, then xa(t) can be exactly recovered from its sample values using the interpolation function

() () ()

…

(1.1)

Thus xa(t) may be expressed as

() ∑∞= −∞ ( )(/) … where xa(/ ) = x (nT ) ≡ x[n] are the samples of x (t).

(1.2)

When the sampling of xa(t) is performed at the minimum sampling rate F s = 2B, the reconstruction formula   −  ∞ becomes:  (1.3) = −∞   −   …  a

 () ∑

a

( )

The sampling rate FN = 2Fmax = 2BHz is called the Nyquist rate.


DSP Lab Manual Consider a continuous-time signal, with amplitude ‘A ’ and frequency ‘F’ Hz. xa(t) = A cos(Ωt) = A cos(2πFt)

…

(1.4)

If Fs is the sampling frequency then F s= 1/T samples per second. T is the sampling interval in seconds. Page | 6 xa[nT]=Acos[ΩnT]

…

(1.5)

= Acos[2πFnT]

…

(1.6)

=Acos[2π

 n]

…

(1.7)

= Acos[2πfn]

…

(1.8)

…

(1.9)

≡ []

= Acos[ωn]

Comparing eqn1.7 and eqn1.8 it can be seen that f =

[  ], and fis called discrete time

frequency with units cycles/sample. 2πf = ω is called discrete-time angular frequency with units radian/sample; ω = ΩT describes the relationship between the D.T. angular frequency and C.T angular frequency(rad/sec) ALGORIHM: Begin Read Amplitude ‘A’ and Frequency ‘F’, t 0,Δt, t 1 define t from t 0 to t 1 in steps of Δt x 1= A1 cos( 2πF 1t ); x 2 = A2 cos( 2πF 2t ) x = x 1 + x 2 choose F s based on max( F 1, F 2 ) define n from n0 = (t 0 /T s ) to n1 = ( t 1 /T s ) in steps of 1 x n1 = A1cos[ 2πF 1nT ] x n2 = A2cos[ 2πF 2nT ] x n = x n1 + x n2

 

compute

() ()

forn0to n1in steps of 1

  ∑  (/)(  /) )  ∑ end End


)   2(−   (2) 2(−  ) 

DSP Lab Manual PROGRAM:

samp_thrm.m

clc clear all close all Page | 7

%Clear Command Window % Clear all the variables % Close all Figure Windows

%

Enter the Amplitude values

A1 = input(‘Amplitude of signal x1(t) = '); A2 = input(‘Amplitude of signal x2(t) = '); %

Enter the frequency values

F1 = input(‘Frequency of signal x1(t) = '); F2 = input(‘Frequency of signal x2(t) = '); %

Define Beginning and End Time

t0 = input(' Beginning Time = '); t1 = input(‘ End Time = '); %

Define the Time Resolution

Delta_t= input('Resolution = ); ’

%

Define x-axis as Time axis

t = t0: Delta_t : t1; % Define Analog Signals x1 = A1 * cos( 2 * pi * F1 * t ); x2 = A2 * cos( 2 * pi * F2 * t ); %

Generate the single composite Analog Signal xa(t)

xa = x1 + x2; %

Sum of the Amplitudes for plotting of xa(t)

A = A1 + A2; %

To plot x1(t), x2(t) and xa(t) signals

figure( 1 ) %Clear Command Window subplot( 3, 1, 1 ) plot( t, x1, 'r' ) axis( [ t0t1 -A1 A1 ] ) ylabel( 'x1(t)' ) title ( ' Continuous Time Signals ' ) subplot( 3, 1, 2 ) plot( t, x2, 'g' ) axis( [ t0 t1 -A2 A2 ] )


DSP Lab Manual ylabel( 'x2(t)' )

Page | 8

subplot( 3, 1, 3 ) plot( t, xa, 'b' ) axis( [ t0 t1 -A A] ) ylabel( 'x3(t)' ) % Choose the Sampling Frequency satisfying Nyquist s Sampling Theorem ’

Fs = input('Enter the Sampling Frequency in Hz = ' ); %

Calculate the Sampling Interval

Ts = 1/Fs; %

Calculate the number of Samples

n0 = t0/Ts n1 = t1/Ts; n = n0:n1; nT = n * Ts; xn1 = A1 * cos( 2 * pi * F1 * nT ); xn2 = A2 * cos( 2 * pi * F2 * nT ); xn = xn1 + xn2; figure(2) subplot( 3, 1, 1 ) stem( n, xn1, 'r' ) axis( [ t0/Ts, t1/Ts, -A1, A1 ] ) ylabel( 'x1[n]' ) title ( ' Discrete Time Signals ' ) subplot( 3, 1, 2 ) stem( n, xn2, 'g' ) axis( [ t0/Ts, t1/Ts, -A2, A2 ] ) ylabel( 'x2[n]' ) subplot( 3, 1, 3 ) stem( n, xn, 'b' ) axis( [ t0/Ts, t1/Ts, -A, A ] ) ylabel( 'x [n]' ) % Find the Maximum Frequency for reconstruction Fmax = max( [ F1 F2 ] ); Theta = zeros( length(nT), length(t) ); tn = nT' * ones( 1, length( t ) ); for row = 1 : length(nT) Theta( row, : ) = t( 1, : ) - tn( row, : ); end


DSP Lab Manual %

Generation of Interpolation Function

Interpolate = sinc( ( 2 * pi * Fmax ) .* Theta ); xr = xn * Interpolate; Page | 9

figure(3) subplot( 2, 1, 1 ) plot( t, xa, 'b' ) axis( [ t0, t1, -A, A ] ) xlabel( ' t--> ' ) ylabel(' xa(t) ' ) title ( ' Original Continous Time Signals ' ) subplot( 2, 1, 2 ) plot( t, xr, 'b' ) xlabel( ' t--> ' ) ylabel(' xr(t) ' ) title ( ' Reconstructed Continous Time Signals ' )

>> Amplitude of signal x1(t) = 1.2 Amplitude of signal x2(t) = 0.8 Frequency in Hz of signal x1(t) = 4 Frequency in Hz of signal x2(t) = 8 Beginning Time = -1 End Time = 1 Resolution =0.0001 Enter the Sampling Frequency in Hz = 64

Figure 1.2Screenshot of the output of Nyquist Sampling Theorem Verification program for perfect reconstruction.


DSP Lab Manual EXPT No: 2 a)

DATE:

LINEAR CONVOLUTION OF TWO SEQUENCES AIM: To execute MATLAB m-file program for finding the response of a given Linear Page | 10 Time Invariant system by linear convolution summation with: (i) (ii) (iii)

Causal sequence Non Causal sequence and To verify commutative, associative and distributive properties.

THEORY: For discrete time systems with input signal x(n) and finite duration unit impulse response h(n), linear convolution summation is defined by

()  ∑=  [][]  ∑= [][]

…

(3.1)

Eqn 3.1 represents the response of a causal system. For a non-causal system the response is given by the convolution summation defined by the eqn3.2  ∑= − [][]  ∑=−  [][] …

y(n) =

(3.2)

If Nx is the length of x[n] and Nh is the length of h[n] then the length of y[n] is N y given by … (3.3) N y = Nx + Nh - 1 Graphical Computation of Linear Convolution Summation:

h(k) resulting in h(-k) • Shifting of h(-k) resulting in h(n-k) • Reflection of

x(k) and h(n-k) • Summation of the product sequence x(k) h(n-k) results in the convolution Sample value for y(n) • Element wise multiplication of the sequences

Example: x(n) = {1, 2, 3, 1}Nx = 4 h(n) = {1, 1, 1} Nh = 3 N y =Nx + Nh - 1= 4 + 3 – 1 = 6


DSP Lab Manual

Page | 11


DSP Lab Manual

Page | 12

y[n] = { 1, 3, 6, 6, 4, 1 }

ALGORIHM: (i) Convolution - Causal

Begin Read x n and hn Compute yn = ∑ x k hn,k , ∀n length of x n, hn, yn display x n, hn, yn End


DSP Lab Manual

PROGRAM:(i)conv_causal.m

Page | 13

clc clear all close all % Enter the samples of Input Sequence x[n] and Unit Impulse Response h[n] x = input('Enter the Discrete Input Signal Sequence x[n] = '); h = input('Enter the Impulse Response Sequence h[n] = '); Nx = length( x ); Nh = length( h ); nx = 0 : 1 : Nx-1; nh = 0 : 1 : Nh-1; % Convolving Input x[n] and Impulse Response h[n] % y = conv( h, x ); Ny = length( y ); ny = 0 : 1 : Ny-1; disp(' The Convolution of x[n] and h[n] is y[n] = ' ) disp( y ) disp(' The number of samples in y[n] is Ny = ' ) disp(Ny ) figure( 1 ) subplot( 3, 1, 1 ) stem(nx, x, 'r', 'filled' ) xlabel(' n--> ' ) ylabel(' x[n] ' ) title(' Input Sequence x[n] ' ) subplot( 3, 1, 2 ) stem(nh, h, 'r', 'filled' ) xlabel(' n--> ' ) ylabel(' h[n] ' ) title(' Impulse Response h[n] ' ) subplot( 3, 1, 3 ) stem(ny, y, 'r', 'filled' ) xlabel(' n--> ' ) ylabel(' y[n] ' )


DSP Lab Manual title(' Output Sequence y[n] ' )

AFTER THE EXECUTION: Page | 14 >>

Enter the Discrete Input Signal Sequence x[n] = [ 1, 2, 3, 4 ] Enter the Impulse Response Sequence h[n] = [ 3 1 2 1 ] The Convolution of x[n] and h[n] is y[n] = 3 7 13 20 12 11 4 The number of samples in y[n] is Ny = 7 >>

Fig 3.1 Screenshot of output of the conv_causal.m program

ALGORIHM: (ii) Convolution - Non Causal Begin Read x n,n x, hn andnh Compute yn = ∑ x khn,k , ∀n Min( ny)= Min( n x)+Min( nh) Max( ny)= Max( n x)+Max( nh) ny = Min( ny):Max( ny) displayx n, hn, yn End


DSP Lab Manual PROGRAM: (ii)conv_noncuasal.m clc clear all close all Page | 15

% Enter the values of x[n] and h[n] x = input(' Enter the input signal x[n] = ' ); nx = input( ' Enter the index values of x[n] = ' ); h = input(' Enter the Unit Impulse Response h[n] = ' ); nh = input( ' Enter the index values of h[n] = ' ); % Convolution % y = conv( h, x ); % Index Sequence of y[n] % ny = ( min(nx)+min(nh) ) : 1 : ( max(nx)+max(nh) ); disp(' The Output Sequence = ' ) disp( y ) disp(ny ) figure(1) subplot( 3, 1, 1 ) stem ( nx, x, 'filled' ) axis( [ min(nx) max(nx) min(x) max(x) ] ) xlabel(' n--> ') ylabel(' x[n] ') title(' Input Signal x[n] ' ) subplot( 3, 1, 2 ) stem ( nh, h, 'filled' ) axis( [ min(nh) max(nh) min(h) max(h) ] ) xlabel(' n--> ') ylabel(' h[n] ') title(' Input Signal h[n] ' ) subplot( 3, 1, 3 ) stem ( ny, y, 'filled' ) axis( [ min(ny) max(ny) min(y) max(y) ] ) xlabel(' n--> ') ylabel(' y[n] ') title(' Input Signal y[n] ' ) AFTER THE EXECUTION: >> Enter the input signal x[n] = [ 3, 4, 5, 2, 1 ] Enter the index values of x[n] = -2:2 Enter the Unit Impulse Response h[n] = [ 1, 3, 2, 1] Enter the index values of h[n] = -1:2 The Output Sequence = 3 13 23 28 21 12 4 1 -3 -2 -1 0 1 2 3 4


DSP Lab Manual >>

Page | 16

Fig 3.2 Screenshot of output of the conv_noncausal.m program

Commutative, Associative & Distributive properties of Linear Convolution clc clear all close all x = [3, 1, 2, 1]; h1 = [ 1, 2, 2, 1 ]; h2 = [ 2, 4, 2, 4 ]; %Commutative Property y1 = conv(x, h1) y2 = conv(h1, x) if y1 == y2 disp('Commutative Property of convolution is proved'); else disp('error in the computation') end


DSP Lab Manual %

Associative Property

y3 = conv(x, conv(h1, h2)) y4 = conv(conv(x, h1), h2) Page | 17

if y3 == y4 disp('Associative Property of convolution is proved'); else disp('error in the computation') end % Distributive Property y5 = conv(x, h1+h2) y6 = conv(x, h1) + conv(x, h2) if y5 == y6 disp('Distributive Property of convolution is proved'); else disp('error in the computation') end

>> y1 = 3

7

10

10

7

4

1

y2 = 3 7

10

10

7

4

1

Commutative Property of convolution is proved y3 = 6

26

54

86 102

96

72

40

18

4

26

54

86 102

96

72

40

18

4

y4 = 6

Associative Property of convolution is proved y5 = 9

21 24

34 19

14

5

21 24

34 19

14

5

y6 = 9


DSP Lab Manual Distributive Property of convolution is proved EXPT No: 2 b)

DATE:

CIRCULAR CONVOLUTION OF TWO SEQUENCES Page | 18

AIM: To execute an m-file program to compute the circular convolution of two given finite-duration discrete-time sequences and to verify the commutative, associative and distributive properties. THEORY: Circular convolution of two finite duration discrete time sequences x[n] and h[n] is given by

()  ∑=  [][(  )]  ∑= [] [(  )]

…

(4.1)

IfNx is the length of x[n] and Nh is the length of h[n] then the length of y[n] is N y given by N y = Max(Nx,Nh ) … (4.2) The computation of y[n] is based on the equalization of the lengths of x[n] and h[n] by appending appropriate zeros with the following logic: if Nx>Nh hn = [ hn, 0, 0, 0, …, 0 ] (Nx -Nh) zeros else xn = [ xn, 0, 0, 0, …, 0] (Nh- Nx ) zeros end

Circular convolution performed in this manner has no practical significance. HoweverLinear convolution result can be obtained from circular convolution using the following steps: The length of y[n] is N y = Max(Nx,Nh).

   

Append ( Nh - 1 ) zeros to x[n]andthus resulting length is (Nx+ Nh - 1 ) Append ( Nx - 1 ) zeros to h[n] andthus resulting length is (Nx+ Nh - 1 ) Perform circular convolution y[n] = x[n] (*) h[n] Resulting y[n] is same as linear convolution between x[n] & h[n]

Linear convolution is the response of an LTI system. Example: Given x1(n) = {1, 1, 2, 1} and x2(n) = {1, 2, 3, 4}. Arrange x1(n) and x2(n) in circular fashion as shown below


DSP Lab Manual

Page | 19

To get x2(-m), rotate x2(m) by 4 samples in clockwise direction. x3(0) = x1(m) x2(-m) = x1(0) x2(0) + x1(1) x2(3) + x1(2) x2(2) + x1(3) x2(1) = 1 + 4 + 6 +2 = 13

x3(0) = 13

Keep x1(m) constant and rotate x2(-m) once to compute further values. To get x3(1) rotate x2(-m) by one sample in anti-clockwise direction x2(1-m)


DSP Lab Manual x3(1) = x1(m) x2(1-m) = x1(0) x2(1) + x1(1) x2(0) + x1(2) x2(3) + x1(3) x2(2) = 2 + 1 + 8 + 3 = 14 x3(1) = 14 Page | 20

To get x3(2) rotate x2(1-m) by one sample in anti-clockwise direction x2(2-m)

x3(2) = x1(m) x2(2-m) = x1(0) x2(2) + x1(1) x2(1) + x1(2) x2(0) + x1(3) x2(3) = 3 + 2 + 2+ 4 = 11 x3(2) = 11 To get x3(3) rotate x2(2-m) by one sample in anti-clockwise direction x2(3-m)

x3(3) = x1(m) x2(3-m) = x1(0) x2(3) + x1(1) x2(2) + x1(2) x2(1) + x1(3) x2(0) = 4 + 3 + 4 + 1 = 12 x3(3) = 12 The convoluted signal is, x3(n) = {13, 14, 11, 12} ALGORIHM: Circular Convolution Begin


DSP Lab Manual Read x n and hn Compute length of x n, hn if (Nx>Nh) hn= [hn, 0,0,…0] else x n= [ x n, 0,0,…0] yn = ∑ x k h((n,k)) , ∀n displayx n, hn, yn

Page | 21

End PROGRAM:circ_conv.m clc clear all close all x = input(' Enter the samples of x[n] = ' ); Nx = length( x ); h = input(' Enter the samples of h[n] = ' ); Nh = length( h ); if ( Nx>Nh ) h = [ h zeros( Nx - Nh ) ]; else x = [ x zeros( Nh - Nx ) ]; end N = max(Nx, Nh ); temp = [ h(1) fliplr( h( 2:N ) ) ]; y(1) = sum( x .* temp ); for k = 2 : N temp = circshift( temp', k-1 )'; y( k ) = sum( x .* temp ); end disp(' Circular convolution of x[n] and h[n] is y[n] ' ) disp( y )

X1 = fft(x); H1 = fft(h); Y1 = X1 .* H1; y1 = ifft(Y1); disp(' Circular convolution of x[n] and h[n] in frequency domain is y1[n] ' ) disp( y1 )

X2 = fft(x, Nx+Nh-1); H2 = fft(h, Nx+Nh-1);


DSP Lab Manual Y2 = X2.*H2; y2= ifft(Y2, Nx+Nh-1); disp(' Linear convolution of x[n] and h[n] in frequency domain is y2[n] ' ) disp( y2 ) Page | 22

figure( 1 ) subplot( 3, 1, 1 ) stem( 0 : Nx-1, x, 'r', 'filled' ) title(' Discrete Time Sequence x[n] ' ) xlabel(' n--> ' ) ylabel(' x[n] ' ) subplot( 3, 1, 2 ) stem( 0 : Nh-1, h, 'r', 'filled' ) title(' Discrete Time Sequence h[n] ' ) xlabel(' n--> ' ) ylabel(' h[n] ' ) subplot( 3, 1, 3 ) stem( 0 : N-1, y, 'r', 'filled' ) title(' Discrete Time Sequence y[n] ' ) xlabel(' n--> ' ) ylabel(' y[n] ' ) AFTER THE EXECUTION:

>> Enter the samples of x[n] = [ 1, 2, 3, 4] Enter the samples of h[n] = [ 1 ,2, 2, 1 ] Circular convolution of x[n] and h[n] is y[n] 17 15 13 15 Circular convolution of x[n] and h[n] in frequency domain is y1[n] 17 15 13 15 Linear convolution of x[n] and h[n] in frequency domain is y2[n] 1 4 9 15 16 11 4


DSP Lab Manual

Page | 23

Fig 4.1 Screenshot of output of the circ_conv.m program

Commutative, Associative & Distributive properties of Convolution clc clear all close all x = [3, 1, 2, 1]; h1 = [ 1, 2, 2, 1 ]; h2 = [ 2, 4, 2, 4 ]; %Commutative Property X1 = fft(x); H1 = fft(h1); Y1x = X1 .* H1; Y1h = H1 .* X1; y1x = ifft(Y1x) y1h = ifft(Y1h) if y1x == y1h disp('Commutative Property proved'); else disp('Error in compuatation'); end


DSP Lab Manual %Associative Property H2 = fft(h2); Yxh = X1 .* ( H1 .* H2 ); Yhx = (X1 .* H1) .* H2; Page | 24 yxh = ifft( Yxh ) yhx = ifft(Yhx) if yxh == yhx disp('Associative Property proved'); else disp('Error in compuatation'); end %Distributive Property Y11 = X1 .* ( H1 + H2); Y12 = X1 .* H1 + X1 .* H2; y11 = ifft(Y11) y12 = ifft(Y12) if y11 == y12 disp('Distributive Property proved'); else disp('Error in compuatation'); end >>y1x = 10

11

11

10

y1h = 10 11 11 10 Commutative Property proved yxh = 126 126 126 126 yhx = 126 126 126 126 Associative Property proved y11 = 28 35 29 34 y12 = 28

35

29

34

Distributive Property proved >>


DSP Lab Manual EXPT No: 3 a)

DATE:

AUTOCORRELATIONOF A SEQUENCE AND VERIFICATION OF ITS PROPERTIES

Page | 25

AIM: To execute an m-file program to find the autocorrelation of a given sequence and to verify its properties. THEORY:The autocorrelation of a sequence x[n] is given by,

 [  ]  ∑∞=−∞ [][   ]  ∑∞=−∞ [  ] []

(5.1)

…

Where the integer parameter l is called the lag indicator. Eqn. 5.1 indicates that the sequence x[n] is shifted by l samples with respect to itself. The operation defined by eqn. 5.1 is identical to convolution except for the time folding. It may be recalled that in convolution one of the sequence is time reversed and then the operations of shifting multiplications and additions are carried out. However in correlation except for the folding he rest of the three operations are carried out. PROPERTIES OF AUTO CORRELATION:

(i)

[  ]attainsmaximum value when the lag parameter l = 0.  []  ∑∞=−∞ [][ ]  ∑∞=−∞ |[]|   … denotes the energy of the signal x[n]. Where 

(ii)

The auto correlation results in a symmetric sequence … i. e.,  []   [  ]

(5.2)

(5.3)

Auto correlation finds application in digital communication. Autocorrelation is a mathematical tool for finding repeating patterns, such as the presence of a periodic signal which has been buried under noise, or identifying the missing fundamental frequency in a signal implied by its harmonic frequencies. It is used frequently in signal processing for analysing functions or series of values, such as time domain signals. Informally, it is the similarity between observations as a function of the time separation between them. Example: x[n] = { 1, 2, 3, 4 } rxx[ 0] => 1 2 3 4 rxx[ 1] => 1 2 3 4 1 2 3 4 1 2 3 4 _________________________ _______________________ 1 + 4 + 9 + 16 = 30 2 + 6 + 12 =20 ________________________ ______________________


DSP Lab Manual rxx[ 2] => 1

2

3

4

1

2

rxx[ 3] =>1 3

1

2

2

3

3

4

4

rxx[ -3] => 1

2

1 3

2

2

4

4

_________________________ 4 _________________________ ALGORIHM: Autocorrelation Begin Read x n Compute

∀

End PROGRAM:Auto_Corr.m clc clear all close all x = input(' Enter the samples of the sequence x[n] = ' ); %

1

2

3

4

3

3 + 8 = 11

length of x n r xx = ∑ x n x (n,l) , E x = ∑ x n, x n E = r xx (0) displayx n, r xx

4

4

________________________

20

3

3

_______________________

1

=

2

_______________________ 4

rxx[ -2] =>

__________________________ 2 + 6 + 12

4 1

___________________________ rxx[ -1] => 1

3

4

________________________ 3 + 8 = 11

Page | 26

2

Compute Auto-Correlation

rxx = xcorr( x );


DSP Lab Manual disp(' The Autocorrelation Sequence rxx[ l ] = ' ) disp(rxx ) N = length( x ); nx = N - 1; Page | 27 l = -nx : 1 : nx; L = length( l ); lx = L - 1; figure( 1 ) subplot ( 2, 1, 1 ) stem ( 0 : nx , x , 'r', 'filled' ) xlabel(' n-- > ' ) ylabel(' x[n] ' ) title ( ' Input D.T. Sequeence x[n] ' ); subplot( 2, 1, 2 ) stem( 0 : lx, rxx, 'r', 'filled' ) xlabel(' l --> ' ) ylabel(' rxx[l] ' ) title(' Autocorrelation Sequence ' ) %

Autocorrelation using Convolution Function

rxx_c = conv( x, fliplr( x ) ); disp(' The Autocorrelation computed using Convolution ' ) disp(rxx_c ) disp(' The Autocorrelation computed using Convolution ' ) disp(rxx_c ) % Autocorrelation with Zero Lag Quantifies Energy in x[n] Ex = sum( x .^ 2 ); disp(' The Energy of the Signal ' ) disp( Ex ); disp(' The Energy of the Signal from autocorrelation' ) disp( rxx( round( L/2 ) ) ) AFTER THE EXECUTION: >> Enter the samples of the sequence x[n] = [ 1, 3, 2, 4 ] The Autocorrelation Sequence rxx[ l ] = 4.0000 14.0000 17.0000 30.0000 17.0000 14.0000 4.0000

The Autocorrelation computed using Convolution 4 14 17 30 17 14 4 The Energy of the Signal 30 The Energy of the Signal from autocorrelation 30 >>


DSP Lab Manual

Page | 28

Fig 5.1 Screenshot of output of the Auto_Corr.m program

EXPT No: 3 b)

DATE:

CROSS CORRELATION OF TWO GIVEN SEQUENCES AND VERIFICATION OF ITS PROPERTIES AIM: To execute an m-file program to find the cross-correlation of two given sequences and to verify its properties. THEORY:Correlation is a mathematical operation that closely resemble convolution. The objective in computing the correlation between two different signals is to measure the degree to which the signals are similar and thus to extract some useful information. The cross-correlation of two real-valued, finite energy, sequences x[n] and y[n] is given by, …  [  ]  ∑∞=−∞ [][   ] [  ]  ∑∞=−∞ [] [   ]… (6.2)   ,±,±,…

(6.1)

Where the integer parameter l is called the lag indicator. Properties of cross-correlation:

) [  ] ≤ √ 

…

(6.3)


DSP Lab Manual )  []  [  ]

(6.4) Where Ex and E y are energies of the signals x[n] and y[n] respectively. The application of correlation finds in RADAR, SONAR and so on.

Page | 29

…

EXAMPLE: Compute the cross correlation of the two sequences:x[n] = { 1, 3, 2, 4 } and y[n] = { 1, 2, 4, 3 }. Solution:

ALGORIHM: Cross-correlation Begin Read x n, yn Compute length of x n , yn r xy = ∑ x n y (n,l) , ∀ ryx = ∑ yn x (n,l) , ∀ displayx n,yn,r xy, ryx End PROGRAM:Cross_Corr.m clc clear all close all x = input(' Enter the samples of the sequence x[n] = ' ); y = input(' Enter the samples of the sequence y[n] = ' ); rxy = xcorr( x, y ); ryx = xcorr( y, x );


DSP Lab Manual % Verify the failure of commutative property with rxy and r yx % disp(' The Cross-Correlation between x[n] and y[n] ' ) disp(rxy ) Page | 30

disp(' The Cross-Correlation between y[n] and x[n] ' ) disp(ryx ) figure( 1 ) subplot( 2, 1, 1 ) stem ( x, 'filled') xlabel(' n--> ' ) ylabel(' x[ n ] ' ) title(' First Sequence x[n] ' )

subplot( 2, 1, 2 ) stem ( y, 'filled') xlabel(' n--> ' ) ylabel(' y[ n ] ' ) title(' Second Sequence y[n] ' ) figure(2 ) subplot( 2, 1, 1 ) stem ( rxy, 'filled') xlabel(' l--> ' ) ylabel(' rxy[ l ] ' ) title(' The Cross-Correlation Sequence rxy[ l ] ' )

subplot( 2, 1, 2 ) stem ( ryx, 'filled' ) xlabel(' l--> ' ) ylabel(' ryx[ l ] ' ) title(' The Cross-Correlation Sequence ryx[ l ] ' )

AFTER THE EXECUTION:

>>

Enter the samples of the sequence x[n] = [ 1, 3, 2, 4 ] Enter the samples of the sequence y[n] = [ 1, 2, 4, 3 ] The Cross-Correlation between x[n] and y[n] 3.0000 13.0000 20.0000 27.0000 23.0000 10.0000 4.0000 The Cross-Correlation between y[n] and x[n] 4.0000 10.0000 23.0000 27.0000 20.0000 13.0000 3.0000 >>


DSP Lab Manual

Page | 31

Fig 6.1 Screenshot of the output of Cross_Corr.m program


DSP Lab Manual EXPT No: 4

DATE:

LINEAR DIFFERENCE EQUATION Page | 32

AIM: To execute an m-file program to find the solution of a given difference equation(i) without initial conditions and (ii) with initial conditions. THEORY:Thegeneral form linear constant coefficient difference equation is given − (7.1) []  ∑− … by: = [   ]  ∑=  [   ]

− []  ∑− =  [ ] 

Applying z-transform to eqn. (7.1)

− ∑− =  []

…

(7.2)

Rearranging eqn(7.2)

 ∑ ()     []  ()  +∑  

(7.3)

H(z) is called transfer function of the system. Linear Constant Coefficient Difference Equations(LCCDEs) model LTI systems. Equation (7.1) governs a recursive or Infinite Impulse Response(IIR) filter. If  ≡  ∀ except for  ,&  =1, then − − (7.4) = 

[]  ∑



…

Equation(7.4) is identical to convolution summation with { bk } being the filter coefficients of a Finite Impulse Response (FIR) Filter. Equation 7.1 models a recursive IIR system. EXAMPLE:

ALGORIHM: LINEAR DIFFERENCE EQUATION WITHOUT INITIAL CONDITIONS(INITIALLY RELAXED SYSTEM) Begin Read ak, bk Compute

,, and h fromrational H(z) using partial fraction expansion.     − ⋯    −  ⋯    − n

hn= IZT { H z } define input signal to be unit step input x n = un solv         displayr,p,k, h, yn End



∑



∑   




DSP Lab Manual

PROGRAM: LCCDE1.m % With ZERO initial conditions OR system initially relaxed. Page | 33

clc clear all close all b = input(' Enter the Numerator Coefficients = ' ); a = input(' Enter the Denomenator Coefficients = ' ); [ r, p, k ] = residuez( b, a ); disp(' Residues = ' ) disp( r ) disp(' Poles = ' ) disp( p ) disp(' Constants = ' ) disp( k ) disp(' ' ) [ h, n ] = impz( b, a ); % Impulse Response of the System disp(' Impulse Response of the System h[n] = ' ) disp( h' ) disp(' ' ) figure( 1 ) subplot( 2, 1, 1 ) stem( n, h ) xlabel(' n--> ' ) ylabel(' h[n] ' ) title(' Impulse Response of the System ' ) % Exciting the System with Unit Step Signal x = ones( 1, 12 ); y = filter( b, a, x ); disp(' The Response of the System y[n] for Step input =' ) disp( y ) disp(' ' ) % Alternatively disp(' Alternatively the Response of the System y[n] for Step input =' ) y1 = filter( h, 1, x ); disp( y1 ) subplot( 2, 1, 2 ) stem( 0:length( y ) - 1, y ) xlabel(' n--> ' ) ylabel(' y[n] ' ) title(' Response of the System ' )


DSP Lab Manual AFTER THE EXECUTION:

>> Page | 34

Enter the Numerator Coefficients = [ 1, 1/4, 1/8 ] Enter the Denomenator Coefficients = [ 1, 1/3, 1/6 ] Residues = 0.1250 + 0.0559i 0.1250 - 0.0559i Poles = -0.1667 + 0.3727i -0.1667 - 0.3727i Constants = 0.7500 Impulse Response of the System h[n] = Columns 1 through 7 1.0000 -0.0833 -0.0139 0.0185 -0.0039 -0.0018 0.0012 Columns 8 through 11 -0.0001 -0.0002 0.0001 0.0000 The Response of the System y[n] for Step input = Columns 1 through 7 1.0000 0.9167 0.9028 0.9213 0.9174 0.9156 0.9169 Columns 8 through 12 0.9168 0.9166 0.9167 0.9167 0.9167 Alternatively the Response of the System y[n] for Step input = Columns 1 through 7 1.0000 0.9167 0.9028 0.9213 0.9174 0.9156 0.9169 Columns 8 through 12 0.9168 0.9166 0.9167 0.9167 0.9167


DSP Lab Manual

Page | 35

Fig 7.1 Screenshot of output of the LCCDE1.m program

ALGORIHM: CONDITIONS

LINEAR

DIFFERENCE

EQUATIONWITH

Begin Read ak, bk Compute hn= IZT { H z } define input signal to be unit step input x n = un  solv            for

e  

∑

∑

display h, yn End PROGRAM: LCCDE2.m clc clear all close all b = input(' Enter the Numerator Constant Coefficients = ' ); a = input(' Enter the Denomenator Constant Coefficients = ' ); ic = input( ' Enter the initial conditions y[-1], y[-2] : ' );


INITIAL

DSP Lab Manual % Impulse Response [ h n ] = impz( b, a );

Page | 36

disp(' Impulse Response ' ) disp( h' ) %

Unit Step Signal as Input

x = ones( 1, 12 ); % Zero Input Response zi = filtic( b, a, ic ); % Total Response = Zero Input Response + Zero State Response y = filter( b, a, x, zi ); disp(' Response of the System for Step Signal ' ) disp( y ) figure( 1 ) subplot( 2, 1, 1 ) stem( n, h ) xlabel(' n--> ' ) ylabel(' h[n] ' ) title(' Impulse Responnse ' ) N = length( y ); ny = 0:N - 1; subplot( 2, 1, 2 ) stem(ny, y ) xlabel(' n--> ' ) ylabel(' y[n] ' ) title(' Response of the System ' ) AFTER THE EXECUTION:

>> Enter the Numerator Constant Coefficients = [ 1, 1/3, 1/6 ] Enter the Denomenator Constant Coefficients = [ 1, 1/4, 1/8 ] Enter the initial conditions y[-1], y[-2] : [ 1, 0 ] Impulse Response Columns 1 through 7

1.0000 0.0833 0.0208 -0.0156 0.0013 0.0016 -0.0006 Columns 8 through 9 -0.0001 0.0001


Dsp Lab -15ecl57 Part_i

Recommend Documents