Default domain for type 2 transform

[sijieh23]

Hi,

I'm using type 2 transform for interpolation on a nonuniform grid. Throughout the document, it says the default domain for the query points is $x_q\in[-\pi, \pi)$. However, I quickly tried it, and it seems like $x_q\in[0, 2\pi)$ is the default one instead of $x_q\in[-\pi, \pi)$.

I'm attaching my code here:

import cupy as cp
import cupyx.scipy.fft as fft
import cufinufft as nuft

n = 1024
x = cp.linspace(0, 2 * pi, n, endpoint=False)
k = fft.fftfreq(n, d=1 / n)

r = 0.05
G = cp.exp(-(r * k)**2 / 4) # filter kernel to suppress high-frequency component

plan = nuft.Plan(2, (n,), n_trans=1, eps=1e-12, isign=+1, dtype=cp.complex128, modeord=1)

u = cp.random.randn(n) * 5
u_hat = fft.fft(u, norm='forward')
u_hat *= G

nq = 10
xq = cp.sort(cp.random.uniform(0, TWOPI - 1e-6, nq))
# xq = cp.sort(cp.random.uniform(-pi, pi - 1e-6, nq))
plan.setpts(xq)
uq = plan.execute(u_hat).real

I'm attaching the results here. The top panel is for $x_q\in[0, 2\pi)$, while the bottom one is for $x_q\in[-\pi, \pi)$. Clearly, the bottom one does not agree with the expected result.

So I'm wondering if I'm missing something here, or if it's just a bug in the document?

Thanks!

[lu1and10]

Hi @sijieh23,

I guess cupyx.scipy.fft implements the DFT basis exp(-i*2*pi*k*n/N), it assumes that the samples live on the index grid n = 0, 1, ..., N-1, the angle grid is 2*pi*n/N(origin at 0).

If your query points are at x_q = -pi + 2*pi*n/N (origin at -pi), the constant shift of -pi introduces a mode wise phase ramp.
To evaluate at the query points with a type-2 NUFFT (which computes sum_k u_hat[k] * exp(+ikx)), you can either:

1. Rephase u_hat by $e^{ik\pi}$ and keep query points in [-pi, pi), add the following after u_hat *= G

k = np.arange(n)          
phase = np.exp(1j*np.pi*k)
u_hat = u_hat * phase

or simply:

2. Keep the u_hat as is and evaluate at shifted points x_q + pi (i.e., move your query points to [0, 2*pi)).

plan.setpts(xq+cp.pi)