12.4 FFT in Two or More Dimensions 521 An alternative way of implementing this algorithm is to form an auxiliary function by copying the even elements of fi into the first N/2 locations,and the odd elements into the next N/2 elements in reverse order.However,it is not easy to implement the alternative algorithm without a temporary storage array and we prefer the above in-place algorithm. Finally,we mention that there exist fast cosine transforms for small N that do not rely on an auxiliary function or use an FFT routine.Instead,they carry out the transform directly,often coded in hardware for fixed N of small dimension [1]. CITED REFERENCES AND FURTHER READING: Brigham,E.O.1974,The Fast Fourier Transform(Englewood Cliffs,NJ:Prentice-Hall),$10-10. Sorensen,H.V..Jones,D.L.,Heideman,M.T.,and Burris,C.S.1987,/EEE Transactions on Acoustics,Speech,and Signal Processing,vol.ASSP-35,pp.849-863. RECIPES I Hou,H.S.1987,IEEE Transactions on Acoustics,Speech,and Signal Processing,vol.ASSP-35, pp.1455-1461 [see for additional references]. Hockney,R.W.1971,in Methods in Computational Physics,vol.9(New York:Academic Press). > 令 Press. Temperton,C.1980.Journal of Computational Physics,vol.34,pp.314-329. Clarke,R.J.1985,Transform Coding of Images,(Reading,MA:Addison-Wesley). Gonzalez,R.C.,and Wintz,P.1987,Digital /mage Processing,(Reading.MA:Addison-Wesley). 9 Chen,W.,Smith,C.H.,and Fralick,S.C.1977,IEEE Transactions on Communications,vol.COM- 25,Pp.1004-1009.[1] 5会 OF SCIENTIFIC( 12.4 FFT in Two or More Dimensions COMPUTING (ISBN 188810920 Given a complex function h(,k2)defined over the two-dimensional grid 0≤k1≤Ni-l,0≤2≤N2-l,we can define its two-dimensional discrete Recipes 10.621 Fourier transform as a complex function H(n1,n2),defined over the same grid, Numerical Recipes 43106 (outside H(n1,n2) exp(2rik2n2/N2)exp(2rikin1/N1)h(k1:k2) k2=0k1=0 North (12.4.1) By pulling the"subscripts 2"exponential outside of the sum overk1,or by reversing the order of summation and pulling the"subscripts 1"outside of the sum over 2, we can see instantly that the two-dimensional FFT can be computed by taking one- dimensional FFTs sequentially on each index of the original function.Symbolically, H(n1,n2)=FFT-on-index-1(FFT-on-index-2 [h(k1,k2)]) (12.4.2) FFT-on-index-2(FFT-on-index-1 [h(1,k2)])12.4 FFT in Two or More Dimensions 521 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). } } } An alternative way of implementing this algorithm is to form an auxiliary function by copying the even elements of fj into the first N/2 locations, and the odd elements into the next N/2 elements in reverse order. However, it is not easy to implement the alternative algorithm without a temporary storage array and we prefer the above in-place algorithm. Finally, we mention that there exist fast cosine transforms for small N that do not rely on an auxiliary function or use an FFT routine. Instead, they carry out the transform directly, often coded in hardware for fixed N of small dimension [1]. CITED REFERENCES AND FURTHER READING: Brigham, E.O. 1974, The Fast Fourier Transform (Englewood Cliffs, NJ: Prentice-Hall), §10–10. Sorensen, H.V., Jones, D.L., Heideman, M.T., and Burris, C.S. 1987, IEEE Transactions on Acoustics, Speech, and Signal Processing, vol. ASSP-35, pp. 849–863. Hou, H.S. 1987, IEEE Transactions on Acoustics, Speech, and Signal Processing, vol. ASSP-35, pp. 1455–1461 [see for additional references]. Hockney, R.W. 1971, in Methods in Computational Physics, vol. 9 (New York: Academic Press). Temperton, C. 1980, Journal of Computational Physics, vol. 34, pp. 314–329. Clarke, R.J. 1985, Transform Coding of Images, (Reading, MA: Addison-Wesley). Gonzalez, R.C., and Wintz, P. 1987, Digital Image Processing, (Reading, MA: Addison-Wesley). Chen, W., Smith, C.H., and Fralick, S.C. 1977, IEEE Transactions on Communications, vol. COM- 25, pp. 1004–1009. [1] 12.4 FFT in Two or More Dimensions Given a complex function h(k1, k2) defined over the two-dimensional grid 0 ≤ k1 ≤ N1 − 1, 0 ≤ k2 ≤ N2 − 1, we can define its two-dimensional discrete Fourier transform as a complex function H(n1, n2), defined over the same grid, H(n1, n2) ≡ N
2−1 k2=0 N
1−1 k1=0 exp(2πik2n2/N2) exp(2πik1n1/N1) h(k1, k2) (12.4.1) By pulling the “subscripts 2” exponential outside of the sum over k 1, or by reversing the order of summation and pulling the “subscripts 1” outside of the sum over k 2, we can see instantly that the two-dimensional FFT can be computed by taking one￾dimensional FFTs sequentially on each index of the original function. Symbolically, H(n1, n2) = FFT-on-index-1(FFT-on-index-2[h(k1, k2)]) = FFT-on-index-2(FFT-on-index-1[h(k1, k2)]) (12.4.2)