Merge branch 'master' into sampler_rework

examples/injection_examples/basic_tutorial.py

@@ -9,15 +9,21 @@ from __future__ import division, print_function
 import tupak
 import numpy as np
 
-
+# Set the duration and sampling frequency of the data segment that we're going to inject the signal into
 time_duration = 4.
 sampling_frequency = 2048.
+
+# Specify the output directory and the name of the simulation.
 outdir = 'outdir'
 label = 'basic_tutorial'
 tupak.utils.setup_logger(outdir=outdir, label=label)
 
-np.random.seed(170809)
+# Set up a random seed for result reproducibility.  This is optional!
+np.random.seed(170801)
 
+# 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.5, tilt_2=1.0, phi_12=1.7, phi_jl=0.3,
                             luminosity_distance=4000., iota=0.4, psi=2.659, phase=1.3, geocent_time=1126259642.413,
                             waveform_approximant='IMRPhenomPv2', reference_frequency=50., ra=1.375, dec=-1.2108)
@@ -29,24 +35,33 @@ waveform_generator = tupak.waveform_generator.WaveformGenerator(time_duration=ti
                                                                 parameters=injection_parameters)
 hf_signal = waveform_generator.frequency_domain_strain()
 
-# Set up interferometers.
+# Set up interferometers.  In this case we'll use three interferometers (LIGO-Hanford (H1), LIGO-Livingston (L1),
+# and Virgo (V1)).  These default to their design sensitivity
 IFOs = [tupak.detector.get_interferometer_with_fake_noise_and_injection(
     name, injection_polarizations=hf_signal, injection_parameters=injection_parameters, time_duration=time_duration,
     sampling_frequency=sampling_frequency, outdir=outdir) for name in ['H1', 'L1', 'V1']]
 
-# Set up prior
+# Set up prior, which is a dictionary
 priors = dict()
-# These parameters will not be sampled
-for key in ['a_1', 'a_2', 'tilt_1', 'tilt_2', 'phi_12', 'phi_jl', 'phase', 'psi', 'iota', 'ra', 'dec', 'geocent_time']:
+# By default we will sample all terms in the signal models.  However, this will take a long time for the calculation,
+# so for this example we will set almost all of the priors to be equall to their injected values.  This implies the
+# prior is a delta function at the true, injected value.  In reality, the sampler implementation is smart enough to
+# not sample any parameter that has a delta-function prior.
+for key in ['a_1', 'a_2', 'tilt_1', 'tilt_2', 'phi_12', 'phi_jl', 'phase', 'psi', 'ra', 'dec', 'geocent_time']:
     priors[key] = injection_parameters[key]
+
+# The above list does *not* include mass_1, mass_2, iota and luminosity_distance, which means those are the parameters
+# that will be included in the sampler.  If we do nothing, then the default priors get used.
 priors['luminosity_distance'] = tupak.prior.create_default_prior(name='luminosity_distance')
 
-# Initialise Likelihood
+# Initialise the likelihood by passing in the interferometer data (IFOs) and the waveoform generator
 likelihood = tupak.likelihood.Likelihood(interferometers=IFOs, waveform_generator=waveform_generator)
 
-# Run sampler
+# Run sampler.  In this case we're going to use the `dynesty` sampler
 result = tupak.sampler.run_sampler(likelihood=likelihood, priors=priors, sampler='dynesty', npoints=1000,
                                    injection_parameters=injection_parameters, outdir=outdir, label=label)
+
+# make some plots of the outputs
 result.plot_corner()
 result.plot_walks()
 result.plot_distributions()

tupak/result.py

@@ -17,13 +17,26 @@ def result_file_name(outdir, label):
     return '{}/{}_result.h5'.format(outdir, label)
 
 
-def read_in_result(outdir, label):
-    """ Read in a saved .h5 data file """
-    filename = result_file_name(outdir, label)
+def read_in_result(outdir=None, label=None, filename=None):
+    """ Read in a saved .h5 data file
+
+    Parameters
+    ----------
+    outdir, label: str
+        If given, use the default naming convention for saved results file
+    filename: str
+        If given, try to load from this filename
+
+    Returns:
+    result: tupak.result.Result instance
+
+    """
+    if filename is None:
+        filename = result_file_name(outdir, label)
     if os.path.isfile(filename):
         return Result(deepdish.io.load(filename))
     else:
-        return None
+        raise ValueError("No information given to load file")
 
 
 class Result(dict):

tupak/sampler.py

@@ -193,7 +193,12 @@ class Sampler(object):
             logging.debug("Command line argument clean given, forcing rerun")
             self.cached_result = None
             return
-        self.cached_result = read_in_result(self.outdir, self.label)
+
+        try:
+            self.cached_result = read_in_result(self.outdir, self.label)
+        except ValueError:
+            self.cached_result = None
+
         if utils.command_line_args.use_cached:
             logging.debug("Command line argument cached given, no cache check performed")
             return