Unable to inject aligned spin waveforms for BBH parameters
Hello @gregory.ashton, @moritz.huebner, @colm.talbot,
It seems that if I try to use IMRPhenomD waveform instead bilby.gw.conversion.generate_all_bbh_parameters function fails to convert the waveform parameters properly, although it works for IMRPhenomPv2 waveform. @sebastian-khan tracked it down to the fact that the conversion function in Bilby has some numerical round off errors in the spin_1x, spin_1y, spin_2x, spin_2y (which should be zero in all the cases). This fails to parse the lalsim_SimInspiralChooseFDWaveform function in LalSimulation.
- The git-hash of the Bilby version:
16:19 bilby INFO : Running bilby version: 0.3.5: (CLEAN) fa5630d 2019-02-04 16:05:57 +1100- I am copying the error message here for your reference:
---------------------------------------------------------------------------
RuntimeError Traceback (most recent call last)
~/from_holodeck/playground_bilby/issue_wfgen/standard_15d_cbc_tutorial_IMRpD.py in <module>
95 injection_parameters=injection_parameters, outdir=outdir,
96 label=label, maxmcmc=2000,
---> 97 conversion_function=bilby.gw.conversion.generate_all_bbh_parameters)
98
99 # Make a corner plot.
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/core/sampler/__init__.py in run_sampler(likelihood, priors, label, outdir, sampler, use_ratio, injection_parameters, conversion_function, plot, default_priors_file, clean, meta_data, save, result_class, **kwargs)
135 injection_parameters=injection_parameters, meta_data=meta_data,
136 use_ratio=use_ratio, plot=plot, result_class=result_class,
--> 137 **kwargs)
138 else:
139 print(implemented_samplers)
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/core/sampler/base_sampler.py in __init__(self, likelihood, priors, outdir, label, use_ratio, plot, skip_import_verification, injection_parameters, meta_data, result_class, **kwargs)
106 self.__fixed_parameter_keys = []
107 self._initialise_parameters()
--> 108 self._verify_parameters()
109 self._verify_use_ratio()
110 self.kwargs = kwargs
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/core/sampler/base_sampler.py in _verify_parameters(self)
245 try:
246 t1 = datetime.datetime.now()
--> 247 self.likelihood.log_likelihood()
248 self._log_likelihood_eval_time = (
249 datetime.datetime.now() - t1).total_seconds()
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/likelihood.py in log_likelihood(self)
228
229 def log_likelihood(self):
--> 230 return self.log_likelihood_ratio() + self.noise_log_likelihood()
231
232 @property
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/likelihood.py in log_likelihood_ratio(self)
159 def log_likelihood_ratio(self):
160 waveform_polarizations =\
--> 161 self.waveform_generator.frequency_domain_strain(self.parameters)
162
163 if waveform_polarizations is None:
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/waveform_generator.py in frequency_domain_strain(self, parameters)
112 transformation_function=utils.nfft,
113 transformed_model=self.time_domain_source_model,
--> 114 transformed_model_data_points=self.time_array)
115
116 def time_domain_strain(self, parameters=None):
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/waveform_generator.py in _calculate_strain(self, model, model_data_points, transformation_function, transformed_model, transformed_model_data_points, parameters)
149 self.parameters = parameters
150 if model is not None:
--> 151 model_strain = self._strain_from_model(model_data_points, model)
152 elif transformed_model is not None:
153 model_strain = self._strain_from_transformed_model(transformed_model_data_points, transformed_model,
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/waveform_generator.py in _strain_from_model(self, model_data_points, model)
158
159 def _strain_from_model(self, model_data_points, model):
--> 160 return model(model_data_points, **self.parameters)
161
162 def _strain_from_transformed_model(self, transformed_model_data_points, transformed_model, transformation_function):
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/source.py in 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)
101 longitude_ascending_nodes, eccentricity, mean_per_ano, delta_frequency,
102 minimum_frequency, maximum_frequency, reference_frequency,
--> 103 waveform_dictionary, approximant)
104
105 h_plus = hplus.data.data
~/anaconda3/lib/python3.7/site-packages/bilby-0.3.5-py3.7.egg/bilby/gw/utils.py in lalsim_SimInspiralChooseFDWaveform(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)
768 longitude_ascending_nodes, eccentricity, mean_per_ano, delta_frequency,
769 minimum_frequency, maximum_frequency, reference_frequency,
--> 770 waveform_dictionary, approximant)
771
772
RuntimeError: Invalid argument- I am also attaching the code (which is only with minimal changes from the default 15D BBH example) that was used to run Bilby:
#!/usr/bin/env python
"""
Tutorial to demonstrate running parameter estimation on a full 15 parameter
space for an injected cbc signal. This is the standard injection analysis script
one can modify for the study of injected CBC events.
"""
from __future__ import division, print_function
import numpy as np
import bilby
# Set the duration and sampling frequency of the data segment that we're
# going to inject the signal into
duration = 4.
sampling_frequency = 2048.
# Specify the output directory and the name of the simulation.
outdir = 'outdir_IMRPhenomD'
label = 'aligned-spin_test-example'
bilby.core.utils.setup_logger(outdir=outdir, label=label)
# Set up a random seed for result reproducibility. This is optional!
np.random.seed(88170235)
# We are going to inject a binary black hole waveform. We first establish a
# dictionary of parameters that includes all of the different waveform
# parameters, including masses of the two black holes (mass_1, mass_2),
# spins of both black holes (a, tilt, phi), etc.
injection_parameters = dict(
mass_1=36., mass_2=29., a_1=0.4, a_2=0.3, tilt_1=0.0, tilt_2=0.0,
phi_12=0.0, phi_jl=0.0, luminosity_distance=2000., iota=0.4, psi=2.659,
phase=1.3, geocent_time=1126259642.413, ra=1.375, dec=-1.2108)
# Fixed arguments passed into the source model
waveform_arguments = dict(waveform_approximant='IMRPhenomD',
reference_frequency=50., minimum_frequency=20.)
# Create the waveform_generator using a LAL BinaryBlackHole source function
# the generator will convert all the parameters
waveform_generator = bilby.gw.WaveformGenerator(
duration=duration, sampling_frequency=sampling_frequency,
frequency_domain_source_model=bilby.gw.source.lal_binary_black_hole,
parameter_conversion=bilby.gw.conversion.convert_to_lal_binary_black_hole_parameters,
waveform_arguments=waveform_arguments)
# Set up interferometers. In this case we'll use two interferometers
# (LIGO-Hanford (H1), LIGO-Livingston (L1). These default to their design
# sensitivity
ifos = bilby.gw.detector.InterferometerList(['H1', 'L1'])
ifos.set_strain_data_from_power_spectral_densities(
sampling_frequency=sampling_frequency, duration=duration,
start_time=injection_parameters['geocent_time'] - 3)
ifos.inject_signal(waveform_generator=waveform_generator,
parameters=injection_parameters)
# For this analysis, we implemenet the standard BBH priors defined, except for
# the definition of the time prior, which is defined as uniform about the
# injected value.
# Furthermore, we decide to sample in chirp mass and mass ratio, due to the
# preferred shape for the associated posterior distributions.
priors = bilby.gw.prior.BBHPriorDict()
priors.pop('mass_1')
priors.pop('mass_2')
priors['chirp_mass'] = bilby.prior.Uniform(
name='chirp_mass', latex_label='$M$', minimum=10.0, maximum=100.0,
unit='$M_{\\odot}$')
priors['mass_ratio'] = bilby.prior.Uniform(
name='mass_ratio', latex_label='$q$', minimum=0.5, maximum=1.0)
priors['geocent_time'] = bilby.core.prior.Uniform(
minimum=injection_parameters['geocent_time'] - 0.1,
maximum=injection_parameters['geocent_time'] + 0.1,
name='geocent_time', latex_label='$t_c$', unit='$s$')
# Initialise the likelihood by passing in the interferometer data (ifos) and
# the waveoform generator, as well the priors.
# The explicit time, distance, and phase marginalizations are turned on to
# improve convergence, and the parameters are recovered by the conversion
# function.
likelihood = bilby.gw.GravitationalWaveTransient(
interferometers=ifos, waveform_generator=waveform_generator, priors=priors,
distance_marginalization=False, phase_marginalization=False, time_marginalization=False)
# Run sampler. In this case we're going to use the `cpnest` sampler
# Note that the maxmcmc parameter is increased so that between each iteration of
# the nested sampler approach, the walkers will move further using an mcmc
# approach, searching the full parameter space.
# The conversion function will determine the distance, phase and coalescence
# time posteriors in post processing.
result = bilby.run_sampler(
likelihood=likelihood, priors=priors, sampler='cpnest', npoints=2000,
injection_parameters=injection_parameters, outdir=outdir,
label=label, maxmcmc=2000,
conversion_function=bilby.gw.conversion.generate_all_bbh_parameters)
# Make a corner plot.
result.plot_corner()