nbodykit.algorithms.convpower.
ConvolvedFFTPower
(source, poles, Nmesh=None, kmin=0.0, dk=None, use_fkp_weights=False, P0_FKP=None)[source]¶Bases: object
Algorithm to compute power spectrum multipoles using FFTs for a data survey with non-trivial geometry.
Due to the geometry, the estimator computes the true power spectrum convolved with the window function (FFT of the geometry).
This estimator implemented in this class is described in detail in Hand et al. 2017 (arxiv:1704.02357). It uses the spherical harmonic addition theorem such that only \(2\ell+1\) FFTs are required to compute each multipole. This differs from the implementation in Bianchi et al. and Scoccimarro et al., which requires \((\ell+1)(\ell+2)/2\) FFTs.
Results are computed when the object is inititalized, and the result is
stored in the poles
attribute. Important meta-data computed
during algorithm execution is stored in the attrs
dict. See the
documenation of run()
.
Note
A full tutorial on the class is available in the documentation here.
Parameters: |
|
---|
References
Methods
load (output[, comm]) |
Load a saved ConvolvedFFTPower result, which has been saved to disk with ConvolvedFFTPower.save() . |
run () |
Compute the power spectrum multipoles. |
save (output) |
Save the ConvolvedFFTPower result to disk. |
to_pkmu (mu_edges, max_ell) |
Invert the measured multipoles \(P_\ell(k)\) into power spectrum wedges, \(P(k,\mu)\). |
load
(output, comm=None)[source]¶Load a saved ConvolvedFFTPower result, which has been saved to
disk with ConvolvedFFTPower.save()
.
The current MPI communicator is automatically used
if the comm
keyword is None
logger
= <logging.Logger object>¶run
()[source]¶Compute the power spectrum multipoles. This function does not return anything, but adds several attributes (see below).
edges
¶array_like – the edges of the wavenumber bins
poles
¶BinnedStatistic
– a BinnedStatistic object that behaves similar to a structured array, with
fancy slicing and re-indexing; it holds the measured multipole
results, as well as the number of modes (modes
) and average
wavenumbers values in each bin (k
)
attrs
¶dict – dictionary holding input parameters and several important quantites computed during execution:
data.W
to randoms.W
data.shotnoise
+ randoms.shotnoise
; this should be subtracted from
the monopole.For further details on the meta-data, see the documentation.
save
(output)[source]¶Save the ConvolvedFFTPower result to disk.
The format is currently json.
Parameters: | output (str) – the name of the file to dump the JSON results to |
---|
to_pkmu
(mu_edges, max_ell)[source]¶Invert the measured multipoles \(P_\ell(k)\) into power spectrum wedges, \(P(k,\mu)\).
Parameters: |
|
---|---|
Returns: | pkmu – a data set holding the \(P(k,\mu)\) wedges |
Return type: |
nbodykit.algorithms.convpower.
get_real_Ylm
(l, m)[source]¶Return a function that computes the real spherical harmonic of order (l,m)
Parameters: | |
---|---|
Returns: | Ylm – a function that takes 4 arguments: (xhat, yhat, zhat) unit-normalized Cartesian coordinates and returns the specified Ylm |
Return type: |
References