what is spectral analysis of signals

cyclic spectral analysis of rolling-element bearing signals facts and fictions and spectral analysis of random signals
JuliyaMadenta Profile Pic
Published Date:15-07-2017
Your Website URL(Optional)
\sm2" i i 2004/2/22 page i i i SPECTRAL ANALYSIS OF SIGNALS Petre Stoica and Randolph Moses PRENTICE HALL, Upper Saddle River, New Jersey 07458 i i i i\sm2" i i 2004/2/22 page 1 i i C H A P T E R 1 Basic Concepts 1.1 INTRODUCTION The essence of the spectral estimation problem is captured by the following informal formulation. From a nite record of a stationary data sequence, estimate how (1.1.1) the total power is distributed over frequency. Spectral analysis nds applications in many diverse elds. In vibration monitoring, the spectral content of measured signals give information on the wear and other characteristics of mechanical parts under study. In economics, meteorology, astron- omy and several other elds, the spectral analysis may reveal \hidden periodicities" in the studied data, which are to be associated with cyclic behavior or recurring processes. In speech analysis, spectral models of voice signals are useful in better understanding the speech production process, and in addition can be used for both speech synthesis (or compression) and speech recognition. In radar and sonar systems, the spectral contents of the received signals provide information on the location of the sources (or targets) situated in the eld of view. In medicine, spectral analysis of various signals measured from a patient, such as electrocardio- gram (ECG) or electroencephalogram (EEG) signals, can provide useful material for diagnosis. In seismology, the spectral analysis of the signals recorded prior to and during a seismic event (such as a volcano eruption or an earthquake) gives useful information on the ground movement associated with such events and may help in predicting them. Seismic spectral estimation is also used to predict sub- surface geologic structure in gas and oil exploration. In control systems, there is a resurging interest in spectral analysis methods as a means of characterizing the dynamical behavior of a given system, and ultimately synthesizing a controller for that system. The previous and other applications of spectral analysis are reviewed in Kay 1988; Marple 1987; Bloomfield 1976; Bracewell 1986; Haykin 1991; Haykin 1995; Hayes III 1996; Koopmans 1974; Priestley 1981; Perci- val and Walden 1993; Porat 1994; Scharf 1991; Therrien 1992; Proakis, Rader, Ling, and Nikias 1992. The textbook Marple 1987 also contains a wellwritten historical perspective on spectral estimation which is worth read- ing. Many of the classical articles on spectral analysis, both applicationdriven and theoretical, are reprinted in Childers 1978; Kesler 1986; these excellent collections of reprints are well worth consulting. There are two broad approaches to spectral analysis. One of these derives its basic idea directly from de nition (1.1.1): the studied signal is applied to a band- pass lter with a narrow bandwidth, which is swept through the frequency band of 1 i i i i\sm2" i i 2004/2/22 page 2 i i 2 Chapter 1 Basic Concepts interest, and the lter output power divided by the lter bandwidth is used as a measure of the spectral content of the input to the lter. This is essentially what the classical (or nonparametric) methods of spectral analysis do. These methods are described in Chapters 2 and 5 of this text (the fact that the methods of Chapter 2 can be given the above lter bank interpretation is made clear in Chapter 5). The second approach to spectral estimation, called the parametric approach, is to postu- late a model for the data, which provides a means of parameterizing the spectrum, and to thereby reduce the spectral estimation problem to that of estimating the parameters in the assumed model. The parametric approach to spectral analysis is treated in Chapters 3, 4 and 6. Parametric methods may o er more accurate spectral estimates than the nonparametric ones in the cases where the data indeed satisfy the model assumed by the former methods. However, in the more likely case that the data do not satisfy the assumed models, the nonparametric meth- ods may outperform the parametric ones owing to the sensitivity of the latter to model misspeci cations. This observation has motivated renewed interest in the nonparametric approach to spectral estimation. Many realworld signals can be characterized as being random (from the ob- server's viewpoint). Brie y speaking, this means that the variation of such a signal outside the observed interval cannot be determined exactly but only speci ed in statistical terms of averages. In this text, we will be concerned with estimating the spectral characteristics of random signals. In spite of this fact, we nd it useful to start the discussion by considering the spectral analysis of deterministic signals (which we do in the rst section of this chapter). Throughout this work, we consider discrete signals (or data sequences). Such signals are most commonly obtained by the temporal or spatial sampling of a continuous (in time or space) signal. The main motivation for focusing on discrete signals lies in the fact that spectral analy- sis is most often performed by a digital computer or by digital circuitry. Chapters 2 to 5 of this text deal with discretetime signals, while Chapter 6 considers the case of discretespace data sequences. In the interest of notational simplicity, the discretetime variablet, as used in this text, is assumed to be measured in units of sampling interval. A similar conven- tion is adopted for spatial signals, whenever the sampling is uniform. Accordingly, the units of frequency are cycles per sampling interval. The signals dealt with in the text are complexvalued. Complexvalued data may appear in signal processing and spectral estimation applications, for instance, as a result of a \complex demodulation" process (this is explained in detail in Chapter 6). It should be noted that the treatment of complexvalued signals is not always more general or more dicult than the analysis of corresponding realvalued signals. A typical example which illustrates this claim is the case of sinusoidal signals considered in Chapter 4. A realvalued sinusoidal signal, cos(t +'), can be rewritten as a linear combination of two complexvalued sinusoidal signals, i( t+' ) i( t+' ) 1 1 2 2 e + e , whose parameters are constrained as follows: = 1 2 1 p = =2, ' = ' = ' and = = . Here i = 1. The fact 2 1 2 1 2 that we need to consider two constrained complex sine waves to treat the case of one unconstrained real sine wave shows that the realvalued case of sinusoidal signals can actually be considered to be more complicated than the complexvalued case Fortunately, it appears that the latter case is encountered more frequently i i i i\sm2" i i 2004/2/22 page 3 i i Section 1.2 Energy Spectral Density of Deterministic Signals 3 in applications, where often both the inphase and quadrature components of the studied signal are available. (For more details and explanations on this aspect, see Chapter 6's introductory section.) 1.2 ENERGY SPECTRAL DENSITY OF DETERMINISTIC SIGNALS Let fy(t);t = 0; 1; 2;:::g denote a deterministic discretetime data sequence. Most commonly, fy(t)g is obtained by sampling a continuoustime signal. For notational convenience, the time index t is expressed in units of sampling interval; that is, y(t) = y (t T ), where y () is the continuous time signal and T is the c s c s sampling time interval. Assume that fy(t)g has nite energy, which means that 1 X 2 jy(t)j 1 (1.2.1) t=1 Then, under some additional regularity conditions, the sequence fy(t)g possesses a discretetime Fourier transform (DTFT) de ned as 1 X it Y () = y(t)e (DTFT) (1.2.2) t=1 i In this text we use the symbol Y (), in lieu of the more cumbersome Y (e ), to denote the DTFT. This notational convention is commented on a bit later, following equation (1.4.6). The corresponding inverse DTFT is then Z  1 it y(t) = Y ()e d (Inverse DTFT) (1.2.3) 2  which can be veri ed by substituting (1.2.3) into (1.2.2). The (angular) frequency is measured in radians per sampling interval. The conversion from to the physical frequency variable  ==T rad/sec can be done in a straightforward manner, as s described in Exercise 1.1. Let 2 (1.2.4) S() = jY ()j (Energy Spectral Density) A straightforward calculation gives Z Z 1 1   X X 1 1  i(ts) S()d = y(t)y (s)e d 2 2   t=1s=1  Z  1 1  X X 1  i(ts) = y(t)y (s) e d 2  t=1s=1 1 X 2 = jy(t)j (1.2.5) t=1 i i i i\sm2" i i 2004/2/22 page 4 i i 4 Chapter 1 Basic Concepts R  1 i(ts) To obtain the last equality in (1.2.5) we have used the fact that e d = 2    (the Kronecker delta). The symbol () will be used in this text to denote the t;s complexconjugate of a scalar variable or the conjugate transpose of a vector or matrix. Equation (1.2.5) can be restated as Z 1  X 1 2 jy(t)j = S()d (1.2.6) 2  t=1 This equality is called Parseval's theorem. It shows that S() represents the dis- tribution of sequence energy as a function of frequency. For this reason, S() is called the energy spectral density. The previous interpretation ofS() also comes up in the following way. Equa- tion (1.2.3) represents the sequence fy(t)g as a weighted \sum" (actually, an inte- 1 it 1 p p gral) of orthonormal sequences f e g ( 2 ;), with weighting Y (). 2 2 1 p Hence, jY ()j \measures" the \length" of the projection of fy(t)g on each of 2 1 p these basis sequences. In loose terms, therefore, jY ()j shows how much (or 2 how little) of the sequence fy(t)g can be \explained" by the orthonormal sequence 1 it p f e g for some given value of . 2 De ne 1 X  (k) = y(t)y (t k) (1.2.7) t=1 It is readily veri ed that 1 1 1 X X X ik  it i(tk) (k)e = y(t)y (t k)e e t=1 k=1 k=1 " "  1 1 X X it is = y(t)e y(s)e t=1 s=1 = S() (1.2.8) which shows that S() can be obtained as the DTFT of the \autocorrelation" (1.2.7) of the niteenergy sequence fy(t)g. The above de nitions can be extended in a rather straightforward manner to the case of random signals treated throughout the remaining text. In fact, the only purpose for discussing the deterministic case in this section was to provide some motivation for the analogous de nitions in the random case. As such, the discussion in this section has been kept brief. More insights into the meaning and properties of the previous de nitions are provided by the detailed treatment of the random case in the following sections. 1.3 POWER SPECTRAL DENSITY OF RANDOM SIGNALS Most of the signals encountered in applications are such that their variation in the future cannot be known exactly. It is only possible to make probabilistic statements about that variation. The mathematical device to describe such a signal is that i i i i\sm2" i i 2004/2/22 page 5 i i Section 1.3 Power Spectral Density of Random Signals 5 of a random sequence which consists of an ensemble of possible realizations, each of which has some associated probability of occurrence. Of course, from the whole ensemble of realizations, the experimenter can usually observe only one realization of the signal, and then it might be thought that the deterministic de nitions of the previous section could be carried over unchanged to the present case. However, this is not possible because the realizations of a random signal, viewed as discretetime sequences, do not have nite energy, and hence do not possess DTFTs. A random signal usually has nite average power and, therefore, can be characterized by an average power spectral density. For simplicity reasons, in what follows we will use the name power spectral density (PSD) for that quantity. The discretetime signal fy(t);t = 0; 1; 2;:::g is assumed to be a sequence of random variables with zero mean: E fy(t)g = 0 for all t (1.3.1) Hereafter,E fg denotes the expectation operator (which averages over the ensemble of realizations). The autocovariance sequence (ACS) or covariance function of y(t) is de ned as  r(k) =E fy(t)y (t k)g (1.3.2) and it is assumed to depend only on the lag between the two samples averaged. The two assumptions (1.3.1) and (1.3.2) imply that fy(t)g is a secondorder stationary sequence. When it is required to distinguish between the autocovariance sequences of several signals, a lower index will be used to indicate the signal associated with a given covariance lag, such as r (k). y The autocovariance sequence r(k) enjoys some simple but useful properties:  r(k) =r (k) (1.3.3) and r(0)  jr(k)j for all k (1.3.4) The equality (1.3.3) directly follows from de nition (1.3.2) and the stationarity assumption, while (1.3.4) is a consequence of the fact that the covariance matrix of fy(t)g, de ned as follows 2 3   r(0) r (1) ::: r (m 1) 6 . 7 . . . 6 7 . r(1) r(0) . 6 7 R = m 6 7 . . . . . .  4 5 . . . r (1) r(m 1) ::: r(1) r(0) 8 9 2 3  y (t 1) 6 7 . = . 6 7 . 6 7 = E y(t 1):::y(t m) (1.3.5) 6 7 . . 4 5 . : ;  y (t m) i i i i\sm2" i i 2004/2/22 page 6 i i 6 Chapter 1 Basic Concepts is positive semide nite for all m. Recall that a Hermitian matrix M is positive  semide nite if a Ma  0 for every vector a (see Section A.5 for details). Since 82 3 9  y (t 1) 6 7 = . . 6 7   . a R a = a E y(t 1):::y(t m) a 6 7 m 4 5 : ;  y (t m)   2 = E fz (t)z(t)g =E jz(t)j  0 (1.3.6) where z(t) = y(t 1):::y(t m)a we see that R is indeed positive semide nite for every m. Hence, (1.3.4) fol- m lows from the properties of positive semide nite matrices (see De nition D11 in Appendix A and Exercise 1.5). 1.3.1 First De nition of Power Spectral Density The PSD is de ned as the DTFT of the covariance sequence: 1 X ik () = r(k)e (Power Spectral Density) (1.3.7) k=1 Note that the previous de nition (1.3.7) of () is similar to the de nition (1.2.8) in the deterministic case. The inverse transform, which recovers fr(k)g from given (), is Z  1 ik r(k) = ()e d (1.3.8) 2  We readily verify that Z  Z  1   X 1 1 ik i(kp) ()e d = r(p) e d =r(k) 2 2   p=1 which proves that (1.3.8) is the inverse transform for (1.3.7). Note that to obtain the rst equality above, the order of integration and summation has been inverted, which is possible under weak conditions (such as under the requirement that () is square integrable; see Chapter 4 in Priestley 1981 for a detailed discussion on this aspect). From (1.3.8), we obtain Z  1 r(0) = ()d (1.3.9) 2   2 Sincer(0) =E jy(t)j measures the (average) power of fy(t)g, the equality (1.3.9) shows that() can indeed be named PSD, as it represents the distribution of the i i i i\sm2" i i 2004/2/22 page 7 i i Section 1.3 Power Spectral Density of Random Signals 7 (average) signal power over frequencies. Put another way, it follows from (1.3.9) that()d=2 is the in nitesimal power in the band (d=2,+d=2), and the total power in the signal is obtained by integrating these in nitesimal contributions. Additional motivation for calling () a PSD is provided by the second de nition of (), given next, which resembles the usual de nition (1.2.2), (1.2.4) in the deterministic case. 1.3.2 Second De nition of Power Spectral Density The second de nition of () is: 8 9 2 N = X 1 it () = lim E y(t)e (1.3.10) N1 : ; N t=1 This de nition is equivalent to (1.3.7) under the mild assumption that the covari- ance sequence fr(k)g decays suciently rapidly, so that N X 1 lim jkjjr(k)j = 0 (1.3.11) N1N k=N The equivalence of (1.3.7) and (1.3.10) can be veri ed as follows: 8 9 2 N N N = X XX 1 1 it  i(ts) lim E y(t)e = lim E fy(t)y (s)ge N1 : ; N1 N N t=1 t=1s=1 N1 X 1 i = lim (N jj)r()e N1N =(N1) 1 N1 X X 1 i i = r()e lim jjr()e N1 N =1 =(N1) = () The second equality is proven in Exercise 1.6, and we used (1.3.11) in the last equality. The above de nition of () resembles the de nition (1.2.4) of energy spec- tral density in the deterministic case. The main di erence between (1.2.4) and (1.3.10) consists of the appearance of the expectation operator in (1.3.10) and the normalization by 1=N; the fact that the \discretetime" variable in (1.3.10) runs over positive integers only is just for convenience and does not constitute an es- sential di erence, compared to (1.2.2). In spite of these di erences, the analogy between the deterministic formula (1.2.4) and (1.3.10) provides further motivation for calling () a PSD. The alternative de nition (1.3.10) will also be quite useful when discussing the problem of estimating the PSD by nonparametric techniques in Chapters 2 and 5. i i i i\sm2" i i 2004/2/22 page 8 i i 8 Chapter 1 Basic Concepts We can see from either of these de nitions that () is a periodic function, with the period equal to 2. Hence, () is completely described by its variation in the interval 2 ; (radians per sampling interval) (1.3.12) Alternatively, the PSD can be viewed as a function of the frequency f = (cycles per sampling interval) (1.3.13) 2 which, according to (1.3.12), can be considered to take values in the interval f 2 1=2; 1=2 (1.3.14) We will generally write the PSD as a function of whenever possible, since this will simplify the notation. As already mentioned, the discretetime sequence fy(t)g is most commonly derived by sampling a continuoustime signal. To avoid aliasing e ects which might be incurred by the sampling process, the continuoustime signal should be (at least, approximately) bandlimited in the frequency domain. To ensure this, it may be necessary to lowpass lter the continuoustime signal before sampling. Let F denote the largest (\signi cant") frequency component in the spectrum of the 0 (possibly ltered) continuous signal, and letF be the sampling frequency. Then it s follows from Shannon's sampling theorem that the continuoustime signal can be \exactly" reconstructed from its samples fy(t)g, provided that F  2F (1.3.15) s 0 In particular, \no" frequency aliasing will occur when (1.3.15) holds (see, e.g., Oppenheim and Schafer 1989). Since the frequency variable, F, associated with the continuoustime signal, is related to f by the equation F =f F (1.3.16) s it follows that the interval of F corresponding to (1.3.14) is   F F s s F 2 ; (cycles/sec) (1.3.17) 2 2 which is quite natural in view of (1.3.15). 1.4 PROPERTIES OF POWER SPECTRAL DENSITIES Since() is a power density, it should be realvalued and nonnegative. That this is indeed the case is readily seen from de nition (1.3.10) of (). Hence, ()  0 for all (1.4.1) i i i i\sm2" i i 2004/2/22 page 9 i i Section 1.4 Properties of Power Spectral Densities 9 From (1.3.3) and (1.3.7), we obtain 1 X ik () =r(0) + 2 Refr(k)e g k=1 where Refg denotes the real part of the bracketed quantity. If y(t), and hence r(k), is real valued then it follows that 1 X () =r(0) + 2 r(k) cos(k) (1.4.2) k=1 which shows that () is an even function in such a case. In the case of complex valued signals, however, () is not necessarily symmetric about the = 0 axis. Thus: For realvalued signals: () =(); 2 ; (1.4.3) For complexvalued signals: in general () = 6 (); 2 ; Remark: The reader might wonder why we did not de ne the ACS as  c(k) =E fy(t)y (t +k)g Comparing with the ACS fr(k)g used in this text, as de ned in (1.3.2), we obtain c(k) =r(k). Consequently, the PSD associated with fc(k)g is related to the PSD corresponding to fr(k)g (see (1.3.7)) via: 1 1 X X ik ik () , c(k)e = r(k)e =() k=1 k=1 It may seem arbitrary as to which de nition of the ACS (and corresponding de ni- tion of PSD) we choose. In fact, from a mathematical standpoint we can use either de nition of the ACS, but the ACS de nition r(k) is preferred from a practical standpoint as we now explain. First, we should stress that the PSD describes the spectral content of the ACS, as seen from equation (1.3.7). The PSD () is sometimes perceived as showing the (in nitesimal) power at frequency in the signal itself, but that is not strictly true. If the PSD represented the power in the signal itself, then we should have had () = () because the signal's spectral content should not depend on the ACS de nition. However, as shown above, in the general complex case () = () = 6 (), which means that the signal power interpretation of the PSD is not (always) correct. Indeed, the PSD () \measures" the power at frequency in the signal's ACS. i i i i\sm2" i i 2004/2/22 page 10 i i 10 Chapter 1 Basic Concepts e(t) y(t) - - H(z) 2  ()  () = jH()j  () e y e Figure 1.1. Relationship between the PSDs of the input and output of a linear system. On the other hand, our motivation for considering spectral analysis is to characterize the average power at frequency in the signal, as given by the sec- ond de nition of the PSD in equation (1.3.10). If c(k) is used as the ACS, its corresponding second de nition of the PSD is 8 9 2 N = X 1 +it () = lim E y(t)e N1 : ; N t=1 which is the average power of y(t) at frequency . Clearly, the second PSD de nition corresponding to r(k) aligns with this average power motivation, while the one forc(k) does not; it is for this reason that we use the de nitionr(k) for the ACS.  Next, we present a useful result which concerns the transfer of PSD through an asymptotically stable linear system. Let 1 X k H(z) = h z (1.4.4) k k=1 1 denote an asymptotically stable linear timeinvariant system. The symbol z 1 denotes the unit delay operator de ned by z y(t) = y(t 1). Also, let e(t) be the stationary input to the system andy(t) the corresponding output, as shown in Figure 1.1. Then fy(t)g and fe(t)g are related via the convolution sum 1 X y(t) =H(z)e(t) = h e(t k) (1.4.5) k k=1 The transfer function of this lter is 1 X ik H() = h e (1.4.6) k k=1 Throughout the text, we will follow the convention of writing H(z) for the con- volution operator of a linear system and its corresponding Z-transform, and H() i i i i\sm2" i i 2004/2/22 page 11 i i Section 1.4 Properties of Power Spectral Densities 11 for its transfer function. We obtain the transfer function H() from H(z) by the i substitution z =e : H() =H(z) i z=e i While we recognize the slight abuse of notation in writing H() instead of H(e ) and in usingz as both an operator and a complex variable, we prefer the simplicity of notation it a ords. From (1.4.5), we obtain 1 1 X X   r (k) = h h E fe(t p)e (t m k)g y p m p=1m=1 1 1 X X  = h h r (m +k p) (1.4.7) p e m p=1m=1 Inserting (1.4.7) in (1.3.7) gives 1 1 1 X X X  i(k+mp) im ip  () = h h r (m +k p)e e e y p e m k=1p=1m=1 " " " 1 1 1 X X X ip  im i = h e h e r ()e p e m p=1 m=1 =1 2 = jH()j  () (1.4.8) e From (1.4.8), we get the following important formula 2 (1.4.9)  () = jH()j  () y e that will be much used in the following chapters. Finally, we derive a property that will be of use in Chapter 5. Let the signals y(t) and x(t) be related by i t 0 y(t) =e x(t) (1.4.10) for some given . Then, it holds that 0  () = ( ) (1.4.11) y x 0 i t 0 In other words, multiplication by e of a temporal sequence shifts its spectral density by the angular frequency . Owing to this interpretation, the process 0 of constructing y(t) as in (1.4.10) is called complex (de)modulation. The proof of (1.4.11) is immediate: since (1.4.10) implies that i k 0 r (k) =e r (k) (1.4.12) y x we obtain 1 X i(0)k  () = r (k)e = ( ) (1.4.13) y x x 0 k=1 which is the desired result. i i i i\sm2" i i 2004/2/22 page 12 i i 12 Chapter 1 Basic Concepts 1.5 THE SPECTRAL ESTIMATION PROBLEM The spectral estimation problem can now be stated more formally as follows. From a nitelength record fy(1);:::;y(N)g of a secondorder (1.5.1) stationary random process, determine an estimate () of its power spectral density (), for 2 ; It would, of course, be desirable that () is as close to () as possible. As we shall see, the main limitation on the quality of most PSD estimates is due to the quite small number of data samples usually available for processing. Note that N will be used throughout this text to denote the number of points of the available data sequence. In some applications, N is small since the cost of obtaining large amounts of data is prohibitive. Most commonly, the value of N is limited by the fact that the signal under study can be considered secondorder stationary only over short observation intervals. As already mentioned in the introductory part of this chapter, there are two main approaches to the PSD estimation problem. The nonparametric approach, discussed in Chapters 2 and 5, proceeds to estimate the PSD by relying essentially only on the basic de nitions (1.3.7) and (1.3.10) and on some properties that di- rectly follow from these de nitions. In particular, these methods do not make any assumption on the functional form of(). This is in contrast with the parametric approach, discussed in Chapters 3, 4 and 6. That approach makes assumptions on the signal under study, which lead to a parameterized functional form of the PSD, and then proceeds by estimating the parameters in the PSD model. The paramet- ric approach can thus be used only when there is enough information about the studied signal, that allows formulation of a model. Otherwise the nonparametric approach should be used. Interestingly enough, the nonparametric methods are close competitors to the parametric ones, even when the model form assumed by the latter is a reasonable description of reality. 1.6 COMPLEMENTS 1.6.1 Coherency Spectrum Let  () yu C () = (1.6.1) yu 1=2  () () yy uu denote the socalled complexcoherency of the stationary signals y(t) and u(t). In the de nition above,  () is the crossspectrum of the two signals, and  () yu yy and () are their respective PSDs. (We implicitly assume in (1.6.1) that () uu yy and  () are strictly positive for all .) Also, let uu 1 X (t) =y(t) h u(t k) (1.6.2) k k=1 i i i i\sm2" i i 2004/2/22 page 13 i i Section 1.6 Complements 13 denote the residues of the least squares problem in Exercise 1.11. Hence, fh g in k equation (1.6.2) satisfy 1 X ik h e ,H() = ()= (): k yu uu k=1 In what follows, we will show that Z   1 2 2 E j(t)j = (1 jC ()j ) ()d (1.6.3) yu yy 2  where jC ()j is the socalled coherency spectrum. We will deduce from (1.6.3) yu that the coherency spectrum shows the extent to which y(t) and u(t) are linearly related to one another, hence providing a motivation for the name given to jC ()j. yu We will also show that jC ()j  1 with equality, for all values, if and only if yu y(t) and u(t) are related as in equation (1.6.2) with (t)  0 (in the mean square sense). Finally, we will show that jC ()j is invariant to linear ltering ofu(t) and yu y(t) (possibly by di erent lters); that is, if u =g u and y =f y where f and g are linear lters, then jC ()j = jC ()j. y u yu P 1 Let x(t) = h u(t k). It can be shown that u(t k) and (t) are k k=1 uncorrelated with one another for allk. (The reader is required to verify this claim; see also Exercise 1.11). Hence x(t) and (t) are also uncorrelated with each other. As y(t) =(t) +x(t); (1.6.4) it then follows that  () = () + () (1.6.5) yy  xx 2 By using the fact that  () = jH()j  (), we can write xx uu Z   1 2 E j(t)j =  ()d  2  Z    1  () uu 2 = 1 jH()j  ()d yy 2  () yy  Z    2 1 j ()j yu = 1  ()d yy 2  () () uu yy  Z    1 2 = 1 jC ()j  ()d yu yy 2  which is (1.6.3). Since the left-hand side in (1.6.3) is nonnegative and the PSD function () yy is arbitrary, we must have jC ()j  1 for all . It can also be seen from (1.6.3) yu that the closer jC ()j is to one, the smaller the residual variance. In particular, if yu jC ()j  1 then(t)  0 (in the mean square sense) and hencey(t) andu(t) must yu be linearly related as in (1.7.11). Owing to the previous interpretation, C () is yu sometimes called the correlation coecient in the frequency domain. i i i i\sm2" i i 2004/2/22 page 14 i i 14 Chapter 1 Basic Concepts Next, consider the ltered signals 1 X y (t) = f y(t k) k k=1 and 1 X u (t) = g u(t k) k k=1 where the lters ff g and fg g are assumed to be stable. As k k  r (p) , E fy (t)u (t p)g y u 1 1 X X   = f g E fy(t k)u (t j p)g k j k=1 j=1 1 1 X X  = f g r (j +p k); k yu j k=1 j=1 it follows that 1 1 1 X X X ik  ij i(j+pk)  () = f e g e r (j +p k)e y u k j yu p=1 k=1 j=1 0 1  1 1 1 X X X ik ij is A = f e g e r (s)e k j yu k=1 j=1 s=1  = F()G () () yu Hence jF()j jG()j j ()j yu jC ()j = = jC ()j y u yu 1=2 1=2 jF()j ()jG()j () yy uu which is the desired result. Observe that the latter result is similar to the invariance of the modulus of the correlation coecient in the time domain, jr (k)j yu ; 1=2 r (0)r (0) yy uu to a scaling of the two signals: y (t) =f y(t) and u (t) =g u(t). 1.7 EXERCISES Exercise 1.1: Scaling of the Frequency Axis In this text, the time variable t has been expressed in units of the sampling interval T (say). Consequently, the frequency is measured in cycles per sampling s interval. Assume we want the frequency units to be expressed in radians per sec- ond or in Hertz (Hz = cycles per second). Then we have to introduce the scaled frequency variables  ==T  2 =T ; =T rad/sec (1.7.1) s s s i i i i\sm2" i i 2004/2/22 page 15 i i Section 1.7 Exercises 15  and f = =2 (in Hz). It might be thought that the PSD in the new frequency variable is obtained by inserting = T into (), but this is wrong. Show that s the PSD, as expressed in units of power per Hz, is in fact given by: 1 X i T k  s ( ) =T ( T ) ,T r(k)e ; j j =T s s s s (1.7.2) k=1 (See Marple 1987 for more details on this scaling aspect.) Exercise 1.2: TimeFrequency Distributions Let y(t) denote a discretetime signal, and let Y () be its discretetime Fourier transform. As explained in Section 1.2, Y () shows how the energy in 1 the whole sequence fy(t)g is distributed over frequency. t=1 Assume that we want to determine how the energy of the signal is distributed in time and frequency. IfD(t;) is a function that characterizes the timefrequency distribution, then it should satisfy the socalled marginal properties: 1 X 2 D(t;) = jY ()j (1.7.3) t=1 and Z  1 2 D(t;)d = jy(t)j (1.7.4) 2  Use intuitive arguments to explain why the previous conditions are desirable prop- erties of a timefrequency distribution. Next, show that the socalled Rihaczek distribution,  it (1.7.5) D(t;) =y(t)Y ()e satis es conditions (1.7.3) and (1.7.4). (For treatments of the timefrequency dis- tributions, the reader is referred to Therrien 1992 and Cohen 1995). Exercise 1.3: Two Useful ZTransform Properties P 1 k (a) Let h be an absolutely summable sequence, and let H(z) = h z k k k=1 be its Ztransform. Find the Ztransforms of the following two sequences: (i) h k P 1  (ii) g = h h k m m=1 mk 1 n  (b) Show that if z is a zero of A(z) = 1 +a z +    +a z , then (1=z ) is a i 1 n i       zero of A (1=z ) (where A (1=z ) = A(1=z ) ). Exercise 1.4: A Simple ACS Example Let y(t) be the output of a linear system as in Figure 1.1 with lter H(z) = 1 1 (1 +b z )=(1 +a z ), and whose input is zero mean white noise with variance 1 1 2 2  (the ACS of such an input is   ). k;0 i i i i\sm2" i i 2004/2/22 page 16 i i 16 Chapter 1 Basic Concepts 2 (a) Find r(k) and () analytically in terms of a , b , and  . 1 1  (b) Verify that r(k) = r (k), and that jr(k)j  r(0). Also verify that when a 1 and b are real, r(k) can be written as a function of jkj. 1 Exercise 1.5: Alternative Proof that jr(k)j  r(0) We stated in the text that (1.3.4) follows from (1.3.6). Provide a proof of that statement. Also, nd an alternative, simple proof of (1.3.4) by using (1.3.8). Exercise 1.6: A Double Summation Formula A result often used in the study of discretetime random signals is the follow- ing summation formula: N N N1 X X X f(t s) = (N jj)f() (1.7.6) t=1 s=1 =N+1 where f() is an arbitrary function. Provide a proof of the above formula. Exercise 1.7: Is a Truncated Autocovariance Sequence (ACS) a Valid ACS? P 1 1 ik Suppose that fr(k)g is a valid ACS; thus, r(k)e  0 for all k=1 k=1 . Is it possible that for some integer p the partial (or truncated) sum p X ik r(k)e k=p is negative for some ? Justify your answer. Exercise 1.8: When Is a Sequence an Autocovariance Sequence? 1 We showed in Section 1.3 that if fr(k)g is an ACS, then R  0 for m k=1 m = 0; 1; 2;:::. We also implied that the rst de nition of the PSD in (1.3.7) satis es ()  0 for all ; however, we did not prove this by using (1.3.7) solely. Show that 1 X ik () = r(k)e  0 for all k=1 if and only if  a R a  0 for every m and for every vector a m Exercise 1.9: Spectral Density of the Sum of Two Correlated Signals Lety(t) be the output to the system shown in Figure 1.2. AssumeH (z) and 1 H (z) are linear, asymptotically stable systems. The inputse (t) ande (t) are each 2 1 2 zero mean white noise, with      2   e (t)    1 1 2   1 E e (s) e (s) =  t;s 1 2 2 e (t)    2 1 2 2 i i i i\sm2" i i 2004/2/22 page 17 i i Section 1.7 Exercises 17 e (t) x (t) 1 1 - H (z) 1 J J J y(t) m - +  e (t) x (t) 2 2 - H (z) 2 Figure 1.2. The system considered in Exercise 1.9. (a) Find the PSD of y(t). (b) Show that for  = 0,  () = () + (). y x x 1 2 2 2 2 2 2 (c) Show that for  = 1 and  = = ,  () = jH () H ()j . y 1 2 1 2 Exercise 1.10: Least Squares Spectral Approximation 1 Assume we are given an ACS fr(k)g or, equivalently, a PSD function k=1 () as in equation (1.3.7). We wish to nd a niteimpulse response (FIR) lter i im as in Figure 1.1, where H() =h +h e +::: +h e , whose input e(t) is 0 1 m zero mean, unit variance white noise, and such that the output sequence y(t) has T PSD  () \close to" (). Speci cally, we wish to nd h = h :::h so that y 0 m the approximation error Z  1 2  = ()  () d (1.7.7) y 2  is minimum. (a) Show that  is a quartic (fourthorder) function of h, and thus no simple closedform solution h exists to minimize (1.7.7). (b) We attempt to reparameterize the minimization problem as follows. We note that r (k)  0 for jkjm; thus, y m X ik  () = r (k)e (1.7.8) y y k=m  Equation (1.7.8) and the fact that r (k) = r (k) mean that  () is a y y y T function of g = r (0):::r (m) . Show that the minimization problem in y y (1.7.7) is quadratic ing; it thus admits a closedform solution. Show that the vector g that minimizes  in equation (1.7.7) gives ( r(k); jkj m r (k) = (1.7.9) y 0; otherwise i i i i\sm2" i i 2004/2/22 page 18 i i 18 Chapter 1 Basic Concepts (c) Can you identify any problems with the \solution" (1.7.9)? Exercise 1.11: Linear Filtering and the CrossSpectrum For two stationary signals y(t) and u(t), with (cross)covariance sequence  r (k) =E fy(t)u (t k)g, the crossspectrum is de ned as: yu 1 X ik  () = r (k)e (1.7.10) yu yu k=1 Let y(t) be the output of a linear lter with input u(t), 1 X y(t) = h u(t k) (1.7.11) k k=1 Show that the input PSD,  (), the lter transfer function uu 1 X ik H() = h e k k=1 and  () are related through the socalled WienerHopf equation: yu  () =H() () (1.7.12) yu uu Next, consider the following least squares (LS) problem, 8 9 2 1 = X min E y(t) h u(t k) (1.7.13) k fh g : ; k k=1 where nowy(t) andu(t) are no longer necessarily related through equation (1.7.11). Show that the lter minimizing the above LS criterion is still given by the Wiener Hopf equation, by minimizing the expectation in (1.7.13) with respect to the real and imaginary parts of h (assume that  () 0 for all ). k uu COMPUTER EXERCISES Exercise C1.12: Computer Generation of Autocovariance Sequences Autocovariance sequences are twosided sequences. In this exercise we develop computer techniques for generating twosided ACSs. Let y(t) be the output of the linear system in Figure 1.1 with lter H(z) = 1 1 (1 +b z )=(1 +a z ), and whose input is zero mean white noise with variance 1 1 2  . i i i i\sm2" i i 2004/2/22 page 19 i i Section 1.7 Exercises 19 2 (a) Find r(k) analytically in terms of a , b , and  (see also Exercise 1.4). 1 1 (b) Plot r(k) for 20 k  20 and for various values of a and b . Notice that 1 1 the tails of r(k) decay at a rate dictated by ja j. 1 2 (c) When a ' b and  = 1, then r(k) '  . Verify this for a = 0:95, 1 1 k;0 1 b = 0:9, and for a = 0:75, b = 0:7. 1 1 1 (d) A quick way to generate (approximately) r(k) on the computer is to use the 2  fact that r(k) =  h(k) h (k) where h(k) is the impulse response of the lter in Figure 1.1 (see equation (1.4.7)) and  denotes convolution. Consider the case where 1 m 1 +b z +    +b z 1 m H(z) = : 1 n 1 +a z +    +a z 1 n 2 Write a Matlab function genacs.m whose inputs areM, ,a andb, wherea andb are the vectors of denominator and numerator coecients, respectively, and whose output is a vector of ACS coecients from 0 toM. Your function M should make use of the Matlab function filter to generate fh g , and k k=0 2  conv to computer(k) = h(k)h (k) using the truncated impulse response sequence. 2 (e) Test your function using  = 1, a = 0:9 and b = 0:8. Try M = 20 and 1 1 M = 150; why is the result more accurate for larger M? Suggest a \rule of thumb" about a good choice of M in relation to the poles of the lter. The above method is a \quick and simple" way to compute an approximation to the ACS, but it is sometimes not very accurate because the impulse response is 2 truncated. Methods for computing the exact ACS from  , a and b are discussed in Exercise 3.2 and also in Kinkel, Perl, Scharf, and Stubberud 1979; Demeure and Mullis 1989. Exercise C1.13: DTFT Computations using TwoSided Sequences In this exercise we consider the DTFT of twosided sequences (including au- tocovariance sequences), and in doing so illustrate some basic properties of autoco- variance sequences. (a) We rst consider how to use the DTFT to determine () from r(k) on a computer. We are given an ACS: ( Mjkj ; jkj M M r(k) = (1.7.14) 0; otherwise Generater(k) forM = 10. Now, in Matlab form a vectorx of lengthL = 256 as: x = r(0);r(1);:::;r(M); 0:::; 0;r(M);:::;r(1) Verify that xf=fft(x) gives ( ) for = 2k=L. (Note that the elements k k of xf should be nonnegative and real.). Explain why this particular choice of x is needed, citing appropriate circular shift and zero padding properties of the DTFT. i i i i

Advise: Why You Wasting Money in Costly SEO Tools, Use World's Best Free SEO Tool Ubersuggest.