SIAM REVIEW c  2004 Society for Industrial and Applied Mathematics Vol. 46, No. 3, pp. 443–454 Accelerating the Nonuniform Fast Fourier Transform ∗ Leslie Greengard † June-Yub Lee ‡ Abstract. The nonequispaced Fourier transform arises in a variety of application areas, from medical imaging to radio astronomy to the numerical solution of partial differential equations. In a typical problem, one is given an irregular sampling of N data in the frequency domain and one is interested in reconstructing the corresponding function in the physical domain. When the sampling is uniform, the fast Fourier transform (FFT) allows this calculation to be computed in O(N log N ) operations rather than O(N 2 ) operations. Unfortunately, when the sampling is nonuniform, the FFT does not apply. Over the last few years, a number of algorithms have been developed to overcome this limitation and are often referred to as nonuniform FFTs (NUFFTs). These rely on a mixture of interpolation and the judicious use of the FFT on an oversampled grid [A. Dutt and V. Rokhlin, SIAM J. Sci. Comput., 14 (1993), pp. 1368–1383]. In this paper, we observe that one of the standard interpolation or “gridding” schemes, based on Gaussians, can be accelerated by a significant factor without precomputation and storage of the interpolation weights. This is of particular value in two- and three- dimensional settings, saving either 10 d N in storage in d dimensions or a factor of about 5–10 in CPU time (independent of dimension). Key words. nonuniform fast Fourier transform, fast gridding, FFT, image reconstruction AMS subject classifications. 42A38, 44A35, 65T50, 65R10 DOI. 10.1137/S003614450343200X 1. Introduction. In this note, we describe an extremely simple and efficient implementation of the nonuniform fast Fourier transform (NUFFT). There are a host of applications of such algorithms, and we refer the reader to the references [2, 6, 8, 11, 13, 14, 17] for examples. We restrict our attention here to one: function (or image) reconstruction from Fourier data as discussed in [6, 8, 11, 14]. Let us begin, however, with a more precise description of the computational task. In two dimensions, we define the nonuniform discrete Fourier transform of types 1 and 2 according to the formulae F (k 1 ,k 2 )= 1 N N−1  j=0 f j e −i(k1,k2)·xj , (1) ∗ Received by the editors July 23, 2003; accepted for publication (in revised form) December 1, 2003; published electronically July 30, 2004. http://www.siam.org/journals/sirev/46-3/43200.html † Courant Institute of Mathematical Sciences, New York University, New York, NY 10012 (greengard@cims.nyu.edu). The work of this author was supported by the Applied Mathematical Sciences Program of the U.S. Department of Energy under contract DEFGO288ER25053. ‡ Department of Mathematics, Ewha Womans University, Seoul, 120-750, Korea (jylee@math. ewha.ac.kr). The work of this author was supported by the Korea Research Foundation under grant 2002-015-CP0044. 443