pynkowski.theory.temperature
1import numpy as np 2import scipy.stats 3 4from .utils import get_μ, define_mu, define_us_for_V 5 6norm = scipy.stats.norm() 7 8 9 10def V0_th_s0(u): 11 """Compute the expected value of the normalised first MF v0 at threshold u for a Gaussian isotropic scalar field normalised for its standard deviation. 12 13 Parameters 14 ---------- 15 u : np.array 16 Thresholds at which the first MF is computed. 17 18 Returns 19 ------- 20 v0 : np.array 21 The expected value of the normalised first MF at thresholds u. 22 23 """ 24 return (1-norm.cdf(u)) 25 26def V1_th_s0(u, μ): 27 """Compute the expected value of the normalised second MF v1 at threshold u for a Gaussian isotropic scalar field normalised for its standard deviation. 28 29 Parameters 30 ---------- 31 u : np.array 32 Thresholds at which the second MF is computed. 33 34 μ : float 35 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 36 37 Returns 38 ------- 39 v1 : np.array 40 The expected value of the normalised second MF at thresholds u. 41 42 """ 43 return np.sqrt(μ) * np.exp(-u**2./2.) / 8. 44 45def V2_th_s0(u, μ): 46 """Compute the expected value of the normalised third MF v2 at threshold u for a Gaussian isotropic scalar field normalised for its standard deviation. 47 48 Parameters 49 ---------- 50 u : np.array 51 Thresholds at which the third MF is computed. 52 53 μ : float 54 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 55 56 Returns 57 ------- 58 v2 : np.array 59 The expected value of the normalised third MF at thresholds u. 60 61 """ 62 return μ * u * np.exp(-u**2. /2.) / np.sqrt(2.*np.pi)**3. 63 64 65 66class TheoryTemperature(): 67 """Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field normalised for its standard deviation defined on the sphere 68 like the temperature anisotropies of the CMB. 69 70 Parameters 71 ---------- 72 us : np.array, optional 73 The thresholds at which the theoretical MFs will be computed. 74 If not given, a range between -5σ and 5σ with steps of 0.1σ is considered, 75 with σ=1 the expected standard deviation of the field. 76 77 cls : np.array, optional 78 The angular power spectrum of the Gaussian isotropic scalar field. 79 Default : None 80 81 μ : float, optional 82 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 83 If both μ and cls are given, an error will be raised. 84 If only cls is given, μ will be computed from input cls. 85 If neither μ nor cls are given, μ is assumed to be 1. 86 Default : None 87 88 average_bin : bool, optional 89 If True, the results of V1 and V2 are the average on each bin, to be compared with binned computations on maps. 90 If False, the results are the evaluation on the center of each bin. 91 The value is always exactly computed for V0, as the computation on maps does not imply binning. 92 Defaul : True 93 94 edges : bool, optional 95 If False, the given 'us' is considered as an array of uniformly distributed thresholds. 96 If True, input 'us' is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 97 In this last case, the thresholds are the central value of the given bins. 98 Neglected if 'us' is not given. 99 Default : False. 100 101 Attributes 102 ---------- 103 us : np.array 104 The thresholds at which the theoretical MFs are computed. 105 106 μ : float 107 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 108 109 """ 110 def __init__(self, us=None, cls=None, μ=None, average_bin=True, edges=False): 111 """Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field defined on the sphere. 112 113 """ 114 if (us is None): 115 Δu = 0.1 116 self.us = np.arange(-5+Δu/2., 5.+Δu/2., Δu) 117 self.dus = Δu*np.ones(self.us.shape[0]) 118 else: 119 us = np.array(us) 120 if us.shape == (1,): 121 self.us = us 122 self.dus = 0. 123 else: 124 if edges: 125 self.dus = (us[1:]-us[:-1]) 126 self.us = (us[1:]+us[:-1])/2. 127 else: 128 self.us = us 129 self.dus = (us[1]-us[0])*np.ones(us.shape[0]) 130 131 132 self.μ = define_mu(cls,μ) 133 if not average_bin: 134 self.dus = 0.*self.dus 135 136 137 def V0(self): 138 """Compute the expected values of the normalised first MF v0 at the different thresholds us. 139 $$\mathbb{E}[ v_{0} ] =1 -\Phi (u)$$ 140 where $\Phi$ is the cumulative normal distribution. 141 """ 142 return (V0_th_s0(self.us)) 143 144 def V1(self): 145 """Compute the expected values of the normalised second MF v1 at the different thresholds us. 146 147 $$\mathbb{E}[ v_{1}(u) ] = {1 \over 8} \exp{(- {u^2 \over 2})} \mu^{1/2}$$ 148 """ 149 us_ = define_us_for_V(self.us,self.dus) 150 v1_ = V1_th_s0(us_, self.μ) 151 152 return np.mean(v1_, axis=1) 153 154 def V2(self): 155 """Compute the expected values of the normalised third MF v2 at the different thresholds us. 156 157 $$\mathbb{E}[v_{2}(u)] = {2 \over \sqrt{(2\pi)^{3}}} \exp(-{u^{2} \over 2}) u$$ 158 """ 159 us_ = define_us_for_V(self.us,self.dus) 160 v2_ = V2_th_s0(us_, self.μ) 161 162 return np.mean(v2_, axis=1) 163 164 165 166__all__ = ["TheoryTemperature"] 167 168__docformat__ = "numpy"
67class TheoryTemperature(): 68 """Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field normalised for its standard deviation defined on the sphere 69 like the temperature anisotropies of the CMB. 70 71 Parameters 72 ---------- 73 us : np.array, optional 74 The thresholds at which the theoretical MFs will be computed. 75 If not given, a range between -5σ and 5σ with steps of 0.1σ is considered, 76 with σ=1 the expected standard deviation of the field. 77 78 cls : np.array, optional 79 The angular power spectrum of the Gaussian isotropic scalar field. 80 Default : None 81 82 μ : float, optional 83 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 84 If both μ and cls are given, an error will be raised. 85 If only cls is given, μ will be computed from input cls. 86 If neither μ nor cls are given, μ is assumed to be 1. 87 Default : None 88 89 average_bin : bool, optional 90 If True, the results of V1 and V2 are the average on each bin, to be compared with binned computations on maps. 91 If False, the results are the evaluation on the center of each bin. 92 The value is always exactly computed for V0, as the computation on maps does not imply binning. 93 Defaul : True 94 95 edges : bool, optional 96 If False, the given 'us' is considered as an array of uniformly distributed thresholds. 97 If True, input 'us' is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. 98 In this last case, the thresholds are the central value of the given bins. 99 Neglected if 'us' is not given. 100 Default : False. 101 102 Attributes 103 ---------- 104 us : np.array 105 The thresholds at which the theoretical MFs are computed. 106 107 μ : float 108 The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. 109 110 """ 111 def __init__(self, us=None, cls=None, μ=None, average_bin=True, edges=False): 112 """Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field defined on the sphere. 113 114 """ 115 if (us is None): 116 Δu = 0.1 117 self.us = np.arange(-5+Δu/2., 5.+Δu/2., Δu) 118 self.dus = Δu*np.ones(self.us.shape[0]) 119 else: 120 us = np.array(us) 121 if us.shape == (1,): 122 self.us = us 123 self.dus = 0. 124 else: 125 if edges: 126 self.dus = (us[1:]-us[:-1]) 127 self.us = (us[1:]+us[:-1])/2. 128 else: 129 self.us = us 130 self.dus = (us[1]-us[0])*np.ones(us.shape[0]) 131 132 133 self.μ = define_mu(cls,μ) 134 if not average_bin: 135 self.dus = 0.*self.dus 136 137 138 def V0(self): 139 """Compute the expected values of the normalised first MF v0 at the different thresholds us. 140 $$\mathbb{E}[ v_{0} ] =1 -\Phi (u)$$ 141 where $\Phi$ is the cumulative normal distribution. 142 """ 143 return (V0_th_s0(self.us)) 144 145 def V1(self): 146 """Compute the expected values of the normalised second MF v1 at the different thresholds us. 147 148 $$\mathbb{E}[ v_{1}(u) ] = {1 \over 8} \exp{(- {u^2 \over 2})} \mu^{1/2}$$ 149 """ 150 us_ = define_us_for_V(self.us,self.dus) 151 v1_ = V1_th_s0(us_, self.μ) 152 153 return np.mean(v1_, axis=1) 154 155 def V2(self): 156 """Compute the expected values of the normalised third MF v2 at the different thresholds us. 157 158 $$\mathbb{E}[v_{2}(u)] = {2 \over \sqrt{(2\pi)^{3}}} \exp(-{u^{2} \over 2}) u$$ 159 """ 160 us_ = define_us_for_V(self.us,self.dus) 161 v2_ = V2_th_s0(us_, self.μ) 162 163 return np.mean(v2_, axis=1)
Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field normalised for its standard deviation defined on the sphere like the temperature anisotropies of the CMB.
Parameters
- us (np.array, optional): The thresholds at which the theoretical MFs will be computed. If not given, a range between -5σ and 5σ with steps of 0.1σ is considered, with σ=1 the expected standard deviation of the field.
- cls (np.array, optional): The angular power spectrum of the Gaussian isotropic scalar field. Default : None
- μ (float, optional): The derivative of the covariance function at the origin for the Gaussian isotropic scalar field. If both μ and cls are given, an error will be raised. If only cls is given, μ will be computed from input cls. If neither μ nor cls are given, μ is assumed to be 1. Default : None
- average_bin (bool, optional): If True, the results of V1 and V2 are the average on each bin, to be compared with binned computations on maps. If False, the results are the evaluation on the center of each bin. The value is always exactly computed for V0, as the computation on maps does not imply binning. Defaul : True
- edges (bool, optional): If False, the given 'us' is considered as an array of uniformly distributed thresholds. If True, input 'us' is considered as a monotonically increasing array of bin edges, including the rightmost edge, allowing for non-uniform distributions of thresholds. In this last case, the thresholds are the central value of the given bins. Neglected if 'us' is not given. Default : False.
Attributes
- us (np.array): The thresholds at which the theoretical MFs are computed.
- μ (float): The derivative of the covariance function at the origin for the Gaussian isotropic scalar field.
111 def __init__(self, us=None, cls=None, μ=None, average_bin=True, edges=False): 112 """Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field defined on the sphere. 113 114 """ 115 if (us is None): 116 Δu = 0.1 117 self.us = np.arange(-5+Δu/2., 5.+Δu/2., Δu) 118 self.dus = Δu*np.ones(self.us.shape[0]) 119 else: 120 us = np.array(us) 121 if us.shape == (1,): 122 self.us = us 123 self.dus = 0. 124 else: 125 if edges: 126 self.dus = (us[1:]-us[:-1]) 127 self.us = (us[1:]+us[:-1])/2. 128 else: 129 self.us = us 130 self.dus = (us[1]-us[0])*np.ones(us.shape[0]) 131 132 133 self.μ = define_mu(cls,μ) 134 if not average_bin: 135 self.dus = 0.*self.dus
Class to compute the expected values of Minkowski functionals (MFs) for a Gaussian isotropic scalar field defined on the sphere.
138 def V0(self): 139 """Compute the expected values of the normalised first MF v0 at the different thresholds us. 140 $$\mathbb{E}[ v_{0} ] =1 -\Phi (u)$$ 141 where $\Phi$ is the cumulative normal distribution. 142 """ 143 return (V0_th_s0(self.us))
Compute the expected values of the normalised first MF v0 at the different thresholds us. $$\mathbb{E}[ v_{0} ] =1 -\Phi (u)$$ where $\Phi$ is the cumulative normal distribution.
145 def V1(self): 146 """Compute the expected values of the normalised second MF v1 at the different thresholds us. 147 148 $$\mathbb{E}[ v_{1}(u) ] = {1 \over 8} \exp{(- {u^2 \over 2})} \mu^{1/2}$$ 149 """ 150 us_ = define_us_for_V(self.us,self.dus) 151 v1_ = V1_th_s0(us_, self.μ) 152 153 return np.mean(v1_, axis=1)
Compute the expected values of the normalised second MF v1 at the different thresholds us.
$$\mathbb{E}[ v_{1}(u) ] = {1 \over 8} \exp{(- {u^2 \over 2})} \mu^{1/2}$$
155 def V2(self): 156 """Compute the expected values of the normalised third MF v2 at the different thresholds us. 157 158 $$\mathbb{E}[v_{2}(u)] = {2 \over \sqrt{(2\pi)^{3}}} \exp(-{u^{2} \over 2}) u$$ 159 """ 160 us_ = define_us_for_V(self.us,self.dus) 161 v2_ = V2_th_s0(us_, self.μ) 162 163 return np.mean(v2_, axis=1)
Compute the expected values of the normalised third MF v2 at the different thresholds us.
$$\mathbb{E}[v_{2}(u)] = {2 \over \sqrt{(2\pi)^{3}}} \exp(-{u^{2} \over 2}) u$$