`
`
`
`
`
`Micro Motion
`
`Micro Motion 1046
`
`1
`
`

`

`DIGITAL SIGNAL PROCESSING
`Principles, Algorithms, and Applications
`
`Second Edition
`
`2
`
`

`

`DIGITAL SIGNAL
`
`
`
`PROCESSING
`Principles, Algorithms,
`
`and Applications
`
`Second Edition
`
`John G. Proakis
`
`Northeastern University
`
`Dimitris G. Manolakis
`
`Worcester Polytechnic Institute
`
`Macmillan Publishing Company
`NEW YORK
`
`Maxwell Macmillan Canada
`TORONTO
`
`Maxwell Macmillan International
`NEWYORK OXFORD SINGAPORE SYDNEY
`
`
`
`3
`
`

`

`Editor(s): John Griffin
`Production Supervisor: Ron Harris
`Production Managcr: Sandra E. Moore
`Text Designer: Eileen Burke
`
`
`
`Copyright © l992 by Macmillan Publishing Company. a division of Macmillan. Inc.
`Printed in the United States of America
`
`All rights reserved. No part of this book may be reproduced or
`transmitted in any form or by any means. electronic or mechanical.
`including photocopying. recording. or any information storage and
`retrieval system. without permission in writing from the publisher.
`
`Earlier edition entitled Introduction to Digital Signal Processing.
`copyright © I988 by Macmillan Publishing Company.
`
`Macmillan Publishing Company is part
`of the Maxwell Communication Group of Companies.
`
`Macmillan Publishing Company
`Roe/Third Avenue. New York. New York l0022
`
`Maxwell Macmillan Canada. Inc.
`l200 Eglinton Avenue East
`Suite 200
`
`Don Mills. Ontario M3C 3Nl
`
`Library of Congress Cataloging in Publication Data
`
`Rev. ed. of: Introduction to digital signal processing. © I988,
`Includes bibliographical references and index.
`ISBN 0—02-3968l5—X
`I. Manolakis. Dimitris
`l. Signal processing—~Digital techniques.
`II. Proakis. John G. Introduction to digital signal processing.
`6.
`III. Title.
`TK5 l 02.5.P6773
`I992
`62l.382'2~—dc20
`'
`
`91— l5536
`CIP
`
`Printing:l2345678
`
`Year:234567890l
`
`4
`
`

`

`- PREFACE
`
`
`
`X.
`
`This book has resulted from our teaching of undergraduate and graduate—level courses
`in digital signal processing over the past several years. In this book we present the
`fundamentals of discrete-time signals, systems, and modern digital processing algo—
`rithms and applications for students in electrical engineering, computer engineering,
`and computer science. As written, the book is suitable for either a one-semester or a
`two-semester undergraduate—level course in discrete systems and digital signal proc-
`essing. It is also intended for use in a one—semester first—year graduate-level course in
`digital sig'r‘ial processing.
`A balanced coverage is provided of both theory and practical applications. Included
`are software implementations of digital filters and digital signal processing algorithms.
`Appropriately designed programs in FORTRAN are provided to aid the student in
`implementing and exercising each algorithm immediately after exposition. A large
`number of well-designed problems and a number of computer experiments are also
`provided to help the student in mastering the subject matter. A solutions manual is
`available for the benefit of the instructor and can be obtained from the publisher.
`It is assumed that the student in electrical and computer engineering has had under—
`graduate courses in advanced calculus (including ordinary differential equations), and
`linear systems for continuous—time signals, including an introduction to the Laplace
`transform. Although the Fourier series and Fourier transforms of periodic and ape—
`riodic signals are described in Chapter 3, we expect that many students may have had
`this material in a prior course.
`In Chapter 1 we describe the operations involved in the analog-to-digital conversion
`of an analog signal. The process of sampling a sinusoid is described in some detail
`and the problem of aliasing is explained. Signal quantization and digital—to—analog
`conversion are also described in general terms in Chapter 1, but the analysis is pre-
`sented in subsequent chapters.
`Chapter 2 is devoted entirely to the characterization and analysis of liner time—
`invariant (shift-invariant) discrete-time systems and discrete-time signals in the time
`domain. The convolution sum is derived and systems are categorized according to the
`duration of their impulse response as finite-duration impulse response (FIR) and
`infinite—duration impulse response (IIR). Linear time—invariant systems characterized
`by difference equations are presented and the solution of difference equations with
`initial conditions is obtained. The chapter concludes with a treatment of discrete—time
`correlation.
`
`Chapter 3 treats the analysis of signals and systems in the frequency domain. Fouri-
`er series and the Fourier transform are presented for both continuous—time and discrete—
`time signals. Linear time—invariant (LTI) discrete systems are characterized in the
`frequency domain by their frequency response function and their response to periodic
`
`5
`
`V
`
`
`
`5
`
`

`

`vi Preface
`
`and aperiodic signals is determined. The frequency response characteristics of ideal
`frequency selection filters are also described.
`The z-transform is introduced in Chapter 4. Both the bilateral and the unilateral
`z-transform are presented, and methods for determining the inverse z-transform are
`described. Use of the z-transform in the analysis of linear time—invariant systems is
`illustrated, and the important properties of systems, such as causality and stability, are
`related to z-domain characteristics.
`
`Chapter 5 is concerned with the analysis and design of discrete-time systems in the
`frequency domain. LTI systems are characterized either by their z-transform or by the
`Fourier transform, and their characteristics are related to the location of the poles and
`zeros of their system function. A number of important types of digital systems are
`described, including resonators, notch filters, comb filters, all—pass filters, and oscil-
`lators. The design of some simple FIR and IIR filters is also considered. In addition,
`the student is introduced to the concepts of minimum—phase, mixed—phase, and max-
`imum-phase systems and the problem of deconvolution.
`Chapter 6 is focused on the sampling of signals in the time and frequency domains.
`In this chapter we establish the sampling theorem for bandlimited continuous—time
`signals, both lowpass and bandpass. Also described in some detail are A/D and D/A
`conversion techniques, including oversampling A/D and D/A converters. The last
`major topic in this chapter is a discussion of sampling of the spectrum of a discrete-
`time signal. This discussion leads to the discrete Fourier transform (DFT) and its use
`in frequency analysis.
`Chapter 7 treats the realization of IIR and FIR systems. This treatment includes
`direct-form, cascade, parallel, lattice, and lattice-ladder realizations. The chapter in—
`cludes a treatment of state-space analysis and structures for discrete—time systems, and
`examines quantization effects in a digital implementation of FIR and IIR systems.
`Techniques for design of digital FIR and IIR filters are presented in Chapter 8. The
`design techniques include both direct design methods in discrete time and techniques
`that involve the conversion of analog filters into digital filters by various transfor-
`mations. Also treated in this chapter is the design of FIR and IIR filters by least-
`squares methods.
`“"-
`Chapter 9 treats the computation of the DFT. Radix—2, radix-4, and split-radix fast
`Fourier transform (FFT) algorithms are derived and applications of the FFT algorithms
`to the computation of convolution and correlation are also described. The Goertzel
`algorithm and the chirp-z transform are introduced as two methods by which the DFT
`can be computed via linear filtering.
`Chapter 10 provides an in-depth treatment of sampling-rate conversion and its
`applications to multirate digital signal processing. In addition to describing decimation
`and interpolation by integer factors, a method of sampling-rate conversion by an ar-
`bitrary factor is presented. Several applications to multirate signal processing are pre—
`sented, including the implementation of digital filters, subband coding of speech sig—
`nals, transmultiplexing, and oversampling A/D and D/A converters.
`Linear prediction and optimum linear (Wiener) filters are topics treated in Chapter
`11. Also included in this chapter are descriptions of the Levinson-Durbin algorithm
`and Schur algorithm for solving the normal equations, as well as AR lattice and ARMA
`lattice-ladder, filters.
`
`Power spectrum estimation is the main topic of Chapter 12. Our coverage of this
`subject includes a description of nonparametric and model-based (parametric) methods
`for power spectrum estimation. Also described are eigen-decomposition—based meth—
`ods, including MUSIC and ESPRIT.
`At Northeastern University, we have used the first six chapters of this book for a
`
`
`
`one—se
`
`For th
`It is a
`
`juncti(
`A (
`discre
`
`then 1:
`introd
`filter (
`In :
`
`ters pi
`may r
`9, fol]
`we hi;
`and pr
`Thl
`
`who 1
`We V
`for as
`H. Le
`
`sugge
`appre
`the m
`We a;
`radix
`Proak
`
`Ackl
`
`We i
`
`gestit
`Engii
`Wins
`State
`
`(D613
`East
`
`Com;
`ence
`
`B. M
`2500
`
`
`
`6
`
`

`

`
`
`i of ideal
`
`tnilateral
`form are
`{stems is
`
`Iiiity. are
`
`ns in the
`
`Jr by the
`oles and
`terns are
`1d oscil-
`1ddilion,
`nd max-
`
`lomains.
`)us-timc
`
`ind D/A
`Fhe last
`iiscrete-
`1 its use
`
`.ncludes
`
`lpter in—
`ms, and
`terns.
`r 8. The
`
`hniques
`ransfor-
`
`y least-
`
`.dix fast
`orithms
`ioerlzcl
`he DPT
`
`and its
`imalion
`I an ar—
`
`are pre-
`:ch sig-
`
`fhapter
`gorithm
`ARMA
`
`of this
`uethods
`l meth-
`
`kfora
`
`
`
`Preface vii
`
`one—semester (junior~level) course in discrete systems and digital signal processing.
`For this purpose we have introduced some simple filter design methods in Chapter 5.
`It is also possible to include the radix-2 decimation—in—time FFT algorithm in con-
`junction with the discussion in Chapter 6 on computation of the BET.
`A one-semester senior-level course for students who have had prior exposure to
`discrete systems may use the material in Chapters 1 through 4 for a quick review and
`then proceed to cover Chapters 5 through 9. In such a course the students may be
`introduced to the we of software packages that are now readily available for digital
`filter design.
`In a first-year graduate-level course in digital signal processing, the first five chap-
`ters provide the student with a good review of discrete-time systems. The instructor
`may move quickly through most of this material and then cover Chapters 6 through
`9, followed by either Chapters‘lO and 11 or by Chapters 11 and 12. In such a course
`we highly-recommend the use of available software packages for digital filter design
`and power spectrum estimation.
`The. authors are indebted to a number of graduate students and faculty colleagues
`who provided valuable suggestions and assistance in preparation of the manuscript.
`We wish to thank Mr. J. Lin for preparing the solutions manual, Dr. A. L. Kok
`for assisting in preparation of the figures, and Drs. J. Deller, V. Ingle, C. Keller,
`H. Lev-Ari. L. Merakos, P. Monticciolo, and M. Schetzen for their comments and
`suggestions on various drafts of the manuscript. We especially wish to express Our
`appreciation to Professors C. L. Nikias, H. J. Trussell, and S. G. Wilson, who reviewed
`the manuscript and provided many important suggestions that have been incorporated.
`We are also indebted to Dr. R. Price for informing us of the different types of split-
`radix FFT algorithms and related suggestions. Finally. we wish to thank Ms. Gloria
`Proakis for typing the entire manuscript.
`
`John G. Proakis
`
`Dimitris G. Manolakis
`
`Acknowledgments
`
`We would like to titanic the following reviewers, who made many valuable sug-
`gestions: Professor Michael Zoltowski (Dept. of Electrical Engineering; Electrical
`Engineering Bldg; Purdue University; West Lafayette,
`IN 47907); Professor
`Winser E. Alexander (Dept. of Electrical and Computer Engineering; North Carolina
`State University; Box 7911; Raleigh, NC 27695-7911); Professor John R. Deller
`(Dept. of Electrical Engineering; 260 Engineering Bldg; Michigan State University;
`East Lansing, MI 48824-1226); Professor Yorarn Bresler (Dept. of Electrical and
`Computer Engineering; University of Illinois#Urbana, Champlain; Coordinated Sci—
`ence Laboratory; 1101 West Springfield Avenue; Urbana, IL 61301); Professor Wasfy
`B. Mikhael (Dept. of Electrical Engineering; University of Central Florida; PO. Box
`25000; Orlando, FL 32316).
`
`7
`
`

`

`
`
`E12
`
`
`
`POWER SPECTRUM ESTIMATION
`
`In this chapter we are concerned with the estimation of the spectral CITHF'cl'CiL‘i‘lsii‘a\ of
`signals that are characterized as random processes. Many of the phenomena am occur
`in nature are best characterized statistically in terms of averages. For exan‘ag‘xte. me»
`teorological phenomena such as the fluctuations in air temperature and pressure are
`best characterized statistically as random processes. Thermal noise voltages '
`.
`
`in resistors and electronic devices are additional examples of physical signalu ,
`:1‘ 1‘ .
`well modeled as random processes.
`Due to the random fluctuations in such signals, we must adopt a statistical \iCVx»
`point, which deals with the average characteristics of random signals. In panicular.
`the autocorrelation function of a random process is the appropriate statisticaé art—rage
`that we will be concerned with for characterizing random signals in the time domain.
`and the Fourier transform of the autocorrelation function, which yields the power
`density spectrum, provides the transformation from the time domain to the frequent}
`domain.
`
`Power spectrum estimation methods have a relatively long history. For a histcrical
`perspective, the reader is referred to the paper by Robinson (1982) and the book by
`Marple (1987). Our treatment of this subject covers the classical power spectrum
`estimation methods based on the periodogram originally introduced by Scittrstcr
`(1898), and the modern model»based or parametric methods that originated with the
`work of Yule (1927). These methods were subsequently developed and applied by
`Walker (l93l), Bartlett (l948), Parzen (1957), Blackman and Tukey (1958;. Burg
`(1967), and others. We also describe the method of Capon (l969) and methods based
`on eigenanalysis of the data correlation matrix.
`
`12.1 ESTIMATION OF SPECTRA FROM FINITE-
`DURATION OBSERVATIONS OF SIGNALS
`
`The basic problem that we consider in this chapter is the estimation of the power
`density spectrum of a signal from observation of the signal over a finite time human}
`As we will see. the finite record length of the data sequence is a major limitation Ht?
`the quality of the power spectrum estimate. When dealing with signals that are statis—
`tically stationary, the longer the data record, the better the estimate that can be 2:»-
`tracted from the data. On the other hand, if the signal statistics are nonstationztr}
`‘
`
`cannot select an arbitrarily long data record to estimate the spectrum. In such a c ,
`the length of the data record that we select is determined by the rapidity of the rim:
`variations in the signal statistics. Ultimately. our goal is to select as short a data record
`864
`
`8
`
`8
`
`

`

`12.1 / Estimation of Spectra from Finite—Duration Observations at Signals
`
`865
`
`as possible that will allow us to resolve spectral cl’iaracteristics of different signal
`components, contained in the data record. that have closely spaced spectra.
`One of the problems that we encounter with classical power spectrum estimation
`methods based on a tinite~length data record is distortion of the spectrum that we are
`attempting to estimate. This problem occurs in both the computation of the spectrum
`for a deterministic signal and the estimation of the power spectrum of a random signal
`Since it
`is easier to observe the effect of the finite length ot‘ the data record en a
`deterministic signal, we treat this case lirst. Thereafter. we consider oniy random
`signals and the estimation of their power spectra
`
`12.1.1 Computation of the Energy Density Spectrum
`
`Let us consider the computation ot‘ the spectrum of a deterministic signal from a finite
`sequence of data. The sequence .\'(/1) is usually the result of sampling a continnons~
`time signal .\‘(,(I) at some uniform sampling rate R. Our objective is to obtain an
`estimate of the true spectrum from a tinite~duiation sequence Alli).
`Recall that if,\‘(t) is a finite-energy signal, that is,
`
`E 2 f T i_\‘u(l‘)iz dr <J 7:
`
`then its Fourier transform exists and is given as
`
`XuiFl : f v\'|l([)(f"_l2’77/‘t (It
`
`From Parseval’s theorem we have
`
`E 2 f 7
`
`l-\1,(f)l3 dz : f
`
`«:4
`
`EXHtFig'3 d1”
`
`.
`
`(Elli)
`
`The quantity [Xutfll2 represents the distribution of signal energy as Lt timction of
`frequency, and hence it is called the energy density spectrum of the signal. that is.
`
`SHtF) : txutrritl
`
`t t 2. l .2;
`
`as described in Chapter 3. Thus the total energy in the signai is simpty the integral of
`S, \tF) over all F lie. the total area under 5‘ \tFil.
`It is also interesting to note that 8‘ \(F) may be viewed as the Fourier transform of
`another function. R_\,_\(T). called the atItor'nrrylutin/iflirtation of the tinite~energy signal
`\,,tt). defined as
`
`RHm : J
`
`‘
`
`.rfftririr + 7) dz
`
`indeed, it easily follows that
`
`l— Rmtm 13“"? in : S‘HtFi : tier/iii
`
`tiltfit
`
`tutti
`
`stl that Rut?) and 5‘ \tF) are a Fourier transform pair.
`Now suppose that we compute the energy density spectrum of the sigtml tuifl from
`it“: samples taken at the rate R samples per second. To ensure that there is no spectrai
`:tiiasing resulting from the sampling process. the signal is assumed to he pretittcred.
`w that. for practical purposes. its bandwidth is limited to 8 hertz. Then the sampling
`tierntency FX is selected such that E 1?» 38.
`
`
`
`9
`
`

`

`866 Ch. 12 / Power Spectrum Estimation
`
`The sampled version of .ra(r) is a sequence .r(n), — 00 < n < 00, which has a Fourier
`transform (voltage spectrum)
`ac
`
`‘
`
`X(w) = 2 x(n>e~fw"
`
`or, equivalently,
`
`00
`
`X(f) = E x(n)e”jz“f”
`
`(12.1.5)
`
`Recall that X(f) may be expressed in terms of the voltage spectrum of the analog
`signal ,1'a(t) as
`
`FX<f> = F5 2 XG(F - kFS)
`
`S
`
`(12.1.6)
`
`where f = F/F3 is the normalized frequency variable.
`In the absence of aliasing, within the fundamental range IFI S FS/Z, we have
`
`XG) = FSXa(F),
`
`S
`
`[Fl 5 F‘s/2
`
`(12.1.7)
`
`Hence the voltage spectrum of the sampled signal is identical to the voltage spectrum
`of the analog signal. As a consequence, the energy density spectrum of the sampled
`signal is
`2
`
`S .
`
`F
`F
`. — = X —
`
`
`
`(F5)
`A)(F5)
`We may proceed further by noting that the autocorrelation of the sampled signal,
`which is defined as
`cc
`
`
`
`= F3 IXa(F)|2
`
`(12.1.8)
`
`rim-(k) = 2 x*(n)x(n ‘+ k)
`n=~0¢
`
`(12.1.9)
`
`has the Fourier transform (Wiener—Khintchine theorem)
`
`51mm = 2 r,.,.(k)e~12rkf
`
`(12.1.10)
`
`Hence the energy density spectrum may be obtained by Fourier transforming the
`autocorrelation of the sequence {x(n)}.
`The relations above lead us to distinguish between two distinct methods for come
`puting the energy density spectrum of a signal xa(t) from its samples x(n). One is the
`direct method, which involves computing the Fourier transform of {x(n)}, and then
`
`SKA-(f) = |X(f)|2x
`
`
`
`2 x(n)e ”12“”
`11:: ~96
`
`2
`
`
`
`(12.l.ll)
`
`The second approach is called the indirect method because it requires two steps. First,
`the autocorrelation ij\._\.(k) is computed from x(n) and then the Fourier tranSform of the
`autocortelation is computed as in (12.1.10) to obtain the energy density SpCCII‘UUL
`
`10
`
`10
`
`

`

`12.1 / Estimation of Spectra from FiniteDuration Observations of Signals 867
`In practice. however. only the finite-duration sequence .r(n), 0 _<. n _<. N — l. is
`available for computing the spectrum of the signal. In effect. limiting the duration of
`the sequence .r(n) to N points is equivalent to multiplying .r(n) by a rectangular win-
`dow. Thus we have
`
`0 E II S N — l
`.\‘(n).
`_
`otherwise
`0
`H”) 2 “WM!” :
`From our discussion of FIR filter design based on the use of windows to limit the
`duration of the impulse response, we recall that multiplication of two sequences is
`equivalent
`to convolution of their voltage spectra. Consequently,
`the frequency-
`domain relation corresponding to (l2. H2) is
`
`W“) = W) Wm
`(l2.l.l3)
`[V2
`4/3 X(or)W(f — or) dd
`Recall from our discussion in Section 6.4.5 that convolution of the window function
`W(f) with X(f) smooths the spectrum X(f). provided that the spectrum W(f) is rela-
`tively narrow compared to X(f). But this condition implies that the window Mn) be
`sufficiently long (i.e.. N must be sufficiently large) such that W(f) is narrow compared
`to X(f). Even if W(f) is narrow compared to X(f). the convolution of X(f) with the
`sidelobes of W(f) results in sidelobe energy in X(f). in frequency bands where the
`true signal spectrum X(f) = 0. This sidelobe energy is called lea/(age. The following
`example illustrates the leakage problem.
`
`EXAMPLE 12.1.1
`
`
`A signal with (voltage) spectrum
`
`otherwise
`X(f) _ {0.
`is convolved with the rectangular window of length N = 6i. Determine the spectrum
`ofX(f) given by (l2.l.l3).
`
`if! S 0-1
`
`_
`
`1’
`
`Solution: The spectral characteristic W(f} for the length N 2 6i rectangular window
`is illustrated in Fig. 8.2(b). Note that the width ofthe main lobe of the window function
`is At» = 4w/6l or Af = 2/6l. which is narrow compared to X(f).
`The convolution of X(f) with W(f) is illustrated in Fig. l2. i. We note that energy
`has leaked into the frequency band 0.l < if! _<. 0.5. where X(f) = l). A part of this
`is due to the width of the main lobe in W(f). which causes a broading or smearing of
`.X’tf) outside the range if! _<. 0.l. However. the sidelobe energy in X(f) is due to the
`presence of the sidelobes of W(f). which are convolved with X(f). The smearing of
`Kit") for if! > 0.l and the sidelobes in the range 0.l S if! S 0.5 constitute the leakage.
`Just as in the case of FIR filter design. we can reduce sidelobe leakage by selecting
`windows that have low sidelobes. This implies that the windows have a smooth time-
`domain cutoff instead of the abrupt cutoff in the rectangular window. Although such
`Window functions reduce sidelobe leakage. they result in an increase in smoothing or
`broadening of the spectral characteristic X(f). For example. the use of a Blackman
`
`11
`
`11
`
`

`

`
`
`
`
`
`
`868 Ch. 12 / Power Spectrum Estimation
`A
`
`10l—
`0
`
`LA.00"l'TlTr
`
` 0L—9
`
`
`
`
`
`Energydensity(dB)
`
`
`_100LJJJA_LJIJJ.LULL[L1L1JLJJ_LJHLJ_[LLJLIJJJJ_LJJJ_1J_LLJ_LJ_>f
`0
`.1
`.2
`.3
`.4
`.5,
`Frequency (cycles/sample)
`
`FIGURE 12.1 Spectrum obtained by convolving an M = 61 rectangular window with
`the ideal lowpass spectrum in Example 12.1.1.
`
`window of length N = 61 in Example 12.1.1 results in the spectral characteristic X(f)
`shown in Fig. 12.2. The sidelobe leakage has certame been reduced, but the spectral
`width has been increased by about 50%.
`The broadening of the spectrum being estimated due to windowing is particularly
`a problem when we wish to resolve signals with closely spaced frequency components.
`For example, the signal with spectral characteristic X(f) = Xl(f) + X2(f), as shown
`in Fig. 12.3, cannot be resolved as two separate signals unless the width of the window
`function is significantly narrower than the frequency separation Af. Thus we observe
`that using smooth time—domain windows reduces leakage at the expense of a decrease
`in frequency resolution.
`It is clear from the discussion above that the energy density spectrum of the win—
`dowed sequence {x(n)} is an approximation of the desired spectrum of the sequence
`
`
`
`
`
`Energydensity(dB)
`
`
`
`0
`
`.l
`
`.2
`
`.3
`
`.4
`
`.5
`
`Frequency (cycles/sample)
`FIGURE 12.2 Spectrum obtained by convolving an M = 61 Blackman window with
`the ideal lowpass spectrum in Example 12.1.1.
`
`12
`
`12
`
`

`

`12.1 / Estimation of Spectra from Finite~Duration Observations of Signals 869
`le)
`
`Xilf)
`
`X20’)
`
`0
`
`‘1 M 1._
`
`*f
`
`0.5
`
`FIGURE 12.3 Two narrowband signal spectra.
`
`{x(n)}. The spectral density obtained from {2201)} is
`
`Sad) = |X(f)|2 =
`
`N—l
`2 fink—1'2““
`n=O
`
`2
`
`
`
`
`
`(12.1.14)
`
`The spectrum given by (12.1.14) can be computed numerically at a set of N fre—
`quency points by means of the DFT. Thus
`N—l
`
`,
`
`m) = 2 awe—WW
`n=O
`
`(12.1.15)
`
`Then
`
`and hence
`
`-
`
`|X(k)|2 = Sfi(f)lf=k/N = S4171)
`
`k
`
`k
`
`
`512.2(3) = Z i(n)e‘-’2"""/N
`
`N—l
`n=0
`
`2
`
`
`
`(12.116)
`
`(12.1.17)
`
`which is a distorted version of the true spectrum Sxx(k/N).
`
`12.1.2 Estimation of the Autocorrelation and Power Spectrum
`of Random Signals: The Periodogram
`
`The finite—energy signals considered in the preceding section possess a Fourier trans—
`form and were characterized in the spectral domain by their energy density spectrum.
`On the other hand, the important class of signals characterized as stationary random
`processes do not have finite energy and hence do not possess a Fourier transform.
`Such signals have finite average power and hence are characterized by a power density
`Spectrum. lf x(t) is a stationary random process, its autocorrelation function is
`
`WAT) = E[x*(t)x(t + 7)]
`
`(12.1.18)
`
`where E H denotes the statistical average. Then via the Wiener~Khintchine theorem,
`the power density spectrum of the stationary random process is the Fourier transform
`of the autocorrelation function, that is,
`
`00
`
`I‘Xx(F) = f w Yxx(T)€ ‘12” dz
`
`(12.1.19)
`
`13
`
`13
`
`

`

`870 Ch. 12 / Power Spectrum Estimation
`
`In practice, we deal with a single realization of the random process from which We
`estimate the power spectrum of the process. We do not know the true autocorrelation
`function vmfi), and as a consequence, we cannot compute the Fourier transform in
`(12.1.19) to obtain I’M(F). On the other hand, from a single realization of the random
`process we can compute the time—average autocorrelation function
`
`T0
`1
`R.T=— x*th+Tdt
`..(>
`27b
`_TO
`(> <
`>
`
`(
`
`1213‘
`-0)
`
`where 2T0 is the observation interval. If the stationary random process is ergodic in
`the first and second moments (mean and autocorrelation function), then
`
`vxxfi) = lim R“,(T)
`T —<)oo
`To
`O
`1
`=lim—f x*txt+Tdt
`“6&20 _% (>(
`>
`
`(12.1.21)
`
`This relation justifies the use of the time—average autocorrelation Rn.("r) as an es-
`timate of the statistical autocorrelation function vxx("r). Furthermore, the Fourier trans—
`form of Rxx("r) provides an estimate PXX(F) of the power density spectrum, that is,
`TO
`I
`RXX(T)€
`‘jZ‘TF" d1"
`To
`1
`To
`U x*(t)x(t + T) 422-12?” d’T
`fig —70
`‘70
`To
`
`2T0
`I—TO X( )6
`
`Pxx(F)
`
`_TO
`
`.
`
`(12.1.22)
`
`= _
`
`
`
`2
`
`t —j21'rFI dt
`
`The actual power density spectrum is the expected value of PXX(F) in the limit as
`To —> 00, that is,
`’”
`
`FXX(F) = lim E[PXX(F)]
`T ~—>oo
`
`'2 F
`To
`1
`o
`
`=limE—f xte‘l"’dt
`TO—éoc
`2T0
`_TU ( )
`
`2
`
`
`
`(12.1.23)
`
`From (12.1.20) and (12.1.22) we again note the two possible approaches to c0m~
`puting PXX(F), the direct method as given by (12.1.22) or the indirect method, in which
`we obtain RH(T) first and then compute the Fourier transform.
`We shall consider the estimation of the power density spectrum from samples of a
`single realization of the random process. In particular, we assume that xa(t) is sampled
`at a rate F5 > 23, where B is the highest frequency contained in the power density
`spectrum of the random process. Thus we obtain a finite—duration sequence x(n),
`0 S n S N — 1, by sampling xa(t). From these samples we may compute the time
`average autocorrelation sequence
`1 N—Iml—l
`
`r;x(m) = N
`2 x*(n)x(n + m), m = 0,1,...,N — 1
`— m n:
`[HO
`Z x*(n)x(n + m), m = —1, —2,..., 1 — N
`N — lml n=|m1
`
`(12.1.24)
`
`r,’.,.(m) =
`
`and then compute the Fourier transform
`
`14
`
`14
`
`

`

`12.1 / Estimation of Spectra from Finite-Duration Observations of Signals
`N—l
`
`871
`
`P410”) = E r.L..—(m)e“/2"f’"
`m= ——N+l
`
`"
`
`(12.1.25)
`
`The normalization factor N —-
`value
`
`[ml in (12.1.24) results in an estimate with mean
`
`E[44171)]
`
`
`1
`N‘imi
`= m0")
`
`N—Iml—l
`2 E[x*(n)x(n + m)]
`"=0
`
`(12.1.26)
`
`where yxx(m) is the true (statistical) autocorrelation sequence of x(n). Hence r,’m(m) is
`an unbiased estimate of the autocorrelation function yxx(m). The variance of the es—
`timate r,’CX(m) is approximately
`N
`GO
`t
`[Iv-[ml]2 ”:23.
`['ymm)!2 + 'yfivm - m)'yxx(n + m)]
`-~——-——
`which is a result given by Jenkins and Watts (1968). Clearly,
`
`var [réc(m)] %
`
`(12.1.27)
`
`lim var [r;x(m)] = 0
`Nam
`
`(12.1.28)
`
`provided that
`
`Do
`
`2 Mn)? < oo
`n=—m
`
`Since E[r)'cx(m)] = yxx(m) and the variance of the estimate converges to zero as
`N -—> 00, the estimate r;X(m) is said to be consistent.
`For large values of the lag parameter m, the estimate ryéx(m) given by ( 12.1.24) has
`a large variance, especially as m approaches N. This is due to the fact that fewer data
`points enter into the estimate for large lags. As an alternative to (12.1.24) we may use
`the estimate
`
`1 N—Iml—l
`rxx(m) = N E x*(n)x(n + m),
`1 ~50
`'
`rxx(m) = N
`x*(n)x(n + m),
`
`n=lml
`
`0 S m S N —-
`
`1
`
`m = —-l, —-2, ..,1 —- N
`
`(12.1.29)
`
`which has a bias of lmlyxx(m)/N, since its mean value is
`
`1 N“ 1m! “ 1
`E[rxx(m)] = N 20 E[x*(n)x(rt + m)]
`N _ imi
`imi
`= T MW) = < - 7V")VXX(M)
`However, this estimate has a smaller variance, given approximately as
`oo
`
`(12.1.30)
`
`var[rxx(m)] ~ _. 2 [Film-(””2 + FY34” _. mnyxUl + ’71)]
`N ,1: ——00
`
`(120-131)
`
`We observe that rxx(m) is asymptotically unbiased, that is,
`
`lim E[r_”(m)] = yxx(m)
`Nam
`
`(12.1.32)
`
`15
`
`15
`
`

`

`872 Ch. 12 / Power Spectrum Estimation
`
`and its variance converges to zero as N ——> 00. Therefore, the estimate 1““..(m) is also a
`consistent estimate of 'yxx(m).
`We shall use the estimate (”(m) given by (12.1.29) in our treatment of power
`spectrum estimation. The corresponding estimate of the power density spectrum is
`N—l
`
`P,.,.(f) =
`
`2
`m=—(N— 1)
`
`rm(m)e_j2"f’"
`
`(12.1.33)
`
`If we substitute for r“.(m) from (12.1.29) into (12.1.33), the estimate P“.(f) may also
`be expressed as
`
`1
`
`P“.(f) = N
`
`N—i
`
`2 X(I’I)€ —j2-rrfn
`n=0
`
`2
`
`
`
`1
`= N 1X(f)l2
`
`(12.1.34)
`
`where X(f) is the Fourier transform of the sample sequence x(n). This well—known
`form of the power density spectrum estimate is called the periodogram. It was origi—
`nally introduced by Schuster (1898) to detect and measure “hidden periodicities” in
`data.
`
`From (12.1.33), the average value of the periodograrn estimate Pxx(f) is
`N— 1
`N— 1
`
`E[P,\x(f)] : E[ 2
`ElPX..(f)l
`E
`
`m= —(N— 1)
`
`=—(N—1)
`N—l
`
`rxx(m)€—j2flfm:l : 3 2
`(I ~ %>vx(m)e"2“fm
`
`m=—(N—I)
`
`E[7'xx(m)]e_j2-rrfm
`
`(12.1.35)
`
`The interpretation that we give to (12.1.35) is that the mean of the estimated spectrum
`is the Fourier transform of the windowed autocorrelation function
`
`m
`
`fixgm) 2 <1 —— U>‘3/J,.,.(171)
`
`N
`
`(12.1.36)
`
`where the window function is the (triangular) Bartlett window. Hence the mean of the
`estimated spectrum is
`
`ElPx1-(f)l
`
`2 1...<m)e‘-"2"fm
`"1: —OO
`1/2
`
`{—1/2 I‘”((x)WB(f -— (x) dOL
`
`(12.1.37)
`
`where WB(f) is the spectral characteristic of the Bartlett window. The relation (1 1.1.37)
`illustrates that the mean of the estimated spectrum is the convolution of the true power
`density spectrum THU") with the Fourier transform WB(f) of the Bartlett window.
`Consequently, the mean of the estimated spectrum is a smoothed version of the true
`spectrum and suffers from the same spectral leakage problems which are due to the
`finite number of data points.
`We observe that the estimated spectrum is asymptotically unbiased, that is,
`N— 1
`m
`
`Um E[ E
`
`N—>3C
`
`m=—(N—1)
`
`r...<m)ewe] = 2 v...<m>eW!“ = F....<f>
`
`m=—:>c
`
`However, in general, the variance of the estimate PXX(f) does not decay to zero as
`N ——> 00. For example, when the data sequence is a Gaussian random process, the
`variance is easily shown to be (see Problem 12.4)
`
`16
`
`16
`
`

`

`12.1 / Estimation of Spectra from Finite—Duration Observations of Signals 873
`2
`Txx(f)[l + (N sin 2m,
`(
`-
`-
`2 W 12 1 38
`
`)
`
`VarlPXX(f)l
`
`z:
`
`which, in the limit as N —~> 00, becomes
`
`lim var[PXX(f)] = FEXQ‘)
`New
`
`(12.1.39)
`
`Hence we conclude that the periodogram is not a consistent estimate offhe true power
`density spectrum (i.e., it does not converge to the true power density spectrum).
`In summary, the estimated autocorrelation rxx(m) is a consistent estimate of the true
`autocorrelation function yxx(m). HOWever, its Fourier transform Pxx(f), the periodo—
`gram, is not a consistent estimate of the true power density spectrum. We observed
`that PXX(f) is an asymptotically unbiased estimate of I‘xx(f), but for a finite—duration
`sequence, the mean value of Pxx(f) contains a bias, which from (12.1.37) is evident
`as a distortion of the true power density spectrum. Thus the estimated spectrum suffers
`from the smoothing effects and the leakage embodied in the Bartlett window. The
`smoothing and leakage ultimately limit our ability to resolve closely spaced spectra.
`The problems of leakage and frequency resolution that we have described above,
`as well as the problem that the periodogram is not a consistent estimate of the power
`spectrum, provide the motivation for the power spectrum estimation methods de-
`scribed in Sections 12.2, 12.3, and 12.4. The methods described in Section 12.2 are
`classical nonparametric methods, which make no assumptions about the data se—
`quence.