Class MeshField#
Defined in File field.hpp
Class Documentation#
-
class MeshField#
Discretely sampled field on a mesh grid from particle catalogues.
Public Functions
-
MeshField(trv::ParameterSet ¶ms)#
Construct the mesh field.
- Parameters
params – Parameter set.
-
~MeshField()#
Destruct the mesh field.
-
void initialise_density_field()#
(Re-)initialise the complex field (and its shadow) on mesh.
This is an explicit method to reset values of trv::MeshField::field (and its interlaced counterpart) to zeros.
-
void finalise_density_field()#
Finalise the complex field (and its shadow) on mesh.
This is an explicit method to free the resources occupied by trv::MeshField::field and may be called outside the class destructor.
-
const fftw_complex &operator[](int gid)#
Return mesh field grid cell value.
- Parameters
gid – Grid index.
- Returns
Field value.
-
void assign_weighted_field_to_mesh(ParticleCatalogue &particles, fftw_complex *weights)#
Assign a weighted field to a mesh by interpolation scheme.
- Parameters
particles – Particle catalogue.
weights – Weight field.
-
void compute_unweighted_field(ParticleCatalogue &particles)#
Compute the unweighted field.
This is is the number density field
\[ n(\vec{x}) = \sum_i \delta^{(\mathrm{D})}(\vec{x} - \vec{x}_i) \,. \]- Parameters
particles – Particle catalogue.
-
void compute_unweighted_field_fluctuations_insitu(ParticleCatalogue &particles)#
Compute the unweighted field fluctuations.
This is typically used for the number density fluctuations in a simulation snapshot in a periodic box, with a global field value calculated from the particle number and box volume subtracted to compute the fluctuations, i.e.
\[ \delta{n}(\vec{x}) = \sum_i \delta^{(\mathrm{D})}(\vec{x} - \vec{x}_i) - N / V \,, \]where \( N \) is the total number of particles and \( V \) the box volume.- Parameters
particles – Particle catalogue.
-
void compute_ylm_wgtd_field(ParticleCatalogue &particles_data, ParticleCatalogue &particles_rand, LineOfSight *los_data, LineOfSight *los_rand, double alpha, int ell, int m)#
Compute the weighted field (fluctuations) further weighted by the reduced spherical harmonics.
For a field (or its fluctuations) \( f \), this is
\[ f_{LM}(\vec{x}) = {\sum_i}{\vphantom{\sum}}' y_{LM}^*(\hat{\vec{x}}) w(\vec{x}) \delta^{(\mathrm{D})}(\vec{x} - \vec{x}_i) \,, \]where if a pair of catalogues are provided,\[ {\sum_i}{\vphantom{\sum}}' = \sum_{i \in \mathrm{data}} - \alpha \sum_{i \in \mathrm{rand}} \]with \( f \equiv \delta{n} \), and otherwise\[ {\sum_i}{\vphantom{\sum}}' = \sum_{i \in \mathrm{data\ or\ rand}} \]with \( f \equiv n \).See also
Eq. (34) in Sugiyama et al. (2019) [1803.02132].
- Parameters
particles_data – (Data-source) particle catalogue.
particles_rand – (Random-source) particle catalogue.
los_data – (Data-source) particle lines of sight.
los_rand – (Random-source) particle lines of sight.
alpha – Alpha contrast.
ell – Degree of the spherical harmonic.
m – Order of the spherical harmonic.
-
void compute_ylm_wgtd_field(ParticleCatalogue &particles, LineOfSight *los, double alpha, int ell, int m)#
Compute the weighted field further weighted by the reduced spherical harmonics.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
- Parameters
particles – Particle catalogue.
los – Particle lines of sight.
alpha – Alpha contrast.
ell – Degree of the spherical harmonic.
m – Order of the spherical harmonic.
-
void compute_ylm_wgtd_quad_field(ParticleCatalogue &particles_data, ParticleCatalogue &particles_rand, LineOfSight *los_data, LineOfSight *los_rand, double alpha, int ell, int m)#
Compute the quadratic weighted field (fluctuations) further weighted by the reduced spherical harmonics.
For a quadratic field (of fluctuations) \( f \), this is
\[ f_{LM}(\vec{x}) = {\sum_i}{\vphantom{\sum}}' y_{LM}^*(\hat{\vec{x}}) w(\vec{x})^2 \delta^{(\mathrm{D})}(\vec{x} - \vec{x}_i) \,, \]where if a pair of catalogues are provided,\[ {\sum_i}{\vphantom{\sum}}' = \sum_{i \in \mathrm{data}} + \alpha^2 \sum_{i \in \mathrm{rand}} \,, \]and otherwise\[ {\sum_i}{\vphantom{\sum}}' = \sum_{i \in \mathrm{data\ or\ rand}} \,. \]See also
Eq. (46) in Sugiyama et al. (2019) [1803.02132].
- Parameters
particles_data – (Data-source) particle catalogue.
particles_rand – (Random-source) particle catalogue.
los_data – (Data-source) particle lines of sight.
los_rand – (Random-source) particle lines of sight.
alpha – Alpha contrast.
ell – Degree of the spherical harmonic.
m – Order of the spherical harmonic.
-
void compute_ylm_wgtd_quad_field(ParticleCatalogue &particles, LineOfSight *los, double alpha, int ell, int m)#
Compute the quadratic weighted field (fluctuations) further weighted by the reduced spherical harmonics.
This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
- Parameters
particles – Particle catalogue.
los – Particle lines of sight.
alpha – Alpha contrast.
ell – Degree of the spherical harmonic.
m – Order of the spherical harmonic.
-
void fourier_transform()#
Fourier transform the field.
If
trv::MeshField.params.interlace
is set to “true”, interlacing is performed where a phase factor is multiplied into the ‘shadow’ complex field before the average of the complex field and its shadow is taken.
-
void inv_fourier_transform()#
Inverse Fourier transform the (FFT-transformed) field.
-
void apply_wide_angle_pow_law_kernel()#
Apply wide-angle correction kernel in configuration space.
This multiplies the field by a power law \( r^{- i - j} \) at order \( (i, j) \) where \( r \) is the grid cell radial distance.
-
void apply_assignment_compensation()#
Apply compensation needed after Fourier transform for assignment schemes.
-
void inv_fourier_transform_ylm_wgtd_field_band_limited(MeshField &field_fourier, std::vector<std::complex<double>> &ylm, double k_band, double dk_band, double &k_eff, int &nmodes)#
Inverse Fourier transform a field \( f \) weighted by the reduced spherical harmonics restricted to a wavenumber band.
This method computes the quantity
\[ F_{LM}(\vec{x}; k) = \frac{(2\pi)^3}{4\pi k^2} \int \frac{\mathrm{d}^3\,\vec{k}'}{(2\pi)^3} \mathrm{e}^{\mathrm{i} \vec{k}' \cdot \vec{x}} \delta^{(\mathrm{D})}(k' - k) y_{LM}(\hat{\vec{k}}) f(\vec{k}) \,. \]See also
Eq. (42) in Sugiyama et al. (2019) [1803.02132].
- Parameters
field_fourier – [in] A Fourier-space field.
ylm – [in] Reduced spherical harmonic on a mesh.
k_band – [in] Band wavenumber.
dk_band – [in] Band wavenumber width.
k_eff – [out] Effective band wavenumber.
nmodes – [out] Number of wavevector modes in band.
-
void inv_fourier_transform_sjl_ylm_wgtd_field(MeshField &field_fourier, std::vector<std::complex<double>> &ylm, trvm::SphericalBesselCalculator &sjl, double r)#
Inverse Fourier transform a field \( f \) weighted by the spherical Bessel function and reduced spherical harmonics.
This method computes the quantity
\[ F_{LM}(\vec{x}; r) = \mathrm{i}^L \int \frac{\mathrm{d}^3\,\vec{k}}{(2\pi)^3} \mathrm{e}^{\mathrm{i} \vec{k} \cdot \vec{x}} j_L(kr) y_{LM}(\hat{\vec{k}}) f(\vec{k}) \,. \]See also
Eq. (49) in Sugiyama et al. (2019) [1803.02132].
- Parameters
field_fourier – A Fourier-space field.
ylm – Reduced spherical harmonic on a mesh.
sjl – Spherical Bessel function interpolator.
r – Separation in configuration space.
-
double calc_grid_based_powlaw_norm(ParticleCatalogue &particles, int order)#
Calculate the normalisation factor \( 1/I_N \) for N-point (for two-point statistics), where
\[ I_N = \int\mathrm{d}^3\,\vec{x} w(\vec{x})^N \bar{n}(\vec{x})^N \,. \]- Parameters
particles – (Typically random-source) particle catalogue.
order – Order \( N \) of the N-point statistics.
- Returns
norm_factor Normalisation factor.
-
MeshField(trv::ParameterSet ¶ms)#