-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathsrd_redshift_distributions.py
152 lines (134 loc) · 6.21 KB
/
srd_redshift_distributions.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
# Niko Sarcevic
# github/nikosarcevic
# ----------
import numpy
from numpy import exp
import pandas
from scipy.integrate import simpson
import yaml
class SRD:
"""
Generate the LSST DESC type redshift distributions
for lens and source sample for year 1 and year 10.
See the LSST DESC Science Requirements Document (SRD)
https://arxiv.org/abs/1809.01669. The model used here
is the Smail type redshift distribution. This class
reads the parameters automatically from a yaml file
included in this repository (lsst_desc_parameters.yaml).
...
Attributes
----------
galaxy_sample: string
galaxy sample for which the redshift distribution will be
calculated. Accepted values are "source_galaxies" and
"lens_galaxies".
forecast_year: string
year that corresponds to the SRD forecast. Accepted values
are "1" and "10"
"""
def __init__(self,
galaxy_sample={},
forecast_year={}):
supported_galaxy_samples = {"lens_sample", "source_sample"}
if galaxy_sample in supported_galaxy_samples:
self.galaxy_sample = galaxy_sample
else:
raise ValueError(f"galaxy_sample must be one of {supported_galaxy_samples}.")
supported_forecast_years = {"1", "10"}
if forecast_year in supported_forecast_years:
self.forecast_year = forecast_year
else:
raise ValueError(f"forecast_year must be one of {supported_forecast_years}.")
# Read in the LSST DESC redshift distribution parameters
with open("parameters/lsst_desc_parameters.yaml", "r") as f:
lsst_desc_parameters = yaml.load(f, Loader=yaml.FullLoader)
self.srd_parameters = lsst_desc_parameters[self.galaxy_sample][self.forecast_year]
self.pivot_redshift = self.srd_parameters["z_0"]
self.alpha = self.srd_parameters["alpha"]
self.beta = self.srd_parameters["beta"]
self.z_start = self.srd_parameters["z_start"]
self.z_stop = self.srd_parameters["z_stop"]
self.z_resolution = self.srd_parameters["z_resolution"]
# Default LSST DESC redshift interval
self.redshift_range = numpy.linspace(self.z_start, self.z_stop, self.z_resolution)
def smail_type_distribution(self,
redshift_range,
pivot_redshift=None,
alpha=None,
beta=None):
"""
Generate the LSST DESC SRD parametric redshift distribution (Smail-type).
For details check LSST DESC SRD paper https://arxiv.org/abs/1809.01669, equation 5.
The redshift distribution parametrisation is dN / dz = z ^ beta * exp[- (z / z0) ^ alpha],
where z is redshift, z0 is pivot redshift, and alpha and beta are power law indices.
----------
Arguments:
redshift_range: array
redshift range
pivot_redshift: float
pivot redshift
alpha: float
power law index in the exponent
beta: float
power law index in the prefactor
Returns:
redshift_distribution: array
A Smail-type redshift distribution over a range of redshifts.
"""
if not pivot_redshift:
pivot_redshift = self.pivot_redshift
if not alpha:
alpha = self.alpha
if not beta:
beta = self.beta
redshift_distribution = [z ** beta * exp(-(z / pivot_redshift) ** alpha) for z in redshift_range]
return redshift_distribution
def get_redshift_distribution(self,
redshift_range=None,
normalised=True,
save_file=True):
"""
Generate the LSST DESC type redshift distribution
for lens and source sample for year 1 and year 10.
See the LSST DESC Science Requirements Document (SRD)
https://arxiv.org/abs/1809.01669, eq. 5. The model is
the Smail type redshift distribution of the form
dN / dz = z ^ beta * exp[- (z / z0) ^ alpha] where
z is the redshift, z0 is the pivot redshift, and
beta and alpha are power law parameters. LSST DESC
has a set of z0, beta, and alpha parameters for
lens and source galaxy sample for year 1 and year 10
science requirements. The value of the parameters can
be found in the SRD paper. The parameters are automatically
read from a yaml file included in this repository
(lsst_desc_parameters.yaml).
----------
Arguments:
redshift_range: array
an array of redshifts over which the redshift distribution
will be defined. If not specified, the SRD default will
be used (redshift interval 0.01 < z < 4.).
normalised: bool
normalise the redshift distribution (defaults to True)
save_file: bool
option to save the output as a .csv (defaults to True).
Saves the redshift range and the corresponding redshift
distribution to a .csv file (with the redshift and the
redshift distribution columns).
Returns:
srd_redshift_distribution (array):
an LSST DESC SRD redshift distribution of a galaxy sample
for a chosen forecast year.
"""
if type(redshift_range) is not numpy.ndarray:
redshift_range = self.redshift_range
redshift_distribution = self.smail_type_distribution(redshift_range)
if normalised:
normalisation = numpy.array(simpson(redshift_distribution, redshift_range))
redshift_distribution = numpy.array(redshift_distribution / normalisation)
if save_file:
dndz_df = pandas.DataFrame({"z": redshift_range, "dndz": redshift_distribution})
dndz_df.to_csv(f"./srd_{self.galaxy_sample}_dndz_year{self.forecast_year}.csv",
index=False)
return redshift_distribution