548 Chapter 13.Fourier and Spectral Applications solution where (f)is a real function.Differentiating with respect to and setting the result equal to zero gives Φ(f)= 1S(f)2 Is(f)2+IN()2 (13.3.6) This is the formula for the optimal filter (f). Notice that equation (13.3.6)involves S,the smeared signal,and N,the noise. The two of these add up to be C,the measured signal.Equation(13.3.6)does not contain U,the"true"signal.This makes for an important simplification:The optimal filter can be determined independently of the determination of the deconvolution function that relates S and U. To determine the optimal filter from equation(13.3.6)we need some way of separately estimating S and N.There is no way to do this from the measured signal C alone without some other information,or some assumption or guess.Luckily,the extra information is often easy to obtain.For example,we can sample a long stretch of data c(t)and plot its power spectral density using equations(12.0.14),(12.1.8),and (12.1.5).This quantity is proportional to the sum 2 S+N,so we have |S(f)2+lN(f)2≈P.(f)=IC(f)20≤f<fe (13.3.7) (More sophisticated methods of estimating the power spectral density will be discussed in 813.4 and 813.7.but the estimation above is almost always good enough 06263 for the optimal filter problem.)The resulting plot (see Figure 13.3.1)will often immediately show the spectral signature of a signal sticking up above a continuous noise spectrum.The noise spectrum may be flat,or tilted,or smoothly varying;it doesn't matter,as long as we can guess a reasonable hypothesis as to what it is. Draw a smooth curve through the noise spectrum.extrapolating it into the region dominated by the signal as well.Now draw a smooth curve through the signal plus noise power.The difference between these two curves is your smooth"model"of the 喻 signal power.The quotient of your model of signal power to your model of signal Numerica 10621 plus noise power is the optimal filter (f).[Extend it to negative values of f by the formula (-f)=(f).]Notice that (f)will be close to unity where the noise is negligible,and close to zero where the noise is dominant.That is how it does its job!The intermediate dependence given by equation(13.3.6)just turns out to be the optimal way of going in between these two extremes. Because the optimal filter results from a minimization problem,the quality of the results obtained by optimal filtering differs from the true optimum by an amount that is secondorder in the precision to which the optimal filter is determined.In other words,even a fairly crudely determined optimal filter(sloppy,say,at the 10 percent level)can give excellent results when it is applied to data.That is why the separation of the measured signal C into signal and noise components S and N can usefully be done"by eye"from a crude plot of power spectral density.All of this may give you thoughts about iterating the procedure we have just described.For example,after designing a filter with response (f)and using it to make a respectable guess at the signal U(f)=(f)C(f)/R(f),you might turn about and regard U(f)as a fresh new signal which you could improve even further with the same filtering technique.548 Chapter 13. Fourier and Spectral Applications Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). solution where Φ(f) is a real function. Differentiating with respect to Φ, and setting the result equal to zero gives Φ(f) = |S(f)| 2 |S(f)| 2 + |N(f)| 2 (13.3.6) This is the formula for the optimal filter Φ(f). Notice that equation (13.3.6) involves S, the smeared signal, and N, the noise. The two of these add up to be C, the measured signal. Equation (13.3.6) does not containU, the “true” signal. This makes for an important simplification: The optimal filter can be determined independently of the determination of the deconvolution function that relates S and U. To determine the optimal filter from equation (13.3.6) we need some way of separately estimating |S| 2 and |N| 2 . There is no way to do this from the measured signal C alone without some other information, or some assumption or guess. Luckily, the extra information is often easy to obtain. For example, we can sample a long stretch of data c(t) and plot its power spectral density using equations (12.0.14), (12.1.8), and (12.1.5). This quantity is proportional to the sum |S| 2 + |N| 2 , so we have |S(f)| 2 + |N(f)| 2 ≈ Pc(f) = |C(f)| 2 0 ≤ f<fc (13.3.7) (More sophisticated methods of estimating the power spectral density will be discussed in §13.4 and §13.7, but the estimation above is almost always good enough for the optimal filter problem.) The resulting plot (see Figure 13.3.1) will often immediately show the spectral signature of a signal sticking up above a continuous noise spectrum. The noise spectrum may be flat, or tilted, or smoothly varying; it doesn’t matter, as long as we can guess a reasonable hypothesis as to what it is. Draw a smooth curve through the noise spectrum, extrapolating it into the region dominated by the signal as well. Now draw a smooth curve through the signal plus noise power. The difference between these two curves is your smooth “model” of the signal power. The quotient of your model of signal power to your model of signal plus noise power is the optimal filter Φ(f). [Extend it to negative values of f by the formula Φ(−f) = Φ(f).] Notice that Φ(f) will be close to unity where the noise is negligible, and close to zero where the noise is dominant. That is how it does its job! The intermediate dependence given by equation (13.3.6) just turns out to be the optimal way of going in between these two extremes. Because the optimal filter results from a minimization problem, the quality of the results obtained by optimal filtering differs from the true optimum by an amount that is second orderin the precision to which the optimal filter is determined. In other words, even a fairly crudely determined optimal filter (sloppy, say, at the 10 percent level) can give excellent results when it is applied to data. That is why the separation of the measured signal C into signal and noise components S and N can usefully be done “by eye” from a crude plot of power spectral density. All of this may give you thoughts about iterating the procedure we have just described. For example, after designing a filter with response Φ(f) and using it to make a respectable guess at the signal U(f) = Φ(f)C(f)/R(f), you might turn about and regard U(f) as a fresh new signal which you could improve even further with the same filtering technique.