Generating waveforms

Through pesummary’s SamplesDict class, we can easily generate a waveform based on a specific set of posterior samples either in the time domain or the frequency domain. This is through the td_waveform and fd_waveform methods. Below we show an example using the publically available GW190814 posterior samples.

First let us plot a waveform in the time domain,

from pesummary.gw.fetch import fetch_open_samples
import matplotlib.pyplot as plt

# First we download and read the publically available posterior samples
f = fetch_open_samples("GW190814", unpack=True, path="GW190814.h5")
samples = f.samples_dict
EOB = samples["C01:SEOBNRv4PHM"]

# Next we generate the plus and cross polarizations in the time domain for the
# maximum likelihood sample
approximant = "SEOBNRv4PHM"
delta_t = 1. / 4096
f_low = 20.
wvfs = EOB.maxL_td_waveform(approximant, delta_t, f_low, f_ref=f_low)

# Alternatively, we may generate the plus and cross polarizations in the time
# domain for a specific sample with the ind kwarg
ind = 100
wvfs = EOB.td_waveform(approximant, delta_t, f_low, f_ref=f_low, ind=ind)

# It is often useful to know what the strain is at a given gravitational wave
# detector. This involves calculating the antenna response function. If we
# wanted to generate the maximum likelihood strain projected onto the LIGO
# Livingtson detector, we may simply run
ht = EOB.maxL_td_waveform(approximant, delta_t, f_low, f_ref=f_low, project="L1")

# In all cases, a `gwpy.timeseries.TimeSeries` object is returned. This can
# therefore easily be plotted with
fig = plt.figure()
plt.plot(ht.times, ht)
plt.savefig("GW190814_SEOBNRv4PHM_td.png")

In the above example only the maximum likelihood waveform is plotted. Sometimes it is useful to know the uncertainty on this waveform. We can calculate and plot the 1 sigma and 2 sigma symmetric confidence intervals of this waveform by taking advantage of the level kwarg,

from pesummary.gw.fetch import fetch_open_samples
import matplotlib.pyplot as plt

f = fetch_open_samples("GW190814", unpack=True, path="GW190814.h5")
samples = f.samples_dict
EOB = samples["C01:SEOBNRv4PHM"]
_ = EOB.downsample(1000)
approximant = "SEOBNRv4PHM"
delta_t = 1. / 4096
f_low = 20.
ht, upper, lower, bound_times = EOB.maxL_td_waveform(
    approximant, delta_t, f_low, f_ref=f_low, project="L1", level=[0.68, 0.95],
    multi_process=4
)
fig = plt.figure()
plt.plot(ht.times, ht)
plt.fill_between(bound_times, upper[0], lower[0], color='r', alpha=0.4)
plt.fill_between(bound_times, upper[1], lower[1], color='r', alpha=0.2)
plt.xlim(1249852257.01 - 0.1, 1249852257.01 + 0.04)
plt.savefig("GW190814_SEOBNRv4PHM_td_with_uncertainty.png")

Here we have chosen to downsample the posterior samples to 1000 samples (_ = EOB.downsample(1000)) and used 4 CPUs (multi_process=4) to speed up waveform generation.

We can also generate a waveform in the frequency domain with,

from pesummary.gw.fetch import fetch_open_samples
import matplotlib.pyplot as plt

# First we download and read the publically available posterior samples
f = fetch_open_samples(
    "GW190814", catalog="GWTC-2", unpack=True, path="GW190814.h5"
)
samples = f.samples_dict
Phenom = samples["C01:IMRPhenomPv3HM"]

# Next we generate the plus and cross polarizations in the frequency domain for
# the maximum likelihood sample
approximant = "IMRPhenomPv3HM"
delta_f = 1. / 256
f_low = 20.
f_high = 1024.
wvfs = Phenom.maxL_fd_waveform(approximant, delta_f, f_low, f_high, f_ref=f_low)

# Alternatively, we may generate the plus and cross polarizations in the
# frequency domain for a specific sample with the ind kwarg
ind = 100
wvfs = Phenom.fd_waveform(approximant, delta_f, f_low, f_high, f_ref=f_low, ind=ind)

# It is often useful to know what the strain is at a given gravitational wave
# detector. This involves calculating the antenna response function. If we
# wanted to generate the maximum likelihood strain projected onto the LIGO
# Livingtson detector, we may simply run
ht = Phenom.maxL_fd_waveform(approximant, delta_f, f_low, f_high, f_ref=f_low, project="L1")

# In all cases, a `gwpy.frequencyseries.FrequencySeries` object is returned.
# This can therefore easily be plotted with
fig = plt.figure()
plt.plot(ht.frequencies, ht)
plt.savefig("GW190814_IMRPhenomPv3HM_fd.png")

For more details about the waveform generator in pesummary see,

pesummary.gw.waveform.antenna_response(samples, ifo)

pesummary.gw.waveform.fd_waveform(samples, approximant, delta_f, f_low, f_high, f_ref=20.0, project=None, ind=0, longAscNodes=0.0, eccentricity=0.0, LAL_parameters=None, mode_array=None, pycbc=False, flen=None, flags={}, debug=False, **kwargs)

Generate a gravitational wave in the frequency domain

Parameters:

• approximant (str) – name of the approximant to use when generating the waveform

• delta_f (float) – spacing between frequency samples

• f_low (float) – frequency to start evaluating the waveform

• f_high (float) – frequency to stop evaluating the waveform

• f_ref (float, optional) – reference frequency

• project (str, optional) – name of the detector to project the waveform onto. If None, the plus and cross polarizations are returned. Default None

• ind (int, optional) – index for the sample you wish to plot

• longAscNodes (float, optional) – longitude of ascending nodes, degenerate with the polarization angle. Default 0.

• eccentricity (float, optional) – eccentricity at reference frequency. Default 0.

• LAL_parameters (LALDict, optional) – LAL dictionary containing accessory parameters. Default None

• mode_array (2d list) – 2d list of modes you wish to include in waveform. Must be of the form [[l1, m1], [l2, m2]]

• pycbc (Bool, optional) – return a the waveform as a pycbc.frequencyseries.FrequencySeries object

• flen (int) – Length of the frequency series in samples. Default is None. Only used when pycbc=True

• debug (bool, optional) – if True, return the model object used to generate the waveform. Only used when the waveform approximant is not in LALSimulation

• flags (dict, optional) – waveform specific flags to add to LAL dictionary

• **kwargs (dict, optional) – all kwargs passed to _calculate_hp_hc_fd

pesummary.gw.waveform.td_waveform(samples, approximant, delta_t, f_low, f_ref=20.0, project=None, ind=0, longAscNodes=0.0, eccentricity=0.0, LAL_parameters=None, mode_array=None, pycbc=False, level=None, multi_process=1, flags={}, debug=False, **kwargs)

Generate a gravitational wave in the time domain

Parameters:

• approximant (str) – name of the approximant to use when generating the waveform

• delta_t (float) – spacing between frequency samples

• f_low (float) – frequency to start evaluating the waveform

• f_ref (float, optional) – reference frequency

• project (str, optional) – name of the detector to project the waveform onto. If None, the plus and cross polarizations are returned. Default None

• ind (int, optional) – index for the sample you wish to plot

• longAscNodes (float, optional) – longitude of ascending nodes, degenerate with the polarization angle. Default 0.

• eccentricity (float, optional) – eccentricity at reference frequency. Default 0.

• LAL_parameters (LALDict, optional) – LAL dictionary containing accessory parameters. Default None

• mode_array (2d list) – 2d list of modes you wish to include in waveform. Must be of the form [[l1, m1], [l2, m2]]

• pycbc (Bool, optional) – return a the waveform as a pycbc.timeseries.TimeSeries object

• level (list, optional) – the symmetric confidence interval of the time domain waveform. Level must be greater than 0 and less than 1

• multi_process (int, optional) – number of cores to run on when generating waveforms. Only used when level is not None,

• flags (dict, optional) – waveform specific flags to add to LAL dictionary. Default {}

• kwargs (dict, optional) – all kwargs passed to _td_waveform