inspiral_range¶
-
gwpy.astro.inspiral_range(psd: FrequencySeries, snr: float =
8
, mass1: float =1.4
, mass2: float =1.4
, fmin: float | None =None
, fmax: float | None =None
, horizon: bool =False
, **kwargs) Quantity [source]¶ Calculate the cosmology-corrected inspiral sensitive distance.
This method returns the distance (in megaparsecs) to which a compact binary inspiral with the given component masses would be detectable given the instrumental PSD. The calculation is defined in Belczynski et. al (2014):
https://dx.doi.org/10.1088/0004-637x/789/2/120
- Parameters:¶
- psd
FrequencySeries
The instrumental power-spectral-density data.
- snr
float
, optional The signal-to-noise ratio for which to calculate range.
- mass1
float
,Quantity
, optional The mass of the first binary component (
float
assumed in solar masses).- mass2
float
,Quantity
, optional The mass of the second binary component (
float
assumed in solar masses).- fmin
float
, optional The lower frequency cut-off of the integral, default:
psd.df
.- fmax
float
, optional The maximum frequency limit of the integral, defaults to the rest-frame innermost stable circular orbit (ISCO) frequency.
- horizon
bool
, optional If
True
, return the maximal ‘horizon’ luminosity distance, otherwise return the angle-averaged comoving distance.- **kwargs
Additional keyword arguments to
CBCWaveform
.
- psd
- Returns:¶
- range
Quantity
The calculated inspiral range [Mpc].
- range
See also
sensemon_range
For the method based on LIGO-T030276, also known as LIGO SenseMonitor.
inspiral-range
The package which does heavy lifting for waveform simulation and cosmology calculations.
Examples
Grab some data for LIGO-Livingston around GW150914 and generate a PSD:
>>> from gwpy.timeseries import TimeSeries >>> hoft = TimeSeries.fetch_open_data("H1", 1126259446, 1126259478) >>> hoff = hoft.psd(fftlength=4)
Now, we can calculate the
inspiral_range()
:>>> from gwpy.astro import inspiral_range >>> r = inspiral_range(hoff, fmin: float = 30) >>> print(r) 70.4612102889 Mpc