c13-9 (Numerical Recipes in C), страница 3

However, the endpointcorrections are equally important in obtaining accurate values of integrals. Narasimhanand Karthikeyan [2] have given a formula that is algebraically equivalent to our trapezoidalformula. However, their formula requires the evaluation of two FFTs, which is unnecessary.The basic idea used here goes back at least to Filon [3] in 1928 (before the FFT!). He usedSimpson’s rule (quadratic interpolation).

Since this interpolation is not left-right symmetric,two Fourier transforms are required. An alternative algorithm for equation (13.9.14) has beengiven by Lyness in [4]; for related references, see [5]. To our knowledge, the cubic-orderformulas derived here have not previously appeared in the literature.Calculating Fourier transforms when the range of integration is (−∞, ∞) can be tricky.If the function falls off reasonably quickly at infinity, you can split the integral at a largeenough value of t. For example, the integration to + ∞ can be writtenZ ∞Z bZ ∞eiωt h(t) dt =eiωt h(t) dt +eiωt h(t) dtaaZbbe=aiωth(b)eiωbh0 (b)eiωbh(t) dt −−···+iω(iω)2(13.9.17)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.

Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).}13.10 Wavelet Transforms591The splitting point b must be chosen large enough that the remaining integral over (b, ∞) issmall. Successive terms in its asymptotic expansion are found by integrating by parts. Theintegral over (a, b) can be done using dftint. You keep as many terms in the asymptoticexpansion as you can easily compute. See [6] for some examples of this idea.

Morepowerful methods, which work well for long-tailed functions but which do not use the FFT,are described in [7-9].Narasimhan, M.S. and Karthikeyan, M. 1984, IEEE Transactions on Antennas & Propagation,vol. 32, pp. 404–408. [2]Filon, L.N.G. 1928, Proceedings of the Royal Society of Edinburgh, vol. 49, pp. 38–47. [3]Giunta, G. and Murli, A. 1987, ACM Transactions on Mathematical Software, vol. 13, pp. 97–107.

[4]Lyness, J.N. 1987, in Numerical Integration, P. Keast and G. Fairweather, eds. (Dordrecht:Reidel). [5]Pantis, G. 1975, Journal of Computational Physics, vol. 17, pp. 229–233. [6]Blakemore, M., Evans, G.A., and Hyslop, J. 1976, Journal of Computational Physics, vol.

22,pp. 352–376. [7]Lyness, J.N., and Kaper, T.J. 1987, SIAM Journal on Scientific and Statistical Computing, vol. 8,pp. 1005–1011. [8]Thakkar, A.J., and Smith, V.H. 1975, Computer Physics Communications, vol. 10, pp. 73–79. [9]13.10 Wavelet TransformsLike the fast Fourier transform (FFT), the discrete wavelet transform (DWT) isa fast, linear operation that operates on a data vector whose length is an integer powerof two, transforming it into a numerically different vector of the same length. Alsolike the FFT, the wavelet transform is invertible and in fact orthogonal — the inversetransform, when viewed as a big matrix, is simply the transpose of the transform.Both FFT and DWT, therefore, can be viewed as a rotation in function space, fromthe input space (or time) domain, where the basis functions are the unit vectors ei ,or Dirac delta functions in the continuum limit, to a different domain.

For the FFT,this new domain has basis functions that are the familiar sines and cosines. In thewavelet domain, the basis functions are somewhat more complicated and have thefanciful names “mother functions” and “wavelets.”Of course there are an infinity of possible bases for function space, almost all ofthem uninteresting! What makes the wavelet basis interesting is that, unlike sines andcosines, individual wavelet functions are quite localized in space; simultaneously,like sines and cosines, individual wavelet functions are quite localized in frequencyor (more precisely) characteristic scale. As we will see below, the particular kindof dual localization achieved by wavelets renders large classes of functions andoperators sparse, or sparse to some high accuracy, when transformed into the waveletdomain.

Analogously with the Fourier domain, where a class of computations, likeconvolutions, become computationally fast, there is a large class of computationsSample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited.

To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).CITED REFERENCES AND FURTHER READING:Stoer, J., and Bulirsch, R. 1980, Introduction to Numerical Analysis (New York: Springer-Verlag),p. 88. [1].