Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

functionality to specify which cls to output + other (param/option name) updates for clarity #88

Merged
merged 3 commits into from
Nov 30, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
66 changes: 44 additions & 22 deletions deepcmbsim/camb_power_spectrum.py
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ def get_noise(self):
provides noise for the TT power spectrum and the polarization power spectra;
shape is (2, max_l_use)
"""
if self.UserParams['noise_type'] == 'white':
if self.UserParams['noise_type'] == 'detector-white':
t_noise = noise.detector_white_noise(self.UserParams['noise_uKarcmin'], self.UserParams['beamfwhm_arcmin'], self.max_l_use,
TT=True)
eb_noise = noise.detector_white_noise(self.UserParams['noise_uKarcmin'], self.UserParams['beamfwhm_arcmin'], self.max_l_use,
Expand All @@ -72,7 +72,7 @@ def get_noise(self):
elif self.UserParams['noise_type'] is None:
return np.zeros((2, self.max_l_use))
else:
print("only detector white noise is currently implemented, via `noise_type = 'white'` in `user_config.yaml`")
print("only detector white noise is currently implemented, via `noise_type = 'detector-white'` in `user_config.yaml`")
return np.zeros((2, self.max_l_use))

def get_cls(self, save_to_dict=None, user_params=True):
Expand All @@ -98,29 +98,51 @@ def get_cls(self, save_to_dict=None, user_params=True):
# main calculation: https://camb.readthedocs.io/en/latest/camb.html#camb.get_results
results = camb.get_results(self.CAMBparams)

outdict = { 'l': range(self.max_l_use + 1) }
if self.UserParams['cls_to_output'] == 'all':
cls_needed = ['clTT', 'clEE', 'clBB', 'clTE', 'clPP', 'clPT', 'clPE']
else:
cls_needed = self.UserParams['cls_to_output']

# https://camb.readthedocs.io/en/latest/results.html#camb.results.CAMBdata.get_total_cls
tt, ee, bb, te = results.get_total_cls(raw_cl=self.normalize_cls, CMB_unit=self.TT_units)[:self.max_l_use + 1].T
if self.UserParams['noise_type'] is not None:
_noise = self.get_noise()
tt += _noise[0]
ee += _noise[1]
bb += _noise[1]
te += _noise[1]
if ('clTT' in cls_needed) or ('clEE' in cls_needed) or ('clEB' in cls_needed) or ('clTE' in cls_needed):
# need to run things to get one/some/all of tt, ee, bb, te
tt, ee, bb, te = results.get_total_cls(raw_cl=self.normalize_cls, CMB_unit=self.TT_units)[:self.max_l_use + 1].T
# add noise
if self.UserParams['noise_type'] is not None:
_noise = self.get_noise()
tt += _noise[0]
ee += _noise[1]
bb += _noise[1]
te += _noise[1]
# now add to outdict
for key in ['clTT', 'clEE', 'clBB', 'clTE']:
if key in cls_needed:
if key == 'clTT':
outdict[key] = tt
elif key == 'clEE':
outdict[key] = ee
elif key == 'clBB':
outdict[key] = bb
elif key == 'clTE':
outdict[key] = te
else:
raise ValueError('somethings wrong.')

# https://camb.readthedocs.io/en/latest/results.html#camb.results.CAMBdata.get_lens_potential_cls
pp, pt, pe = results.get_lens_potential_cls(raw_cl=self.normalize_cls)[:self.max_l_use + 1].T
lvals = range(self.max_l_use + 1)

outdict = {
'l': lvals,
'clTT': tt,
'clEE': ee,
'clBB': bb,
'clTE': te,
'clPP': pp,
'clPT': pt,
'clPE': pe
}
if ('clPP' in cls_needed) or ('clPT' in cls_needed) or ('clPE' in cls_needed):
pp, pt, pe = results.get_lens_potential_cls(raw_cl=self.normalize_cls)[:self.max_l_use + 1].T
# now add to outdict
for key in ['clPP', 'clPT', 'clPE']:
if key in cls_needed:
if key == 'clPP':
outdict[key] = pp
elif key == 'clPT':
outdict[key] = pt
elif key == 'clPE':
outdict[key] = pe
else:
raise ValueError('somethings wrong.')

if bool(self.UserParams["verbose"]):
time_end = dt.now()
Expand Down
26 changes: 13 additions & 13 deletions deepcmbsim/cl_plotting.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,8 +15,8 @@ class flatmap:
number of pixels on each side of a flat map
degrees : float
number of degrees on each side of a flat map
seed : int, default -1
random seed to provide `namaster` (potential bug in `namaster`, does not always have the
namaster_seed : int, default -1
random namaster_seed to provide `namaster` (potential bug in `namaster`, does not always have the
behavior that is described at https://namaster.readthedocs.io/en/latest/pymaster.html
cl_dict : dict, optional
dictionary of power spectra. Necessary to use the primary method of this class.
Expand All @@ -40,7 +40,7 @@ class flatmap:
function for plotting a realization of a flat-sky map with the given number of pixels
representing the given size of sky patch based on the provided power spectra
"""
def __init__(self, pixels, degrees, seed = -1, cl_dict = None):
def __init__(self, pixels, degrees, namaster_seed = -1, cl_dict = None):
"""
initialize a map-making class
Parameters
Expand All @@ -49,32 +49,32 @@ def __init__(self, pixels, degrees, seed = -1, cl_dict = None):
number of pixels on each side of a flat map
degrees : float
number of degrees on each side of a flat map
seed : int, default -1
random seed to provide `namaster` (potential bug in `namaster`, does not always have the
namaster_seed : int, default -1
random namaster_seed to provide `namaster` (potential bug in `namaster`, does not always have the
behavior that is described at https://namaster.readthedocs.io/en/latest/pymaster.html
cl_dict : dict, optional
dictionary of power spectra. Necessary to use the primary method of this class.
If not specified, can still use the hidden mapping function.
"""
self.pixels, self.degrees = pixels, degrees # pixels on each side, degrees on each side
self.seed, self.cl_dict = seed, cl_dict # random seed if applicable, dictionary of CLs if you want to include them
self.namaster_seed, self.cl_dict = namaster_seed, cl_dict # random namaster_seed if applicable, dictionary of CLs if you want to include them
self.reso = self.degrees / self.pixels
self.nx, self.lx, self.lx_rad, self.ny, self.ly, self.ly_rad = [int(self.pixels), int(self.degrees), self.degrees*np.pi/180]*2
self.ticks, self.labels = np.linspace(0, self.pixels, self.degrees+1), np.arange(0, self.degrees+1, dtype=float)

def _flatmap(self, cl_array, spin_array = [0], seed = None):
seed = seed if seed is not None else self.seed
return nmt.synfast_flat(self.nx, self.ny, self.lx_rad, self.ly_rad, np.array([x for x in cl_array]), spin_array, seed = seed)
def _flatmap(self, cl_array, spin_array = [0], namaster_seed = None):
namaster_seed = namaster_seed if namaster_seed is not None else self.namaster_seed
return nmt.synfast_flat(self.nx, self.ny, self.lx_rad, self.ly_rad, np.array([x for x in cl_array]), spin_array, namaster_seed = namaster_seed)

def flatmap(self, maps_out, seed = None):
def flatmap(self, maps_out, namaster_seed = None):
"""

Parameters
----------
maps_out : ["T", "E", "B", "P", "TT", "EE", "BB", "TE", "PP", "PT", "PE", "clTT", "clEE", "clBB", "clTE", "clPP", "clPT", "clPE", "TEB", "TQU"]
map(s) that you would like to produce
seed : int, optional
random seed to provide `namaster` (potential bug in `namaster`, does not always have the
namaster_seed : int, optional
random namaster_seed to provide `namaster` (potential bug in `namaster`, does not always have the
behavior that is described at https://namaster.readthedocs.io/en/latest/pymaster.html
Returns
-------
Expand All @@ -98,6 +98,6 @@ def flatmap(self, maps_out, seed = None):
else:
print("not a valid map specification")
return None
return self._flatmap(cl_arr, spin_arr, seed=seed)
return self._flatmap(cl_arr, spin_arr, namaster_seed=namaster_seed)
else:
print("if you don't want to restrict to a `cl_dict` dictionary, use `self._flatmap` instead")
10 changes: 7 additions & 3 deletions deepcmbsim/settings/user_config.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -3,14 +3,18 @@ FORCAMB: # these simply overwrite any values in BASECAMBPARAMS
ITERABLES: # these provide iterables that can overwrite their corresponding values in BASECAMBPARAMS in loop settings
InitPower.r : [0.01, 0.1, 3] # for nested CAMBparams attributes, use the dot structure
Alens : [.8, 1.2, 3]
seed: 0
namaster_seed: 0 # seed for the map realization using namaster
verbose: 1
normalize_cls: False #raw_cl – return Cl rather than l*(l+1)*Cl/2π (Cl alone is not conventional)
TT_units: muK #return uK**2 units for the TT, EE, BB, TE Cls
outfile_dir: "outfiles"
noise_type: "white"
noise_type: "detector-white" # only option for now; adds detector white noise to the output; see noise.detector_white_noise for more details.
noise_uKarcmin: 5 # noise level in uK*arcmin
beamfwhm_arcmin: 3 # size of beam in arcmin
extra_l: 300
max_l_use: 10000 # max_l_use will differ from max_l and max_l_tensor by "extra_l" because
# according to the CAMB documentation errors affect the last "100 or so" multipoles
# according to the CAMB documentation errors affect the last "100 or so" multipoles
# ---
# decide on what to output; either 'all' or a subset list of
# ['clTT', 'clEE', 'clBB', 'clTE', 'clPP', 'clPT', 'clPE']
cls_to_output: 'all'