Skip to content
Snippets Groups Projects
Commit ce56c551 authored by Colm Talbot's avatar Colm Talbot
Browse files

set up a base lal FD waveform source function

parent 59869044
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 138 additions and 158 deletions
bilby/gw/source.py
+ 138
− 158
View file @ ce56c551
......@@ -19,8 +19,8 @@ except ImportError:
def lal_binary_black_hole(
frequency_array, mass_1, mass_2, luminosity_distance, a_1, tilt_1, phi_12, a_2, tilt_2, phi_jl,
iota, phase, **kwargs):
frequency_array, mass_1, mass_2, luminosity_distance, a_1, tilt_1,
phi_12, a_2, tilt_2, phi_jl, iota, phase, **kwargs):
""" A Binary Black Hole waveform model using lalsimulation
Parameters
......@@ -56,67 +56,82 @@ def lal_binary_black_hole(
-------
dict: A dictionary with the plus and cross polarisation strain modes
"""
waveform_kwargs = dict(
waveform_approximant='IMRPhenomPv2', reference_frequency=50.0,
minimum_frequency=20.0, maximum_frequency=frequency_array[-1])
waveform_kwargs.update(kwargs)
waveform_approximant = waveform_kwargs['waveform_approximant']
reference_frequency = waveform_kwargs['reference_frequency']
minimum_frequency = waveform_kwargs['minimum_frequency']
maximum_frequency = waveform_kwargs['maximum_frequency']
delta_frequency = frequency_array[1] - frequency_array[0]
frequency_bounds = ((frequency_array >= minimum_frequency) *
(frequency_array <= maximum_frequency))
if mass_2 > mass_1:
return None
luminosity_distance = luminosity_distance * 1e6 * utils.parsec
mass_1 = mass_1 * utils.solar_mass
mass_2 = mass_2 * utils.solar_mass
return _base_lal_cbc_waveform(
frequency_array=frequency_array, mass_1=mass_1, mass_2=mass_2,
luminosity_distance=luminosity_distance, iota=iota, phase=phase,
a_1=a_1, a_2=a_2, tilt_1=tilt_1, tilt_2=tilt_2, phi_12=phi_12,
phi_jl=phi_jl, **waveform_kwargs)
if tilt_1 == 0 and tilt_2 == 0:
spin_1x = 0
spin_1y = 0
spin_1z = a_1
spin_2x = 0
spin_2y = 0
spin_2z = a_2
else:
iota, spin_1x, spin_1y, spin_1z, spin_2x, spin_2y, spin_2z = (
lalsim_SimInspiralTransformPrecessingNewInitialConditions(
iota, phi_jl, tilt_1, tilt_2, phi_12, a_1, a_2, mass_1,
mass_2, reference_frequency, phase))
longitude_ascending_nodes = 0.0
eccentricity = 0.0
mean_per_ano = 0.0
waveform_dictionary = None
approximant = lalsim_GetApproximantFromString(waveform_approximant)
hplus, hcross = lalsim_SimInspiralFD(
mass_1, mass_2, spin_1x, spin_1y, spin_1z, spin_2x, spin_2y,
spin_2z, luminosity_distance, iota, phase,
longitude_ascending_nodes, eccentricity, mean_per_ano, delta_frequency,
minimum_frequency, maximum_frequency, reference_frequency,
waveform_dictionary, approximant)
def lal_binary_neutron_star(
frequency_array, mass_1, mass_2, luminosity_distance, chi_1, chi_2,
iota, phase, lambda_1, lambda_2, **kwargs):
""" A Binary Neutron Star waveform model using lalsimulation
h_plus = hplus.data.data
h_cross = hcross.data.data
Parameters
----------
frequency_array: array_like
The frequencies at which we want to calculate the strain
mass_1: float
The mass of the heavier object in solar masses
mass_2: float
The mass of the lighter object in solar masses
luminosity_distance: float
The luminosity distance in megaparsec
chi_1: float
Dimensionless aligned spin
chi_2: float
Dimensionless aligned spin
iota: float
Orbital inclination
phase: float
The phase at coalescence
ra: float
The right ascension of the binary
dec: float
The declination of the object
geocent_time: float
The time at coalescence
psi: float
Orbital polarisation
lambda_1: float
Dimensionless tidal deformability of mass_1
lambda_2: float
Dimensionless tidal deformability of mass_2
h_plus = h_plus[:len(frequency_array)] * frequency_bounds
h_cross = h_cross[:len(frequency_array)] * frequency_bounds
kwargs: dict
Optional keyword arguments
return {'plus': h_plus, 'cross': h_cross}
Returns
-------
dict: A dictionary with the plus and cross polarisation strain modes
"""
waveform_kwargs = dict(
waveform_approximant='TaylorF2', reference_frequency=50.0,
minimum_frequency=20.0, maximum_frequency=frequency_array[-1])
a_1 = abs(chi_1)
a_2 = abs(chi_2)
tilt_1 = np.arccos(np.sign(chi_1))
tilt_2 = np.arccos(np.sign(chi_2))
phi_12 = 0.0
phi_jl = 0.0
waveform_kwargs.update(kwargs)
return _base_lal_cbc_waveform(
frequency_array=frequency_array, mass_1=mass_1, mass_2=mass_2,
luminosity_distance=luminosity_distance, iota=iota, phase=phase,
a_1=a_1, a_2=a_2, tilt_1=tilt_1, tilt_2=tilt_2, phi_12=phi_12,
phi_jl=phi_jl, lambda_1=lambda_1, lambda_2=lambda_2, **waveform_kwargs)
def lal_eccentric_binary_black_hole_no_spins(
frequency_array, mass_1, mass_2, eccentricity, luminosity_distance, iota, phase, **kwargs):
""" Eccentric binary black hole waveform model using lalsimulation (EccentricFD)
frequency_array, mass_1, mass_2, eccentricity, luminosity_distance,
iota, phase, **kwargs):
"""
Eccentric binary black hole waveform model using lalsimulation (EccentricFD)
Parameters
----------
......@@ -141,11 +156,62 @@ def lal_eccentric_binary_black_hole_no_spins(
-------
dict: A dictionary with the plus and cross polarisation strain modes
"""
waveform_kwargs = dict(
waveform_approximant='EccentricFD', reference_frequency=10.0,
minimum_frequency=10.0, maximum_frequency=frequency_array[-1])
waveform_kwargs.update(kwargs)
return _base_lal_cbc_waveform(
frequency_array=frequency_array, mass_1=mass_1, mass_2=mass_2,
luminosity_distance=luminosity_distance, iota=iota, phase=phase,
eccentricity=eccentricity, **waveform_kwargs)
def _base_lal_cbc_waveform(
frequency_array, mass_1, mass_2, luminosity_distance, iota, phase,
a_1=0.0, a_2=0.0, tilt_1=0.0, tilt_2=0.0, phi_12=0.0, phi_jl=0.0,
lambda_1=0.0, lambda_2=0.0, eccentricity=0.0, **waveform_kwargs):
""" Generate a cbc waveform model using lalsimulation
Parameters
----------
frequency_array: array_like
The frequencies at which we want to calculate the strain
mass_1: float
The mass of the heavier object in solar masses
mass_2: float
The mass of the lighter object in solar masses
luminosity_distance: float
The luminosity distance in megaparsec
a_1: float
Dimensionless primary spin magnitude
tilt_1: float
Primary tilt angle
phi_12: float
Azimuthal angle between the component spins
a_2: float
Dimensionless secondary spin magnitude
tilt_2: float
Secondary tilt angle
phi_jl: float
Azimuthal angle between the total and orbital angular momenta
iota: float
Orbital inclination
FIXME: this should be theta_jn
phase: float
The phase at coalescence
eccentricity: float
Binary eccentricity
lambda_1: float
Tidal deformability of the more massive object
lambda_2: float
Tidal deformability of the less massive object
kwargs: dict
Optional keyword arguments
Returns
-------
dict: A dictionary with the plus and cross polarisation strain modes
"""
waveform_approximant = waveform_kwargs['waveform_approximant']
reference_frequency = waveform_kwargs['reference_frequency']
minimum_frequency = waveform_kwargs['minimum_frequency']
......@@ -162,17 +228,25 @@ def lal_eccentric_binary_black_hole_no_spins(
mass_1 = mass_1 * utils.solar_mass
mass_2 = mass_2 * utils.solar_mass
spin_1x = 0.0
spin_1y = 0.0
spin_1z = 0.0
spin_2x = 0.0
spin_2y = 0.0
spin_2z = 0.0
if tilt_1 == 0 and tilt_2 == 0:
spin_1x = 0
spin_1y = 0
spin_1z = a_1
spin_2x = 0
spin_2y = 0
spin_2z = a_2
else:
iota, spin_1x, spin_1y, spin_1z, spin_2x, spin_2y, spin_2z = (
lalsim_SimInspiralTransformPrecessingNewInitialConditions(
iota, phi_jl, tilt_1, tilt_2, phi_12, a_1, a_2, mass_1,
mass_2, reference_frequency, phase))
longitude_ascending_nodes = 0.0
mean_per_ano = 0.0
waveform_dictionary = None
waveform_dictionary = lal.CreateDict()
lalsim_SimInspiralWaveformParamsInsertTidalLambda1(waveform_dictionary, lambda_1)
lalsim_SimInspiralWaveformParamsInsertTidalLambda2(waveform_dictionary, lambda_2)
approximant = lalsim_GetApproximantFromString(waveform_approximant)
......@@ -183,8 +257,11 @@ def lal_eccentric_binary_black_hole_no_spins(
minimum_frequency, maximum_frequency, reference_frequency,
waveform_dictionary, approximant)
h_plus = hplus.data.data * frequency_bounds
h_cross = hcross.data.data * frequency_bounds
h_plus = hplus.data.data
h_cross = hcross.data.data
h_plus = h_plus[:len(frequency_array)] * frequency_bounds
h_cross = h_cross[:len(frequency_array)] * frequency_bounds
return {'plus': h_plus, 'cross': h_cross}
......@@ -248,103 +325,6 @@ def supernova_pca_model(
return {'plus': h_plus, 'cross': h_cross}
def lal_binary_neutron_star(
frequency_array, mass_1, mass_2, luminosity_distance, chi_1, chi_2,
iota, phase, lambda_1, lambda_2, **kwargs):
""" A Binary Neutron Star waveform model using lalsimulation
Parameters
----------
frequency_array: array_like
The frequencies at which we want to calculate the strain
mass_1: float
The mass of the heavier object in solar masses
mass_2: float
The mass of the lighter object in solar masses
luminosity_distance: float
The luminosity distance in megaparsec
chi_1: float
Dimensionless aligned spin
chi_2: float
Dimensionless aligned spin
iota: float
Orbital inclination
phase: float
The phase at coalescence
ra: float
The right ascension of the binary
dec: float
The declination of the object
geocent_time: float
The time at coalescence
psi: float
Orbital polarisation
lambda_1: float
Dimensionless tidal deformability of mass_1
lambda_2: float
Dimensionless tidal deformability of mass_2
kwargs: dict
Optional keyword arguments
Returns
-------
dict: A dictionary with the plus and cross polarisation strain modes
"""
waveform_kwargs = dict(
waveform_approximant='TaylorF2', reference_frequency=50.0,
minimum_frequency=20.0, maximum_frequency=frequency_array[-1])
waveform_kwargs.update(kwargs)
waveform_approximant = waveform_kwargs['waveform_approximant']
reference_frequency = waveform_kwargs['reference_frequency']
minimum_frequency = waveform_kwargs['minimum_frequency']
maximum_frequency = waveform_kwargs['maximum_frequency']
delta_frequency = frequency_array[1] - frequency_array[0]
frequency_bounds = ((frequency_array >= minimum_frequency) *
(frequency_array <= maximum_frequency))
if mass_2 > mass_1:
return None
luminosity_distance = luminosity_distance * 1e6 * utils.parsec
mass_1 = mass_1 * utils.solar_mass
mass_2 = mass_2 * utils.solar_mass
spin_1x = 0
spin_1y = 0
spin_1z = chi_1
spin_2x = 0
spin_2y = 0
spin_2z = chi_2
longitude_ascending_nodes = 0.0
eccentricity = 0.0
mean_per_ano = 0.0
waveform_dictionary = lal.CreateDict()
lalsim_SimInspiralWaveformParamsInsertTidalLambda1(waveform_dictionary, lambda_1)
lalsim_SimInspiralWaveformParamsInsertTidalLambda2(waveform_dictionary, lambda_2)
approximant = lalsim_GetApproximantFromString(waveform_approximant)
hplus, hcross = lalsim_SimInspiralFD(
mass_1, mass_2, spin_1x, spin_1y, spin_1z, spin_2x, spin_2y,
spin_2z, luminosity_distance, iota, phase,
longitude_ascending_nodes, eccentricity, mean_per_ano, delta_frequency,
minimum_frequency, maximum_frequency, reference_frequency,
waveform_dictionary, approximant)
h_plus = hplus.data.data * frequency_bounds
h_cross = hcross.data.data * frequency_bounds
h_plus = h_plus[:len(frequency_array)]
h_cross = h_cross[:len(frequency_array)]
return {'plus': h_plus, 'cross': h_cross}
def roq(frequency_array, mass_1, mass_2, luminosity_distance, a_1, tilt_1,
phi_12, a_2, tilt_2, phi_jl, iota, phase, **waveform_arguments):
"""
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment