Mathematical definitions of transforms¶
We use notation with a general space dimensionality \(d\), which will be 1, 2, or 3, in our library. The arbitrary (ie nonuniform) points in space are denoted \(\mathbf{x}_j \in \mathbb{R}^d\), \(j=1,\dots,M\). We will see that for type 1 and type 2, without loss of generality one could restrict to the periodic box \([-\pi,\pi)^d\). For type 1 and type 3, each such NU point carries a given associated strength \(c_j\in\mathbb{C}\). Type 1 and type 2 involve the Fourier “modes” (Fourier series coefficients) with integer frequency indices lying in the Cartesian product set
where
For instance, \(K_{10} = \{-5,-4,\dots,4\}\), whereas \(K_{11} = \{-5,-4,\dots,5\}\), and so \(K_{10,11} = \{(-5,-5),(-4,-5),\dots,(4,-5),(-5,-4),(-4,-4),\dots,(3,5),(4,5)\}\). Note that the ordering in the last case is with the first index “fast”, second “slow”; this matches the storage ordering in the library interface. Thus, in the 1D case \(K\) is an interval containing \(N_1\) integer indices, in 2D it is a list of \(N_1N_2\) index pairs (which may be thought of as a rectangle of frequencies), and in 3D it is a list of \(N_1N_2N_3\) index triplets (which may be thought of as a cuboid).
Then the type 1 (nonuniform to uniform, aka “adjoint”) NUFFT evaluates
This can be viewed as evaluating a set of Fourier series coefficients due to sources with strengths \(c_j\) at the arbitrary locations \(\mathbf{x}_j\). Either sign of the imaginary unit in the exponential can be chosen in the interface. Note that our normalization differs from that of references [DR,GL].
The type 2 (U to NU, aka “forward”) NUFFT evaluates
This is the adjoint of the type 1, ie the evaluation of a given Fourier series at a set of arbitrary points. Both type 1 and type 2 transforms are invariant under translations of the NU points by multiples of \(2\pi\), thus one could require that all NU points live in the origin-centered box \([-\pi,\pi)^d\). In fact, the NU points may be any real numbers: each coordinate is folded periodically into \([-\pi,\pi)^d\), so round-off error grows with \(|\mathbf{x}_j|\).
Finally, the type 3 (NU to NU) transform does not have restrictions on the NU points, and there is no periodicity. Let \(\mathbf{x}_j\in\mathbb{R}^d\), \(j=1,\ldots,M\), be NU locations, with strengths \(c_j \in \mathbb{C}\), and let \(\mathbf{s}_k\), \(k=1,\ldots,N\) be NU frequencies. Then the type 3 transform evaluates:
For all three transforms, the computational effort scales like the product of the space-bandwidth products (real-space width times frequency-space width) in each dimension. For type 1 and type 2 this means near-linear scaling in the total number of modes \(N := N_1\dots N_d\). However, be warned that for type 3 this means that, even if \(N\) and \(M\) are small, if the product of the tightest intervals enclosing the coordinates of \(\mathbf{x}_j\) and \(\mathbf{s}_k\) is large, the algorithm will be inefficient. For such NU points, a direct sum should be used instead.
We emphasise that the NUFFT tasks that this library performs should not be confused with either the discrete Fourier transform (DFT), the (continuous) Fourier transform (although it may be used to approximate this via a quadrature rule), or the inverse NUFFT (the iterative solution of the linear system arising from nonuniform Fourier sampling, as in, eg, MRI). It is also important to know that, for NU points, the type 1 is not the inverse of the type 2. See the references for clarification.