astropaint package

Submodules

astropaint.paint_bucket module

library for simulating semi-analytic mock maps of CMB secondary anisotropies

class astropaint.paint_bucket.Canvas(catalog, nside, mode='healpy', R_times=1, inclusive=False)

Bases: object

healpy or flat-sky canvas with the location of the halos to paint the signal on

Cl
class Disc(catalog, nside, R_times, inclusive)

Bases: object

R_times
catalog
center_D_a
center_ang
center_index
center_vec
gen_cent2pix_hat(halo_list='All')
gen_cent2pix_mpc(halo_list='All')
gen_cent2pix_mpc_vec(halo_list='All')
gen_cent2pix_rad(halo_list='All')
gen_center_ang(halo_list='All', snap2pixel=False)

generate the angle (theta, phi) of the halo centers if snap2pixel is True, the angular coordinate of the halo center pixel will be returned

gen_center_ipix(halo_list='All')

generate ipix of the halo centers

gen_center_vec(halo_list='All')
gen_pixel_ang(halo_list='All')
gen_pixel_index(halo_list='All')
gen_pixel_vec(halo_list='All')

generate the unit vectors pointing to the disc pixels

Returns:
Return type:None
inclusive
nside
pix2cent_mpc
pix2cent_rad
pix2cent_vec
pixel_ang
pixel_index
Dl
R_max
add_cmb(Cl='LCDM', mode='TT', lmax=None, inplace=True, weight=1, overwrite=False, *args, **kwargs)

Add CMB to the pixels If Cl array is not provided, a Lambda-CDM power spectrum is loaded from disc

Parameters:
  • Cl – Input CMB Power Spectrum The default keyword ‘LCDM’ loads a CAMB generated power spectrum from disc
  • mode – ‘TT’ for Temperature Note: ‘EE’ and ‘BB’ for polarization are currently unavailable
  • lmax – Maximum ell mode to include in the CMB power spectrum
  • inplace – If True, the result will be added to canvas.pixels Otherwise the generated CMB map will be returned as output
  • weight – The multiplicative factor for the generated CMB pixels
  • args*args to be passed to healpy.synfast(*args)
  • kwargs*kwargs to be passed to healpy.synfast(**kwargs)
Returns:

  • None or np.ndarray
  • The generated CMB map has the same NSIDE as the canvas.pixels map

add_noise(Nl='Planck', mode='TT', frequency=[217], lmax=None, sigma_n=None, fwhm_b=None, apply_beam=False, inplace=True, weight=1, *args, **kwargs)

Add Noise to the pixels. The Nl can be either provided directly as an array, or as a keyword using the name of an experiment. The list of all available experiments can be found in astropaint/lib/noise.yml.

Parameters:
  • Nl – Input Noise Power Spectrum. - If Nl is an array, it will be used as the noise power spectrum Nl. - If Nl is a string (e.g. ‘Planck’, ‘SO’, or ‘S4’) the noise configuration for the selected frequency channel will be read from the noise.yml file. Arbitrary noise configurations can be added to this file and then called here (see lib.misc.get_experiment_Nl for details). - If Nl is set to ‘white’, then sigma_n and fwhm_b will be used to construct a white noise power spectrum on the fly (see lib.misc.get_custom_Nl for details).
  • mode – ‘TT’ for Temperature Note: ‘EE’ and ‘BB’ for polarization are currently unavailable
  • [GHz] (frequency) – Frequency channels at which the noise will be calculated. If a list is provided, the final noise will be the combination of all channels.
  • lmax – Maximum ell mode to include in the noise power spectrum
  • sigma_n

    For Nl=’white’, this will be used as the noise level sigma [uK-arcmin], i.e.

    w_inverse = arcmin2rad(sigma_n) ** 2 Nl = w_inverse (for each ell)

  • fwhm_b

    if apply_beam=True, this will be used as the beam Full-Width-Half-Maximum (FWHM), i.e.

    fwhm = arcmin2rad(np.asarray(fwhm)) sigma_theta = fwhm ** 2 / 8 / np.log(2) B2l = np.exp(-L * (L + 1) * sigma_theta

  • apply_beam – If True, the Noise power spectrum will be divided by the Beam spectrum B2l
  • inplace – If True, the result will be added to canvas.pixels Otherwise the generated noise map will be returned as output
  • weight – The multiplicative factor for the generated CMB pixels
  • args*args to be passed to healpy.synfast(*args)
  • kwargs*kwargs to be passed to healpy.synfast(**kwargs)
Returns:

  • None or np.ndarray
  • The generated Noise map has the same NSIDE as the canvas.pixels map

alm
almxfl(fl, inplace=True)

wrapper for healpy.almxfl multiplies the alms by the array fl (lmin=0, to lmax=3*nside+1)

beam_smooth(fwhm_b=0.0, sigma_b=None, *args, **kwargs)

Smoothes canvas.pixels with a gaussian beam using healpy.sphtfunc.smoothing

Parameters:
  • fwhm_b
  • sigma_b
  • args
  • kwargs
Returns:

Return type:

None

catalog
clean()

Clean the canvas and set all pixels to zero

Returns:
Return type:None
cmap
cutouts(halo_list='all', lon_range=[-1, 1], lat_range=None, xpix=200, ypix=None, apply_func=None, func_kwargs={})

Generate cutouts of angular size lon_range x lat_range around halo center with xpix & ypix pixels on each side.

*This method uses Healpy’s projector.CartesianProj class to perform the cartesian projection.

Parameters:
  • halo_list – index of halos to consider (e.g. [1,2,5,10]). goes through all the halos in the catalog when set to “all”.
  • lon_range – range of longitutes to cut around the halo center in degrees. e.g. [-1,1] cuts out 1 degree on each side of the halo. same as lon_range in healpy
  • lat_range – range of longitutes to cut around the halo center in degrees. by default (None) it is set equal to lon_range. same as lat_range in healpy
  • xpix – number of pixels on the x axis same as xpix in healpy
  • ypix – number of pixels on the y axis by default (None) it is set equal to xrange same as ypix in healpy
  • apply_func – function to apply to the patch after cutting it out. THe first argument of the function must be the input patch.
func_kwargs: dictionary or list of dictionaries

dictionary of keyword arguments to pass to apply_func

Example usage:

cutouts(apply_func=linear_transform,
slope=2, intercept=1)

if the func_kwargs are scalars, they will be the same for all the halos in apply_func. If they are arrays or lists, their lengths must be the same as the catalog size.

Returns:
Return type:generator
ell
get_Cl(save_alm=True, lmax=None)

find the power spectrum of the map (.pixels)

get_alm_from_pixels()

find the alm coefficients of the map

get_pixels_from_alm()

get the map from the alm coefficients

instantiate_discs()

instantiate the discs attribute using the Disc class Useful when the disc generators are exhausted and need to be reset

lmax
load_map_from_file(filename=None, prefix=None, suffix=None, inplace=True)

save the healpy map to file

Parameters:
  • filename (str) – custom file name; overrides the prefix and suffix and default file name
  • prefix (str) – prefix string to be added to the beginning of the default file name
  • suffix (str) – suffix string to be added to the end of default file name
  • inplace (bool) – if True, canvas.pixels will be loaded with the map from file
mask_alm(lmin=0, lmax=None, inplace=True)

only keep alms in the range between lmin and lmax (inclusive)

npix
nside
pixels
remove_cmb()

remove cmb from the pixels

remove_noise()

Remove canvas.noise from the pixels

save_Cl_to_file(prefix=None, suffix=None, filename=None)

save the map power spectrum to file

Parameters:
  • prefix (str) – prefix string to be added to the beginning of the default file name
  • suffix (str) – suffix string to be added to the end of default file name
  • filename (str) – custom file name; overrides the prefix and suffix and default file name
save_map_to_file(filename=None, prefix=None, suffix=None, overwrite=True)

save the healpy map to file

Parameters:
  • filename (str) – custom file name; overrides the prefix and suffix and default file name
  • prefix (str) – prefix string to be added to the beginning of the default file name
  • suffix (str) – suffix string to be added to the end of default file name
show_discs(projection='mollweide', graticule=False, *args, **kwargs)
show_halo_centers(projection='mollweide', graticule=True, marker='o', color=None, s=None, *args, **kwargs)
show_map(projection='mollweide', apply_func=None, *args, **kwargs)
stack_cutouts(halo_list='all', lon_range=[-1, 1], lat_range=None, xpix=200, ypix=None, inplace=True, parallel=False, apply_func=None, func_kwargs={})

Stack cutouts of angular size lon_range x lat_range around halo center with xpix & ypix pixels on each side. apply_func is applied to each cutout before stacking

*This method uses Healpy’s projector.CartesianProj class to perform the cartesian projection.

Parameters:
  • halo_list – index of halos to consider (e.g. [1,2,5,10]). goes through all the halos in the catalog when set to “all”.
  • lon_range – range of longitutes to cut around the halo center in degrees. e.g. [-1,1] cuts out 1 degree on each side of the halo. same as lon_range in healpy
  • lat_range – range of longitutes to cut around the halo center in degrees. by default (None) it is set equal to lon_range. same as lat_range in healpy
  • xpix – number of pixels on the x axis same as xpix in healpy
  • ypix – number of pixels on the y axis by default (None) it is set equal to xrange same as ypix in healpy
  • inplace – if True, the result is saved in canvas.stack. Otherwise the stack is returned as outpu.
  • with_ray – if True, stacking is be performed in parallel using ray.
  • apply_func

    function to apply to the patch after cutting it out. THe first argument of the function must be the input patch.

    Example:

    def linear_transform(patch, slope, intercept):
    return slope * patch + intercept
  • func_kwargs

    keyword arguments to pass to apply_func

    Example usage:

    stack_cutouts(halo_list=np.arange(10),
    apply_func=linear_transform, inplace=True, slope=2, intercept=1)

    This will apply the function linear_transform to the first 10 halos and stacks them together in canvas.stack.

    if the func_kwargs are scalars, they will be the same for all the halos in apply_func. If they are arrays or lists, their lengths must be the same as the catalog size.

Returns:

  • np.ndarray
  • or None (if inplace=True)

ud_grade(nside_out, inplace=True, **kwargs)

“wrapper for healpy.ud_grade() function changes the nside of canvas.pixels to nside_out

Parameters:
  • nside_out (int) – nside of the new map
  • inplace (bool) – if True, change canvas.pixels, otherwise return the new map
  • kwargs (kwargs or dict) – kwargs to be passeed to hp.ud_grade()
class astropaint.paint_bucket.Catalog(data=None, calculate_redshifts=False, default_redshift=0)

Bases: object

halo catalog containing halo masses, locations, velocities, and redshifts

x, y, z: [Mpc] v_x, v_y, v_z: [km/s] M_200c: [M_sun]

build_dataframe(calculate_redshifts=False, default_redshift=0)
cut_D_a(D_min=0.0, D_max=inf, inplace=True)

Cut the catalog according the the given angular diameter distance range

Parameters:
  • [Mpc] (D_max) – minimum halo angular diameter distance to keep
  • [Mpc] – maximum halo angular diameter distance to keep
Returns:

  • None
  • catalog.data will only contain halos with angular diameter distance D_a in the range
  • D_min < D_a < D_max

cut_D_c(D_min=0.0, D_max=inf, inplace=True)

Cut the catalog according the the given comoving distance range

Parameters:
  • [Mpc] (D_max) – minimum halo comoving distance to keep
  • [Mpc] – maximum halo comoving distance to keep
Returns:

  • None
  • catalog.data will only contain halos with comoving distance D_c in the range D_min < D_c <
  • D_max

cut_M_200c(mass_min=0.0, mass_max=inf, inplace=True)

Cut the catalog according the the given mass range

Parameters:
  • [M_sun] (mass_max) – minimum halo mass to keep
  • [M_sun] – maximum halo mass to keep
Returns:

  • None
  • catalog.data will only contain halos with mass M in the range mass_min < M < mass_max

cut_R_200c(R_min=0.0, R_max=inf, inplace=True)

Cut the catalog according the the given radius range

Parameters:
  • [Mpc] (R_max) – minimum halo radius to keep
  • [Mpc] – maximum halo radius to keep
Returns:

  • None
  • catalog.data will only contain halos with radius R in the range R_min < R < R_max

cut_R_ang_200c(R_ang_min=0.0, R_ang_max=inf, inplace=True)

Cut the catalog according the the given angular radius range

Parameters:
  • [arcmin] (R_ang_max) – minimum halo angular radius to keep
  • [arcmin] – maximum halo angular radius to keep
Returns:

  • None
  • catalog.data will only contain halos with angular radius R_ang_200c in the range
  • R_ang_min < R_ang_200c < R_ang_max

cut_lon_lat(lon_range=[0, 360], lat_range=[-90, 90], inplace=True)

Cut the catalog according the the given longitude and latitude range

Parameters:
  • [deg] (lat_range) – range of longitutes to keep
  • [deg] – rane of latitudes to keep
Returns:

  • None
  • catalog.data will only contain halos with longitutes in the range lon_range and
  • latitudes in the range lat_range

cut_mask(mask, threshold=0.5, inplace=True)

cut the catalog according to the input mask halos outside of the mask (mask<threshold) will be discarded

cut_redshift(redshift_min=0.0, redshift_max=inf, inplace=True)

Cut the catalog according the the given redshift range

Parameters:
  • redshift_min – minimum halo redshift to keep
  • redshift_max – maximum halo redshift to keep
Returns:

  • None
  • catalog.data will only contain halos with redshift in the range
  • redshift_min < redshift < redshift_max

cut_theta_phi(theta_range=[0, 3.141592653589793], phi_range=[0, 6.283185307179586], inplace=True)

Cut the catalog according the the given longitude and latitude range

Parameters:
  • [rad] (phi_range) – range of longitutes to keep
  • [rad] – rane of latitudes to keep
Returns:

  • None
  • catalog.data will only contain halos with theta in the range theta_range and
  • phi in the range phi_range

data
generate_random_box(box_size=50, v_max=100, mass_min=100000000000000.0, mass_max=1000000000000000.0, n_tot=50000, put_on_shell=False, inplace=True)
generate_random_shell(shell_radius=50, v_max=100, mass_min=100000000000000.0, mass_max=1000000000000000.0, n_tot=50000, inplace=True)
generate_test_box(configuration=['all'], distance=100, mass=1000000000000000.0, inplace=True)
load_from_csv(sample_name='MICE')

load sample data using the name of dataset

move_to_box_center()

move the observer from (0,0,0) to the center of the box (Lx/2, Ly/2, Lz/2) to make coordinates symmetric *Not recommended for light-cone catalogs

replicate(mode='rotate')

Replicate an octant to get a whole-sky box

Parameters:mode
save_to_csv(sample_name)

load sample data using the name of dataset

class astropaint.paint_bucket.Painter(template)

Bases: object

Painter object sprays a signal over the canvas using a template

calculate_template(R, catalog, halo_list=[0], **template_kwargs)

calculate the 1D profile of the template as a function of R [Mpc]

plot_template(R, catalog, halo_list=[0], axis=None, **template_kwargs)

plot the 1D profile of the template as a function of R [Mpc]

spray(canvas, distance_units='Mpc', n_cpus=1, cache=False, lazy=False, lazy_grid=None, **template_kwargs)

#TODO: add example

Parameters:
  • canvas
  • distance_units
  • with_ray
  • cache
  • lazy
  • template_kwargs
template

Module contents