diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index dd5fbb4e7f103401579b787c47c54c54bdd14714..b88a77cdbcf2dfc1106d3725cd8638c6aebb4d00 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -22,11 +22,11 @@ exitcode-jessie:
- pip install coverage
- pip install coverage-badge
- coverage erase
- - coverage run -a test/prior_tests.py
- - coverage run -a test/tests.py
- - coverage run -a test/waveform_generator_tests.py
- - coverage run -a test/noise_realisation_tests.py
- - coverage html --include=tupak/*
+ - coverage run --include=tupak/* -a test/prior_tests.py
+ - coverage run --include=tupak/* -a test/tests.py
+ - coverage run --include=tupak/* -a test/waveform_generator_tests.py
+ - coverage run --include=tupak/* -a test/noise_realisation_tests.py
+ - coverage html
- coverage-badge -o coverage_badge.svg -f
artifacts:
paths:
diff --git a/README.md b/README.md
index 08f2f969933964909670540a79f82d24f6a831d5..3fd3060dbd9daab3d4282c8f41685a4eb22970b2 100644
--- a/README.md
+++ b/README.md
@@ -1,5 +1,5 @@
[](https://git.ligo.org/Monash/tupak/commits/master)
-[](
+[](
https://monash.docs.ligo.org/tupak/)
# Tupak
diff --git a/examples/other_examples/alternative_likelihoods.py b/examples/other_examples/alternative_likelihoods.py
index 1b6c0b9f120bd32d1d8bd32884683fa57c5b5238..329f9f339053ad67bc598a8cb83be0e169ab6f2a 100644
--- a/examples/other_examples/alternative_likelihoods.py
+++ b/examples/other_examples/alternative_likelihoods.py
@@ -20,7 +20,7 @@ outdir = 'outdir'
# use of the `tupak` waveform_generator to make the signal (more on this later)
# But, one could make this work without the waveform generator.
-class GaussianGravitationalWaveTransient():
+class GaussianLikelihood():
def __init__(self, x, y, waveform_generator):
self.x = x
self.y = y
@@ -76,7 +76,7 @@ waveform_generator = tupak.waveform_generator.WaveformGenerator(time_duration=ti
# Now lets instantiate a version of out GravitationalWaveTransient, giving it the time, data
# and waveform_generator
-likelihood = GaussianGravitationalWaveTransient(time, data, waveform_generator)
+likelihood = GaussianLikelihood(time, data, waveform_generator)
# From hereon, the syntax is exactly equivalent to other tupak examples
# We make a prior
diff --git a/test/prior_tests.py b/test/prior_tests.py
index bbed21f609404280932b7038c49a0e014e2cc3a4..d44f2fe075e63c962a76a9641fa30f09e9ea72ba 100644
--- a/test/prior_tests.py
+++ b/test/prior_tests.py
@@ -198,6 +198,7 @@ class TestFillPrior(unittest.TestCase):
def setUp(self):
self.likelihood = Mock()
self.likelihood.parameters = dict(a=0, b=0, c=0, d=0, asdf=0, ra=1)
+ self.likelihood.non_standard_sampling_parameter_keys = dict(t=8)
self.priors = dict(a=1, b=1.1, c='string', d=tupak.prior.Uniform(0, 1))
self.priors = tupak.prior.fill_priors(self.priors, self.likelihood)
diff --git a/test/tests.py b/test/tests.py
index 07f09b9bb5c57b4abfd372b0dc0e24334b8e9430..fabe4559b6d68cbdd980de3c611d21ce16c7da8f 100644
--- a/test/tests.py
+++ b/test/tests.py
@@ -57,8 +57,10 @@ class Test(unittest.TestCase):
priors['luminosity_distance'] = tupak.prior.Uniform(
name='luminosity_distance', minimum=dL - 10, maximum=dL + 10)
- result = tupak.sampler.run_sampler(likelihood, priors, sampler='nestle')
- self.assertAlmostEqual(np.mean(result.samples), dL, delta=np.std(result.samples))
+ result = tupak.sampler.run_sampler(
+ likelihood, priors, sampler='nestle', verbose=False, npoints=100)
+ self.assertAlmostEqual(np.mean(result.samples), dL,
+ delta=np.std(result.samples))
if __name__ == '__main__':
diff --git a/tupak/likelihood.py b/tupak/likelihood.py
index dd496935391f813b3ee6ab5429707dc724fb8e02..43b9937ff4af4c28c4a7538807842fe4fed2dfd9 100644
--- a/tupak/likelihood.py
+++ b/tupak/likelihood.py
@@ -11,7 +11,20 @@ import tupak
import logging
-class GravitationalWaveTransient(object):
+class Likelihood(object):
+ """ Empty likelihood class to be subclassed by other likelihoods """
+
+ def log_likelihood():
+ return np.nan
+
+ def noise_log_likelihood():
+ return np.nan
+
+ def log_likelihood_ratio(self):
+ return self.log_likelihood() - self.noise_log_likelihood()
+
+
+class GravitationalWaveTransient(Likelihood):
""" A gravitational-wave transient likelihood object
This is the usual likelihood object to use for transient gravitational
@@ -114,7 +127,7 @@ class GravitationalWaveTransient(object):
return log_l.real
def log_likelihood(self):
- return self.log_likelihood() + self.noise_log_likelihood()
+ return self.log_likelihood_ratio() + self.noise_log_likelihood()
def setup_distance_marginalization(self):
if 'luminosity_distance' not in self.prior.keys():
@@ -136,6 +149,62 @@ class GravitationalWaveTransient(object):
bounds_error=False, fill_value=-np.inf)
+class BasicGravitationalWaveTransient(Likelihood):
+ """ A basic gravitaitonal wave transient likelihood
+
+ The simplest frequency-domain gravitational wave transient likelihood. Does
+ not include distance/phase marginalization.
+
+ Parameters
+ ----------
+ interferometers: list
+ A list of `tupak.detector.Interferometer` instances - contains the
+ detector data and power spectral densities
+ waveform_generator: `tupak.waveform_generator.WaveformGenerator`
+ An object which computes the frequency-domain strain of the signal,
+ given some set of parameters
+
+ Returns
+ -------
+ Likelihood: `tupak.likelihood.Likelihood`
+ A likehood object, able to compute the likelihood of the data given
+ some model parameters
+
+ """
+ def __init__(self, interferometers, waveform_generator):
+ self.interferometers = interferometers
+ self.waveform_generator = waveform_generator
+
+ def noise_log_likelihood(self):
+ log_l = 0
+ for interferometer in self.interferometers:
+ log_l -= 2. / self.waveform_generator.time_duration * np.sum(
+ abs(interferometer.data) ** 2 /
+ interferometer.power_spectral_density_array)
+ return log_l.real
+
+ def log_likelihood(self):
+ log_l = 0
+ waveform_polarizations = self.waveform_generator.frequency_domain_strain()
+ if waveform_polarizations is None:
+ return np.nan_to_num(-np.inf)
+ for interferometer in self.interferometers:
+ log_l += self.log_likelihood_interferometer(
+ waveform_polarizations, interferometer)
+ return log_l.real
+
+ def log_likelihood_interferometer(self, waveform_polarizations,
+ interferometer):
+ signal_ifo = interferometer.get_detector_response(
+ waveform_polarizations, self.waveform_generator.parameters)
+
+ log_l = - 2. / self.waveform_generator.time_duration * np.vdot(
+ interferometer.data - signal_ifo,
+ (interferometer.data - signal_ifo)
+ / interferometer.power_spectral_density_array)
+ return log_l.real
+
+
def get_binary_black_hole_likelihood(interferometers):
""" A rapper to quickly set up a likelihood for BBH parameter estimation
@@ -157,3 +226,4 @@ def get_binary_black_hole_likelihood(interferometers):
parameters={'waveform_approximant': 'IMRPhenomPv2', 'reference_frequency': 50})
likelihood = tupak.likelihood.GravitationalWaveTransient(interferometers, waveform_generator)
return likelihood
+
diff --git a/tupak/prior.py b/tupak/prior.py
index a1c8665acae05f56b7605fbce4357b9cf188da2f..4546075d44ff6958d69629b2e552546c21511ab4 100644
--- a/tupak/prior.py
+++ b/tupak/prior.py
@@ -433,7 +433,7 @@ def parse_keys_to_parameters(keys):
return parameters
-def fill_priors(prior, likelihood, parameters=None):
+def fill_priors(prior, likelihood):
"""
Fill dictionary of priors based on required parameters of likelihood
@@ -446,9 +446,9 @@ def fill_priors(prior, likelihood, parameters=None):
dictionary of prior objects and floats
likelihood: tupak.likelihood.GravitationalWaveTransient instance
Used to infer the set of parameters to fill the prior with
- parameters: list
- list of parameters to be sampled in, this can override the default
- priors for the waveform generator
+
+ Note: if `likelihood` has `non_standard_sampling_parameter_keys`, then this
+ will set-up default priors for those as well.
Returns
-------
@@ -470,8 +470,8 @@ def fill_priors(prior, likelihood, parameters=None):
missing_keys = set(likelihood.parameters) - set(prior.keys())
- if parameters is not None:
- for parameter in parameters:
+ if getattr(likelihood, 'non_standard_sampling_parameter_keys', None) is not None:
+ for parameter in likelihood.non_standard_sampling_parameter_keys:
prior[parameter] = create_default_prior(parameter)
for missing_key in missing_keys:
diff --git a/tupak/result.py b/tupak/result.py
index 77bd42bb4c84080e1b22b36832023321841ee214..3ef71a34cebf5d37ba05652b026c0922aa3726e7 100644
--- a/tupak/result.py
+++ b/tupak/result.py
@@ -67,6 +67,9 @@ class Result(dict):
else:
return ''
+ def get_result_dictionary(self):
+ return dict(self)
+
def save_to_file(self, outdir, label):
""" Writes the Result to a deepdish h5 file """
file_name = result_file_name(outdir, label)
@@ -80,7 +83,7 @@ class Result(dict):
logging.info("Saving result to {}".format(file_name))
try:
- deepdish.io.save(file_name, self)
+ deepdish.io.save(file_name, self.get_result_dictionary())
except Exception as e:
logging.error(
"\n\n Saving the data has failed with the following message:\n {} \n\n"
diff --git a/tupak/sampler.py b/tupak/sampler.py
index 805cd87bf10b8e695f66a1d3001af36fb4e60c97..c6e9565b9574ec9ff0b86ddfc4b88d241373cbcf 100644
--- a/tupak/sampler.py
+++ b/tupak/sampler.py
@@ -245,8 +245,9 @@ class Nestle(Sampler):
def run_sampler(self):
nestle = self.external_sampler
self.external_sampler_function = nestle.sample
- if self.kwargs.get('verbose', True):
- self.kwargs['callback'] = nestle.print_progress
+ if 'verbose' in self.kwargs:
+ if self.kwargs['verbose']:
+ self.kwargs['callback'] = nestle.print_progress
self.kwargs.pop('verbose')
self.verify_kwargs_against_external_sampler_function()
@@ -444,7 +445,7 @@ def run_sampler(likelihood, priors=None, label='label', outdir='outdir',
if priors is None:
priors = dict()
- priors = fill_priors(priors, likelihood, parameters=likelihood.non_standard_sampling_parameter_keys)
+ priors = fill_priors(priors, likelihood)
tupak.prior.write_priors_to_file(priors, outdir)
if implemented_samplers.__contains__(sampler.title()):
@@ -468,7 +469,7 @@ def run_sampler(likelihood, priors=None, label='label', outdir='outdir',
result.injection_parameters = injection_parameters
if conversion_function is not None:
conversion_function(result.injection_parameters)
- result.fixed_parameter_keys = [key for key in priors if isinstance(key, prior.DeltaFunction)]
+ result.fixed_parameter_keys = sampler.fixed_parameter_keys
# result.prior = prior # Removed as this breaks the saving of the data
result.samples_to_data_frame(likelihood=likelihood, priors=priors, conversion_function=conversion_function)
result.kwargs = sampler.kwargs
diff --git a/tupak/waveform_generator.py b/tupak/waveform_generator.py
index 30667d84050c34500034fbd736ceb6c49fce8f88..7ffb7c84e394514addd5d49eb7325c6b63a1538b 100644
--- a/tupak/waveform_generator.py
+++ b/tupak/waveform_generator.py
@@ -89,7 +89,7 @@ class WaveformGenerator(object):
if self.parameter_conversion is not None:
for key in added_keys:
self.parameters.pop(key)
- return model_time_series[key]
+ return model_time_series
@property
def frequency_array(self):