You are not logged in to this journal. Log in

This paper is concerned with the evaluation of the partial Fourier transform $u_x=\sum_{|k|<c(x)}e^{2\pi\imath xk/N}f_k$ in one and two dimensions. The motivating application is the wave extrapolation procedure in reflection seismology. As the summation in the frequency variable depends on the location $x$, the standard fast Fourier transform does not apply here. Direct summation requires quadratic complexity, which is extremely expensive for large values of $N$. The main idea of our approach is to decompose the summation domain in the $(x,|k|)$ space hierarchically into dyadic squares or cubes. The computation associated with each square or cube is accelerated by the fractional Fourier transform in one dimension and by the butterfly algorithm for the sparse Fourier transform in two dimensions. The resulting algorithm in one dimension has an $O(N\log^2N)$ complexity for an input of size $N$ and is numerically exact. The algorithm in two dimensions takes $O(N^2\log^2N)$ steps for an input of size $N\times N$ and is also highly accurate. Numerical examples are included to demonstrate the efficiency of these algorithms.