Fix bug

?

examples/gw_examples/data_examples/GW150914.prior

+mass_ratio = Uniform(name='mass_ratio', minimum=0.125, maximum=1, boundary='reflective')
+chirp_mass = Uniform(name='chirp_mass', minimum=25, maximum=31, unit='$M_{\odot}$', boundary='reflective')
+mass_1 = Constraint(name='mass_1', minimum=10, maximum=80)
+mass_2 = Constraint(name='mass_2', minimum=10, maximum=80)
+a_1 = Uniform(name='a_1', minimum=0, maximum=0.99, boundary='reflective')
+a_2 = Uniform(name='a_2', minimum=0, maximum=0.99, boundary='reflective')
+tilt_1 = Sine(name='tilt_1', boundary='reflective')
+tilt_2 = Sine(name='tilt_2', boundary='reflective')
+phi_12 = Uniform(name='phi_12', minimum=0, maximum=2 * np.pi, boundary='periodic')
+phi_jl = Uniform(name='phi_jl', minimum=0, maximum=2 * np.pi, boundary='periodic')
+luminosity_distance = PowerLaw(alpha=2, name='luminosity_distance', minimum=50, maximum=2000, unit='Mpc', latex_label='$d_L$')
+dec =  Cosine(name='dec', boundary='reflective')
+ra =  Uniform(name='ra', minimum=0, maximum=2 * np.pi, boundary='periodic')
+theta_jn =  Sine(name='theta_jn', boundary='reflective')
+psi =  Uniform(name='psi', minimum=0, maximum=np.pi, boundary='periodic')
+phase =  Uniform(name='phase', minimum=0, maximum=2 * np.pi, boundary='periodic')

examples/open_data_examples/GW150914.py → examples/gw_examples/data_examples/GW150914.py

 #!/usr/bin/env python
 """
-Tutorial to demonstrate running parameter estimation on GW150914 using open
-data.
+Tutorial to demonstrate running parameter estimation on GW150914
 
 This example estimates all 15 parameters of the binary black hole system using
-commonly used prior distributions.  This will take a few hours to run.
+commonly used prior distributions. This will take several hours to run. The
+data is obtained using gwpy, see [1] for information on how to access data on
+the LIGO Data Grid instead.
+
+[1] https://gwpy.github.io/docs/stable/timeseries/remote-access.html
 """
 from __future__ import division, print_function
 import bilby
+from gwpy.timeseries import TimeSeries
 
+logger = bilby.core.utils.logger
 outdir = 'outdir'
 label = 'GW150914'
-time_of_event = bilby.gw.utils.get_event_time(label)
-
-# This sets up logging output to understand what bilby is doing
-bilby.core.utils.setup_logger(outdir=outdir, label=label)
-
-# Here we import the detector data. This step downloads data from the
-# LIGO/Virgo open data archives. The data is saved to an `Interferometer`
-# object (here called `H1` and `L1`). A Power Spectral Density (PSD) estimate
-# is also generated and saved to the same object. To check the data and PSD
-# makes sense, for each detector a plot is created in the `outdir` called
-# H1_frequency_domain_data.png and LI_frequency_domain_data.png. The two
-# objects are then placed into a list.
-# For GW170817, 170608 and 170814 add the following line to select clean data.
-# kwargs = {"tag": 'CLN'}
-interferometers = bilby.gw.detector.get_event_data(label)
+
+# Data set up
+trigger_time = 1126259462
+
+roll_off = 0.4  # Roll off duration of tukey window in seconds, default is 0.4s
+duration = 4  # Analysis segment duration
+post_trigger_duration = 2  # Time between trigger time and end of segment
+end_time = trigger_time + post_trigger_duration
+start_time = end_time - duration
+
+psd_duration = 32 * duration
+psd_start_time = start_time - psd_duration
+psd_end_time = start_time
+
+# We now use gwpy to obtain analysis and psd data and create the ifo_list
+ifo_list = bilby.gw.detector.InterferometerList([])
+for det in ["H1", "L1"]:
+    logger.info("Downloading analysis data for ifo {}".format(det))
+    ifo = bilby.gw.detector.get_empty_interferometer(det)
+    data = TimeSeries.fetch_open_data(det, start_time, end_time)
+    ifo.strain_data.set_from_gwpy_timeseries(data)
+
+    logger.info("Downloading psd data for ifo {}".format(det))
+    psd_data = TimeSeries.fetch_open_data(det, psd_start_time, psd_end_time)
+    psd_alpha = 2 * roll_off / duration
+    psd = psd_data.psd(
+        fftlength=duration,
+        overlap=0,
+        window=("tukey", psd_alpha),
+        method="median"
+    )
+    ifo.power_spectral_density = bilby.gw.detector.PowerSpectralDensity(
+        frequency_array=psd.frequencies.value, psd_array=psd.value)
+    ifo_list.append(ifo)
+
+logger.info("Saving data plots to {}".format(outdir))
+bilby.core.utils.check_directory_exists_and_if_not_mkdir(outdir)
+ifo_list.plot_data(outdir=outdir, label=label)
 
 # We now define the prior.
-# We have defined our prior distribution in a file packaged with BILBY.
+# We have defined our prior distribution in a local file, GW150914.prior
 # The prior is printed to the terminal at run-time.
 # You can overwrite this using the syntax below in the file,
 # or choose a fixed value by just providing a float value as the prior.
-prior = bilby.gw.prior.BBHPriorDict(filename='GW150914.prior')
+priors = bilby.gw.prior.BBHPriorDict(filename='GW150914.prior')
 
-# In this step we define a `waveform_generator`. This is out object which
+# In this step we define a `waveform_generator`. This is the object which
 # creates the frequency-domain strain. In this instance, we are using the
 # `lal_binary_black_hole model` source model. We also pass other parameters:
-# the waveform approximant and reference frequency.
+# the waveform approximant and reference frequency and a parameter conversion
+# which allows us to sample in chirp mass and ratio rather than component mass
 waveform_generator = bilby.gw.WaveformGenerator(
     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_approximant': 'IMRPhenomPv2',
                         'reference_frequency': 50})
 
 # In this step, we define the likelihood. Here we use the standard likelihood
 # function, passing it the data and the waveform generator.
 likelihood = bilby.gw.likelihood.GravitationalWaveTransient(
-    interferometers, waveform_generator)
+    ifo_list, waveform_generator, priors=priors, time_marginalization=True,
+    phase_marginalization=True, distance_marginalization=True)
 
 # Finally, we run the sampler. This function takes the likelihood and prior
 # along with some options for how to do the sampling and how to save the data
-result = bilby.run_sampler(likelihood, prior, sampler='dynesty',
-                           outdir=outdir, label=label)
+result = bilby.run_sampler(
+    likelihood, priors, sampler='dynesty', outdir=outdir, label=label,
+    nlive=1000, walks=100, n_check_point=10000, check_point_plot=True,
+    conversion_function=bilby.gw.conversion.generate_all_bbh_parameters)
 result.plot_corner()