Digital Signal Processing_S. Salivahanan, A. Vallavaraj and C. Gnanapriya

Scilab Textbook Companion for Digital Signal Processing by S. Salivahanan, A. Vallavaraj And C. Gnanapriya1 Created by Priya Sahani B TECH EXTC Electrical Engineering V.J.T.I College Teacher Rizwn Ahmed Cross-Checked by Lavitha Pareira July 17, 2017

1

Funded by a grant from the National Mission on Education through ICT, http://spoken-tutorial.org/NMEICT-Intro. This Textbook Companion and Scilab codes written in it can be downloaded from the ”Textbook Companion Project” section at the website http://scilab.in

Book Description Title: Digital Signal Processing Author: S. Salivahanan, A. Vallavaraj And C. Gnanapriya Publisher: Tata McGraw - Hill, New Delhi Edition: 1 Year: 2008 ISBN: 978-0-07-463996-2

1

Book Description Title: Digital Signal Processing Author: S. Salivahanan, A. Vallavaraj And C. Gnanapriya Publisher: Tata McGraw - Hill, New Delhi Edition: 1 Year: 2008 ISBN: 978-0-07-463996-2

1

Scilab numbering policy used in this document and the relation to the above book. Exa Example (Solved example) Eqn Equation (Particular equation of the above book) AP Appendix to Example(Scilab Code that is an Appednix to a particular

Example of the above book) For example, Exa 3.51 means solved example 3.51 of this book. Sec 2.3 means a scilab code whose theory is explained in Section 2.3 of the book.

2

Contents List of Scilab Codes

4

1 Classifications of signals and systems

5

2 Fourier Analysis of Preiodic and Aperiodic Continuous Time Signals and Systems

7

3 Applications of Laplace Transform to System Analysis

16

4 Z Transforms

18

5 Linear Time Invariant Systems

22

6 Discrete and Fast Fourier Transforms

24

7 Finite Impulse Response Filters

41

8 Infinite Impulse Response Filters

43

9 Realisation of Digital Linear Systems

51

10 Effects of Finite Word Length in Digital Filters

53

11 Multirate Digital Signal Processing

57

3

12 Spectral Estimation

63

4

List of Scilab Codes Exa 1.2.a Exa 1.2.b Exa 1.2.c Exa 1.2.d Exa 2.1 Exa 2.2 Exa 2.3 Exa 2.4 Exa 2.5 Exa 2.6 Exa 2.8 Exa 3.10 Exa 3.11 Exa 3.12 Exa 4.2 Exa 4.4 Exa 4.13 Exa 4.14 Exa 4.16 Exa 4.19 Exa 5.20 Exa 5.21 Exa 6.1 Exa 6.2 Exa 6.3 Exa 6.4 Exa 6.5

Rectangular wave . . . . . . . . . . . . . . . Rectangular wave . . . . . . . . . . . . . . . Cosine wave . . . . . . . . . . . . . . . . . . Ramp wave . . . . . . . . . . . . . . . . . . Fourier Series of Periodic Square Wave . . . Fourier Series of Periodic Rectangular Wave Fourier Series of Periodic Half Wave Rectified Sine Wave . . . . . . . . . . . . . . . . . . . Fourier Series of Periodic Triangular Wave . Fourier Series of Periodic Rectangular Pulse Fourier Series of Square Wave . . . . . . . . Complex fourier series representation . . . . Poles and Zeros . . . . . . . . . . . . . . . . Poles and zeros . . . . . . . . . . . . . . . . Poles and zeros . . . . . . . . . . . . . . . . Z transform . . . . . . . . . . . . . . . . . . Z transform . . . . . . . . . . . . . . . . . . Convolution . . . . . . . . . . . . . . . . . . Convolution . . . . . . . . . . . . . . . . . . Cross correlation . . . . . . . . . . . . . . . System response . . . . . . . . . . . . . . . System response . . . . . . . . . . . . . . . Poles and zeros . . . . . . . . . . . . . . . . Linear and Circular convolution . . . . . . . FIR filter response . . . . . . . . . . . . . . Convolution . . . . . . . . . . . . . . . . . . Convolution . . . . . . . . . . . . . . . . . . Convolution . . . . . . . . . . . . . . . . . . 5

5 5 6 6 7 8 9 10 11 12 14 16 16 17 18 19 19 20 20 21 22 22 24 25 26 26 27

Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa Exa

6.6 6.8 6.9 6.10 6.11 6.12 6.13 6.14 6.15 6.16 6.17 6.18 6.19 6.20 6.21 6.22 6.23 6.24 6.25 6.26 6.27 6.34 6.35 6.36 6.37 7.3 7.4

Exa 8.1 Exa 8.2 Exa 8.3 Exa 8.4 Exa 8.5 Exa 8.6

Convolution . . . . . . . . . . . . . . . . . . DFT . . . . . . . . . . . . . . . . . . . . . . DFT . . . . . . . . . . . . . . . . . . . . . . DFT . . . . . . . . . . . . . . . . . . . . . . DFT . . . . . . . . . . . . . . . . . . . . . . DFT . . . . . . . . . . . . . . . . . . . . . . Inverse DFT . . . . . . . . . . . . . . . . . . Inverse DFT . . . . . . . . . . . . . . . . . . DIT FFT . . . . . . . . . . . . . . . . . . . DIT FFT . . . . . . . . . . . . . . . . . . . DIT FFT . . . . . . . . . . . . . . . . . . . DIT FFT . . . . . . . . . . . . . . . . . . . DIF FFT . . . . . . . . . . . . . . . . . . . DIF FFT . . . . . . . . . . . . . . . . . . . DIF FFT . . . . . . . . . . . . . . . . . . . DIF FFT . . . . . . . . . . . . . . . . . . . IFFT . . . . . . . . . . . . . . . . . . . . . . IFFT . . . . . . . . . . . . . . . . . . . . . . IFFT . . . . . . . . . . . . . . . . . . . . . . IFFT . . . . . . . . . . . . . . . . . . . . . . IFFT . . . . . . . . . . . . . . . . . . . . . . Overlap Add Convolution . . . . . . . . . . Overlap Save Convolution . . . . . . . . . . Cross Correlation . . . . . . . . . . . . . . . Circular Correlation . . . . . . . . . . . . . Low pass filter using fourier series method . Low pass filter using Type 1 frequency sampling technique . . . . . . . . . . . . . . . . IIR filter Design byBackward Difference For Derivative method . . . . . . . . . . . . . . IIR filter Design byBackward Difference For Derivative method . . . . . . . . . . . . . . IIR filter Design byBackward Difference For Derivative method . . . . . . . . . . . . . . IIR filter Design by Impulse Invariant method IIR filter Design by Impulse Invariant method IIR filter Design by Impulse Invariant method

6

28 29 29 30 30 31 31 32 32 32 33 33 33 34 34 34 35 35 36 36 36 37 38 38 39 41 42 43 43 44 44 45 45

Exa 8.7 Exa 8.8 Exa 8.9 Exa 8.10 Exa 8.11 Exa 8.12 Exa 8.14 Exa 8.15 Exa 9.4 Exa 9.5.a Exa 9.5.b Exa 10.2 Exa 10.3 Exa 10.4 Exa 10.5 Exa 11.1 Exa 11.2 Exa 11.4 Exa 11.5 Exa 11.6 Exa 12.2 Exa 12.4

IIR filter Design by Bilinear Transformation method . . . . . . . . . . . . . . . . . . . . IIR filter Design by Bilinear Transformation method . . . . . . . . . . . . . . . . . . . . IIR filter Design by Bilinear Transformation method . . . . . . . . . . . . . . . . . . . . IIR filter Design by Bilinear Transformation method . . . . . . . . . . . . . . . . . . . . Butterworth Filter using Impulse Invariant transformation . . . . . . . . . . . . . . . . Butterworth Filter using Bilinear transformation . . . . . . . . . . . . . . . . . . . . . . Filter transformation . . . . . . . . . . . . . Filter transformation . . . . . . . . . . . . . Cascade Realisation . . . . . . . . . . . . . Parallel Realisation . . . . . . . . . . . . . . Parallel Realisation . . . . . . . . . . . . . . Output Quantisation Noise . . . . . . . . . Deadband Interval . . . . . . . . . . . . . . Deadband Interval . . . . . . . . . . . . . . Output Quantisation Noise for Cascade realisation . . . . . . . . . . . . . . . . . . . . . Time Decimation . . . . . . . . . . . . . . . Interpolation . . . . . . . . . . . . . . . . . Polyphase Decomposition . . . . . . . . . . Decimator implementation . . . . . . . . . . Decimator implementation . . . . . . . . . . Power Spectrum . . . . . . . . . . . . . . . Frequency resolution . . . . . . . . . . . . .

7

46 46 46 47 47 48 50 50 51 51 52 53 54 54 55 57 57 58 59 60 63 64

Chapter 1 Classifications of signals and systems

Scilab code Exa 1.2.a Rectangular wave

1 2 3 4 5 6

/ / Example 1 . 2 ( a )

clc ; clear ; t=-5:0.01:5; x = 1 * ( abs (2*t +3) <0.5); plot ( t , x ) ; t i t l e ( ’ x ( t )=re c t (2 t +3) ’ ) ;

Scilab code Exa 1.2.b Rectangular wave

1 2 3 4 5 6

/ / Exam ple 1 . 2 ( b )

clc ; clear ; t=-5:0.01:5; x = 2 * ( abs (t -1/4) <0.5); plot ( t , x ) ; t i t l e ( ’ x ( t )=2 ∗ r e c t ( t − 1/4) ’ ) ;

8

Scilab code Exa 1.2.c Cosine wave

1 2 3 4 5 6 7

/ / E xa mp le 1 . 2 ( c )

clc ; clear ; pi=22/7; t=-5:0.01:5; x = cos ( 2 * p i * t - 5 0 * p i ) ; plot ( t , x ) ; t i t l e ( ’ x ( t )=co s (2 ∗ p i ∗ t −∗ p i ) ’ ) ;

Scilab code Exa 1.2.d Ramp wave

1 2 3 4 5 6 7 8

/ / Exam ple 1 . 2 ( d )

clc ; clear ; t=-5:0.01:5; x=-0.5*(t-4); plot ( t , x ) ; t i t l e ( ’ x ( t )=r ( − 0.5 t+2) ’ ) ; z o o m _r e ct ( [ - 5 0 5 5 ]) ;

9

Chapter 2 Fourier Analysis of Preiodic and Aperiodic Continuous Time Signals and Systems

Scilab code Exa 2.1 Fourier Series of Periodic Square Wave

1 2 3 4 5 6 7

/ / E xa mp le 2 . 1 clc ; clear ; close ; A=1;T=2; w0=2*%pi/T;

/ / C a l c u l a t i o n o f t r i g n o m e t r i c f o u r i e r s e r i e s co − efficients 8 a 0 = A / T * ( i n t e g r a t e( ’ −1 ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’ +1 ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ −1 ’ , ’ t ’ , T / 4 , T / 2 ) ) ; 9 for n = 1 : 1 0 ; 10 a ( 1 , n ) = 2 * A / T * ( i n t e g r a t e( ’ − c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’+co s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ − c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , T / 4 , T / 2 ) ) ; 11 b ( 1 , n ) = 2 * A / T * ( i n t e g r a t e( ’ − s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’ +s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ − s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , T / 4 , T / 2 ) ) ;

10

12 end 13 14 / / D i s p l a y i n g f o u r i e r c o e f f i c i e n t s 15 disp (T , ’ f u n da m e nt a l p e r i o d T= ’ ,A , ’ Assu mpti on : A m p l i t u d e A= ’ ) ; 16 disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; 17 disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ; 18 19 x = [ - A * ones ( 1 , 25 ) A *ones ( 1 , 50 ) - A *ones ( 1 , 2 5 ) ] //

F un ct io n f o r p l o t i n g p ur po se 20 21 22 23 24 25 26

t=-T/2:0.01*T:T/2-0.01; subplot ( 3 1 1 ) ; plot ( t , x ) ; t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; subplot ( 3 1 2 ) ; plot2d3 ( a ) ; t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; subplot ( 3 1 3 ) ; plot2d3 ( b ) ; t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.2 Fourier Series of Periodic Rectangular Wave

1 2 3 4 5 6 7


// C a lc u la t io n o f t r ig n o me t ri c f o u r i e r efficients 8 a 0 = A / T * i n t e g r a t e( ’ 1 ’ , ’ t ’ , - T / 4 , T / 4 ) ; 9 10 11 12 13 14

s e r i e s co −

for n = 1 : 1 0 ; a ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ cos (n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) ; b ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) ; end

// Di sp la yi ng f o u r i e r c o e f f i c i e n t s 11

15 disp (T , ’ f u n da m e nt a l p e r i o d T= ’ ,A , ’ Assu mpti on : A m p l i t u d e A= ’ ) ; 16 disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; 17 disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ; 18 19 x =[ zeros ( 1 , 25 ) A *ones (1,50) zeros ( 1 , 2 5 ) ] ; 20 t = - T / 2 : 0 . 0 1 * T : T / 2 - 0 . 0 1 ; 21 subplot ( 3 1 1 ) ; plot ( t , x ) ; 22 t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; 23 subplot ( 3 1 2 ) ; plot2d3 ( a ) ; 24 t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; 25 subplot ( 3 1 3 ) ; plot2d3 ( b ) ; 26 t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.3 Fourier Series of Periodic Half Wave Rectified Sine Wave

1 2 3 4 5 6 7


/ / C a l c u l a t i o n o f t r i g n o m e t r i c f o u r i e r s e r i e s co − efficients 8 a 0 = A / T * i n t e g r a t e( ’ si n (w0 ∗ t ) ’ , ’ t ’ , 0 , T / 2 ) ; 9 for n = 1 : 1 0 ; 10 a ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ si n (w0 ∗ t ) ∗ c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ ,0,T/2); 11 b ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ si n (w0 ∗ t ) ∗ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ ,0,T/2); 12 end 13 14 / / D i s p l a y i n g f o u r i e r c o e f f i c i e n t s 15 disp (T , ’ f u n da m e nt a l p e r i o d T= ’ ,A , ’ Assu mpti on : A m p l i t u d e A= ’ ) ;

12

16 17 18 19 20 21 22 23 24 25 26 27

disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ;

t=0:0.01*T:T/2; x = [ A * sin ( w 0 * t ) zeros ( 1 , 5 0 ) ] ; t=0:0.01*T:T; subplot ( 3 1 1 ) ; plot ( t , x ) ; t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; subplot ( 3 1 2 ) ; plot2d3 ( a ) ; t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; subplot ( 3 1 3 ) ; plot2d3 ( b ) ; t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.4 Fourier Series of Periodic Triangular Wave

1 2 3 4 5 6 7


/ / C a l c u l a t i o n o f t r i g n o m e t r i c f o u r i e r s e r i e s co − efficients 8 a 0 = 4 * A / T * ( i n t e g r a t e( ’ t − 0 .5 ∗ T ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’ t ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ − t +0.5 ∗ T ’ , ’ t ’ ,T/4,T/2)); 9 for n = 1 : 1 0 ; 10 a ( 1 , n ) = 2 * 4 * A / T * ( i n t e g r a t e( ’ ( t − 0. 5 ∗ T) ∗ c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’ t ∗ c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ ( − t +0. 5 ∗ T) ∗ c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , T / 4 , T /2)); 11 b ( 1 , n ) = 2 * 4 * A / T * ( i n t e g r a t e( ’ ( t − 0. 5 ∗ T) ∗ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , - T / 4 ) + i n t e g r a t e( ’ t ∗ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 4 , T / 4 ) + i n t e g r a t e( ’ ( − t +0. 5 ∗ T) ∗ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , T / 4 , T /2));

13

12 end 13 14 / / D i s p l a y i n g f o u r i e r c o e f f i c i e n t s 15 disp (T , ’ f u n da m e nt a l p e r i o d T= ’ ,A , ’ Assu mpti on : A m p l i t u d e A= ’ ) ; 16 disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; 17 disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ; 18 19 t = - T / 2 : 0 . 0 1 * T : T / 2 ; 20 x = [ - 4 * A / T *t ( 1 : 2 5) - 2* A 4 * A / T * t ( 2 6: 7 5 ) - 4* A / T * t (76:101)+2*A]; 21 subplot ( 3 1 1 ) ; plot ( t , x ) ; 22 t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; 23 subplot ( 3 1 2 ) ; plot2d3 ( a ) ; 24 t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; 25 subplot ( 3 1 3 ) ; plot2d3 ( b ) ; 26 t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.5 Fourier Series of Periodic Rectangular Pulse

1 2 3 4 5 6 7

/ / E xa mp le 2 . 5 clc ; clear ; close ; A=1;T=2;d=0.1; w0=2*%pi/T;

// C a lc u la t io n o f t r ig n o me t ri c f o u r i e r efficients 8 a 0 = A / T * i n t e g r a t e( ’ 1 ’ , ’ t ’ , - T / 4 , T / 4 ) ; 9 10 11 12 13 14

s e r i e s co −

for n = 1 : 1 0 ; a ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ cos (n ∗ w0 ∗ t ) ’ , ’ t ’ , - d / 2 , d / 2 ) ; b ( 1 , n ) = 2 * A / T * i n t e g r a t e( ’ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - d / 2 , d / 2 ) ; end

// Di sp la yi ng f o u r i e r c o e f f i c i e n t s 14

15 disp (d , ’ p u l s e w id th d= ’ ,T , ’ f u n da m en t a l p e r i o d T= ’ , A , ’ A s s u m p t i o n : A m p l i t u d e A= ’ ) ; 16 disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; 17 disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ; 18 19 n = round ( 5 0 * d / T ) ; //

V ar ia bl e used f o r p l o t t i ng p u l se s a c cu r a te l y 20 21 22 23 24 25 26 27

x =[ zeros ( 1 ,5 0 - n ) A *ones ( 1 , 2 * n + 1 ) zeros ( 1 , 5 0 - n ) ] t=-T/2:0.01*T:T/2; subplot ( 3 1 1 ) ; plot ( t , x ) ; t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; subplot ( 3 1 2 ) ; plot2d3 ( a ) ; t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; subplot ( 3 1 3 ) ; plot2d3 ( b ) ; t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.6 Fourier Series of Square Wave

1 2 3 4 5 6 7


/ / C a l c u l a t i o n o f t r i g n o m e t r i c f o u r i e r s e r i e s co − efficients 8 a 0 = A / T * ( i n t e g r a t e( ’ −1 ’ , ’ t ’ , - T / 2 , 0 ) + i n t e g r a t e( ’ +1 ’ , ’ t ’ , 0 , T / 2 ) ) ; 9 for n = 1 : 1 0 10 a ( 1 , n ) = 2 * A / T * ( i n t e g r a t e( ’ − c o s ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , 0 ) + i n t e g r a t e( ’+co s (n ∗ w0 ∗ t ) ’ , ’ t ’ , 0 , T / 2 ) ) ; 11 b ( 1 , n ) = 2 * A / T * ( i n t e g r a t e( ’ − s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , - T / 2 , 0 ) + i n t e g r a t e( ’+ s i n ( n ∗ w0 ∗ t ) ’ , ’ t ’ , 0 , T / 2 ) ) ; 12 end 13 a = clean ( a ) ; b = clean ( b ) ; / / F u n ct i o n u s ed t o

15

ro u nd s m a l l e n t i t i e s t o z e r o 14 15 / / C a l c u l a t i o n o f e x p o n e n ti a l f o u r i e r

s e r i e s co −

efficients 16 function y = f ( t ) , y = c o m p l e x (cos ( n * w 0 * t ) , - sin ( n * w 0 * t ) ) , e n d f u n c t i o n; 17 for n = - 1 0 : 1 0 18 c ( 1 , n + 1 1 ) = A / T * ( - 1 * intc ( - T / 2 , 0 , f ) + intc ( 0 , T / 2 , f ) ) ; 19 end 20 c = clean ( c ) ; / / F u n ct i o n u s ed t o

ro un d s ma ll e n t i t i e s t o z er o 21 22 / / C a l c u l a t i o n o f t r i g n o m e t r i c

f o u r i e r s e r i e s co − e f f i c i e n t s f ro m e x p o n e n t i a l f o u r i e s e r i e s coefficients

23 24 25 26 27 28

a01=c(1); a 1 = 2 * real ( c ( 1 2 : 2 1 ) ) ; b 1 = - 2 * imag ( c ( 1 2 : 2 1 ) ) ;

29 30 31 32

33 34 35 36 37 38 39 40 41

// Di sp la yi ng f o u r i e r c o e f f i c i e n t s disp (T , ’ f u n da m e nt a l p e r i o d T= ’ ,A , ’ Assu mpti on : A m p l i t u d e A= ’ ) ; disp ( ’ T ig no me tr ic f o u r i e r s e r i e s co − e f f i c i e n t s : ’ ) ; disp ( a 0 , ’ a 0= ’ ) ; disp (a , ’ a n= ’ ) ; disp (b , ’ b n= ’ ) ; disp ( ’ E xp on en ti al f o u r i e r s e r i e s co − e f f i c i e n t s ’ ) ; disp ( c ( 1 1 ) , ’ c 0= ’ ) ; disp ( c ( 1 2 : 2 1 ) , ’ c n= ’ ) ; disp ( c ( 1 0 : - 1 : 1 ) , ’ c −n= ’ ) ; disp ( ’ T ri g n o m et ri c f o u r i e r s e r i e s co − e f f i c i e n t s fr om ex po n en ti al c o e f f i c i e n t s : ’ ) ; disp ( a 0 1 , ’ a 0= ’ ) ; disp ( a 1 , ’ a n= ’ ) ; disp ( b 1 , ’ b n= ’ ) ; disp ( ’ T he c o − e f f i f c i e n t s o b t a i n e d a r e s am e b y b o t h methods ’ ) x = [ - A * ones ( 1 , 50 ) A *ones ( 1 , 5 1 ) ] ; t=-T/2:0.01*T:T/2; n=-10:10; subplot ( 3 1 1 ) ; plot ( t , x ) ; t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ;

16

42 43 44 45 46 47 48 49 50 51 52

subplot ( 3 1 2 ) ; plot2d3 ( a ) ; t i t l e ( ’ C o e f f i c i e n t s an ’ ) ; x l a b e l ( ’ n ’ ) ; subplot ( 3 1 3 ) ; plot2d3 ( b ) ; t i t l e ( ’ C o e f f i c i e n t s bn ’ ) ; x l a b e l ( ’ n ’ ) ; figure ; subplot ( 3 1 1 ) ; plot ( t , x ) ; t i t l e ( ’ x ( t ) ’ ) ; x l a b e l ( ’ t i m e t ’ ) ; subplot ( 3 1 2 ) ; plot2d3 (n , abs ( c ) ) ; t i t l e ( ’ Ma gn itu de o f C o e f f i c i e n t s | c | ’ ) ; x l a b e l ( ’ n ’ ) ; subplot ( 3 1 3 ) ; plot2d3 (n , atan ( c ) ) ; t i t l e ( ’ Phas e o f C o e f f i c i e n t s / c ’ ) ; x l a b e l ( ’ n ’ ) ;

Scilab code Exa 2.8 Complex fourier series representation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

/ / E xa mp le 2 . 8 clc ; clear ; close ; t = poly (0 , ’ t ’ ) ;

//cn=3/(4+(n ∗ %pi ) ˆ2 ) // Total energy Pt=0.669; Preq=0.999*Pt; // Required energy c0=3/(4+(0*%pi)^2); disp ( c 0 , ’ c0= ’ ) ; P =( abs ( c 0 ) ) ^ 2 ; c=[];n=0; while P < P r e q n=n+1; c(n)=3/(4+(n*%pi)^2); disp ( c ( n ) , ’ c n= ’ ,n , ’ n= ’ ) ; P = P + 2 * ( abs ( c ( n ) ) ) ^ 2 ; end disp ( P t , ’ T o t a l p o we r P t= ’ ) ; disp ( P r e q , ’ 9 9 . 9% o f t o t a l po we r P reqd= ’ ) ; disp (n , ’ To i c l u d e 9 9 .9% o f e ne rg y , we n eed t o r e t a i n n t e rm s w he re n= ’ ) ;

17

18

Chapter 3 Applications of Laplace Transform to System Analysis

Scilab code Exa 3.10 Poles and Zeros

1 2 3 4 5 6 7 8 9 10 11

/ / E x am pl e 3 . 1 0

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; I=3*s/(s+2)/(s+4); disp (I , ’ G iv en T r a n s f e r F u n c t i o n : ’ ) ; z e r o = roots ( numer ( I ) ) ; p o l e = roots ( denom ( I ) ) ; disp ( z e r o , ’ Z er os o f t r a n s f e r f u n c t i o n : disp ( p o l e , ’ P o l es o f t r a n s f e r f u n c t i o n : plzr ( I ) ;

Scilab code Exa 3.11 Poles and zeros

1 / / E x am pl e 3 . 1 1 2

19

’ ); ’ );

3 4 5 6 7 8 9 10 11

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; F=4*(s+1)*(s+3)/(s+2)/(s+4); disp (F , ’ G iv en T r a n s f e r F u n c t i o n : ’ ) ; z e r o = roots ( numer ( F ) ) ; p o l e = roots ( denom ( F ) ) ; disp ( z e r o , ’ Z er os o f t r a n s f e r f u n c t i o n : disp ( p o l e , ’ P o l es o f t r a n s f e r f u n c t i o n : plzr ( F ) ;

’ ); ’ );


1 2 3 4 5 6 7 8 9 10 11

/ / E x am pl e 3 . 1 2

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; F=10*s/(s^2+2*s+2); disp (F , ’ G iv en T r a n s f e r F u n c t i o n : ’ ) ; z e r o = roots ( numer ( F ) ) ; p o l e = roots ( denom ( F ) ) ; disp ( z e r o , ’ Z er os o f t r a n s f e r f u n c t i o n : disp ( p o l e , ’ P o l es o f t r a n s f e r f u n c t i o n : plzr ( F ) ;

20

’ ); ’ );

Chapter 4 Z Transforms

Scilab code Exa 4.2 Z transform

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21

/ / E xa mp le 4 . 2

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x 1 =[3 1 2 5 7 0 1]; n1 = - 3:3; X1=x1*(z^-n1)’; x 2 =[2 4 5 7 0 1 2]; n2 = - 2:4; X2=x2*(z^-n2)’; x3 =[1 2 5 4 0 1]; n3 =0:5; X3=x3*(z^-n3)’; x4 =[0 0 1 2 5 4 0 1]; n4 =0:7; X4=x4*(z^-n4)’; X5=z^0; X6=z^-5; X7=z^5; disp ( X 1 , ’ x1 ( n ) = { 3 ,1 ,2 ,5 ,7 ,0 ,1 } X1( z )= ’ ) ; disp ( X 2 , ’ x2 ( n ) = { 2 ,4 ,5 ,7 ,0 ,1 ,2 } X2( z )= ’ ) ; disp ( X 3 , ’ x3 ( n ) = { 1 ,2 ,5 ,4 ,0 ,1 } X3( z )= ’ ) ; disp ( X 4 , ’ x4 ( n ) = { 0 ,0 ,1 ,2 ,5 ,4 ,0 ,1 } X4( z )= ’ ) ; disp ( X 5 , ’ x 5 ( n )= d e l t a ( n ) X5 ( z )= ’ ) ; disp ( X 6 , ’ x6 ( n )= d e l t a ( n − 5) X6( z )=’ ) ;

21

22 disp ( X 7 , ’ x 7 ( n ) = d e l t a ( n +5 ) X7 ( z )= ’ ) ;

Scilab code Exa 4.4 Z transform

1 2 3 4 5 6 7

/ / E xa mp le 4 . 4 clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x = [1 3 0 0 6 -1]; n = -1:4; X=x*(z^-n)’; disp (X , ’ x ( n ) = { 1 ,3 ,0 ,0 ,6 , − 1 } X( z )= ’ ) ;

Scilab code Exa 4.13 Convolution

1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E x am pl e 4 . 1 3

15

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x 1 = [4 -2 1 ]; n 1 = 0: length ( x 1 ) - 1 ; X1=x1*(z^-n1)’; x2 =[1 1 1 1 1]; n2 =0: length ( x 2 ) - 1 ; X2=x2*(z^-n2)’; X3=X1*X2; l = coeff ( numer ( X 3 ) ) ; x3=l(:,$:-1:1); disp ( X 1 , ’ x 1 ( n ) = { 4 , − 2 ,1 } X1( z )= ’ ) ; disp ( X 2 , ’ x 2 ( n ) = { 1 ,1 ,1 ,1 ,1 } X2( z )= ’ ) ; disp ( X 3 , ’ Z t r a ns f o r m o f c o n vo l u ti o n o f t he two s i g n a l s X3 ( z )= ’ ) ; disp ( x 3 , ’ C on vo lu ti on r e s u l t o f t he two s i g n a l s = ’ )

22


1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E x am pl e 4 . 1 4

15

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x 1 =[ 2 1 0 0 .5 ]; n 1 =0 :length ( x 1 ) - 1 ; X1=x1*(z^-n1)’; x 2 = [2 2 1 1 ]; n 2 =0 :length ( x 2 ) - 1 ; X2=x2*(z^-n2)’; X3=X1*X2; l = coeff ( numer ( X 3 ) ) ; x3=l(:,$:-1:1); disp ( X 1 , ’ x1 ( n ) = { 2 , 1 , 0 , 0 . 5 } X1( z )= ’ ) ; disp ( X 2 , ’ x 2 ( n ) = { 2 ,2 ,1 ,1 } X2( z )= ’ ) ; disp ( X 3 , ’ Z t r a ns f o r m o f c o n vo l u ti o n o f t he two s i g n a l s X3 ( z )= ’ ) ; disp ( x 3 , ’ C on vo lu ti on r e s u l t o f t he two s i g n a l s = ’ )

Scilab code Exa 4.16 Cross correlation

1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E x am pl e 4 . 1 6

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x 1 = [1 2 3 4 ]; n 1 =0 :length ( x 1 ) - 1 ; X1=x1*(z^-n1)’; x 2 = [4 3 2 1 ]; n 2 =0 :length ( x 2 ) - 1 ; X2=x2*(z^-n2)’; X2_=x2*(z^n2)’; X3=X1*X2_; l = coeff ( numer ( X 3 ) ) ; x3=l(:,$:-1:1); disp ( X 1 , ’ x 1 ( n ) = { 4 , − 2 ,1 } X1( z )= ’ ) ; disp ( X 2 , ’ x 2 ( n ) = { 4 , − 2 ,1 } X2( z )= ’ ) ;

23

15 disp ( X 3 , ’ Z t ra ns fo rm o f c r o s s c r r e l a t i o n o f t h e two s i g n a l s X3 ( z )= ’ ) ; 16 disp ( x 3 , ’ C r o s s c o r r e l a t i o n r e s u l t o f t h e two s i g n a l s = ’)

Scilab code Exa 4.19 System response

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

/ / E x am pl e 4 . 1 9

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; h = [ 1 2 3 ]; n 1 = 0: length ( h ) - 1 ; H=h*(z^-n1)’; y =[1 1 2 -1 3]; n2 =0: length ( y ) - 1 ; Y=y*(z^-n2)’; X=Y/H; l = coeff ( numer ( X ) ) ; x=l(:,$:-1:1); disp (H , ’ h ( n ) = { 1 ,2 ,3 } H( z )= ’ ) ; disp (Y , ’ y ( n ) = { 1 ,1 ,2 , − 1 , 3 } Y( z )= ’ ) ; disp (X , ’ Z t r a n s fo r m o f i n p u t s e q ue n c e X( z )= ’ ) ; disp (x , ’ I n pp u t S e qu e nc e = ’ )

24

Chapter 5 Linear Time Invariant Systems

Scilab code Exa 5.20 System response

1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E x am pl e 5 . 2 0

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; x = [ -1 1 0 -1] ;n =0 :length ( x ) - 1 ; X=x*(z^-n)’; H=0.2-0.5*z^-2+0.4*z^-3 Y=H*X; l = coeff ( numer ( Y ) ) ; y=l(:,$:-1:1); disp (X , ’ I n pu t s e q ue n c e x ( n ) = { − 1 , 1 , 0 , − 1 } X( z )= ’ ) ; disp (H , ’ S ys te m T r a n s f e r F u n c t i o n H( z )= ’ ) ; disp (Y , ’ Z t r a n s fo r m o f o u tp ut r e s p o n s e Y( z )= ’ ) ; disp (y , ’ D i g i t a l o u tp ut s e q ue n c e y= ’ )


1 / / E x am pl e 5 . 2 1

25

2 3 4 5 6 7 8 9 10

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; H=(1+z^-1)/(1+3/4*z^-1+1/8*z^-2); p o l e = roots ( numer ( H ) ) ; z e r o = roots ( denom ( H ) ) ; disp (H , ’ S ys te m T r a n s f e r F u n c t i o n H( z )= ’ ) ; disp ( z e r o , ’ S ystem z e r o s a r e a t ’ ) ; disp ( p o l e , ’ S ystem p o l e s a re a t ’ ) ; plzr ( H ) ;

26

Chapter 6 Discrete and Fast Fourier Transforms

Scilab code Exa 6.1 Linear and Circular convolution

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19

/ / E xa mp le 6 . 1 clc ; clear ; close ; x1 =[1 1 2 2]; x2 =[1 2 3 4]; y l e n g t h = length ( x 1 ) ;

/ / C a l c ul a t i o n o f l i n e a r c o n vo l u ti o n z = convol ( x 1 , x 2 ) ;

/ / C a l c ul a t i o n o f c i r c u l a r c o n vo l u ti o n for n = 1 : y l e n g t h y(n)=0; for k=1: ylength , l=n-k+1; if l <= 0 then l=l+ylength; end y(n)=y(n)+(x1(k)*x2(l)); end end

/ / C a l c u l a t i o n o f c i r c u l a r c o n v o l ut i o n u s i ng DFT and 27

IDFT 20 21 22 23 24 25 26 27

X1 = fft ( x 1 , - 1 ) ; X2 = fft ( x 2 , - 1 ) ; Y1=X1.*X2; y1 = fft ( Y 1 , 1 ) ; y1 = clean ( y 1 ) ; disp (z , ’ L i ne a r C on vo l ut io n s e qu en c e i s z ( n ) : ’ ) ; disp (y , ’ C i r c u l a r C on vo lu ti on s eq ue nc e i s y ( n ) : ’ ) ; disp ( y 1 , ’ C i r c u l a r C on vo l ut io n s e qu e nc e c a l c u l a t e d u s i n g DFT−IDFT metho d i s y ( n ) : ’ ) ;

Scilab code Exa 6.2 FIR filter response

1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E xa mp le 6 . 2 clc ; clear ; close ; x = [ 1 2 ]; h = [1 2 4];

/ / C a l c ul a t i o n o f l i n e a r c o n vo l u ti o n

15 16 17 18 19 20

y = convol ( x , h ) ; disp (x , ’ I np ut S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ I m pu l se r e s p n o s e o f FIR f i l t e r h ( n ) : ’ ) ; disp (y , ’ O utput s e qu e nc e i s y ( n ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ I n p u t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( h ) ; t i t l e ( ’ I m p u l s e R e s p o n s e h [ n ] : ’ ) ; y l a b e l ( ’ Amplitu de −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( y ) ; t i t l e ( ’ O u tp u t S e q e n c e y [ n ] = x [ n ] ∗ h [ n ] : ’ ) ; y l a b e l ( ’

28

Amplitude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )


1 2 3 4 5 6 7 8 9 10 11 12 13 14

/ / E xa mp le 6 . 3 clc ; clear ; close ; x = [1 1 1]; h = [1 1 1];


15 16 17 18 19 20

y = convol ( x , h ) ; disp (x , ’ F i r s t S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ S ec on d S eq ue nc e i s h ( n ) : ’ ) ; disp (y , ’ O utput s e qu e nc e i s y ( n ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( h ) ; t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( y ) ; t i t l e ( ’ C o n v o l u t i o n S e q e n c e y [ n ] = x [ n ] ∗ h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )


1 / / E xa mp le 6 . 4 2 3 clc ; clear ; close ;

29

4 5 6 7 8 9 10 11 12 13 14

a=0.5; n=1:50; x = ones ( 1 , 5 0 ) ; h=a^n;

15 16 17 18


19 20 21 22 23 24

for i = 1 : 5 0 y ( 1 , i ) = sum ( h ( 1 : i ) ) ; end disp ( ’ F i r s t S eq ue nc e i s x ( n )=u ( n ) ’ ) ; disp (a , ’ S ec on d S e qu e nc e i s h ( n ) =a ˆ n ∗ u ( n ) wh ere a= ’ ) ; disp (y , ’ O utput s e qu e nc e i s y ( n ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( h ) ; t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( y ) ; t i t l e ( ’ C o n v o l u t i o n S e q e n c e y [ n ] = x [ n ] ∗ h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )


1 2 3 4 5 6 7

/ / E xa mp le 6 . 5 clc ; clear ; close ; x = [1 2 3 ]; x m in = 0 ; nx = x mi n :length (x) +xmin -1; h = [ 1 2 -2 - 1] ; h mi n = - 1; n h =length (h) +hmin -1;

/ / C a l c ul a t i o n o f l i n e a r c o n vo l u ti o n y = convol ( x , h ) ;

30

8 9 10 11 12 13 14 15

y m i n = x m i n + h m i n ; n y = y m i n :length (y) +ymin -1;

16 17 18 19 20 21

disp (x , ’ F i r s t S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ S ec on d S eq ue nc e i s h ( n ) : ’ ) ; disp (y , ’ O utput s e qu e nc e i s y ( n ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( n x , x ) ; t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( n h , h ) ; t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( n y , y ) ; t i t l e ( ’ C o n v o l u t i o n S e q e n c e y [ n ] = x [ n ] ∗ h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

/ / E xa mp le 6 . 6 clc ; clear ; close ; x = [1 1 0 1 1 ]; x mi n = -2; n x= x mi n :length (x) +xmin -1; h = [ 1 -2 -3 4 ]; h m in = - 3 ; nh =length (h) +hmin -1;

/ / C a l c ul a t i o n o f l i n e a r c o n vo l u ti o n y = convol ( x , h ) ; y m i n = x m i n + h m i n ; n y = y m i n :length (y) +ymin -1;

disp (x , ’ F i r s t S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ S ec on d S eq ue nc e i s h ( n ) : ’ ) ; disp (y , ’ O utput s e qu e nc e i s y ( n ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( n x , x ) ;

31

16 t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) 17 subplot ( 3 , 1 , 2 ) ; 18 plot2d3 ( n h , h ) ; 19 t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) 20 subplot ( 3 , 1 , 3 ) ; 21 plot2d3 ( n y , y ) ; 22 t i t l e ( ’ C o n v o l u t i o n S e q e n c e y [ n ] = x [ n ] ∗ h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )

Scilab code Exa 6.8 DFT

1 2 3 4 5 6 7 8 9 10

/ / E xa mp le 6 . 8 clc ; clear ; close ; L=3;A=1/4; x = A * ones ( 1 , L ) ;

/ / C a l c u l a t i o n o f DFT X = fft ( x , - 1 ) ; X = clean ( X ) ; disp (x , ’ Given S eq ue nc e i s x ( n ) : ’ ) ; disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) :


1 2 3 4 5 6 7

/ / E xa mp le 6 . 9 clc ; clear ; close ; L=3;A=1/5; n=-1:1; x = A * ones ( 1 , L ) ;

/ / C a l c u l a t i o n o f DFT 32

’ );

8 9 10 11

X = fft ( x , - 1 ) ; X = clean ( X ) ; disp (x , ’ Given S eq ue nc e i s x ( n ) : ’ ) ; disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) :

’ );


1 2 3 4 5 6 7 8 9 10 11 12

/ / E x am pl e 6 . 1 0 clc ; clear ; close ; x =[1 1 2 2 3 3];

/ / C a l c u l a t i o n o f DFT

13 14 15 16 17 18

X = fft ( x , - 1 ) ; X = clean ( X ) ; disp (x , ’ Given S eq ue nc e i s x ( n ) : ’ ) ; disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) : ’ ) ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ G i v e n S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ; ’ ) ; x l a b e l ( ’ n−−> ; ’ ) ; subplot ( 3 , 1 , 2 ) ; plot2d3 ( abs ( X ) ) ; t i t l e ( ’ M a g n i t u de S p e c t ru m | X ( k ) | ’ ) ; x l a b e l ( ’ k−−> ; ’ ) ; subplot ( 3 , 1 , 3 ) ; plot2d3 ( atan ( X ) ) ; t i t l e ( ’ P ha s e S p ec tr u m / X ( k ) ’ ) ; x l a b e l ( ’ k−−> ; ’ ) ;


1 / / E x am pl e 6 . 1 1

33

2 3 4 5 6 7 8 9 10 11 12

clc ; clear ; close ; N=8;A=1/4; n=0:N-1; x=A^n;

/ / C a l c u l a t i o n o f DFT X = fft ( x , - 1 ) ; X = clean ( X ) ; disp (x , ’ Given S eq ue nc e i s x ( n ) : ’ ) ; disp (N , ’N= ’ ) disp (X , ’ N− p o i nt DFT o f t he S eq ue nc e i s X( k ) :


1 2 3 4 5 6 7 8 9 10 11

/ / E x am pl e 6 . 1 2 clc ; clear ; close ; N=4; n=0:N-1; x = cos ( % p i / 4 * n ) ;

/ / C a l c u l a t i o n o f DFT X = fft ( x , - 1 ) ; X = clean ( X ) ; disp (x , ’ Given S eq ue nc e i s x ( n ) : ’ ) ; disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) :

Scilab code Exa 6.13 Inverse DFT

1 2 3 4 5

/ / E x am pl e 6 . 1 3 clc ; clear ; close ; X = [1 2 3 4];

/ / C a l c u l a t i o n o f IDFT x = fft ( X , 1 ) ;

34

’ );

’ );

6 x = clean ( x ) ; 7 disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) : 8 disp (x , ’ S eq ue nc e i s x ( n ) : ’ ) ;

’ );

Scilab code Exa 6.14 Inverse DFT

1 2 3 4 5 6 7 8

/ / E x am pl e 6 . 1 4 clc ; clear ; close ; X = [3 2 + %i 1 2 - %i ] ;

/ / C a l c u l a t i o n o f IDFT x = fft ( X , 1 ) ; x = clean ( x ) ; disp (X , ’DFT o f t h e S eq ue nc e i s X( k ) : disp (x , ’ S eq ue nc e i s x ( n ) : ’ ) ;

Scilab code Exa 6.15 DIT FFT

1 2 3 4 5 6 7

/ / E x am pl e 6 . 1 5 clc ; clear ; x =[1 2 3 4 4 3 2 1 ]; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


1 / / E x am pl e 6 . 1 6 2 3 clc ; clear ;

35

’ );

4 5 6 7

x =[0 1 2 3 4 5 6 7]; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


1 2 3 4 5 6 7 8

/ / E x am pl e 6 . 1 7

clc ; clear ; n=0:7; x=2^n; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


1 2 3 4 5 6 7

/ / E x am pl e 6 . 1 8 clc ; clear ; x = [0 1 2 3]; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;

Scilab code Exa 6.19 DIF FFT

1 / / E x am pl e 6 . 1 9 2

36

3 4 5 6 7

clc ; clear ; x =[1 2 3 4 4 3 2 1 ]; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


1 2 3 4 5 6 7 8

/ / E x am pl e 6 . 2 0

clc ; clear ; n=0:7; x=2^n; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


1 2 3 4 5 6 7 8

/ / E x am pl e 6 . 2 1

clc ; clear ; n=0:7; x=n+1; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;


37

1 2 3 4 5 6 7 8

/ / E x am pl e 6 . 2 1 clc ; clear ; n=0:3; x = cos ( n * % p i / 2 ) ; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;

Scilab code Exa 6.23 IFFT

1 2 3 4 5 6 7

/ / E x am pl e 6 . 2 3 clc ; clear ; X = [ 6 - 2+ 2* % i -2 - 2 -2 * %i ] ; x = clean ( i f f t ( X ) ) ; disp (X , ’X ( k ) = ’ ) ; disp (x , ’ x ( n )= ’ ) ;


1 / / E x am pl e 6 . 2 4 2 3 clc ; clear ; 4 X = [2 0 - 5. 82 8 - 2 .4 14 * % i 0 - 0. 17 2 - 0. 41 4* % i 0 - 0 . 17 2 +0 . 4 14 * % i 0 - 5 . 82 8 + 2. 4 14 * % i ] ; 5 x = round ( clean ( i f f t ( X ) ) ) ; 6 disp (X , ’X ( k ) = ’ ) ; 7 disp (x , ’ x ( n )= ’ ) ;

38


1 / / E x am pl e 6 . 2 5 2 3 clc ; clear ; 4 X = [ 2 5 5 4 8 . 6 3+ 1 6 6 .0 5 * % i - 5 1+ 1 02 * % i - 7 8. 6 3 +4 6 . 05 * % i - 85 - 7 8. 6 3 - 4 6. 0 5* % i - 51 - 1 0 2* % i 4 8 .6 3 - 1 6 6 . 05 * % i ] ; 5 x = round ( clean ( i f f t ( X ) ) ) ; 6 disp (X , ’X ( k ) = ’ ) ; 7 disp (x , ’ x ( n )= ’ ) ;


1 / / E x am pl e 6 . 2 6 2 3 clc ; clear ; 4 X = [ 36 - 4 +9 . 65 6 * % i - 4+ 4* % i - 4 +1 . 65 6 * % i - 4 - 4 - 1 .6 56 * % i - 4 -4 * %i - 4 - 9 .6 56 * % i ] ; 5 x = round ( clean ( i f f t ( X ) ) ) ; 6 disp (X , ’X ( k ) = ’ ) ; 7 disp (x , ’ x ( n )= ’ ) ;


1 2 3 4 5 6 7 8 9

/ / E x am pl e 6 . 2 7

clc ; clear ; t=0:0.0025:0.0175; f=50; x = sin ( 2 * % p i * f * t ) ; X = clean ( fft ( x ) ) ; disp (x , ’ x ( n )= ’ ) ; disp (X , ’X ( k ) = ’ ) ;

39

Scilab code Exa 6.34 Overlap Add Convolution

1 2 3 4 5 6

/ / E x am pl e 6 . 3 4 clc ; clear ; close ; h = [2 2 1]; x =[3 0 -2 0 2 1 0 -2 -1 0]; // l e n g th o f i m pu l se M = length ( h ) ;

response / / l e n g t h o f FFT/ IFFT

7 L=2^M;

operation 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27

N=L-M+1; xl = length ( x ) ; // number o f i t e r a t i o n s K = ceil ( x l / N ) ; h = [ h zeros ( 1 , L - M ) ] ; x = [ x x ( 1 : K * N - xl ) ] ; H = fft ( h ) ; y = zeros ( 1 , M - 1 ) ; for k = 0 : K - 1 x k = [ x ( k * N + 1 : ( k + 1 ) * N ) zeros ( 1 , M - 1 ) ] ; Xk = fft ( x k ) ; Yk=H.*Xk; yk=ifft(Yk); yk = clean ( y k ) ; y = [ y ( 1: k * N ) y ( k * N + 1: k * N + M - 1 ) + yk ( 1 : M - 1 ) y k ( M : L ) ]; disp ( k + 1 , ’ S eg m en t = ’ ) ; disp ( x k , ’ xk ( n )= ’ ) ; disp ( y k , ’ yk ( n )= ’ ) ; end y=y(1:xl+M-1); disp (y , ’ O utput S eq ue nc e i s y ( n ) : ’ ) ;

40

Scilab code Exa 6.35 Overlap Save Convolution

1 2 3 4 5 6

/ / E x am pl e 6 . 3 5 clc ; clear ; close ; h = [2 2 1]; x =[3 0 -2 0 2 1 0 -2 -1 0]; M = length ( h ) ; // l e n g th o f i m pu l se

response 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25

L=2^M; / / l e n g t h o f FFT/ IFFT o p e r a t i o n N=L-M+1; xl = length ( x ) ; K = ceil ( x l / N ) ; // number o f i t e r a t i o n s h = [ h zeros ( 1 , L - M ) ] ; x =[ zeros ( 1 ,M - 1) x x ( 1: K * N - xl ) ] ; H = fft ( h ) ; for k = 0 : K - 1 xk=x(k*N+1:(k+1)*N+M-1); Xk = fft ( x k ) ; Yk=H.*Xk; yk=ifft(Yk); yk = clean ( y k ) ; y = [ y k ( 1: k * N ) y k ( M : L ) ]; disp ( k + 1 , ’ S eg m en t = ’ ) ; disp ( x k , ’ xk ( n )= ’ ) ; disp ( y k , ’ yk ( n )= ’ ) ; end disp (y , ’ O utput S eq ue nc e i s y ( n ) : ’ ) ;

Scilab code Exa 6.36 Cross Correlation

1 / / E x am pl e 6 . 3 6 2 3 clc ; clear ; close ; 4 x = [1 0 0 1];

41

5 6 7 8 9 10 11 12

h = [4 3 2 1]; y l e n g t h = length ( x ) + length ( h ) - 1 ; x l e n g t h = length ( x ) ; x =[ zeros (1 , length ( h ) -1) x zeros (1 , length ( h ) - 1 ) ] ; y=0;

// C a l c u l a t i o n o f c r o s s c o r r e l a t i o n for n = 1 : y l e n g t h ; y ( n ) = x * [ zeros ( 1 ,n - 1 ) h zeros (1, ylength -n)] ’;

// t h i s i n s t r u c t i o n p er fo rm s c r o s s c o rr e la t i on o f x & h 13 14 15 16 17 18 19 20 21

end

22 23 24 25 26 27

disp (x , ’ F i r s t S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ S ec on d S eq ue nc e i s h ( n ) : ’ ) ; disp (y , ’ C o r r e l a t i o n S eq ue nc e y [ n ] i s ’ ) ; figure ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( h ) ; t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( y ) ; t i t l e ( ’ C o r r e l a t i o n S e q en c e y [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )

Scilab code Exa 6.37 Circular Correlation

1 / / E x am pl e 6 . 3 7 2 3 clc ; clear ; close ; 4 x = [1 0 0 1];

42

5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27

h = [4 3 2 1]; y l e n g t h = length ( x ) ; y=0;

/ / C a lc u la t io n o f c i r c u l a r

correlation

for n=1: ylength , y(n)=0; for k=1: ylength , l=k-n+1; if l <= 0 then l=l+ylength; end y(n)=y(n)+(x(k)*h(l)); end y(n)=y(n)/4; end

28 29 30 31 32 33

disp (x , ’ F i r s t S eq ue nc e i s x ( n ) : ’ ) ; disp (h , ’ S ec on d S eq ue nc e i s h ( n ) : ’ ) ; disp (y , ’ C o r r e l a t i o n S eq ue nc e y [ n ] i s ’ ) ; figure ; subplot ( 3 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ F i r s t S e q e n c e x [ n ] : ’ ) ; y l a b e l ( ’ Ampli tude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 2 ) ; plot2d3 ( h ) ; t i t l e ( ’ S e c o n d S e q e n c e h [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ ) subplot ( 3 , 1 , 3 ) ; plot2d3 ( y ) ; t i t l e ( ’ C o r r e l a t i o n S e q en c e y [ n ] : ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ n−−> ’ )

43

Chapter 7 Finite Impulse Response Filters

Scilab code Exa 7.3 Low pass filter using fourier series method

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

/ / E xa mp le 7 . 3 clc ; clear ; close ; / / p a s sb a nd f r e q u e n c y fp=2000; F=9600; / / s a m p l in g f r e q u a n c y

/ / C a l c u l a t i o n o f f i l t e r co − e f f i c i e n t s a 0 = 1 / F * i n t e g r a t e( ’ 1 ’ , ’ t ’ , - f p , f p ) ; for n = 1 : 1 0 ; a ( 1 , n ) = 2 / F * i n t e g r a t e( ’ cos (2 ∗ %p i ∗ n ∗ f /F) ’ , ’ f ’ , - f p , f p ) ; end h = [ a( : , $ : - 1: 1) / 2 a 0 a / 2] ;

/ / D i s p l a y i n g f i l t e r c o− e f f i c i e n t s disp (F , ’ S a mp li n g f r e q u e n c y F= ’ , f p , ’ Assu mpti on : P as sb an d f r e q u e n cy f p= ’ ) ; 16 disp ( ’ F i l t e r co − e f f i c i e n t s : ’ ) ; 17 disp ( a 0 , ’ h ( 0 )= ’ ) ; disp ( a / 2 , ’h (n )=h(− n ) = ’ ) ; 18 19 n = - 1 0 : 1 0 ; 20 plot2d3 ( n , h ) ;

44

21 t i t l e ( ’ F i l t e r );

t r a n s f e r f u n c t i o n h ( n ) ’ ) ; x l a b e l ( ’ n−−> ’

Scilab code Exa 7.4 Low pass filter using Type 1 frequency sampling technique

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19

/ / E xa mp le 7 . 4 clc ; clear ; close ; M=7;w=2*%pi/M;

/ / C a l c u l a t i o n o f f i l t e r co − e f f i c i e n t s k = [0 1 6]; for n = 0 : M - 1 h ( n + 1 ) = sum ( exp ( - % i * 3 * w * k ) . * exp ( % i * w * k * n ) ) / M ; end h = clean ( h ) ;

/ / D i s p l a y i n g f i l t e r c o− e f f i c i e n t s disp (M , ’ F i l t e r O rd er M= ’ ) ; disp ( ’ F i l t e r co − e f f i c i e n t s : ’ ) ; disp (h , ’ h ( n )= ’ ) ; plot2d3 ( h ) ; title( ’ F i l t e r );

t r a n s f e r f u n c t i o n h ( n ) ’ ) ; x l a b e l ( ’ n−−> ’

45

Chapter 8 Infinite Impulse Response Filters

Scilab code Exa 8.1 IIR filter Design byBackward Difference For Derivative method

1 2 3 4 5 6 7 8

/ / E xa mp le 8 . 1

9 10

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=1/(s+2); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ; disp ( ’ U si ng Backward d i f f e r e n c e f or mu l a f o r deriva tive : ’) disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 / / E xa mp le 8 . 2 2 clc ; clear ; close ;

46

3 4 5 6 7 8

9 10

s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=1/(s^2+16); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ; disp ( ’ U si ng Backward d i f f e r e n c e f or mu l a f o r deriva tive : ’) disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6 7 8

/ / E xa mp le 8 . 3

9 10

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=1/((s+0.1)^2+9); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ; disp ( ’ U si ng Backward d i f f e r e n c e f or mu l a f o r deriva tive : ’) disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;

Scilab code Exa 8.4 IIR filter Design by Impulse Invariant method

1 2 3 4 5 6 7

/ / E xa mp le 8 . 4 clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=(s+0.2)/((s+0.2)^2+9); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ;

47

8 disp ( ’ U s in g I m p u ls e I n v a r i a n t T e ch n iq u e : ’ ) 9 disp ( H s , ’H( s )= ’ ) ; 10 disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6 7 8 9 10

/ / E xa mp le 8 . 5

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=1/(s+1)/(s+2); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ; disp ( ’ U s in g I m p u ls e I n v a r i a n t T e ch n iq u e : ’ ) disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6 7 8 9 10

/ / E xa mp le 8 . 6

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=1; Hs=1/(s+0.5)/(s^2+0.5*s+2); Hz = horner ( H s , ( 1 - 1 / z ) / T ) ; disp ( ’ U s in g I m p u ls e I n v a r i a n t T e ch n iq u e : ’ ) disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;

48

Scilab code Exa 8.7 IIR filter Design by Bilinear Transformation method

1 2 3 4 5 6 7 8 9 10

/ / E xa mp le 8 . 7

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=0.276; Hs=(s+0.1)/((s+0.1)^2+9); Hz = ss2tf ( cls2dls ( tf2ss ( H s ) , T ) ) ; disp ( ’ U si ng B i l i n e a r T r an s f or m a ti o n : ’ ) ; disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6 7 8 9

/ / E xa mp le 8 . 8

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=0.1; Hs=2/(s+1)/(s+2); Hz = ss2tf ( cls2dls ( tf2ss ( H s ) , T ) ) ; disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6

/ / E xa mp le 8 . 9 clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=0.1; wr=0.25*%pi;

// G iven c ut o f f f r eq u e nc y 49

7 8 9 10 11 12

f c = 2 / T * tan ( w r / 2 ) ; Hs=fc/(s+fc); Hz = ss2tf ( cls2dls ( tf2ss ( H s ) , T ) ) ; disp ( ’ U si ng B i l i n e a r T r an s f or m a ti o n : ’ ) ; disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;


1 2 3 4 5 6 7 8 9 10

/ / E x am pl e 8 . 1 0

clc ; clear ; close ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; T=0.1; Hs=1/(s+1)^2; Hz = ss2tf ( cls2dls ( tf2ss ( H s ) , T ) ) ; disp ( ’ U si ng B i l i n e a r T r an s f or m a ti o n : ’ ) ; disp ( H s , ’H( s )= ’ ) ; disp ( H z , ’H( z )= ’ ) ;

Scilab code Exa 8.11 Butterworth Filter using Impulse Invariant transformation

1 2 3 4 5 6 7 8 9 10 11

/ / E x am pl e 8 . 1 1

clc ; clear ; close ; rp=0.707 rs=0.2 wp=%pi/2; ws=3*%pi/4; T=1; fp=wp/T; fs=ws/T; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ;

// p a ss ba nd r i p p l e // s t op ba nd r i p p l e / / p a s sb a nd f r e q u e n c y / / s t o pb a nd f r e q u e n c y

50

12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32

hs=1;

//Calculat ing the order of f i l t e r n u m = log ( ( r s ^ -2 - 1) / ( r p ^ - 2 - 1) ) ; d e n = 2 * log ( f s / f p ) ; N = ceil ( n u m / d e n ) ;

// C a l c u l a t i o n o f cu t − o f f f r eq u en c y f c = f p / ( r p ^ - 2 - 1) ^ ( 0 . 5/ N ) ;

// Calc ulat ing

f i l t e r respons e

i f m o du lo ( N , 2 ) = = 1 then b = - 2 * sin ( % p i / ( 2 * N ) ) ; hs=hs*fc/(s+fc); end for k = 1 : N / 2 b = 2 * sin ( ( 2 * k - 1 ) * % p i / ( 2 * N ) ) ; hs=hs*fc^2/(s^2+b*fc*s+fc^2); end hs = clean ( h s ) ; s y s = syslin ( ’ c ’ , h s ) ; hz = horner ( ss2tf ( dscr ( s y s , T ) ) , 1 / z ) ;

//

c o n v e r t i n g H( s ) t o H( z ) 33 34 35 36 37 38

// Displ aying


[ h z m , f r ] = frmag ( h z , 2 5 6 ) ; disp ( h z , ’ F i l t e r T ra n s f e r f u n c t i o n : ’ ) ; plot ( f r , h z m ) ; t i t l e ( ’ L o wp as s B u tt e rw o rt h F i l t e r R es po ns e ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ N o rm a li s ed f r e q u e n c y f / f s −−> ’ ) ;

Scilab code Exa 8.12 Butterworth Filter using Bilinear transformation

1 / / E x am pl e 8 . 1 2 2 clc ; clear ; close ;

51

3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31

// p a ss ba nd // s t op ba nd / / p a s sb a nd / / s t o pb a nd

rp=0.9 rs=0.2 wp=%pi/2; ws=3*%pi/4; T=1; f p = 2 / T * tan ( w p / 2 ) ; f s = 2 / T * tan ( w s / 2 ) ; s = poly (0 , ’ s ’ ) ; z = poly (0 , ’ z ’ ) ; hs=1;

ripp le ripp le frequency frequency

//Calculat ing the order of f i l t e r n u m = log ( ( r s ^ -2 - 1) / ( r p ^ - 2 - 1) ) ; d e n = 2 * log ( f s / f p ) ; N = ceil ( n u m / d e n ) ;

// C a l c u l a t i o n o f cu t − o f f f r eq u en c y f c = f p / ( r p ^ - 2 - 1) ^ ( 0 . 5/ N ) ;

// Calc ulat ing


i f m o du lo ( N , 2 ) = = 1 then hs=hs*fc/(s+fc); end for k = 1 : N / 2 b = 2 * sin ( ( 2 * k - 1 ) * % p i / ( 2 * N ) ) ; hs=hs*fc^2/(s^2+b*fc*s+fc^2); end hs = clean ( h s ) ; s y s = syslin ( ’ c ’ , h s ) ; hz = ss2tf ( cls2dls ( tf2ss ( s y s ) , T ) ) ;

// conve rting

H( s ) to H( z ) 32 33 34 35 36 37

// Displ aying


[ h z m , f r ] = frmag ( h z , 2 5 6 ) ; disp ( h z , ’ F i l t e r T ra n s f e r f u n c t i o n : ’ ) ; plot ( f r , h z m ) ; t i t l e ( ’ L o wp as s B u tt e rw o rt h F i l t e r R es po ns e ’ ) ; y l a b e l ( ’ Amplit ude −−> ’ ) ; x l a b e l ( ’ N o rm a li s ed f r e q u e n c y f / f s −−> ’ ) ;

52

Scilab code Exa 8.14 Filter transformation transformation

1 2 3 4 5 6 7 8 9 10

/ / E x am am pl pl e 8 . 1 4

c l c ; clear ; close ; s = poly ( 0 , ’ s ’ ) ; fc=1; / / As Assumed c u t o f f / / G i ve ve n d a t a Q=10;f0=2; Hs=1/(s^2+2*s+1); l=fc*(s^2+f0^2)/(s*f0/Q); H s 1 = horner ( H s , l ) ; disp ( H s , ’ Lo Low p a s s f i l t e r H ( s ) = ’ ) ; B an d p a s s f i l t e r H ( s )= ’ ) ; disp ( H s 1 , ’ Ba

frequency

Scilab code Exa 8.15 Filter transformation transformation

1 2 3 4

/ / E x am am pl pl e 8 . 1 5 c l c ; clear ; close ; s = poly ( 0 , ’ s ’ ) ; fc=1;

low pass

/ / As A s s u m e d c ut u t o f f f r e qu q u e n cy cy o f

filter

/ / As As s u m e d c ut u t o f f f r e q u e n cy cy o f high pass f i l t e r

5 f0=5; 6 7 8 9

Hs=fc/(s+fc); H s 1 = horner ( H s , f c * f 0 / s ) ; )= ’ , f c , ’ Lo Low p a s s f i l t e r w i t h f c = ’ ) ; disp ( H s , ’H( s )= disp ( H s 1 , ’H( s )= ’ , f 0 , ’ H i gh g h p a s s f i l t e r w i th th f c = ’ ) ;

53

Chapter 9 Realisation of Digital Linear Systems

Scilab code Exa 9.4 Cascade Cascade Realisation Realisation

1 2 3 4 5 6 7

/ / E xa x a mp m p le le 9 . 4 c l c ; clear ; close ; z = poly ( 0 , ’ z ’ ) ; Hz=2*(z+2)/(z*(z-0.1)*(z+0.5)*(z+0.4)); H = dscr ( H z , 0 . 1 ) ; st e m F u n c t i o n H ( z ) = ’ ) ; disp ( H z , ’ S y st disp ( H , ’ S ys ys t e m F u un n cctt io i o n f o r c a sc sc a d dee r e a l i s a t i o n Hk ( z )= ’ ) ;

Scilab code Exa 9.5.a Parallel Parallel Realisation Realisation

1 2 3 4 5

/ / E xa xa mp mp le le 9 . 5 . a c l c ; clear ; close ; z = poly ( 0 , ’ z ’ ) ; s = poly ( 0 , ’ s ’ ) ; Hz=3*(2*z^2+5*z+4)/(2*z+1)/(z+2);

54

6 H = pfss ( H z / z ) ; 7 f o r k=1: length ( H ) 8 H ( k ) = clean ( H ( k ) ) ; H 1 ( k ) = z * horner ( H ( k ) , z ) ; 9 10 disp ( H 1 ( k ) , ’ S y s t e m F u n ct c t i on on f o r p a r a l l e l Hk( z )= ’ ) ; 11 e n d 12 disp ( H z , ’ S y st st e m F u n c t i o n H ( z ) = ’ ) ;

reali satio n

Scilab code Exa 9.5.b Parallel Parallel Realisation Realisation

1 2 3 4 5 6 7 8 9 10

/ / E xa xa mp mp le le 9 . 5 . b

11 12

c l c ; clear ; close ; z = poly ( 0 , ’ z ’ ) ; s = poly ( 0 , ’ s ’ ) ; Hz=3*z*(5*z-2)/(z+1/2)/(3*z-1); H = pfss ( H z / z ) ; f o r k=1: length ( H ) H ( k ) = clean ( H ( k ) ) ; H 1 ( k ) = z * horner ( H ( k ) , z ) ; disp ( H 1 ( k ) , ’ S y s t e m F u n ct c t i on on f o r p a r a l l e l Hk( z )= ’ ) ; end disp ( H z , ’ S y st st e m F u n c t i o n H ( z ) = ’ ) ;

55

reali satio n

Chapter 10 Effects of Finite Word Length in Digital Filters

Scilab code Exa 10.2 Output Quantisation Noise

1 2 3 4 5 6 7 8 9 10 11 12 13

/ / E x am pl e 1 0 . 2

clc ; clear ; close ; z = poly (0 , ’ z ’ ) ; H=0.5*z/(z-0.5); B=8; pn=2^(-2*B)/12; X = H * horner ( H , 1 / z ) / z ; r = roots ( denom ( X ) ) ; rl = length ( r ) ; rc = coeff ( denom ( X ) ) q1=[];q2=[]; for n = 1 : r l

/ / N o i s e p o we r

// Loop t o s e p a r a t e p o l e s

ins id e the unit c i r c l e 14 15 16 17 18

if ( abs ( r ( n ) ) < 1 ) then q 1 = [ q1 r ( n ) ] ; else q 2 = [ q2 r ( n ) ] ; end

56

19 20 21 22 23 24 25 26

end P = numer ( X ) / r c ( length ( r c ) ) ; Q1 = poly ( q 1 , ’ z ’ ) ; Q2 = poly ( q 2 , ’ z ’ ) ; / / R es i du e C a l c u l a t i o n I = residu ( P , Q 1 , Q 2 ) ; po=pn*I; / / O ut pu t N o i s e p ow er disp ( p n , ’ I n p u t N o i s e p ow er ’ ) ; disp ( p o , ’ O ut pu t N o i s e p ow er ’ ) ;

Scilab code Exa 10.3 Deadband Interval

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16

/ / E x am pl e 1 0 . 3 clc ; clear ;

/ / y ( n ) = 0 .9 y ( n − 1)+x(n) / / I n p u t x ( n ) =0 n=-1;y=12; / / I n i t i a l C o n d i t i o n y ( − 1)=12 flag=1; while n <8 n=n+1; y = [ y 0 . 9* y ( n + 1 ) ] ; yr = round ( y ) ; end disp (n , ’n= ’ ) ; disp (y , ’ y ( n )− e x a c t ’ ) ; disp ( y r , ’ y ( n )− r o u n d e d ’ ) ; disp ( [ - y r ( n + 2) y r ( n + 2) ] , ’ Deadband i n t e r v a l

Scilab code Exa 10.4 Deadband Interval

1 / / E x am pl e 1 0 . 4 2 3 clc ; clear ;

57

’)

4 5 6 7

/ / y ( n ) = 0 .9 y ( n − 1)+x(n) a=0.9; l = ceil ( 0 . 5 / ( 1 - abs ( a ) ) ) ; disp ([ - l l ], ’ Deadband i n t e r v a l

’)

Scilab code Exa 10.5 Output Quantisation Noise for Cascade realisation

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16

/ / E x am pl e 1 0 . 5

clc ; clear ; close ; x = poly (0 , ’ x ’ ) ; z = poly (0 , ’ z ’ ) ; H1=1/(1-0.9/z); H2=1/(1-0.8/z); H=H1*H2; pn=x/12;

//x=2ˆ−2B

/ / I n p u t N o i s e p ow er

// C a l c u l a t i o n o f o ut pu t n o i s e f o r H1 ( z ) X 1 = H * horner ( H , 1 / z ) / z ; r1 = roots ( denom ( X 1 ) ) ; r c 1 = coeff ( denom ( X 1 ) ) ; q1=[];s1=[]; for n=1: length ( r 1 )

// Loop t o s e p a r a t e

p o l e s i n s i d e t he u n i t c i r c l e 17 18 19 20 21 22 23 24 25 26

if ( abs ( r 1 ( n ) ) < 1 ) then q 1 = [ q1 r 1 ( n ) ]; else s 1 = [ s1 r 1 ( n ) ]; end

end P1 = numer ( X 1 ) / r c 1 ( length ( r c 1 ) ) ; Q1 = poly ( q 1 , ’ z ’ ) ; S1 = poly ( s 1 , ’ z ’ ) ; I1 = abs ( residu ( P 1 , Q 1 , S 1 ) ) ;

Calculation 58

// Resi due

27 28 29 30 31 32 33 34

/ / O ut pu t N o i s e p ow er

po1=pn*I1;

// C a l c u l a t i o n o f o ut pu t n o i s e f o r H2 ( z ) X 2 = H 2 * horner ( H 2 , 1 / z ) / z ; r2 = roots ( denom ( X 2 ) ) ; r c 2 = coeff ( denom ( X 2 ) ) ; q2=[];s2=[]; for n=1: length ( r 2 )

// Loop t o s e p a r a t e

p o l e s i n s i d e t he u n i t c i r c l e 35 36 37 38 39 40 41 42 43 44

if ( abs ( r 2 ( n ) ) < 1 ) then q 2 = [ q2 r 2 ( n ) ]; else s 2 = [ s2 r 2 ( n ) ]; end

end P2 = numer ( X 2 ) / r c 2 ( length ( r c 2 ) ) ; Q2 = poly ( q 2 , ’ z ’ ) ; S2 = poly ( s 2 , ’ z ’ ) ; I2 = abs ( residu ( P 2 , Q 2 , S 2 ) ) ;

// Resi due

Calculation / / O ut pu t N o i s e

45 p o 2 = p n * I 2 ;

power 46 47 48 49 50 51 52

po=po1+po2; disp ( p n , ’ I n p u t N o i s e p ow er ’ ) ; disp ( I 1 , ’ I1= ’ ) ; disp ( I 2 , ’ I2= ’ ) ; disp ( p o 1 , ’ O ut pu t N o i s e p ow er f o r H1 ( z ) ’ ) ; disp ( p o 2 , ’ O ut pu t N o i s e p ow er f o r H2 ( z ) ’ ) ; disp ( p o , ’ T o t a l O ut pu t N o i s e p ow er ’ ) ;

59

Chapter 11 Multirate Digital Signal Processing

Scilab code Exa 11.1 Time Decimation

/ / E x am pl e 1 1 . 1

1 2 3 4 5 6 7

8 9 10 11

clc ; clear ; close ; x = [ 0 : 6 0 : 6] ; y = x ( 1 : 3 : length ( x ) ) ; disp (x , ’ I np ut s i g n a l x ( n )= ’ ) ; disp (y , ’ O utput s i g n a l o f d e c im a t io n p r o c es s by f a c t o r t h re e y ( n ) ’ ) ; subplot ( 2 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ I np ut s i g n a l x ( n ) ’ ) ; subplot ( 2 , 1 , 2 ) ; plot2d3 ( y ) ; t i t l e ( ’ O utput s i g n a l y ( n ) ’ ) ;

Scilab code Exa 11.2 Interpolation

1 / / E x am pl e 1 1 . 2

60

2 3 4 5 6 7 8 9 10 11 12 13 14

clc ; clear ; close ; x=0:5; y=[]; for i=1: length ( x ) y(1,2*i)=x(i); end disp (x , ’ I np ut s i g n a l x ( n )= ’ ) ; disp (y , ’ Output s i g n a l o f i n t e r p o l a t i o n p r o c e s s w i t h f a c t o r two y ( n ) ’ ) ; subplot ( 2 , 1 , 1 ) ; plot2d3 ( x ) ; t i t l e ( ’ I np ut s i g n a l x ( n ) ’ ) ; subplot ( 2 , 1 , 2 ) ; plot2d3 ( y ) ; t i t l e ( ’ O utput s i g n a l y ( n ) ’ ) ;

Scilab code Exa 11.4 Polyphase Decomposition

1 2 3 4 5 6 7 8 9 10 11 12 13

/ / E x am pl e 1 1 . 4

clc ; clear ; z = poly (0 , ’ z ’ ) ; num=1-4*z^-1; den=1+5*z^-1; H=num/den; num1=num*(1-5*z^-1); den1=den*(1-5*z^-1); H1=num1/den1; c = coeff ( numer ( n u m 1 ) ) ; c l e n g t h = length ( c ) ; c = [ c zeros (1 , pmodulo (clength ,2))];

//make

l en gt h o f ’ c ’ m ul ti pl e o f 2 14 c 0 = [ ] ; c 1 = [ ] ; 15 for n=1: ceil ( c l e n g t h / 2 )

//loop to

s e p a r a t e e ve n a nd odd p ow er s o f z 16

c 0 =[ c 0 c ( 2* n - 1) 0 ];

61

17 18 19 20 21 22 23

c 1 =[ c 1 c ( 2* n ) 0 ]; end E0 = poly ( c 0 , ’ z ’ , ’ c o e f f ’ ) / z ^ n / d e n 1 ; E1 = poly ( c 1 , ’ z ’ , ’ c o e f f ’ ) / z ^ ( n - 2 ) / d e n 1 ; disp ( ’ P o l y p h a s e C om p on en t s ’ ) disp ( E 0 , ’ E0( z ) ’ ) ; disp ( E 1 , ’ E1( z ) ’ ) ;

Scilab code Exa 11.5 Decimator implementation

1 / / E x am pl e 1 1 . 5 2 3 clc ; clear ; 4 5 function [ N , R ] = f u n c ( F s , F p , F t , F t i , d p , d s , M ) 6 dF=(Fs-Fp)/Ft;

// N o rm a li s ed t r a n s i t i o n b an dw id th 7

8

N = round ((-20* log10 ( sqroot ( d p * d s ) ) - 1 3 ) / ( 1 4 . 6 * d F ) ) // FIR F i l t e r l e n g th ; R=N*Fti/M;

// Number o f M u l t i p l i c a t i o n s p er s ec on d 9 10 11 12 13 14 15 16 17 18 19 20

endfunction

// S am p li ng r a t e o f i n p ut s i g n a l / / P as sb an d f r e q u e n c y / / S to pb an d f r e q u e n c y / / Pa ss ba nd r i p p l e / / St op ba nd r i p p l e / / D e c im a ti o n F a c t o r / / I n pu t s a mp l in g r a t e / / S i n g l e s t a g e i m p le m en t a ti o n

Ft=20000; Fp=40; Fs=50; dp=0.01; ds=0.002; M=100; Fti=Ft;

[N1,R1]=func(Fs,Fp,Ft,Fti,dp,ds,M);

62

21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40

/ /Two s t a g e i m p l e m en t a t i o n // S ta ge 1 F ( z ) w it h d e ci m at i on f a c t o r 50 / / P as sb an d f r e q u e n c y Fpf=Fp; Fsf=190; / / S to pb an d f r e q u e n c y // P a ss ba nd r i p p l e dpf=0.005; dsf=0.002; // S to pb an d r i p p l e / / D e c im a t io n F a c t or Mf=50; Fti=Ft; // I n pu t s a mp l in g r a t e

[N2f,R2f]=func(Fsf,Fpf,Ft,Fti,dpf,dsf,Mf);

/ / S t ag e 2 G( z ) w it h d e c i ma t i o n f a c t o r 2 / / P as sb an d f r e q u e n c y Fpg=50*Fp; Fsg=50*Fs; / / S to pb an d f r e q u e n c y // P a ss ba nd r i p p l e dpg=0.005; dsg=0.002; // S to pb an d r i p p l e / / D e c im a t io n F a c t or Mg=2; Fti=Ft/50; // I n pu t s a mp l in g r a t e

[N2g,R2g]=func(Fsg,Fpg,Ft,Fti,dpg,dsg,Mg); N2=N2f+50*N2g+2; // Total f i l t e r leng th / / T o t a l Number o f R2=R2f+R2g;

41

42 43 44 45

M u l t i p l i c a t i o n s p er s ec on d disp ( R 1 , ’ Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 1 , ’ FIR f i l t e r l e n g t h = ’ , ’ F or S i n g l e s t a g e implementation : ’ ); disp ( ’ F o r Two s t a g e i m p l e m e nt a t i o n : ’ ) ; disp ( R 2 f , ’ Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 f , ’ FIR f i l t e r l e n g t h = ’ , ’ For F( z ) : ’ ) ; disp ( R 2 g , ’ Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 g , ’ FIR f i l t e r l e n g t h = ’ , ’ For G( z ) : ’ ) ; disp ( R 2 , ’ T ot al Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 , ’ O v e r a l l FIR f i l t e r l e n g t h = ’ ) ;

Scilab code Exa 11.6 Decimator implementation

1 / / E x am pl e 1 1 . 6

63

2 3 clc ; clear ; 4 5 6 function [ N , R ] = f u n c ( F s , F p , F t , F t i , d p , d s , M ) dF=(Fs-Fp)/Ft; 7

// N o rm a li s ed t r a n s i t i o n b an dw id th 8

9

N = round ((-20* log10 ( sqroot ( d p * d s ) ) - 1 3 ) / ( 1 4 . 6 * d F ) ) ; // FIR F i l t e r l e n g th R=N*Fti/M;

// Number o f M u l t i p l i c a t i o n s p er s ec on d 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34

endfunction

// S am p li ng r a t e o f i n p ut s i g n a l / / P as sb an d f r e q u e n c y / / S to pb an d f r e q u e n c y // P a ss ba nd r i p p l e / / St op ba nd r i p p l e / / D e c im a t io n F a c t or / / I n pu t s a mp l in g r a t e / / S i n g l e s t a g e i m p le m en t a ti o n

Ft=10000; Fp=150; Fs=180; dp=0.002; ds=0.001; M=20; Fti=Ft;

[N1,R1]=func(Fs,Fp,Ft,Fti,dp,ds,M);

/ /Two s t a g e i m p l e m en t a t i o n // S ta ge 1 F ( z ) w it h d e ci m at i on f a c t o r 50 Fpf=Fp; / / P as sb an d f r e q u e n c y / / S to pb an d f r e q u e n c y Fsf=720; dpf=0.001; // P a ss ba nd r i p p l e // S to pb an d r i p p l e dsf=0.001; Mf=10; / / D e c im a t io n F a c t or // I n pu t s a mp l in g r a t e Fti=Ft;

[N2f,R2f]=func(Fsf,Fpf,Ft,Fti,dpf,dsf,Mf);

/ / S t ag e 2 G( z ) w it h d e c i ma t i o n f a c t o r 2 Fpg=10*Fp; / / P as sb an d f r e q u e n c y / / S to pb an d f r e q u e n c y Fsg=10*Fs; 64

35 36 37 38 39 40 41

// P a ss ba nd r i p p l e dpg=0.001; dsg=0.001; // S to pb an d r i p p l e / / D e c im a t io n F a c t or Mg=2; Fti=Ft/10; // I n pu t s a mp l in g r a t e [N2g,R2g]=func(Fsg,Fpg,Ft,Fti,dpg,dsg,Mg); N2=N2f+10*N2g+2; // Total f i l t e r leng th / / T o t a l Number o f R2=R2f+R2g; M u l t i p l i c a t i o n s p er s ec on d

42 43 disp ( R 1 , ’ Number o f

44 45 46 47

M u l t i p l i c a t i o n s p er s ec on d = ’ , N 1 , ’ FIR f i l t e r l e n g t h = ’ , ’ F or S i n g l e s t a g e implementation : ’ ); disp ( ’ F o r Two s t a g e i m p l e m e nt a t i o n : ’ ) ; disp ( R 2 f , ’ Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 f , ’ FIR f i l t e r l e n g t h = ’ , ’ For F( z ) : ’ ) ; disp ( R 2 g , ’ Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 g , ’ FIR f i l t e r l e n g t h = ’ , ’ For G( z ) : ’ ) ; disp ( R 2 , ’ T ot al Number o f M u l t i p l i c a t i o n s p er s ec on d = ’ , N 2 , ’ O v e r a l l FIR f i l t e r l e n g t h = ’ ) ;

65

Chapter 12 Spectral Estimation

Scilab code Exa 12.2 Power Spectrum

1 2 3 4 5 6 7 8

/ / E x am pl e 1 2 . 2

9 10 11 12 13 14 15 16 17 18

clc ; clear ; close ; N=8;n=0:N-1; f1=0.6;f2=0.62; x = cos ( 2 * % p i * f 1 * n ) + cos ( 2 * % p i * f 2 * n ) ; L1=8; for k = 0 : L 1 - 1 P 1 ( k + 1 ) = 1 / N * abs ( x * ( cos ( % p i * n * k / L 1 ) - % i * sin ( % p i * n * k/L1))’)^2 end L2=16; for k = 0 : L 2 - 1 P 2 ( k + 1 ) = 1 / N * abs ( x * ( cos ( % p i * n * k / L 2 ) - % i * sin ( % p i * n * k/L2))’)^2; end L3=32; for k = 0 : L 3 - 1 P 3 ( k + 1 ) = 1 / N * abs ( x * ( cos ( % p i * n * k / L 3 ) - % i * sin ( % p i * n * k/L3))’)^2; end subplot ( 3 1 1 ) ;

66

Digital Signal Processing_S. Salivahanan, A. Vallavaraj and C. Gnanapriya

Recommend Documents