Accurate Estimators of Correlation Functions in Fourier Space

TL;DR

This paper addresses aliasing in FFT-based estimates of Fourier-space statistics for cosmological data when small-scale power leaks into grid-based estimates. It proposes an efficient solution by combining interlaced density grids with higher-order mass-assignment kernels (notably PCS) to suppress aliasing and recover accurate amplitudes and phases up to the Nyquist frequency, $k_{ m Nyq}$. The authors also develop a fast bispectrum estimator using FFT-compatible integrals and demonstrate that aliasing can be reduced to well below the 0.01% level across relevant scales, outperforming or matching more costly methods like POWMES. The approach yields substantial practical benefits for analyzing large simulations and surveys, enabling unbiased Fourier-space statistics with reduced computational resources, and a public implementation is provided.

Abstract

Efficient estimators of Fourier-space statistics for large number of objects rely on Fast Fourier Transforms (FFTs), which are affected by aliasing from unresolved small scale modes due to the finite FFT grid. Aliasing takes the form of a sum over images, each of them corresponding to the Fourier content displaced by increasing multiples of the sampling frequency of the grid. These spurious contributions limit the accuracy in the estimation of Fourier-space statistics, and are typically ameliorated by simultaneously increasing grid size and discarding high-frequency modes. This results in inefficient estimates for e.g. the power spectrum when desired systematic biases are well under per-cent level. We show that using interlaced grids removes odd images, which include the dominant contribution to aliasing. In addition, we discuss the choice of interpolation kernel used to define density perturbations on the FFT grid and demonstrate that using higher-order interpolation kernels than the standard Cloud in Cell algorithm results in significant reduction of the remaining images. We show that combining fourth-order interpolation with interlacing gives very accurate Fourier amplitudes and phases of density perturbations. This results in power spectrum and bispectrum estimates that have systematic biases below 0.01% all the way to the Nyquist frequency of the grid, thus maximizing the use of unbiased Fourier coefficients for a given grid size and greatly reducing systematics for applications to large cosmological datasets.

Figures (9)

• Figure 1: Left panel: one-dimensional, real-space window functions $W^{(p)}(x)$ of the four mass assignment schemes described in the text. Different curves show Nearest Grid Point (NGP, dotted, black), Cloud In Cell (CIC, dashed, blue), Triangular Shaped Cloud (TSC, dotdashed, green) and Piecewise Cubic Spline interpolation (PCS, continuous, red). Right panel: the same window functions in Fourier space, $W^{(p)}(k)$.
• Figure 2: Residual factor $F_{\rm res, 1D}(k)$, in one spatial dimension, as a function of the ratio of the wavenumber $k$ to the Nyquist frequency, evaluated for the four interpolation schemes considered. Continuous and dashed curves assume respectively results with and without interlacing, eq. (\ref{['eq:res1Di']}) and eq. (\ref{['eq:res1Dnoi']}), respectively.
• Figure 3: Scatter plots comparing the (log) amplitude of the Fourier-space density, $\ln|\delta({\bf k})|$ obtained by interpolation on a grid and the same quantity obtained by direct summation. Top row assumes the CIC mass assignment scheme, while bottom row the PCS mass assignment scheme. Left-hand column assumes no correction by interlacing, while the right-hand column does. In all panels light shaded (cyan and orange) points correspond to all modes in the grid of linear size $N_G=100$ while dark shaded (blue and red) points are restricted to modes of $|{\bf k}|\le 0.7~{k_{\rm Nyq}}$.
• Figure 4: Same as figure \ref{['fig:amplitudes']}, but comparing the phases of the Fourier-space mode, $\phi=\arctan\,[{\rm Im}\, \delta({\bf k})\, /\, {\rm Re}\, \delta({\bf k})]$.
• Figure 5: Comparison between the power spectrum measured by FFT of a grid interpolated field and the same quantity obtained by direct summation as a function of the wavenumber in units of the Nyquist frequency. The left-hand panel show the absolute value of the relative difference between the cross-power spectrum of the grid-based, window-corrected field with the direct summation density field, $P_{\rm G-DS}\sim \langle\delta^{\rm G}({\bf k})\,\delta(-{\bf k})\rangle$ to the auto-power spectrum $P_{DS}\sim \langle|\delta({\bf k})|^2\rangle$. The right-hand panel shows instead the relative difference between the auto-power spectra $P_{G}\sim \langle|\delta^{\rm G}({\bf k})|^2\rangle$ and $P_{DS}$. Results in the top panels assume a (linear) FFT grid of $N_G=100$ while the bottom panels assume $N_G=256$. All measurements are performed on the $N$-body simulation described in ColombiEtal2009.