API documentation

Core BLIP logic

`run_blip` module

class blip.run_blip.LISA(params, inj)[source]

Bases: LISAdata, Model

Generic class for getting data and setting up the prior space and likelihood.

compute_inj_fdata_spectra()[source]

A function to compute the injected spectra as evaluated at the simulated data frequencies.

This is needed to allow us to accurately plot and compare injected spectra vs. simulated data vs. recovered models.

diag_spectra()[source]: A function to do simple diagnostics. Plot the expected spectra and data.

generate_inj_fdata_responses()[source]

A function to generate the response functions for each injection at the simulated data frequencies.

This is needed to allow us to accurately plot and compare injected spectra vs. simulated data vs. recovered models.

load_blip_injection_data()[source]: Function to load in data created via BLIP’s own Injection routines, with associated Injection file. This is useful because we save true values, frozen spectra, etc. in addition to the raw data stream, so loading autogenerated data in this way allows BLIP’s plotting and post-processing routines to return significantly more informative outputs.

make_data_correlation_matrix()[source]

Uses the generated time-domain data series to construct a data correlation matrix.

Used to be the initialization of the (now defunct) likelihoods.py

makedata()[source]: Just a wrapper function to use the methods the LISAdata class to generate data. Return Frequency domain data.

plot_spectra()[source]: A function to make a plot of the data spectrum. For use with external (non-autogenerated) data, where we cannot calculate the intrinsic components.

blip.run_blip.blip(paramsfile, *, resume)[source]

Run BLIP on a given parameter file with given command-line options.

Parameters:

paramsfile (str) – Path to INI parameter file.
resume (bool) – CLI option to resume from a checkpointed run.

blip.run_blip.main()[source]

blip.run_blip.run_pipeline(parsed_params, resume, pre_sample_hook=None)[source]

Run the Bayesian pipeline.

Parameters:

parsed_params (tuple) – Parameters for the run as output by parse_config().
resume (bool) – Whether to resume a checkpointed run.
pre_sample_hook (Callable, optional) – Function to call on the analysis model before sampling. Can be used to perform arbitrary changes to the model. By default None

`submodel` module

class blip.src.submodel.SubmodelKind(*values)[source]

Bases: Enum

Kind of submodel specifier string.

Most are ‘spectral_spatial’.

There are only two exceptions: - ‘population’, which is shorthand for ‘population_population’; - ‘noise’ and ‘fixednoise’, because there is no spatial model for noise.

NOISE = 1: ‘noise’ or ‘fixednoise’ submodel.

POPULATION = 2: ‘population’ submodel, shorthand for ‘population_population’.

SPECTRAL_SPATIAL = 3

‘powerlaw_fixedgalaxy’.

Type:: General spectral_spatial submodel specifier. Example

class blip.src.submodel.SubmodelSpec(name: str, kind: SubmodelKind, is_injection: bool, spectral: str | None, spatial: str | None, count: str, raw_name: str, truevals: dict, fixedvals: dict, alias: str | None)[source]

Bases: object

Parsed (hence valid) specification of a submodel. Use this in submodel.__init__().

To get a specification from a string, use blip.config.parse_model_spec().

alias: str | None: Raw name of aliased injection submodel. Only used in analysis submodels.

count: str

‘powerlaw_isgwb-2’ -> ‘(2)’.

Type:: Duplicate count enclosed by parenthesis. Example

fixedvals: dict: Fixed values of parameters. Only for analysis (not injection) submodels.

is_injection: bool: True if the submodel is intended to be used as part of an injection. False if it is intended for analysis.

kind: SubmodelKind: Kind of specifier string.

name: str

‘population-1’ -> ‘population’.

Type:: Submodel name without duplicate count. Example

raw_name: str: Complete submodel name. Used as unique id.

spatial: str | None: Spatial model.

spectral: str | None: Spectral model.

truevals: dict: True parameter values. For analysis submodels, this is derived from the aliased injection submodel.

blip.src.submodel.gen_blm_parameters(blmax)[source]

Function to make the blm parameter name strings for all blms of a given lmax, in the correct order.

Parameters:: (int) (blmax)
Returns:: blm_parameters (list of str)
Return type:: Ordered list of blm parameter name strings

class blip.src.submodel.submodel(params, inj, spec: SubmodelSpec, fs, f0, tsegmid, injection=False, suffix='', parallel_response=False)[source]

Bases: fast_geometry, clebschGordan, instrNoise

Modular class that can represent either an injection or an analysis model. Will have different attributes depending on use case.

Includes all information required to generate an injection or a likelihood/prior.

New models (injection or analysis) should be added here.

bdspec_prior(theta)[source]

Prior for spectrum of bulge or disk of MW (bdspec).

Parameters:: theta (array (3,)) – samples from unit cube
Returns:: rescaled samples in order: log(Omega_ref), log(f_cut), and log(f_scale)
Return type:: theta

broken_powerlaw_fixed_a1_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter broken power law spectral model. Fixed low-frequency slope.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

broken_powerlaw_fixed_a1_spectrum(fs, log_omega0, alpha_2, log_fbreak, delta)[source]

Function to calculate a broken power law spectrum, with a fixed low-frequency slope and variable turnover scale.

Parameters:

floats) (fs (array of)
(float) (log_fbreak)
(float)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting broken power law spectrum

broken_powerlaw_fixed_a1delta_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter broken power law spectral model. Fixed low-frequency slope.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

broken_powerlaw_fixed_a1delta_spectrum(fs, log_omega0, alpha_2, log_fbreak)[source]

Function to calculate a broken power law spectrum, with a fixed low-frequency slope and turnover scale.

Parameters:

floats) (fs (array of)
(float) (log_fbreak)
(float)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting broken power law spectrum

broken_powerlaw_prior(theta)[source]

Prior function for a stochastic signal search with a broken power law spectral model.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha_1, log(Omega_ref), alpha_2, and log(f_break).
Return type:: float

broken_powerlaw_spectrum(fs, alpha_1, log_omega0, alpha_2, log_fbreak)[source]

Function to calculate a broken power law spectrum.

Parameters:

floats) (fs (array of)
(float) (log_fbreak)
(float)
(float)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting broken power law spectrum

compute_Omega0_from_Sgw(fs, Sgw)[source]

Wrapper function to generically calculate the associated stochastic gravitational wave dimensionless GW energy density Omega(f): for a spectral model given in terms of the PSD (S_gw)

Parameters:

floats) (fs (array of)
(list) (sgw_args)

Returns:

Sgw (array of floats)

Return type:

the resulting GW PSD

compute_Sgw(fs, omegaf_args)[source]

Wrapper function to generically calculate the associated stochastic gravitational wave PSD (S_gw): for a spectral model given in terms of the dimensionless GW energy density Omega(f)

Parameters:

floats) (fs (array of)
(list) (omegaf_args)

Returns:

Sgw (array of floats)

Return type:

the resulting GW PSD

compute_cov_asgwb(theta)[source]

Computes the covariance matrix contribution from a generic anisotropic stochastic GW signal.

Parameters:: (float) (theta)
Returns:: cov_sgwb (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for an anisotropic SGWB submodel.

compute_cov_fixed(dummy_theta)[source]

Wrapper to allow for “models” that are fixed a priori.

Parameters:: (NoneType) (dummy_theta)
Returns:: cov_fixed (array)
Return type:: The precomputed 3 x 3 x frequency x time covariance matrix for the fixed model.

compute_cov_fixed_asgwb(theta)[source]

Computes the covariance matrix contribution from an anisotropic stochastic GW signal with a known (assumed) sky distribution.

Parameters:: (float) (theta)
Returns:: cov_sgwb (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for an anisotropic SGWB submodel.

compute_cov_fixedspec_parameterized_asgwb(theta)[source]

Computes the covariance matrix contribution from a explicitly parameterized (i.e. not a generic spherical harmonic model), pixel-basis anisotropic stochastic GW signal.

Assumes a fixed spectral model. Only compatible with fixedspec spectral models.

Parameters:: (float) (theta)
Returns:: cov_sgwb (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for an anisotropic SGWB submodel.

compute_cov_isgwb(theta)[source]

Computes the covariance matrix contribution from a generic isotropic stochastic GW signal.

Parameters:: (float) (theta)
Returns:: cov_sgwb (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for an isotropic SGWB submodel.

compute_cov_noise(theta)[source]

Computes the noise covariance for a given draw of log_Np, log_Na

Parameters:: (float) (theta)
Returns:: cov_noise (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for the detector noise submodel.

compute_cov_parameterized_asgwb(theta)[source]

Computes the covariance matrix contribution from a explicitly parameterized (i.e. not a generic spherical harmonic model), pixel-basis anisotropic stochastic GW signal.

Parameters:: (float) (theta)
Returns:: cov_sgwb (array)
Return type:: The corresponding 3 x 3 x frequency x time covariance matrix for an anisotropic SGWB submodel.

compute_skymap_alms(blm_params)[source]

Function to compute the anisotropic skymap a_lms from the blm parameters.

Parameters:: floats) (blm_params (array of complex)
Returns:: alm_vals (array of complex floats)
Return type:: the corresponding alms

compute_summed_pixel_response(pixelmap)[source]

Function to compute the integrated, skymap-convolved anisotropic response for an arbitrary skymap in the pixel basis.

Parameters:: array) (pixelmap (healpy)
Returns:: summ_response_mat (array)
Return type:: the sky-integrated response (3 x 3 x frequency x time)

compute_summed_response(alms)[source]

Function to compute the integrated, skymap-convolved anisotropic response

Parameters:: floats) (alms (array of complex)
Returns:: summ_response_mat (array)
Return type:: the sky/alm-integrated response (3 x 3 x frequency x time)

fixed_model_wrapper_prior(theta)[source]

Wrapper function to allow fixed “models” to function within Model.prior

Parameters:: theta (float) – Shoudl always be None or [].
Returns:: theta – empty list
Return type:: float

fixed_truncated_powerlaw_spectrum(fs)[source]

Function to calculate a tanh-truncated power law spectrum with all parameters fixed.

Parameters:: floats) (fs (array of)
Returns:: spectrum (array of floats)
Return type:: the resulting truncated power law spectrum

fixedpowerlaw_prior(theta)[source]

Prior function for a power law with fixed slope.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha and log(Omega0)
Return type:: float

fixedpowerlaw_spectrum(fs, log_omega0)[source]

Function to calculate a simple power law spectrum, fixed to the alpha=2/3 prediction for the stellar origin binary background.

Parameters:

floats) (fs (array of)
(float) (log_omega0)

Returns:

spectrum (array of floats)

Return type:

the resulting power law spectrum

fixedsky_prior(theta)[source]

Fixed sky prior transform. Just serves as a wrapper for the spectral prior, as no additional foofaraw is necessary.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled for the spectral parameters.
Return type:: float

flatlowpowerlaw_prior(theta)[source]

Prior function for an isotropic stochastic backgound analysis. Imposes an upper constraint on the prior to aid in spectral separation.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha and log(Omega0)
Return type:: float

instr_noise_prior(theta)[source]

Prior function for only instrumental noise

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, omega_ref, Np and Na
Return type:: float

isotropic_prior(theta)[source]

Isotropic prior transform. Just serves as a wrapper for the spectral prior, as no additional foofaraw is necessary.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled for the spectral parameters.
Return type:: float

lmcspec_prior(theta)[source]

Prior function for a stochastic signal search with a 3-parameter truncated power law spectral model.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

lmcspecbpl_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter broken power law spectral model. Tailored for the LMC DWD SGWB spectrum.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

lmcspecbplad_prior(theta)[source]

Prior function for a stochastic signal search with a 3-parameter broken power law spectral model. Tailored for the LMC DWD SGWB spectrum.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

lmcspecfbpl_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter broken power law spectral model. Tailored for the LMC DWD SGWB spectrum. In contrast to lmcspecbpl, this variant fixes the smoothing parameter delta and allows alpha_1 to vary.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), log(f_scale)
Return type:: float

lowpowerlaw_prior(theta)[source]

Prior function for an isotropic stochastic backgound analysis. Imposes an upper constraint on the prior to aid in spectral separation.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha and log(Omega0)
Return type:: float

mask_and_norm_pixel_skymap(skymap)[source]

Function that takes in a modeled astrophysical skymap and:

masks to the pixels where we have computed the response
normalizes the masked map such that the sky integral = 1

Parameters:: array) (skymap (healpy)

mw1parameter_prior(theta)[source]

Hierarchical anisotropic prior transform. Combines a generic spectral prior function with the hierarchical astrophysical prior.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled for both the spectral and spatial parameters.
Return type:: float

mw2parameter_prior(theta)[source]

Hierarchical anisotropic prior transform. Combines a generic spectral prior function with the hierarchical astrophysical prior.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled for both the spectral and spatial parameters.
Return type:: float

mwspec3par_prior(theta)[source]

Prior function for a stochastic signal search with a 3-parameter truncated power law spectral model.

Bounds are astrophysically-motivated and tailored to expectations of the MW foreground.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), and log(f_scale)
Return type:: float

mwspec4par_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter truncated power law spectral model.

Bounds are astrophysically-motivated and tailored to expectations of the MW foreground.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), and log(f_scale)
Return type:: float

mwspec_prior(theta)[source]

Prior function for a stochastic signal search with a 2-parameter truncated power law spectral model.

Bounds are astrophysically-motivated and tailored to expectations of the MW foreground.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), and log(f_cut)
Return type:: float

powerlaw_prior(theta)[source]

Prior function for an isotropic stochastic backgound analysis.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha and log(Omega0)
Return type:: float

powerlaw_spectrum(fs, alpha, log_omega0)[source]

Function to calculate a simple power law spectrum.

Parameters:

floats) (fs (array of)
(float) (log_omega0)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting power law spectrum

process_astro_skymap_injection(skymap)[source]

Function that takes in an astrophysical pixel skymap and:

calculates all associated sph quantities
computes corresponding blm parameter truevals
convolves with response

Parameters:: array) (skymap (healpy)

process_astro_skymap_model(skymap)[source]

Function that takes in an astrophysical pixel skymap and:

calculates all associated sph quantities
convolves with response
sets sample-time response to be the map-convolved

This is intended for use with models that assume a fixed spatial distribution (e.g., fixedgalaxy, hotpixel).

Parameters:: array) (skymap (healpy)

recompute_response(f0=None, tsegmid=None)[source]

Function to recompute the LISA response matrices if needed.

When we save the Injection object, we delete the LISA response of each injection, as to do otherwise takes up egregious amounts of disk space. This allows us to recompute them identically as desired.

Parameters:

(array) (tsegmid)
(array)

Returns:

response_mat (array)

Return type:

The associated response for this submodel.

robson19_foreground_spectrum(fs, logA)[source]

Function to calculate an analytical spectrum for the Galactic foreground of the form given in Robson et al. (2019) (arXiv:1803.01944)

NOTE: this is given in terms of PSD amplitude A, as opposed to the usual units used in BLIP (dimensionless GW energy density)

Parameters:

floats) (fs (array of)
(float) (A)

Returns:

spectrum (array of floats)

Return type:

the resulting analytical foreground spectrum

robson19_foreground_varied_spectrum(fs, alpha, logA, log_fknee)[source]

Function to calculate an analytical spectrum for the Galactic foreground of the form given in Robson et al. (2019) (arXiv:1803.01944)

This version also varies the slope alpha and break (“knee”) frequency f_knee.

NOTE: this is given in terms of PSD amplitude A, as opposed to the usual units used in BLIP (dimensionless GW energy density)

Parameters:

floats) (fs (array of)
(float) (log_fknee)
(float)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting analytical foreground spectrum

robson19foreground_prior(theta)[source]

robson19foregroundvaried_prior(theta)[source]

sobbh_powerlaw_prior(theta)[source]

Prior function for a power law with fixed slope, with astrophysical prior bounds tailored to the expected SOBBH ISGWB amplitude (see, e.g., Babak+2023)

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha and log(Omega0)
Return type:: float

sph_convolve_inj_response_mat(fdata_flag=False)[source]

Function to convolve the sph response matrix with an injected spherical harmonic skymap.

Parameters:: (bool) (fdata_flag)
Return type:: (none)

sph_prior(theta)[source]

Spherical harmonic anisotropic prior transform. Combines a generic spectral prior function with the spherical harmonic priors for the desired lmax.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled for both the spectral and spatial parameters.
Return type:: float

truncated_powerlaw_2par_prior(theta)[source]

Prior function for a stochastic signal search with a 2-parameter truncated power law spectral model.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as log(Omega_ref) and log(f_cut)
Return type:: float

truncated_powerlaw_2par_spectrum(fs, log_omega0, log_fcut)[source]

Function to calculate a tanh-truncated power law spectrum with a set truncation scale and low-f slope.

Parameters:

floats) (fs (array of)
(float) (log_fcut)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting truncated power law spectrum

truncated_powerlaw_3par_prior(theta)[source]

Prior function for a stochastic signal search with a 3-parameter truncated power law spectral model.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), and log(f_cut)
Return type:: float

truncated_powerlaw_3par_spectrum(fs, alpha, log_omega0, log_fcut)[source]

Function to calculate a tanh-truncated power law spectrum with a set truncation scale.

Parameters:

floats) (fs (array of)
(float) (log_fcut)
(float)
(float)

Returns:

spectrum (array of floats)

Return type:

the resulting truncated power law spectrum

truncated_powerlaw_4par_prior(theta)[source]

Prior function for a stochastic signal search with a 4-parameter truncated power law spectral model.

Parameters:: theta (float) – A list or numpy array containing samples from a unit cube.
Returns:: theta – theta with each element rescaled. The elements are interpreted as alpha, log(Omega_ref), log(f_cut), and log(f_scale)
Return type:: float

truncated_powerlaw_4par_spectrum(fs, alpha, log_omega0, log_fcut, log_fscale)[source]

Function to calculate a tanh-truncated power law spectrum.

Parameters:

floats) (fs (array of)
(float) (log_fcut)
(float)
(float)
log_fscale (log of the cutoff scale factor in Hz)

Returns:

spectrum (array of floats)

Return type:

the resulting truncated power law spectrum

truncated_powerlaw_fixedalpha_spectrum(fs, log_omega0, log_fcut, log_fscale)[source]

Function to calculate a tanh-truncated power law spectrum with a set low-f slope.

Parameters:

floats) (fs (array of)
(float) (log_fcut)
(float)
log_fscale (log of the cutoff scale factor in Hz)

Returns:

spectrum (array of floats)

Return type:

the resulting truncated power law spectrum

twothirdspowerlaw_spectrum(fs, log_omega0)[source]

Function to calculate a simple power law spectrum, fixed to the alpha=2/3 prediction for the stellar origin binary background.

Parameters:

floats) (fs (array of)
(float) (log_omega0)

Returns:

spectrum (array of floats)

Return type:

the resulting power law spectrum

wrapper_convolve_inj_response_mat(fdata_flag=False)[source]

A wrapper function for the ISGWB and pixel basis cases, the skymaps are convolved implicitly when calculating the response.

Parameters:: (bool) (fdata_flag)
Return type:: (none)

`models` module

class blip.src.models.Injection(params, inj, fs, f0, tsegmid)[source]

Bases: object

Simulation model. This is a container of submodels used for synthesizing LISA data.

Parameters:

params (dict) – Configuration dictionaries from parse_config.
inj (dict) – Configuration dictionaries from parse_config.
fs (array (nfreqs,)) – Frequency grid in Hz.
f0 (array (nfreqs,)) – Frequency grid scaled by transfer frequency (f0 = fs/2*FSTAR).
tsegmid (array (ntimes,)) – Time grid in seconds.

add_component(name_args)[source]

Wrapper function for the injection component creation process, to allow for parallelization.

Parameters:: (tuple) (name_args)
Returns:: cm (submodel object)
Return type:: Injection component

compute_convolved_spectra(component_name, fs_new=None, channels='11', return_fs=False, imaginary=False)[source]

Wrapper to return the frozen injected detector-convolved GW spectra for the desired channels.

Useful note - these frozen spectra are computed in diag_spectra(), as they are calculated and saved at the analysis frequencies.

Also note that this is meant for plotting purposes only, and includes interpolation/absolute values that are not desirable in a data generation/analysis environment.

Parameters:

(str) (channels)
(array) (fs_new)
(str)
(bool) (imaginary)
(bool)

Returns:

PSD (array) (Power spectral density of the specified channels’ auto/cross-correlation at the desired frequencies.)
fs (array, optional) (The PSD frequencies, if return_fs==True.)

extract_and_save_skymap_data(map_data_path=None)[source]

plot_injected_spectra(component_name, fs_new=None, ax=None, convolved=False, legend=False, channels='11', return_PSD=False, scale='log', flim=None, ymins=None, **plt_kwargs)[source]

Wrapper to plot the injected spectrum component on the specified matplotlib axes (or current axes if unspecified).

Parameters:

(str) (scale)
(array) (fs_new)
axes) (ax (matplotlib)
(bool) (return_PSD)
(bool)
(str)
(bool)
(str)
(tuple) (flim)
(list) (ymins)
(kwargs) (**plt_kwargs)

Returns:

PSD plot on specified axes.
PSD (array, optional) (Power spectral density of the specified channels’ auto/cross-correlation at the desired frequencies.)

plot_skymaps(component_name, save_figures=True, return_mapdata=False, **plt_kwargs)[source]

Function to plot the injected skymaps.

NOTE - will need to be generalized when I add the astro injections

class blip.src.models.Model(params, inj, fs, f0, tsegmid, rmat)[source]

Bases: object

Analysis model. This is a container of submodels with a prior and a likelihood that can be sampled.

Parameters:

params (dict) – Configuration dictionaries from parse_config.
inj (dict) – Configuration dictionaries from parse_config.
fs (array (nfreqs,)) – Frequency grid in Hz.
f0 (array (nfreqs,)) – Frequency grid scaled by transfer frequency (f0 = fs/2*FSTAR).
tsegmid (array (ntimes,)) – Time grid in seconds.
rmat (complex array (nfreqs, ntimes, 3, 3)) – Data correlation matrix.

likelihood(theta)[source]

Unified likelihood function to compare the combined covariance contributions of a generic set of noise/SGWB models to the data.

Parameters:: (list) (theta)
Returns:: loglike (float)
Return type:: resulting joint log likelihood

prior(unit_theta)[source]

Unified prior function to interatively perform prior draws for each submodel in the proper order

Parameters:: (array) (unit_theta)
Returns:: theta (list)
Return type:: transformed prior draws for all submodels in sequence

tree_flatten()[source]

classmethod tree_unflatten(aux_data, children)[source]

blip.src.models.bespoke_inv(A)[source]

compute inverse without division by det; …xv3xc3 input, or array of matrices assumed

Credit to Eelco Hoogendoorn at stackexchange for this piece of wizardy. This is > 3 times faster than numpy’s det and inv methods used in a fully vectorized way as of numpy 1.19.1

https://stackoverflow.com/questions/21828202/fast-inverse-and-transpose-matrix-in-python

LISA instrument and SGWB response

Legacy `geometry` module

class blip.src.geometry.geometry(params: dict, inj: dict, injection: bool)[source]

Bases: sph_geometry

Module containing geometry methods. The methods here include calculation of antenna patterns for a single doppler channel, for the three michelson channels or for the AET TDI channels and calculation of noise power spectra for various channel combinations.

doppler_response(f0, theta, phi, tsegmid, tsegstart)[source]

Calculate antenna pattern/ detector transfer functions for a GW originating in the direction of (theta, phi) for the u doppler channel of an orbiting LISA with satellite position vectors rs1, rs2, rs3. Return the detector response for + and x polarization. Note that f0 is (pi*L*f)/c and is input as an array.

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
theta (phi) – Sky position values.
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

Rplus, Rcross – Plus and cross antenna Patterns for the given sky direction for each time in midpoints.

Return type:

float

isgwb_aet_response(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using A, E and T TDI channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:: f0 (float) – A numpy array of scaled frequencies (see above for def)

Return s

R1, R2 and R3float: Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization.

isgwb_aet_response_parallel(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using A, E and T TDI channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:: f0 (float) – A numpy array of scaled frequencies (see above for def)

Return s

R1, R2 and R3float: Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization.

isgwb_frequency_response_wrapper(ii)[source]

Wrapper function to help with parallelization of the response function calculations.

Arguments

ii (int) : Frequency index

Returns

response_ii : Response matrix in that frequency bin

isgwb_mich_response(f0, tsegmid)[source]

Calculate the Antenna pattern/detector transfer function for an isotropic SGWB using basic michelson channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegstart: A numpy array of segment start times :type tsegstart: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float

Returns:: response_tess – 4D array of covariance matrices for antenna patterns of the three channels, integrated over sky direction and averaged over polarization, across all frequencies and times.
Return type:: float

isgwb_mich_response_parallel(f0, tsegmid)[source]

Parallel version of isgwb_mich_response.

Calculate the Antenna pattern/detector transfer function for an isotropic SGWB using basic michelson channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegstart: A numpy array of segment start times :type tsegstart: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float

Returns:: response_tess – 4D array of covariance matrices for antenna patterns of the three channels, integrated over sky direction and averaged over polarization, across all frequencies and times.
Return type:: float

isgwb_xyz_response(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using X, Y and Z TDI channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array

Parameters:: f0 (float) – A numpy array of scaled frequencies (see above for def)
Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization.
Return type:: float

isgwb_xyz_response_parallel(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using X, Y and Z TDI channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array

Parameters:: f0 (float) – A numpy array of scaled frequencies (see above for def)
Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization.
Return type:: float

lisa_orbits(tsegmid)[source]

Define LISA orbital positions at the midpoint of each time integration segment using analytic MLDC orbits.

Parameters:: tsegmid (array) – A numpy array of the tsegmid for each time integration segment.
Returns:: rs1, rs2, rs3 – Arrays of satellite positions for each segment midpoint in timearray. e.g. rs1[1] is [x1,y1,z1] at t=midpoint[1]=timearray[1]+(segment length)/2.
Return type:: array

pixel_aet_response(f0, tsegmid, skymap_inj)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using A,E,T TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.
Return type:: float

pixel_aet_response_parallel(f0, tsegmid, skymap_inj)[source]

Parallel implementation of pixel_art_response(). Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using A,E,T TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.
Return type:: float

pixel_frequency_response_wrapper(ii)[source]

Wrapper function to help with parallelization of the response function calculations.

Arguments

ii (int) : Frequency index

Returns

response_ii : Response matrix in that frequency bin

pixel_mich_response(f0, tsegmid, skymap_inj)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using basic michelson channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: response_mat – 4D array of covariance matrices for antenna patterns of the three channels, integrated over sky direction and averaged over polarization, across all frequencies and times.
Return type:: float

pixel_mich_response_parallel(f0, tsegmid, skymap_inj)[source]

Parallel version of pixel_mich_response(). Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using basic michelson channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: response_mat – 4D array of covariance matrices for antenna patterns of the three channels, integrated over sky direction and averaged over polarization, across all frequencies and times.
Return type:: float

pixel_xyz_response(f0, tsegmid, skymap_inj)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using X,Y,Z TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.
Return type:: float

pixel_xyz_response_parallel(f0, tsegmid, skymap_inj)[source]

Parallel implementation of pixel_xyz_response(). Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using X,Y,Z TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.
Return type:: float

unconvolved_pixel_aet_response(f0, tsegmid, masked_skymap)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using A,E,T TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param masked_skymap: A pixel map in healpy ordering of GW power on the sky, masked to cover only the areas of interest. :type masked_skymap: healpy pixel map

Returns:: aet_response_mat – 4D array of covariance matrices for antenna patterns of the three channels, averaged over polarization, across all frequencies, times, and sky directions of interest.
Return type:: float

unconvolved_pixel_aet_response_parallel(f0, tsegmid, masked_skymap)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using A,E,T TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param masked_skymap: A pixel map in healpy ordering of GW power on the sky, masked to cover only the areas of interest. :type masked_skymap: healpy pixel map

Returns:: aet_response_mat – 4D array of covariance matrices for antenna patterns of the three channels, averaged over polarization, across all frequencies, times, and sky directions of interest.
Return type:: float

unconvolved_pixel_frequency_response_wrapper(ii)[source]

Wrapper function to help with parallelization of the response function calculations.

Arguments

ii (int) : Frequency index

Returns

response_ii : Response matrix in that frequency bin

unconvolved_pixel_mich_response(f0, tsegmid, masked_skymap)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using basic michelson channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param masked_skymap: A pixel map in healpy ordering of GW power on the sky, masked to cover only the areas of interest. :type masked_skymap: healpy pixel map

Returns:: response_mat – 4D array of covariance matrices for antenna patterns of the three channels, averaged over polarization, across all frequencies, times, and sky directions of interest.
Return type:: float

unconvolved_pixel_mich_response_parallel(f0, tsegmid, masked_skymap)[source]

Parallel version of unconvolved_pixel_mich_response(). Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using basic michelson channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param skymap: A pixel map in healpy ordering of GW power on the sky :type skymap: healpy pixel map

Returns:: response_mat – 4D array of covariance matrices for antenna patterns of the three channels, integrated over sky direction and averaged over polarization, across all frequencies and times.
Return type:: float

unconvolved_pixel_xyz_response(f0, tsegmid, masked_skymap)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using X,Y,Z TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param masked_skymap: A pixel map in healpy ordering of GW power on the sky, masked to cover only the areas of interest. :type masked_skymap: healpy pixel map

Returns:: xyz_response_mat – 4D array of covariance matrices for antenna patterns of the three channels, averaged over polarization, across all frequencies, times, and sky directions of interest.
Return type:: float

unconvolved_pixel_xyz_response_parallel(f0, tsegmid, masked_skymap)[source]

Calculate the Antenna pattern/detector transfer function for a pixel-basis anisotropic SGWB using X,Y,Z TDI channels. Note that we only evaluate the response to sky directions with power in them. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array. :param f0: A numpy array of scaled frequencies (see above for def) :type f0: float :param tsegmid: A numpy array of segment midpoints :type tsegmid: float :param masked_skymap: A pixel map in healpy ordering of GW power on the sky, masked to cover only the areas of interest. :type masked_skymap: healpy pixel map

Returns:: xyz_response_mat – 4D array of covariance matrices for antenna patterns of the three channels, averaged over polarization, across all frequencies, times, and sky directions of interest.
Return type:: float

`fast_geometry` module

class blip.src.fast_geometry.fast_geometry(params, nthreads=1, use_gpu=False)[source]

Bases: sph_geometry

Module containing fast, unified geometry methods. The methods here include calculation of antenna patterns for the Michelson, XYZ, and AET TDI channels. The calculations are multilayered and abstracted so as to avoid (almost) all repeated calculations when considering multiple submodels/injection components.

Parameters:

(dict) (params)
(int) (nthreads)

apply_response_slice_convolutions(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii)[source]: Helper function to perform the sky integrals/masking as needed.

assign_responses_to_submodels()[source]: Method to attach the unique response functions to their (potentially non-unique) respective submodels.

calculate_response_functions(f0, tsegmid, submodels, tdi_lev, plot_flag=False)[source]

Prototype federated function for response function calculations, designed to eliminate all redundant calculations when computing the response functions for multiple submodels.

Parameters:

array) (tsegmid (float)
array)
object) (submodels (models.submodel)
(str) (tdi_lev)
(bool) (plot_flag)

Returns:

response_mat(s) (numpy array/s)

Return type:

the appropriate response function for each submodel, attached to that submodel (the function does not return the response arrays directly)

construct_aet_response_mat(xyz_response_mat)[source]

Calculate the Antenna pattern/detector transfer function for a generic isotropic or anisotropic SGWB using A,E,T TDI channels.

Parameters:: array) (xyz_response_mat (float)
Returns:: aet_response_mat (float array)
Return type:: Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averagedover polarization.

frequency_response_wrapper(ii)[source]

Wrapper function to help with parallelization of the response function calculations.

Parameters:: (int) (ii)
Returns:: response_ii (float array)
Return type:: Response matrix in that frequency bin

get_aet_from_xyz()[source]: Compute the LISA response functions for the AET Time-Delay Interferometry channels from the XYZ channels.

get_geometric_Fplus_Fcross(arm_hat_ij, mhat_prod, nhat_prod)[source]

Wrapper function to reduce repeated calculations for the geometric components of the antenna pattern functions

Parameters:

arrays) (nhat_prod (float)
mhat_prod (outer products of phi hat and theta hat, respectively)
arrays)

Returns:

Fplus_ij, Fcross_ij (float arrays)

Return type:

Geometric components of the plus and cross antenna pattern functions for LISA arms i and j

get_xyz_from_michelson()[source]: Compute the LISA response functions for the XYZ Time-Delay Interferometry channels from the Michelson channels.

isgwb_wrapper(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii, ghost_arg=None)[source]: Wrapper function to take sky integral and return response array slice

pix_convolved_asgwb_wrapper(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii, skymap)[source]: Wrapper function to convolve with pixel-basis skymap and return response array slice

pix_masked_unconvolved_asgwb_wrapper(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii, mask)[source]: Wrapper function which just returns its inputs, as the unconvolved case doesn’t integrate over the sky

pix_unconvolved_asgwb_wrapper(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii, ghost_arg=None)[source]: Wrapper function which just returns its inputs, as the unconvolved case doesn’t integrate over the sky

sph_asgwb_wrapper(F1_ii, F2_ii, F3_ii, F12_ii, F13_ii, F23_ii, ghost_arg=None)[source]: Wrapper function to convolve with Ylms and return response array slice

unpack_wrapper(ii, Rf)[source]

Wrapper function to generically assigen the calculated per-frequency response function slices to their appropriate unique response array.

Parameters:

(int) (ii)
(array) (Rf)

blip.src.fast_geometry.get_model_responses(Model_obj)[source]

Wrapper function to call methods for computing the response function.

This exists so that the methods in question can be easily jettisoned when the Model object in JAX JIT compilation.

The experimental `faster_geometry` module

Functions to compute the LISA response to stochastic backgrounds.

The main procedure is calculate_response_functions(), which is used to interface with the rest of BLIP.

The main objective is to compute the LISA response for a given sky direction, frequency and time. This is done in mich_response_unconvolved().

blip.src.faster_geometry.calculate_response_functions(freqs, times, submodels, params, plot_flag=False)[source]

Compute SGWB responses for each submodel and write it to the submodel attributes.

This procedure avoids redundant calculations by only computing the response for a given pixel once. It also uses a sparse time-frequency grid and linearly interpolates the results in-between.

It mirrors the functionality in fast_geometry.calculate_response_functions().

Parameters:

freqs (array (nfreqs,)) – frequencies.
times (array (ntimes,)) – times.
submodels (list[submodel]) –
The submodels whose response matrices will be computed. They must be in pixel basis and have a skymap attribute, not assumed normalized.

These objects also act as output parameters: after the procedure has run, each of them will receive a response matrix as an attribute, which is an array of shape (3, 3, nfreqs, ntimes).

The attributes used for the output depend on the value of plot_flag: the response matrix will be written to sm.response_mat and sm.inj_response_mat if it is False, or to sm.fdata_response_mat if it is True.
params (dict) – Parameter dictionary from parse_config().
plot_flag (bool, optional) – If True, don’t write the output response matrix to sm.response_mat or sm.inj_response_mat, but to submodel.fdata_response_mat. This is useful for plotting simulation and analysis models together (hence the name). Defaults to False.

Warning

This is experimental, tested only for healpix responses.

blip.src.faster_geometry.compute_orbits(times, use_lisaorbits=False, betaphase=0)[source]

Compute orbit information at specified time array.

Parameters:

times (array 1D) – Times at which the positions of the spacecraft and link vectors should be computed.

Returns:

array 1D – the input times array
array (3, ntimes, 3) – Spacecraft positions in ecliptic cartesian coordinates as an array, where the first dimension specifies the spacecraft.
array (6, ntimes, 3) – Single link unit vectors in lisaorbits order.

blip.src.faster_geometry.get_arm_orientations(t, sc, orbits)[source]

Get unit vectors for left and right arm of a given spacecraft.

The numbering convention is the same as in lisaorbits.

Parameters:

t (float) – time
sc (int) – Spacecraft number (1, 2, 3). Must be known at JAX trace-time.
orbits (tuple) – orbital information from compute_orbits().

Returns:

array (3,) – Unit vector pointing away from the given spacecraft (1 -> 3, 2 -> 1, 3 -> 2).
array (3,) – Unit vector pointing away from the given spacecraft (1 -> 2, 2 -> 3, 3 -> 1).

blip.src.faster_geometry.get_link_vectors(t, orbits)[source]

Look up link unit vectors in orbits, in a way that is jax-traceable.

This function does not perform interpolation. It will return the ‘last known’ unit vectors.

Parameters:

t (float) – time.
orbits (tuple) – orbital information returned by compute_orbits().

Returns:

The unit vectors in 3D (second axis) for each of the single links (first axis) in lisaorbits order.

Return type:

array (6, 3)

blip.src.faster_geometry.get_ortho_basis_ecliptic_3d(lam, beta)[source]

Get right-handed orthonormal basis (n, l, m).

This is the basis in Romano & Cornish 2017 eq (2.4).

Parameters:

lam (float) – ecliptic longitude
beta (float) – ecliptic latitude

Returns:

A tuple of arrays of shape (3,).

Return type:

tuple (n, l, m)

blip.src.faster_geometry.get_sc_positions(t, orbits)[source]

Look up S/C positions in orbits, in a way that is jax-traceable.

This function does not perform interpolation. It will return the ‘last known’ position of the spacecraft.

Parameters:

t (float) – time.
orbits (tuple) – orbital information returned by compute_orbits().

Returns:

The positions of the three spacecraft (first axis) in cartesian coordinates (second axis).

Return type:

array (3, 3)

blip.src.faster_geometry.get_vecs_all_sky(nside)[source]

Compute array of all unit vectors in the sky.

Parameters:: nside (int) – Healpix nside. Should be a power of 2.
Returns:: Unit vectors for each healpix direction.
Return type:: array (npix, 3)

blip.src.faster_geometry.mich_antenna_pattern(t, f, n, polarization: str, channel, orbits)[source]

Compute Michelson (TDI gen 0) antenna pattern.

Checked against Banagiri+21 and Romano & Cornish 2017.

Parameters:

t (float) – time
f (float) – frequency
n (array (3,)) – Unit vector in the direction of the GW source.
polarization (str) – Should be “plus” or “cross”. Must be known at JAX trace time.
channel (int) – Channel index 0, 1, 2. Must be known at JAX trace time.
orbits (tuple) – orbital information returned by compute_orbits().

Returns:

The antenna pattern.

Return type:

complex

blip.src.faster_geometry.mich_detector_tensor(f, u, v, n, r)[source]

Michelson channel detector tensor.

Checked against Banagiri+21 eq (15).

Parameters:

f (float) – frequency in Hz
u (array (3,)) – normalized vector in the direction of the first arm
v (array (3,)) – normalized vector in the direction of the second arm
n (array (3,)) – normalized vector in the direction of the GW source
r (array (3,)) – position of vertex S/C in barycentric ecliptic cartesian coordinates

Returns:

detector tensor

Return type:

complex array (3, 3)

blip.src.faster_geometry.mich_response_unconvolved(t, f, n, orbits)[source]

Unconvolved Michelson (TDI gen 0) sky SGWB response.

Parameters:

t (float) – time
f (float) – frequency
n (array (3,)) – normalized vector in the direction of the GW source.
orbits (tuple) – orbital information returned by compute_orbits().

Returns:

Unconvolved response matrix for the three data channels.

Return type:

complex array (3, 3)

Notes

The quantity computed here is exactly

\[\frac{1}{2} \sum_{A=+,\times}\left(F_I^A(f, \mathbf{n})^* F_J^A(f, \mathbf{n}) \right)\]

where $I$ and $J$ stand for TDI channels, and $F_I^A$ are antenna pattern functions.

This is integrated against the sky map to produce the GW time-frequency correlation matrix.

The convention here agrees with Criswell+25 eq. (6) (but for a complex conjugate), which is a corrected version of Banagiri+21 eq. (18).

blip.src.faster_geometry.timing_transfer_fn(f, costheta)[source]

Timing transfer function for two-way photon propagation.

Checked against Banagiri+21 eq (16) and Cornish & Rubbo 2003 eq (37). Also agrees with Romano & Cornish 2017 eq (5.27) up to a constant 2L/c. This seems due to the conversion between strain and timing measurements, eq (5.4) in the living review.

Parameters:

f (float) – frequency in Hz
costheta (float) – cosine of angle between arm and sky direction

Returns:

the transfer function.

Return type:

complex

`sph_geometry` module

class blip.src.sph_geometry.sph_geometry[source]

Bases: clebschGordan

asgwb_aet_response(f0, tsegmid, set_almax=None)[source]

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using X,Y,Z TDI channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

asgwb_aet_response_parallel(f0, tsegmid, set_almax=None)[source]

Parallel implementation of asgwb_aet_response.

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using X,Y,Z TDI channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

asgwb_frequency_response_wrapper(ii)[source]

Wrapper function to help with parallelization of the response function calculations.

Arguments

ii (int) : Frequency index

Returns

response_ii : Response matrix in that frequency bin

asgwb_mich_response(f0, tsegmid, set_almax=None)[source]

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using michelson channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

asgwb_mich_response_parallel(f0, tsegmid, set_almax=None)[source]

Parallel implementation of asgwb_mich_response.

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using michelson channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

asgwb_xyz_response(f0, tsegmid, set_almax=None)[source]

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using X,Y,Z TDI channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

asgwb_xyz_response_parallel(f0, tsegmid, set_almax=None)[source]

Parallel implementation of asgwb_xyz_response.

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using X,Y,Z TDI channels, and using a spherical harmonic decomposition. Note that the response function to power is integrated over sky direction with the appropriate spherical harmonics, and averaged over polarozation. The angular integral is numerically done by divvying up the sky into a healpix grid.

Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (float) – A numpy array of the time series midpoints.
set_almax (int) – Allows the user to manually set the almax used in the response function calculations. If None, almax will default to the globally set self.almax. Otherwise almax will be the value given.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.

Return type:

float

`makeLISAdata` module

class blip.src.makeLISAdata.LISAdata(params, inj)[source]

Bases: object

Class for lisa data. Includes methods for generation of gaussian instrumental noise, and generation of isotropic stochastic background. Signal models should be added as methods here.

add_sgwb_data(injmodel, key, tbreak=0.0)[source]

Function to simulate the SGWB time series, given a spectrum.

Parameters:

component) (injmodel (Injection)
(int) (key)

process_external_data()[source]: Read external data from params.datafile assuming it is in the file format specified by params.datafileformat.

read_spectrum()[source]

Read an input frequency domain data file. Returns the fourier transform of the data from the three channels and an array of reference frequencyes

Returns:: rA, rE, rT, fdata
Return type:: float

tser2fser(h1, h2, h3, timearray)[source]

Convert time domain data to fourier domain and return ffts. The convention is that the the ffts are divided by the sampling frequency and corrected for windowing. A hann window is applied by default when moving to the fourier domain. The ffts are also normalized so that thier square gives the PSD.

Parameters:

h1 (float) – time series data for the three input channels
h2 (float) – time series data for the three input channels
h3 (float) – time series data for the three input channels
timearray (float) – times corresponding to data in h1, h2, h3

Returns:

r1, r2, r3 (float) – frequency series data for the three input channels
fdata (float) – Reference frequency series
tsegstart (float) – Segmented time array giving segment start points
tsegmid (float) – Segmented time array giving segment midpoints

blip.src.makeLISAdata.cholesky_symm(m)[source]

Cholesky decomposition of response matrix.

Parameters:

m (complex array (3, 3, ...)) –

Response matrix. Assumed to have the form

\[\begin{split}\begin{pmatrix}A & B & B \\ B^* & A & B \\ B^* & B^* & A\end{pmatrix}\end{split}\]

Returns:

Lower-triangular matrix L such that m = L L^H.

Return type:

complex array (3, 3, …)

blip.src.makeLISAdata.get_data_length(dur, seglen, fs)[source]

Compute data array length.

The data is required to have this length, aligned with segment sizes, for use in LISAdata.tser2fser().

Parameters:

dur (float) – total data duration in seconds, as set by the duration parameter.
seglen (float) – segment length set by seglen.
fs (float) – data sampling rate in Hz. Inverse of dt.

Returns:

length of analysis data array.

Return type:

int

blip.src.makeLISAdata.get_data_tf_grid(dur, seglen, fs)[source]

Compute time-frequency grid parameters used for data analysis.

Parameters:

dur (float) – total data duration in seconds, as set by the duration parameter.
seglen (float) – segment length set by seglen.
fs (float) – data sampling rate in Hz. Inverse of dt.

Returns:

int – nsegs: number of (overlapping) segments
int – Nperseg: number of time-domain points in each segment

`instrNoise` module

class blip.src.instrNoise.instrNoise[source]

Bases: object

Class with methods to calcualte instrumental noises

aet_noise_spectrum(freqs, f0, Np=9e-42, Na=3.6e-49)[source]

Calculates A, E, and T channel noise spectra for a stationary lisa. Following the defintions in Adams & Cornish, http://iopscience.iop.org/article/10.1088/0264-9381/18/17/308

Parameters:

freqs (float) – A numpy array of frequencies
(optional) (Na) – Position noise value
(optional) – Acceleration noise level

Returns:

SAA, SEE, STT – Frequencies arrays with the noise PSD for the A, E and T TDI channels

Return type:

float

freqdomain_gaussianData(Sh, freqs, fs=1, dur=100000.0)[source]

Script to generate freq Domain gaussian data of a given spectral density.

Parameters:

Sh ((float)) – A frequency array with the desired power spectral density
freqs ((float)) – An array with corresponding frequencies to Sh
fs ((float)) – SampleRate in Hz
dur ((int)) – Duration in seconds

Returns:

ht (float)
frequency domain gaussian.

fundamental_noise_spectrum(freqs, Np=9e-42, Na=3.6e-49)[source]

Compute position and acceleration noise PSDs in strain units.

Position noise is also known as optical metrology system (OMS) noise, while acceleration noise is test mass (TM) noise. There is one of these per MOSA.

The default values agree with the 2017 LISA proposal noise budget. In the white noise regions corresponding to the coefficients, these are 15 pm/√Hz for OMS and 3 fm/s²/√Hz for TM.

Parameters:

freqs (array) – frequencies in Hz
Np (float, optional) – Coefficient of position noise as equivalent strain, by default 9e-42
Na (float, optional) – Coefficient of acceleration noise as equivalent strain, by default 3.6e-49

Returns:

array – Sp: position noise PSD in strain units
array – Sa: acceleration noise PSD in strain units

gaussianData(Sh, freqs, fs=1, dur=100000.0)[source]

Script for generation time series noise drawn from a gaussian process of a given spectral density. Adapted from gaussian_noise.m from stamp

Parameters:

Sh ((float)) – A frequency array with the desired power spectral density
freqs ((float)) – An array with corresponding frequencies to Sh
fs ((float)) – SampleRate in Hz
dur ((int)) – Duration in seconds

Returns:

ht (float)
Array with time series data of duration, dur with the prescribed spectrum Sh

gen_aet_noise()[source]

Generate interferometric A, E and T channel TDI (time-domain) noise, using freqDomain.fundamental_noise_spectrum

Returns:: h1_noi, h2_noi, h3_noi – Time series data for the three TDI channels
Return type:: float

gen_michelson_noise()[source]

Generate interferometric michelson (time-domain) noise, using freqDomain.fundamental_noise_spectrum

Returns:: h1, h2, h3 – Time series data for the three michelson channels
Return type:: float

gen_noise_cov()[source]

Generate interferometric (time-domain) noise, using a frequency domain covariance spectrum matrix rather than time delays in time domain. ———

h1_noi, h2_noi, h3_noifloat: Time series data for the three TDI channels

gen_xyz_noise()[source]

Generate interferometric A, E and T channel TDI (time-domain) noise, using freqDomain.fundamental_noise_spectrum

Returns:: h1_noi, h2_noi, h3_noi – Time series data for the three TDI channels
Return type:: float

mich_noise_spectrum(freqs, f0, Np=9e-42, Na=3.6e-49)[source]

Calculates michelson channel noise spectra for a stationary lisa. Following the defintions in Adams & Cornish, http://iopscience.iop.org/article/10.1088/0264-9381/18/17/308. We assume that there is no phase noise.

Parameters:

freqs (float) – A numpy array of frequencies
(optional) (Na) – Position noise value
(optional) – Acceleration noise level

Returns:

SAA, SEE, STT – Frequencies arrays with the noise PSD for the A, E and T TDI channels

Return type:

float

xyz_noise_spectrum(freqs, f0, Np=9e-42, Na=3.6e-49)[source]

Calculates X,Y,Z channel noise spectra for a stationary lisa. Following the defintions in Adams & Cornish, http://iopscience.iop.org/article/10.1088/0264-9381/18/17/308

Parameters:

freqs (float) – A numpy array of frequencies
(optional) (Na) – Position noise value
(optional) – Acceleration noise level

Returns:

SAA, SEE, STT – Frequencies arrays with the noise PSD for the A, E and T TDI channels

Return type:

float

`orbitinglisa` module

class blip.src.orbitinglisa.orbitinglisa[source]

Bases: object

Module containing geometry methods accounting for an orbiting satellite. The methods here include a base orbit, a single doppler channel, for the three michelson channels or for the AET TDI channels and calculation of noise power spectra for various channel combinations.

isgwb_oaet_response(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using A, E and T TDI channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarization. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note that f0 is (pi*L*f)/c and is input as an array.

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarizationevaluated at the midpoint of each time segment.

Return type:

arrays

isgwb_omich_response(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions for an orbiting LISA to an isotropic SGWB using basic michelson channels. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarization. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note also that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization, evaluated at the midpoint of each time segment.

Return type:

float

isgwb_oxyz_response(f0, tsegmid)[source]

Calcualte the Antenna pattern/ detector transfer function functions to an isotropic SGWB using X, Y and Z TDI channels for an orbiting LISA. Note that since this is the response to an isotropic background, the response function is integrated over sky direction and averaged over polarozation. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization, evaluated at the midpoint of each time segment.

Return type:

float

lisa_orbits(tsegmid)[source]

Define LISA orbital positions at the midpoint of each time integration segment using analytic MLDC orbits.

Parameters:: tsegmid (array) – A numpy array of the tsegmid for each time integration segment.
Returns:: rs1, rs2, rs3 – Arrays of satellite positions for each segment midpoint in timearray. e.g. rs1[1] is [x1,y1,z1] at t=midpoint[1]=timearray[1]+(segment length)/2.
Return type:: array

orbiting_aet_response(f0, theta, phi, tsegmid)[source]

Calculate Antenna pattern/ detector transfer functions for a GW originating in the direction of (theta, phi) for the A, E and T TDI channels of an orbiting LISA. Return the detector responses for + and x polarization. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
theta (phi) – Sky position values.
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

RAplus, RAcross, REplus, REcross, RTplus, RTcross – Plus and cross antenna Patterns for the given sky direction for the three channels for each time in midpoints.

Return type:

arrays

orbiting_doppler_response(f0, theta, phi, tsegmid)[source]

Calculate antenna pattern/ detector transfer functions for a GW originating in the direction of (theta, phi) for the u doppler channel of an orbiting LISA with satellite position vectors rs1, rs2, rs3. Return the detector response for + and x polarization. Note that f0 is (pi*L*f)/c and is input as an array.

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
theta (phi) – Sky position values.
tsegmid (array) – A numpy array of the midpoints for each time integration segment.
rs1 (array) – Satellite position vectors.
rs2 (array) – Satellite position vectors.
rs3 (array) – Satellite position vectors.

Returns:

Rplus, Rcross – Plus and cross antenna Patterns for the given sky direction for each time in midpoints.

Return type:

float

orbiting_michelson_response(f0, theta, phi, tsegmid)[source]

Calculate Antenna pattern/ detector transfer function for a GW originating in the direction of (theta, phi) at a given time for the three Michelson channels of an orbiting LISA. Return the detector response for + and x polarization. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:

f0 (float) – A numpy array of scaled frequencies (see above for def)
theta (phi) – Sky position values.
rs1 (arrays) – Satellite position vectors.
rs2 (arrays) – Satellite position vectors.
rs3 (arrays) – Satellite position vectors.
tsegmid (array) – A numpy array of the midpoints for each time integration segment.

Returns:

R1plus, R1cross, R2plus, R2cross, R3plus, R3cross – Plus and cross antenna Patterns for the given sky direction for the three channels for each time in midpoints.

Return type:

arrays

tdi_aniso_sph_sgwb_response(f0)[source]

Calculate the Antenna pattern/ detector transfer function functions to acSGWB using A, E and T TDI channels, and using a spherical harmonic decomposition. Note that the response function is integrated over sky direction with the appropriate legandre polynomial, and averaged over polarozation. Finally note that the spherical harmonic coeffcients correspond to strain sky distribution, while the legandre polynomials describe the power sky. The angular integral is a linear and rectangular in the cos(theta) and phi space. Note that f0 is (pi*L*f)/c and is input as an array

Parameters:: f0 (float) – A numpy array of scaled frequencies (see above for def)
Returns:: R1, R2 and R3 – Antenna Patterns for the given sky direction for the three channels, integrated over sky direction and averaged over polarization. The arrays are 2-d, one direction corresponds to frequency and the other to the l coeffcient.
Return type:: float

Astrophysics

`astro` module

class blip.src.astro.Galaxy_Model(nside, grid_spec='interval', grid_res=0.33, gal_rad=16, gal_height=8, max_rh=4, max_zh=2, fix_rh=None)[source]

Bases: object

Class to support parameterized inference of the Galactic white dwarf binary spatial distribution.

mw_mapmaker_1par(zh)[source]

Generate a galactic white dwarf binary foreground modeled after Breivik et al. (2020), consisting of a bulge + disk. The default values are those given in Breivik+20; the model itself traces back to electromagnetic studies of the Galaxy by McMillan et al. (2011), Bissantz and Gerhard (2002), and Binney et al. (1997). zh is the vertical scale height in kpc. The distribution is azimuthally symmetric in the galactocentric frame.

Designed for speed, as it is intended for use during sampling. Relies on pre-computed galaxy grid that is produced as part of Galaxy_Model() initialization.

Returns:: skymap – Healpy GW power skymap of the Milky Way white dwarf binary distribution.
Return type:: float

mw_mapmaker_2par(rh, zh)[source]

Generate a galactic white dwarf binary foreground modeled after Breivik et al. (2020), consisting of a bulge + disk. The default values are those given in Breivik+20; the model itself traces back to electromagnetic studies of the Galaxy by McMillan et al. (2011), Bissantz and Gerhard (2002), and Binney et al. (1997). rh is the radial scale height in kpc, zh is the vertical scale height in kpc. The distribution is azimuthally symmetric in the galactocentric frame.

Designed for speed, as it is intended for use during sampling. Relies on pre-computed galaxy grid that is produced as part of Galaxy_Model() initialization.

Returns:: skymap – Healpy GW power skymap of the Milky Way white dwarf binary distribution.
Return type:: float

class blip.src.astro.Population(params, inj, frange, popdict, seed=None, map_only=False)[source]

Bases: object

Class for handling binary populations. Includes methods for loading, parsing, calculating physical quantities, and preparing data for use in some of the makeLISAdata.py signal simulation methods.

Sgw_wrapper(frange, spoof_arg=None)[source]

This is a wrapper function to allow the population spectrum to play well with some of the generic Injection-handling code.

Evaluated at the injection frequencies.

Sgw_wrapper_true(frange, spoof_arg=None)[source]: This is a wrapper function to allow the population spectrum to play well with some of the generic Injection-handling code. Evaluated at the data frequencies.

classmethod file2map(popfile, nside, df, fmin, fmax, SNR_cut=7, **read_csv_kwargs)[source]

Wrapper function to get a skymap directly from a population catalogue file.

Parameters:

lats (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
longs (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
PSDs (array of floats) – Corresponding catalogue binary PSDs (assumed monochromatic)
t_obs (astropy Quantity, time units) – Observation time.
nside (int) – Healpix nside to use for skymap. Must be power of 2 < 2**32.

Returns:

Healpix skymap of GW power on the sky logskymap (array of floats) : Healpix skymap of log GW power on the sky

Return type:

skymap (array of floats)

classmethod file2spec(popfile, frange, t_obs, SNR_cut=7, plot=False, return_median=False, **read_csv_kwargs)[source]

Wrapper function to calculate the foreground spectrum directly from a population catalogue file.

Parameters:

popfile (str) – ‘/path/to/binary/population/data/file.csv’
frange (1D array of floats) – Frequencies at which to calculate binned PSD
t_obs (astropy Quantity, time units) – Observation time.

Returns:

Resulting PSD of unresolved binary background/foreground for all f in frange

Return type:

fg_PSD (array of floats)

static filter_by_snr(data, SNRs, SNR_cut=7, get_type='unresolved')[source]

Function to filter DWD data by SNR. Can return either unresolved (SNR < SNR_cut) or resolved (SNR > SNR_cut) binaries.

Parameters:

data (1D array of floats) – Binary population data of your choice (or list thereof), corresponding to the given SNRs
SNRs (1D array of floats) – SNR value for each system corresponding to data
SNR_cut (float) – Value of SNR that delineates resolved and unresolved binaries. Default is SNR = 7.
get_type (str) – Whether to return the resolved or unresolved binaries. Default is unresolved.

Returns:

Filtered arrays of frequencies and strains.

Return type:

data_filt

static gen_summed_map(lats, longs, PSDs, nside, return_log=False)[source]

Function to get a skymap from a collection of binary sky coordinates and (monochromatic) PSDs. Note that this function will process all binaries given to it; SNR filtering must be done beforehand.

Parameters:

lats (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
longs (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
PSDs (array of floats) – Corresponding catalogue binary PSDs (assumed monochromatic)
nside (int) – Healpix nside to use for skymap. Must be power of 2 < 2**32.
return_log (bool) – If True, also return the log skymap, with slight zero-buffering.

Returns:

Healpix skymap of GW power on the sky logskymap (array of floats) : Healpix skymap of log GW power on the sky

Return type:

skymap (array of floats)

classmethod gen_summed_spectrum(fs, hs, cos_incs, frange, t_obs, plot=False, saveto=None, return_median=False)[source]

Function to calculate the foreground spectrum arising from a set of monochromatic strains and associated frequencies.

Binaries passed to this function are assumed to all contribute to the foreground.

Parameters:

fs (1D array of floats) – Binary frequencies. Assumed monochromatic.
hs (1D array of floats) – Binary strain amplitudes.
cos_incs (1D array of floats) – Cosine of the binary inclinations
t_obs (astropy Quantity, time units) – Observation time.
frange (1D array of floats) – Frequencies at which to calculate binned PSD

Returns:

Resulting PSD of unresolved binary background/foreground for all f in frange

Return type:

fg_PSD_binned (array of floats)

static get_binary_psd(hs, cos_incs, df)[source]

Function to calculate PSD of catalogue binaries. Assumed monochromatic.

We assume a definition of amplitude such that A = h0, so

h(t) = h+(t) + hx(t),

h+(t) = (1+cos^2(i)) * A * cos(2*omega*t + phi0),

and

hx(t) = 2 * cosi * A * sin(2omega*t + phi0).

The strain power is then

<h(t)^2> = (1+cos^2(i))^2 * A^2 + 4 * cos^2(i) * A^2

which for optimal inclination (i=0, face-on) yields

<h(t)^2> = 8A^2.

The PSD contribution from the monochromatic binary at frequency resolution df = 1/Tobs is then

PSD = (1/df) * < h(t)^2 >

which in the optimal inclination case is

PSD = 8 * (1/df) * A^2.

Note that by combining the + and x contributions prior to convolution with the LISA response functions, we implicitly assume that the overall population produces an unpolarized stochastic signal. This is statistically true for any stochastic signal produced by a population of binaries with uniformly distributed inclinations, but the assumption may break down in some cases.

Also note that there is an alternate definition for the amplitude A, such that A = 2h0. If used here, this will result in an erroneous factor of 4 in the PSD.

Parameters:

hs (1D array of floats) – Binary strains.
cos_incs (1D array of floats) – cosine of binary inclinations
df (astropy Quantity, frequency units) – Binning frequency resolution.

Returns:

Monochromatic PSDs for each binary

Return type:

binary_psds (1D array of floats)

classmethod get_snr(fs, hs, cos_incs, t_obs, noise_PSD='default')[source]

Function to get SNRs of catalogue binaries, given their frequencies/strains, observation time, and a detector PSD. Assumes monochromatic systems.

Parameters:

fs (1D array of floats) – Binary frequencies. Assumed monochromatic.
hs (1D array of floats) – Binary strains.
t_obs (astropy Quantity, time units) – Observation time.
PSD (varies) – If ‘default’, uses the Legwork LISA PSD without confusion noise. If an array, used as the detector PSD at each frequency value in fs

Returns:

Binary SNRs.

Return type:

SNRs (array of floats)

static load_population(popfile, fmin, fmax, coldict={'f': 'f', 'h': 'h', 'lat': 'lat', 'long': 'long'}, unitdict={'f': Unit('Hz'), 'lat': Unit('rad'), 'long': Unit('rad')}, sep=' ', seed=None, **read_csv_kwargs)[source]

Function to load a population file and store relevant data. For now this assumes columns with labels [‘f’,’h’,lat’,’long’]. Assumes LEGWORK definition of binary strain (see https://arxiv.org/abs/2111.08717)

Parameters:

popfile (str) – ‘/path/to/binary/population/data/file.csv’
fmin (float) – Minimum analysis frequency
fmax (float) – Maximum analysis frequency
coldict (dict) – Dictionary explaining which columns correspond to which quantities (allows for users to specify different column names)
unitdict (dict) – Dictionary specifying units of each column.
sep (str) – File delimiter. Overwritten if sep or delimiter is specified in **read_csv_kwargs. Default is ‘ ‘.
**read_csv_kwargs – Optional keyword arguments to be passed to pd.read_csv()

Returns:

Binary frequencies hs (array) : Binary strain amplitudes lats (array) : Binary latitudes in degrees longs (array) : Binary longitudes in degrees

Return type:

fs (array)

omegaf_wrapper(fs, spoof_arg=None)[source]: This is a wrapper function to allow the pupulation spectrum to play well with some of the generic Injection-handling code.

classmethod pop2map(pop, nside, df, fmin, fmax, SNR_cut=7)[source]

Function to get a skymap from a catalogue of binaries.

Parameters:

lats (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
longs (arrays of floats) – Latitudes and longitudes of catalogue binaries. IMPORTANT: Must be given in ecliptic coordinates and units of degrees!
PSDs (array of floats) – Corresponding catalogue binary PSDs (assumed monochromatic)
t_obs (astropy Quantity, time units) – Observation time.
nside (int) – Healpix nside to use for skymap. Must be power of 2 < 2**32.

Returns:

Healpix skymap of GW power on the sky logskymap (array of floats) : Healpix skymap of log GW power on the sky

Return type:

skymap (array of floats)

classmethod pop2spec(pop, frange, t_obs, SNR_cut=7, plot=False, return_median=False, saveto=None)[source]

Function to calculate the foreground spectrum arising from a population catalogue of unresolved DWD binaries.

Parameters:

popfile (str) – ‘/path/to/binary/population/data/file.csv’
frange (1D array of floats) – Frequencies at which to calculate binned PSD
t_obs (astropy Quantity, time units) – Observation time.
SNR_cut (float) – SNR above which a binary will be assumed to be individually resolveable and subtracted. Default SNR=7.
return_median (bool) – If True, also return a running median of the spectrum (Useful for smoothing, plotting). Default False.

Returns:

Resulting PSD of unresolved binary background/foreground for all f in frange

Return type:

fg_PSD (array of floats)

rebin_PSD(fs_new)[source]: Function to correctly interpolate the population spectrum to new frequencies without violating conservation of energy

blip.src.astro.generate_galactic_bulge(nside, rho_c=1, r_cut=2.1, r0=0.075, alpha=1.8, q=0.5)[source]

Generate the bulge component of a galactic white dwarf binary foreground modeled after the Breivik et al. (2020) bulge+disk spatial model. The default values are those given in Breivik+20; the model itself traces back to electromagnetic studies of the Galaxy by McMillan et al. (2011), Bissantz and Gerhard (2002), and Binney et al. (1997). The distribution is azimuthally symmetric in the galactocentric frame.

Parameters:

(int) (nside)
(float) (q)
(float)
(float)
1.8. (alpha (float) Bulge density power law index. Default)
(float)

Returns:

astro_map (float) – Healpy GW power skymap of the Galactic bulge contribution to the foreground.
log_DWD_FG_map (float) – Healpy log GW power skymap. For plotting purposes.

blip.src.astro.generate_galactic_disk(rh, zh, nside)[source]

Generate the disk component of a galactic white dwarf binary foreground modeled after the Breivik et al. (2020) bulge+disk spatial model. The default values are those given in Breivik+20; the model itself traces back to electromagnetic studies of the Galaxy by McMillan et al. (2011), Bissantz and Gerhard (2002), and Binney et al. (1997). The distribution is azimuthally symmetric in the galactocentric frame.

Parameters:

(float) (zh)
(float)
(int) (nside)

Returns:

astro_map (float) – Healpy GW power skymap of the Galactic disk contribution to the foreground.
log_DWD_FG_map (float) – Healpy log GW power skymap. For plotting purposes.

blip.src.astro.generate_galactic_foreground(rh, zh, nside)[source]

Generate a galactic white dwarf binary foreground modeled after Breivik et al. (2020), consisting of a bulge + disk. The default values are those given in Breivik+20; the model itself traces back to electromagnetic studies of the Galaxy by McMillan et al. (2011), Bissantz and Gerhard (2002), and Binney et al. (1997). rh is the radial scale height in kpc, zh is the vertical scale height in kpc. Thin disk has rh=2.9kpc, zh=0.3kpc; Thick disk has rh=3.31kpc, zh=0.9kpc. Defaults to thin disk. The distribution is azimuthally symmetric in the galactocentric frame.

Parameters:

(float) (zh)
(float)
(int) (nside)

Returns:

astro_map (float) – Healpy GW power skymap of the DWD galactic foreground.
log_DWD_FG_map (float) – Healpy log GW power skymap. For plotting purposes.

blip.src.astro.generate_point_source(ang_coord1, ang_coord2, nside, convention='healpy', pad=True)[source]

Generates a point source skymap.

Parameters:

ang_coord1 (float) – angular coordinates of the point source in radians. Either theta, phi or ra, dec (see convention variable)
ang_coord2 (float) – angular coordinates of the point source in radians. Either theta, phi or ra, dec (see convention variable)
nside (int) – Healpy nside (skymap resolution)
convention (str) – Angle specification convention. Can be ‘healpy’ (Healpy polar theta, aziumuthal phi) or ‘radec’ (standard astronomical RA/DEC). Default is theta/phi.
pad (bool) – Whether to allow a small amount of power to artifically bleed into adjacent pixels to avoid numerical error issues later on. Only needed for single-pixel case.

Returns:

astro_map (array of floats)

Return type:

healpy skymap

blip.src.astro.generate_point_sources(coord_list, nside, convention='healpy')[source]

Generates a skymap with a flexible number of point sources.

Parameters:

coord_list (list of tuples) – List of (ang_coord1,ang_coord2) tuples, one tuple per source. Each tuple gives angular coordinates of their respective point sourc as either (theta, phi) or (ra, dec) (see convention variable).
nside (int) – Healpy nside (skymap resolution)
convention (str) – Angle specification convention. Can be ‘healpy’ (Healpy polar theta, aziumuthal phi) or ‘radec’ (standard astronomical RA/DEC). Default is theta/phi.

Returns:

astro_map (array of floats)

Return type:

healpy skymap

blip.src.astro.generate_sdg(nside, ra=80.21496, dec=-69.37772, D=50, r=2.1462, N=2169264)[source]

Generates the stochastic DWD signal from a a generic toy model spherical dwarf galaxy (SDG). Default values are for the LMC.

Parameters:

ra (float) – Right ascension and declination.
dec (float) – Right ascension and declination.
D (float) – Distance to SDG in kpc.
r (float) – radius of SDG in kpc
N (int) – Number of DWD systems in the SDG

Returns:

skymap – Healpy GW power skymap of the stochastic DWD signal.

Return type:

float

blip.src.astro.generate_two_point_source(theta_1, phi_1, theta_2, phi_2, nside)[source]

Generates a two-point-source skymap.

Depreciation note: Keeping until the angular resolution study is finished, then will depreciate in favor of generate_point_sources() (below)/

Parameters:

theta_1 (float) – angular coordinates of the 1st point source in radians
phi_1 (float) – angular coordinates of the 1st point source in radians
theta_2 (float) – angular coordinates of the 2nd point source in radians
phi_2 (float) – angular coordinates of the 2nd point source in radians

Returns:

astro_map (array of floats)

Return type:

healpy skymap

blip.src.astro.skymap_pix2sph(skymap, blmax)[source]

Transform a pixel-basis skymap into the b_lm spherical harmonic basis

Returns:: astro_blms – Spherical harmonic healpy expansion of the galactic foreground
Return type:: float

Sampler interfaces

Interface to `dynesty`

class blip.src.dynesty_engine.dynesty_engine[source]

Bases: object

Class for interfacing with dynesty sampler. This method also contains the priors definition for all models written to work with the dynesty sampler.

classmethod define_engine(lisaobj, params, nlive, nthread, randst, pool=None, resume=False)[source]

load_engine(randst, pool)[source]

static run_engine(engine)[source]

static run_engine_with_checkpointing(engine, parameters, interval, checkpoint_file, step=1000)[source]

Interface to `emcee`

class blip.src.emcee_engine.emcee_engine[source]

Bases: object

Class for interfacing with dynesty sampler. This method also contains the priors definition for all models written to work with the dynesty sampler.

classmethod define_engine(model, nlive, randst, pool=None)[source]

Defines the emcee engine.

Parameters:

model (Model object containing all relevant methods and parameters for evaluating the unified prior and likelihood)
(int) (nlive)
randst (initial random state)

classmethod logpost(theta, prior_transform, loglike)[source]

Wrapper to calculate log-posterior at theta.

Redefining this to place theta on the unit cube, consistent with the other implemented samplers.

static run_engine(engine, model, init_samples, Nburn, Nsamples)[source]

Runs the emcee sampler.

Parameters:

engine (the emcee EnsembleSampler engine)
array) (init_samples Ndim x Nwalkers)
(int) (Nsamples)
(int)

Returns:

post_samples (array)

Return type:

the final posterior samples

Interface to `numpyro`

class blip.src.numpyro_engine.numpyro_engine[source]

Bases: object

Class for interfacing with the Numpyro HMC (NUTS) sampler.

classmethod define_engine(lisaobj, Nburn, Nsamples, Nthreads, prog, seed, gpu=False)[source]

load_engine()[source]

static run_engine(engine, lisaobj, rng_key)[source]

static run_engine_with_checkpointing(engine, lisaModel, rng_key, chain, checkpoint_file, Ntotal, checkpoint_at, resume=False)[source]

Runs the numpyro sampler with sampler state checkpointing.

Parameters:

in] ([fill)
(bool) (resume) – ‘end’ (only saves sampler state at the very end of the run) ‘warmup’ (saves after warmup phase and at end) ‘interval’ (saves after warmup, at end, and after every checkpoint_interval number of samples) Note: Generally not worth checkpointing while sampling for large models/datasets, as the recompliation and GPU off/onloading time will exceed the sampling time.
(bool)

Returns:

post_samples (array)

Return type:

Posterior samples.

blip.src.numpyro_engine.numpyro_model(Model)[source]

Wrapper to translate our unified prior and likelihood to Numpyro-friendly input.

Parameters:

str) (parameters (list of)
(function) (log_likelihood)
(function)

blip.src.numpyro_engine.numpyro_model_sph(Model)[source]

Wrapper to translate our unified prior and likelihood to Numpyro-friendly input.

The _sph version here accounts for the fact that the phase parameter prior should be periodic.

Parameters:

str) (parameters (list of)
(function) (log_likelihood)
(function)

Miscellaneous

`utils` module

blip.src.utils.catch_color_duplicates(Object, color_pool=None, sacred_labels=[])[source]

Function to catch duplicate plotting colors and reassign from a default or user-specified pool of matplotlib colors.

Parameters:

Object (Model or Injection with attached submodels.)
color_pool (List of matplotlib color namestrings; see https://matplotlib.org/stable/gallery/color/named_colors.html)
sacred_labels (List of submodel names whose colors should be treated as inviolate.)

blip.src.utils.catch_duplicates(names)[source]

Function to catch duplicate names so we don’t overwrite keys while building a Model or Injection

Parameters:: str) (names (list of)
Returns:: names (list of str)
Return type:: model or injection submodel names, with duplicates numbered

blip.src.utils.ensure_color_matching(Model, Injection)[source]

Function to ensure linked Model and Injection models share a color in the final posterior fitmaker plot.

(i.e., pairwise matching between submodels and injection components that share a name.)

Parameters:

Model (Model object)
Injection (Injection object)

blip.src.utils.gen_suffixes(names)[source]

Function to generate appropriate parameter suffixes so repeated parameters are clearly linked to their respective submodel configurations.

Parameters:: str) (names (list of)
Returns:: suffixes (list of str)
Return type:: parameter suffixes for each respective model or injection submodel

blip.src.utils.get_robson19_shape_pars_from_tobs(T_obs, shape_pars=['$\\\\alpha_{\\\\rm shape}$', '$\\\\beta$', '$\\\\kappa$', '$\\\\gamma$', '$f_{\\\\rm knee}$'])[source]

Function to return the time-dependent foreground shape parameters of the Robson+19 analytic foreground spectrum, given an observing time in years. The function will round the given duration to the closest observing time in Table 1 of https://arxiv.org/abs/1803.01944.

Parameters:

(float) (T_obs)
str) (shape_pars (list of) – r”$alpha_{rm shape}$”, r”$beta$”, r”$kappa$”, r”$gamma$”, r”$f_{rm knee}$’

Returns:

shapevals (dict)

Return type:

Dictionary of ‘parameter_name’:val for the shape parameters of the Robson 19 analytic foreground model.

blip.src.utils.log_manager(level)[source]

Context manager to clean up bits of the code where we want e.g., healpy to be quieter. Adapted from code by Martin Heinz (https://martinheinz.dev/blog/34)

Parameters:: level (logging level (DEBUG, INFO, WARNING, ERROR))

`clebschGordan` module

class blip.src.clebschGordan.clebschGordan(blmax)[source]

Bases: object

Class with methods for manipulating clebsch-gordon coeffcients.

blm_2_alm(blms_in)[source]: Convert complex blm values to alm complex values. This will contain both -ve m values too in the standard order

blm_params_2_blms(blm_params)[source]: convert blm parameter values where amplitudes and phases are seperate to complex blm values.

blms_2_blm_params(blms)[source]: Convert complex-valued blms to our parameter notation where amplitudes and phases are separate.

calc_beta()[source]: Method to calculate beta array to convert from blm to alm

calc_blm_full(blms_in)[source]

Convert samples in blm space to blm complex values including negetive m vals

Input: blms ordered dictionary Ouput: blms_full, list including blms with negative m vals

idxtoalm(lmax, ii)[source]: index –> (l, m) function which works for negetive indices too

`hierarchical` module

class blip.src.hierarchical.postprocess(rundir, params, inj, parameters)[source]

Bases: LISAdata

blm_decompose(blm_samples)[source]

Function to decompose complex blms into real-valued components in the proper order.

Parameters:: blm_samples (array) – complex blms
Returns:: real-valued components
Return type:: blms (array)

breivik2020_log_likelihood(theta, post_dist)[source]

Likelihood of theta = {rh,zh} given set of BLIP alm posteriors. NOTE: Assumes blm posteriors (as transformed to alms) are normal and only uses their moments to determine the likelihood. Assumes blm covariance is diagonal.

Parameters:

rh (float) – Breivik+2020 radial scale height
zh (float) – Breivik+2020 vertical scale height
post_dist (scipy mv norm) – Scipy multivariate normal created from the blm posteriors

Returns:

log likelihood of theta = {rh,zh}

Return type:

loglike (float)

breivik2020_log_prior(theta, bounds=array([[2., 4.], [0.05, 2.]]))[source]

Prior for the breivik2020 model. Uniform on user-specified bounds in kpc. Default bounds are reasonable for the Milky Way.

Parameters:

theta (array) – [rh,zh], Breivik+2020 radial and vertical scale height
bounds (array) – [[rhmin,rhmax],[zhmin,zhmax]]

Returns:

log prior of theta = {rh,zh}

Return type:

logprior (float)

breivik2020_log_prob(theta, post_dist, bounds=array([[2, 4], [0, 2]]))[source]

Log probability for the Breivik+2020 model. Prior is uniform on user-specified bounds in kpc; default bounds are reasonable for the Milky Way. Likelihood of theta = {rh,zh} given set of BLIP alm posteriors. NOTE: Assumes blm posteriors (as transformed to alms) are normal and only uses their moments to determine the likelihood. Assumes blm covariance is diagonal.

Parameters:

theta (array) – [rh,zh], Breivik+2020 radial and vertical scale height
post_dist (scipy mv norm) – Scipy multivariate normal created from the blm posteriors
bounds (array) – [[rhmin,rhmax],[zhmin,zhmax]]

Returns:

log posterior probability of theta = {rh,zh}

Return type:

logp (float)

breivik2020_mapmaker(rh, zh)[source]

This is a streamlined version of makeLISAdata.generate_galactic_foreground(), sacrificing compactness for speed. Requires initialization via init_breivik2020_grid(), above.

Generate a galactic white dwarf binary foreground modeled after Breivik et al. (2020), consisting of a bulge + disk. rh is the radial scale height in kpc, zh is the vertical scale height in kpc. The distribution is azimuthally symmetric in the galactocentric frame.

Returns:: DWD_FG_map – Healpy GW power skymap of the DWD galactic foreground.
Return type:: float

hierarchical_sampler(model='breivik2020', Nwalkers=50, Nsamples=10000, Nburn=1000, rng=None, Nthread=1)[source]: Function to perform hierarchical sampling to constrain a spatial model of choice.

init_breivik2020_grid(grid_spec='interval', grid_res=0.33, gal_rad=16, gal_height=8)[source]

Function to initialize a grid on which to generate simple parameterized density models of the galactic DWD distribution.

Parameters:

grid_spec (str) – Determines the nature of grid_res (below). Can be ‘interval’ or ‘npoints’. If ‘interval’, grid_res is the dx=dy=dz grid interval in kpc. If ‘npoints’, grid_res is the number of number of points along x and y. (Note that the number of points along z will be scaled to keep dx=dy=dz if gal_rad and gal_height are different.)
grid_res (float) – Grid resolution as defined above. If grid_spec=’npoints’, type must be int.
gal_rad (float) – Max galactic radius of the grid in kpc. Grid will be definded on -gal_rad <= x,y <= +gal_rad.
gal_height (float) – Max galactic height of the grid in kpc. Grid will be definded on -gal_height <= z <= +gal_height.

post2dist(post)[source]

Function to generate a scipy multivariate normal distribution from BLIP blm posteriors. NOTE: Assumes blm posteriors (as transformed to alms) are normally distributed and only uses their moments to determine the distribution. Assumes alm covariance is diagonal.

Parameters:: post (array) – blm posterior samples
Returns:: Scipy multivariate normal created from the blm posteriors
Return type:: dist (scipy mv norm)

samples2alm(post)[source]

Function to convert BLIP blm amplitude+phase posterior samples to the alm basis

Parameters:: post (array) – blm posterior samples
Returns:: alm posterior samples
Return type:: alm_samples (array)

samples2blm(post)[source]

Function to convert BLIP blm amplitude+phase posterior samples to the alm basis

Parameters:: post (array) – blm posterior samples
Returns:: alm posterior samples
Return type:: alm_samples (array)

`config` module

class blip.config.Option(name: str, desc: str, default: str | None = None, required: bool = False)[source]

Bases: object

Configuration option.

default: str | None = None

desc: str

name: str

required: bool = False

blip.config.parse_config(paramsfile: str, resume: bool)[source]

Parse a configuration file paramsfile.

Parameters:

paramsfile – the config filename.
resume – a command-line flag to run_blip. It affects whether the configuration file is valid because it only makes sense to resume checkpointed runs if they have fixed seeds.

Returns:

a tuple of dictionaries (params, inj, misc) containing a valid configuration where all parameters have been converted from str to the appropriate types.

For invalid input configurations: raises TypeError, ValueError, KeyError, AssertionError, or configparse.NoOptionError, FileNotFoundError.

blip.config.parse_model_spec(specs: str, *, is_injection: bool, truevals_all: dict, fixedvals_all: dict | None = None, alias_all: dict | None = None, injection_specs: list[SubmodelSpec] | None = None) → list[SubmodelSpec][source]

Parse a model specifier string into its component submodels.

Parameters:

specs (str) – model specifier string, with submodels separated by “+” and no spaces.
is_injection (bool) – if True, parse this as an injection model. This will ignore the parameters alias_all and injection_specs, since the injections cannot be aliased to anything.
truevals_all (dict) – the truevals dictionary in the config params, specifying injected values. Not checked for correctness (this is for now the job of submodel.__init__).
fixedvals_all (dict | None, optional) – the fixedvals dictionary in the config params, specifying parameters to hold fixed in likelihood evaluation. Ignored if is_injection=True.
alias_all (dict | None, optional) – the alias dictionary in the config params. Ignored if is_injection=True.
injection_specs (list[SubmodelSpec] | None, optional) – the already-parsed specifiers for the injection submodels. Required to correctly parse the analysis model. Ignored if is_injection=True.

Returns:

The specifications to build the submodels.

Return type:

list[SubmodelSpec]

Raises:

ValueError – If the model specifier string contains invalid syntax.

blip.config.preprocess_truevals(tv_raw)[source]

API documentation

Core BLIP logic

run_blip module

submodel module

models module

LISA instrument and SGWB response

Legacy geometry module

Return s

Return s

fast_geometry module

The experimental faster_geometry module

sph_geometry module

makeLISAdata module

instrNoise module

orbitinglisa module

Astrophysics

astro module

Sampler interfaces

Interface to dynesty

Interface to emcee

Interface to numpyro

Miscellaneous

utils module

clebschGordan module

hierarchical module

config module

`run_blip` module

`submodel` module

`models` module

Legacy `geometry` module

`fast_geometry` module

The experimental `faster_geometry` module

`sph_geometry` module

`makeLISAdata` module

`instrNoise` module

`orbitinglisa` module

`astro` module

Interface to `dynesty`

Interface to `emcee`

Interface to `numpyro`

`utils` module

`clebschGordan` module

`hierarchical` module

`config` module