Source code for wale.VarianceCalculator

import numpy as np
from scipy.integrate import simpson
import pyccl as ccl
from .FilterFunctions import top_hat_filter, starlet_filter


[docs] class Variance: """ A class to compute linear and nonlinear variance using power spectrum interpolators and a specific cosmological model. Attributes: cosmo (Cosmology): An instance of a cosmology class providing necessary cosmological functions and parameters. PK_interpolator_linear (Interpolator): An interpolator instance for linear power spectrum calculations. PK_interpolator_nonlinear (Interpolator): An interpolator instance for nonlinear power spectrum calculations. model (str): The name of the cosmological model to be used for variance calculations. """ def __init__(self, cosmo, filter_type, pk): """ Initializes the Variance class with cosmology and parameters for P(k) calculation. Calculates the non-linear power spectrum, including cosmic variance noise if volume is specified. Parameters: cosmo (Cosmology_function): An instance of the cosmology class. z_values (array-like): Redshifts for lensing planes/primary calculations. volume (float, optional): Volume for cosmic variance calculation in (Mpc/h)^3. Defaults to None (no CV noise). delta_A0 (float, optional): Parameter for additional non-Gaussian noise. Defaults to 1.9. """ self.cosmo = cosmo self.filter_type = filter_type print("Variance module initialized...") self.pk = pk
[docs] def linear_sigma2(self, redshift, R1, R2=None, **kwargs): """ Calculates the linear variance σ² for given scales and redshift, considering the specified model adjustments. Parameters: redshift (float): The redshift at which to evaluate the variance. R1 (float): The first scale radius. R2 (float, optional): The second scale radius. Defaults to R1 if not specified. Returns: float: The linear variance σ² at the given scales and redshift. """ if R2 is None: R2 = R1 else: R2 = R2 # pk = self.PK_interpolator_linear.P(redshift, self.cosmo.k_values) pk = ( ccl.linear_matter_power( self.cosmo.cosmoccl, self.cosmo.k, 1 / (1 + redshift) ) * self.cosmo.h**3 ) if self.filter_type == "tophat": w1_2D = top_hat_filter(self.cosmo.k, R1) w2_2D = top_hat_filter(self.cosmo.k, R2) w2 = w1_2D * w2_2D elif self.filter_type == "starlet": w1_2D = starlet_filter(self.cosmo.k, R1) w2_2D = starlet_filter(self.cosmo.k, R2) w2 = w1_2D * w2_2D constant = 1.0 / 2.0 / np.pi integrand = self.cosmo.k_values * pk * w2 * constant return simpson(integrand, x=self.cosmo.k_values)
[docs] def nonlinear_sigma2(self, redshift, R1, R2=None, **kwargs): """ Calculates the nonlinear variance σ² for given scales and redshift, considering the specified model adjustments. Parameters: redshift (float): The redshift at which to evaluate the variance. R1 (float): The first scale radius. R2 (float, optional): The second scale radius. Defaults to R1 if not specified. Returns: float: The nonlinear variance σ² at the given scales and redshift. """ if R2 is None: R2 = R1 else: R2 = R2 # pk = kwargs.get("pk", self.pk[redshift]) pk = kwargs["pk"] if "pk" in kwargs else self.pk[redshift] # pk = self.pk[redshift] k = self.cosmo.k * self.cosmo.h if self.filter_type == "tophat": w1_2D = top_hat_filter(self.cosmo.k, R1) w2_2D = top_hat_filter(self.cosmo.k, R2) w2 = w1_2D * w2_2D elif self.filter_type == "starlet": w1_2D = starlet_filter(self.cosmo.k, R1) w2_2D = starlet_filter(self.cosmo.k, R2) w2 = w1_2D * w2_2D constant = 1.0 / 2.0 / np.pi integrand = k * pk * w2 * constant return simpson(integrand, x=k) / self.cosmo.h
[docs] def get_sig_slice(self, z, R1, R2): """ Calculates the slice variance σ² for the given scales and redshift in the nonlinear regime. Parameters: z (float): The redshift at which to evaluate the slice variance. R1 (float): The first scale radius. R2 (float): The second scale radius. Returns: float: The slice variance σ² at the given scales and redshift. """ if self.filter_type == "tophat": sigslice = ( self.nonlinear_sigma2(z, R1) + self.nonlinear_sigma2(z, R2) - 2.0 * self.nonlinear_sigma2(z, R1, R2) ) return sigslice elif self.filter_type == "starlet": sigslice = ( self.nonlinear_sigma2(z, R1) + self.nonlinear_sigma2(z, R2) - 2.0 * self.nonlinear_sigma2(z, R1, R2) ) return sigslice