FML
|
This namespace deals with computing correlations functions (N-point functions) in real and fourier space. More...
Namespaces | |
namespace | PAIRCOUNTS |
This namespace deals with paircounting and contains general auto and cross paircount methods in 1D, 2D or 3D for simulation boxes and survey data (with randoms). | |
Data Structures | |
class | PolyspectrumBinning |
Class for storing the results from general polyspectrum. More... | |
class | PowerSpectrumBinning |
Class for holding the results after binning up a power-spectrum. More... | |
Typedefs | |
template<int N> | |
using | FFTWGrid = FML::GRID::FFTWGrid< N > |
template<int N> | |
using | MonospectrumBinning = PolyspectrumBinning< N, 2 > |
template<int N> | |
using | BispectrumBinning = PolyspectrumBinning< N, 3 > |
template<int N> | |
using | TrispectrumBinning = PolyspectrumBinning< N, 4 > |
Functions | |
template<int N, class T > | |
void | compute_power_spectrum (int Ngrid, T *part, size_t NumPart, size_t NumPartTotal, PowerSpectrumBinning< N > &pofk, std::string density_assignment_method, bool interlacing) |
Assign particles to grid using density_assignment_method = NGP,CIC,TSC,PCS,... Fourier transform, decolvolve the window function for the assignement above. More... | |
template<int N, class T > | |
void | compute_power_spectrum_direct_summation (int Ngrid, const T *part, size_t NumPart, PowerSpectrumBinning< N > &pofk) |
Brute force (but aliasing free) computation of the power spectrum. More... | |
template<int N> | |
void | bin_up_power_spectrum (const FFTWGrid< N > &fourier_grid, PowerSpectrumBinning< N > &pofk) |
Compute the power-spectrum of a fourier grid. More... | |
template<int N> | |
void | bin_up_deconvolved_power_spectrum (const FFTWGrid< N > &fourier_grid, PowerSpectrumBinning< N > &pofk, std::string density_assignment_method) |
Compute the power-spectrum of a fourier grid. More... | |
template<int N> | |
void | bin_up_cross_power_spectrum (const FFTWGrid< N > &fourier_grid_1, const FFTWGrid< N > &fourier_grid_2, PowerSpectrumBinning< N > &pofk) |
Compute the cross power-spectrum of two fourier grids. More... | |
template<int N> | |
void | compute_power_spectrum_multipoles_fourier (const FFTWGrid< N > &fourier_grid, std::vector< PowerSpectrumBinning< N > > &Pell, std::vector< double > line_of_sight_direction) |
Compute power-spectrum multipoles from a Fourier grid where we have a fixed line_of_sight_direction (typical coordinate axes like \( (0,0,1) \)). More... | |
template<int N, class T > | |
void | compute_power_spectrum_multipoles (int Ngrid, FML::PARTICLE::MPIParticles< T > &part, double velocity_to_displacement, std::vector< PowerSpectrumBinning< N > > &Pell, std::string density_assignment_method, bool interlacing) |
A simple power-spectrum estimator for multipoles in simulations. More... | |
template<int N, class T , int ORDER> | |
void | compute_polyspectrum (int Ngrid, T *part, size_t NumPart, size_t NumPartTotal, PolyspectrumBinning< N, ORDER > &polyofk, std::string density_assignment_method, bool interlacing) |
Computes the polyspectrum \( P(k_1,k_2,\ldots,k_{\rm ORDER}) = \left<\delta(k_1)\cdots\delta(k_{\rm ORDER})\right> \) from particles. More... | |
template<int N, int ORDER> | |
void | compute_polyspectrum (const FFTWGrid< N > &fourier_grid, PolyspectrumBinning< N, ORDER > &polyofk) |
Computes the polyspectrum \( P(k_1,k_2,\ldots,k_{\rm ORDER}) = \left<\delta(k_1)\cdots\delta(k_{\rm ORDER})\right> \) from a fourier grid. More... | |
template<int N> | |
void | compute_monospectrum (const FFTWGrid< N > &fourier_grid, MonospectrumBinning< N > &pofk) |
Computes the monospectrum \( P(k_1,k_2) = \left<\delta(k_1)\delta(k_2)\right> \) from a fourier grid. More... | |
template<int N, class T > | |
void | compute_monospectrum (int Ngrid, T *part, size_t NumPart, size_t NumPartTotal, MonospectrumBinning< N > &pofk, std::string density_assignment_method, bool interlacing) |
Computes the monospectrum \( P(k_1,k_2) = \left<\delta(k_1)\delta(k_2)\right> \) from particles. More... | |
template<int N, class T > | |
void | compute_bispectrum (int Ngrid, T *part, size_t NumPart, size_t NumPartTotal, BispectrumBinning< N > &bofk, std::string density_assignment_method, bool interlacing) |
Computes the bispectrum \( B(k_1,k_2,k_3) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\right> \) from particles. More... | |
template<int N> | |
void | compute_bispectrum (const FFTWGrid< N > &fourier_grid, BispectrumBinning< N > &bofk) |
Computes the bispectrum \( B(k_1,k_2,k_3) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\right> \) from a fourier grid. More... | |
template<int N> | |
void | compute_trispectrum (const FFTWGrid< N > &fourier_grid, TrispectrumBinning< N > &tofk) |
Computes the trispectrum \( T(k_1,k_2,k_3,k_4) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\delta(k_4)\right> \) from a fourier grid. More... | |
template<int N, class T > | |
void | compute_trispectrum (int Ngrid, T *part, size_t NumPart, size_t NumPartTotal, TrispectrumBinning< N > &tofk, std::string density_assignment_method, bool interlacing) |
Computes the truspectrum \( T(k_1,k_2,k_3,k_4) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\delta(k_4)\right> \) from particles. More... | |
template<int N> | |
void | bin_up_cross_power_spectrum (FFTWGrid< N > &fourier_grid_1, FFTWGrid< N > &fourier_grid_2, PowerSpectrumBinning< N > &pofk) |
template<int N, class T > | |
void | compute_power_spectrum_direct_summation (int Ngrid, T *part, size_t NumPart, PowerSpectrumBinning< N > &pofk) |
template<int N, int ORDER> | |
void | compute_polyspectrum_bincount (int Nmesh, PolyspectrumBinning< N, ORDER > &polyofk) |
This method is used by compute_polyspectrum. More... | |
template<int N> | |
void | compute_monospectrum (const FFTWGrid< N > &fourier_grid, PolyspectrumBinning< N, 2 > &pofk) |
template<int N> | |
void | compute_bispectrum (const FFTWGrid< N > &fourier_grid, PolyspectrumBinning< N, 3 > &bofk) |
template<int N> | |
void | compute_trispectrum (const FFTWGrid< N > &fourier_grid, PolyspectrumBinning< N, 4 > &tofk) |
This namespace deals with computing correlations functions (N-point functions) in real and fourier space.
using FML::CORRELATIONFUNCTIONS::BispectrumBinning = typedef PolyspectrumBinning<N, 3> |
Definition at line 61 of file ComputePowerSpectrum.h.
using FML::CORRELATIONFUNCTIONS::FFTWGrid = typedef FML::GRID::FFTWGrid<N> |
Definition at line 39 of file ComputePowerSpectrum.h.
using FML::CORRELATIONFUNCTIONS::MonospectrumBinning = typedef PolyspectrumBinning<N, 2> |
Definition at line 58 of file ComputePowerSpectrum.h.
using FML::CORRELATIONFUNCTIONS::TrispectrumBinning = typedef PolyspectrumBinning<N, 4> |
Definition at line 64 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::bin_up_cross_power_spectrum | ( | const FFTWGrid< N > & | fourier_grid_1, |
const FFTWGrid< N > & | fourier_grid_2, | ||
PowerSpectrumBinning< N > & | pofk | ||
) |
Compute the cross power-spectrum of two fourier grids.
The result has no scales. Get scales by calling pofk.scale(boxsize) which does k \( k \to k/B \) and \( P(k) \to B^N P(k) \) where \( B\) is the boxsize once spectrum has been computed. The method assumes the two grids are fourier transforms of real grids (i.e. \( f(-k) = f^*(k) \)) and we only bin up the real part of \( f_1(k)f_2^*(k) \). The imaginary part is also binned up, but not returned.
N | Dimension of the grid |
[in] | fourier_grid_1 | Grid in fourier space |
[in] | fourier_grid_2 | Grid in fourier space |
[out] | pofk | Binned cross power-spectrum |
void FML::CORRELATIONFUNCTIONS::bin_up_cross_power_spectrum | ( | FFTWGrid< N > & | fourier_grid_1, |
FFTWGrid< N > & | fourier_grid_2, | ||
PowerSpectrumBinning< N > & | pofk | ||
) |
Definition at line 630 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::bin_up_deconvolved_power_spectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PowerSpectrumBinning< N > & | pofk, | ||
std::string | density_assignment_method | ||
) |
Compute the power-spectrum of a fourier grid.
The result has no scales. Get scales by calling pofk.scale(boxsize) which does \( k \to k/B \) and \( P(k) \to B^N P(k) \) where \( B\) is the boxsize once spectrum has been computed. The method assumes the two grids are fourier transforms of real grids (i.e. \( f(-k) = f^*(k) \)). This method divides the power-spectrum by the square of the window function corresponding to the density assignement method for every mode we bin up.
N | Dimension of the grid |
[in] | fourier_grid | Grid in fourier space |
[out] | pofk | Binned power-spectrum |
[in] | density_assignment_method | The density assignment method we use to create the grid |
Definition at line 534 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::bin_up_power_spectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PowerSpectrumBinning< N > & | pofk | ||
) |
Compute the power-spectrum of a fourier grid.
The result has no scales. Get scales by calling pofk.scale(boxsize) which does \( k \to k/B \) and \( P(k) \to B^N P(k) \) where \( B\) is the boxsize once spectrum has been computed. The method assumes the two grids are fourier transforms of real grids (i.e. \( f(-k) = f^*(k) \)).
N | Dimension of the grid |
[in] | fourier_grid | Grid in fourier space |
[out] | pofk | Binned power-spectrum |
Definition at line 585 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_bispectrum | ( | const FFTWGrid< N > & | fourier_grid, |
BispectrumBinning< N > & | bofk | ||
) |
Computes the bispectrum \( B(k_1,k_2,k_3) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\right> \) from a fourier grid.
This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{3} \) integrals. This is just an alias for compute_polyspectrum<N, 3>
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | fourier_grid | A fourier grid. |
[out] | bofk | The binned bispectrum. We required it to be initialized with the number of bins, kmin and kmax. |
void FML::CORRELATIONFUNCTIONS::compute_bispectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PolyspectrumBinning< N, 3 > & | bofk | ||
) |
Definition at line 1270 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_bispectrum | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
size_t | NumPartTotal, | ||
BispectrumBinning< N > & | bofk, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
Computes the bispectrum \( B(k_1,k_2,k_3) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\right> \) from particles.
Note that with interlacing we change the particle positions, but when returning they should be in the same state as when we started. This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{3} \) integrals. This is just an alias for compute_polyspectrum<N, 3>. Shot-noise is subtracted.
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. |
[in] | NumPartTotal | Total number of particles on all tasks. |
[out] | bofk | The binned bispectrum. We required it to be initialized with the number of bins, kmin and kmax. |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) |
[in] | interlacing | Use interlaced grids when computing density field for alias reduction |
Definition at line 1043 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_monospectrum | ( | const FFTWGrid< N > & | fourier_grid, |
MonospectrumBinning< N > & | pofk | ||
) |
Computes the monospectrum \( P(k_1,k_2) = \left<\delta(k_1)\delta(k_2)\right> \) from a fourier grid.
This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{2} \) integrals. This is just an alias for compute_polyspectrum<N, 2>
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | fourier_grid | A fourier grid. |
[out] | pofk | The binned monospectrum. We required it to be initialized with the number of bins, kmin and kmax. |
void FML::CORRELATIONFUNCTIONS::compute_monospectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PolyspectrumBinning< N, 2 > & | pofk | ||
) |
Definition at line 1265 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_monospectrum | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
size_t | NumPartTotal, | ||
MonospectrumBinning< N > & | pofk, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
Computes the monospectrum \( P(k_1,k_2) = \left<\delta(k_1)\delta(k_2)\right> \) from particles.
Note that with interlacing we change the particle positions, but when returning they should be in the same state as when we started. This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{2} \) integrals. This is just an alias for compute_polyspectrum<N, 2>. Shot-noise is subtracted
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. |
[in] | NumPartTotal | Total number of particles on all tasks. |
[out] | pofk | The binned monospectrum. We required it to be initialized with the number of bins, kmin and kmax. |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) |
[in] | interlacing | Use interlaced grids when computing density field for alias reduction |
Definition at line 1029 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_polyspectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PolyspectrumBinning< N, ORDER > & | polyofk | ||
) |
Computes the polyspectrum \( P(k_1,k_2,\ldots,k_{\rm ORDER}) = \left<\delta(k_1)\cdots\delta(k_{\rm ORDER})\right> \) from a fourier grid.
This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{\rm ORDER} \) integrals. If one is to compute many spectra with the same Ngrid and binning then one can precompute N123 in polyofk and set it using polyofk.set_bincount(N123). This speeds up the polyspectrum estimation by a factor of 2 by avoiding half of the fourier transforms.
N | The dimension of the particles. |
ORDER | The order. 2 is the power-spectrum, 3 is the bispectrum, 4 is the trispectrum. |
[in] | fourier_grid | Grid in fourier space |
[out] | polyofk | The binned polyspectrum. |
Definition at line 1280 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_polyspectrum | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
size_t | NumPartTotal, | ||
PolyspectrumBinning< N, ORDER > & | polyofk, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
Computes the polyspectrum \( P(k_1,k_2,\ldots,k_{\rm ORDER}) = \left<\delta(k_1)\cdots\delta(k_{\rm ORDER})\right> \) from particles.
Note that with interlacing we change the particle positions, but when returning they should be in the same state as when we started. This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{\rm ORDER} \) integrals. If one is to compute many spectra with the same Ngrid and binning then one can precompute N123 in polyofk and set it using polyofk.set_bincount(N123). This speeds up the polyspectrum estimation by a factor of 2. Shot-noise term is subtracted if ORDER=2 (power) or ORDER=3 (bi). Doing this for higher order requires all lower order correlators and we currently only bin up the polyspectrum and the power-spectrum.
N | The dimension of the particles. |
ORDER | The order. 2 is the power-spectrum, 3 is the bispectrum, 4 is the trispectrum. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. |
[in] | NumPartTotal | Total number of particles on all tasks. |
[out] | polyofk | The binned polyspectrum. We required it to be initialized with the number of bins, kmin and kmax. |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) |
[in] | interlacing | Use interlaced grids when computing density field for alias reduction |
Definition at line 1070 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_polyspectrum_bincount | ( | int | Nmesh, |
PolyspectrumBinning< N, ORDER > & | polyofk | ||
) |
This method is used by compute_polyspectrum.
It computes the number of generalized triangles of the bins needed to normalize the polyspectra up to symmetry (i.e. we only compute it for \( k_1 \leq k_2 \leq k_3 \) and only for valid triangle configurations \( k_1 + k_2 \geq k_3 \) and then set rest using symmetry). If one is to compute many spectra with the same Nmesh and binning then one can precompute N123 and set it using polyofk.set_bincount(N123) (which sets polyofk.bincount_is_set = true and avoid a call to this function) This speeds up the polyspectrum estimation by a factor of 2.
N | The dimension we work in |
ORDER | The order (mono = 2, bi = 3, tri = 4) |
[in] | Nmesh | The size of the grid we us |
[out] | polyofk | The binning (we compute and store the volumes of each bin, N123, in this binning). |
Definition at line 1147 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_power_spectrum | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
size_t | NumPartTotal, | ||
PowerSpectrumBinning< N > & | pofk, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
Assign particles to grid using density_assignment_method = NGP,CIC,TSC,PCS,... Fourier transform, decolvolve the window function for the assignement above.
With interlacing bin to a grid and an interlaced grid (displaced by dx/2 in all directions) Fourier transform both and add the together to cancel the leading aliasing contributions bin up power-spectrum and subtract shot-noise 1/NumPartTotal. Note that with interlacing we change the particle positions, but when returning they should be in the same state as when we started.
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. |
[in] | NumPartTotal | Total number of particles on all tasks. |
[out] | pofk | The binned power-spectrum. We required it to be initialized with the number of bins, kmin and kmax. |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) |
[in] | interlacing | Use interlaced grids for alias reduction. Twice as expensive, but allows us to use a smaller grid (so its actually faster if used correctly). |
Definition at line 891 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_power_spectrum_direct_summation | ( | int | Ngrid, |
const T * | part, | ||
size_t | NumPart, | ||
PowerSpectrumBinning< N > & | pofk | ||
) |
Brute force (but aliasing free) computation of the power spectrum.
Loop over all grid-cells and all particles and add up contribution and subtracts shot-noise term. Since we need to combine all particles with all cells this is not easiy parallelizable with MPI so we assume all CPUs have exactly the same particles when this is run on more than 1 MPI tasks (so best run just using OpenMP). This method scales as \( \mathcal{O}({\rm Npart}\cdot {\rm Nmesh}^N) \) so will be slow!
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. With MPI assumes all tasks have the same particles. |
[out] | pofk | The binned power-spectrum. We required it to be initialized with the number of bins, kmin and kmax. |
void FML::CORRELATIONFUNCTIONS::compute_power_spectrum_direct_summation | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
PowerSpectrumBinning< N > & | pofk | ||
) |
Definition at line 704 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_power_spectrum_multipoles | ( | int | Ngrid, |
FML::PARTICLE::MPIParticles< T > & | part, | ||
double | velocity_to_displacement, | ||
std::vector< PowerSpectrumBinning< N > > & | Pell, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
A simple power-spectrum estimator for multipoles in simulations.
Displacing particles from realspace to redshift-space using their velocities along each of the coordinate axes. Result is the mean of this. Deconvolving the window-function and subtracting shot-noise (1/NumPartTotal) for monopole.
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Particles in the form of a MPIParticle container |
[in] | velocity_to_displacement | Factor to convert a velocity to a displacement. This is \( \frac{c}{aH(a)B} \) for peculiar and \( \frac{c}{H(a)B} \) (where \( B\) is the boxsize) for comoving velocities At \( z = 0 \) velocity_to_displacement = \( \frac{1}{100 {\rm Boxsize}} \) when Boxsize is in Mpc/h |
[out] | Pell | Vector of power-spectrum binnings. The size of Pell is the maximum ell to compute. All binnings has to have nbins, kmin and kmax set. At the end Pell[ ell ] is a binning of P_ell(k). |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) to use. |
[in] | interlacing | Use interlaced grids for alias reduction when computing the density field |
Definition at line 777 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_power_spectrum_multipoles_fourier | ( | const FFTWGrid< N > & | fourier_grid, |
std::vector< PowerSpectrumBinning< N > > & | Pell, | ||
std::vector< double > | line_of_sight_direction | ||
) |
Compute power-spectrum multipoles from a Fourier grid where we have a fixed line_of_sight_direction (typical coordinate axes like \( (0,0,1) \)).
Pell contains \( P_0,P_1,P_2,\ldots,P_{\ell_{\rm max}} \) where \( \ell_{\rm max} = n-1 \) is the size of Pell
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | fourier_grid | The fourier grid to compute the multipoles from. |
[out] | Pell | Vector of power-spectrum binnings. The size of Pell is the maximum ell to compute. All binnings has to have nbins, kmin and kmax set. At the end Pell[ ell ] is a binning of P_ell(k). |
[in] | line_of_sight_direction | The line of sight direction, e.g. \( (0,0,1) \) for the z-axis. |
Definition at line 425 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_trispectrum | ( | const FFTWGrid< N > & | fourier_grid, |
PolyspectrumBinning< N, 4 > & | tofk | ||
) |
Definition at line 1275 of file ComputePowerSpectrum.h.
void FML::CORRELATIONFUNCTIONS::compute_trispectrum | ( | const FFTWGrid< N > & | fourier_grid, |
TrispectrumBinning< N > & | tofk | ||
) |
Computes the trispectrum \( T(k_1,k_2,k_3,k_4) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\delta(k_4)\right> \) from a fourier grid.
This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{4} \) integrals. This is just an alias for compute_polyspectrum<N, 4>
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | fourier_grid | A fourier grid. |
[out] | tofk | The binned trispectrum. We required it to be initialized with the number of bins, kmin and kmax. |
void FML::CORRELATIONFUNCTIONS::compute_trispectrum | ( | int | Ngrid, |
T * | part, | ||
size_t | NumPart, | ||
size_t | NumPartTotal, | ||
TrispectrumBinning< N > & | tofk, | ||
std::string | density_assignment_method, | ||
bool | interlacing | ||
) |
Computes the truspectrum \( T(k_1,k_2,k_3,k_4) = \left<\delta(k_1)\delta(k_2)\delta(k_3)\delta(k_4)\right> \) from particles.
Note that with interlacing we change the particle positions, but when returning they should be in the same state as when we started. This method allocates nbins FFTWGrids at the same time and performs \( 2{\rm nbins} \) fourier transforms and does \( {\rm nbins}^{4} \) integrals. This is just an alias for compute_polyspectrum<N, 4> The shot-noise term (something like (B(123) + cyc)/n + (P1P2 + cyc)/n^2 + 1/n^3) is not subtracted as it requires the bispectrum which we currently don't bin up together with the trispectrum
N | The dimension of the particles. |
T | The particle class. Must have a get_pos() method. |
[in] | Ngrid | Size of the grid to use. |
[in] | part | Pointer to the first particle. |
[in] | NumPart | Number of particles on the local task. |
[in] | NumPartTotal | Total number of particles on all tasks. |
[out] | tofk | The binned trispectrum. We required it to be initialized with the number of bins, kmin and kmax. |
[in] | density_assignment_method | The density assignment method (NGP, CIC, TSC, PCS or PQS) |
[in] | interlacing | Use interlaced grids when computing density field for alias reduction |
Definition at line 1057 of file ComputePowerSpectrum.h.