change_sampled_parameters.py

#!/usr/bin/env python
"""
Tutorial to demonstrate running parameter estimation sampling in non-standard
parameters for an injected signal.

This example estimates the masses using a uniform prior in chirp mass,
mass ratio and redshift.
The cosmology is according to the Planck 2015 data release.
"""
from __future__ import division, print_function
import tupak
import numpy as np


tupak.core.utils.setup_logger(log_level="info")

duration = 4.
sampling_frequency = 2048.
outdir = 'outdir'

np.random.seed(151226)

injection_parameters = dict(
    total_mass=66., mass_ratio=0.9, 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=2000, iota=0.4, psi=2.659,
    phase=1.3, geocent_time=1126259642.413, ra=1.375, dec=-1.2108)

waveform_arguments = dict(waveform_approximant='IMRPhenomPv2',
                          reference_frequency=50.)

# Create the waveform_generator using a LAL BinaryBlackHole source function
waveform_generator = tupak.gw.waveform_generator.WaveformGenerator(
    sampling_frequency=sampling_frequency, duration=duration,
    frequency_domain_source_model=tupak.gw.source.lal_binary_black_hole,
    parameter_conversion=
        tupak.gw.conversion.convert_to_lal_binary_black_hole_parameters,
    waveform_arguments=waveform_arguments)

# Set up interferometers.
ifos = tupak.gw.detector.InterferometerList(['H1', 'L1', 'V1', 'K1'])
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)

# Set up prior
priors = tupak.gw.prior.BBHPriorSet()
priors.pop('mass_1')
priors.pop('mass_2')
priors.pop('luminosity_distance')
priors['chirp_mass'] = tupak.prior.Uniform(
    name='chirp_mass', latex_label='$m_c$', minimum=13, maximum=45,
    unit='$M_{\\odot}$')
priors['symmetric_mass_ratio'] = tupak.prior.Uniform(
    name='symmetric_mass_ratio', latex_label='q', minimum=0.1, maximum=0.25)
priors['redshift'] = tupak.prior.Uniform(
    name='redshift', latex_label='$z$', minimum=0, maximum=0.5)
# These parameters will not be sampled
for key in ['a_1', 'a_2', 'tilt_1', 'tilt_2', 'phi_12', 'phi_jl', 'psi',
            'ra', 'dec', 'geocent_time', 'phase']:
    priors[key] = injection_parameters[key]
priors.pop('iota')
priors['cos_iota'] = np.cos(injection_parameters['iota'])
print(priors)

# Initialise GravitationalWaveTransient
likelihood = tupak.gw.likelihood.GravitationalWaveTransient(
    interferometers=ifos, waveform_generator=waveform_generator)

# Run sampler
result = tupak.core.sampler.run_sampler(
    likelihood=likelihood, priors=priors, sampler='dynesty', outdir=outdir,
    injection_parameters=injection_parameters, label='DifferentParameters',
    conversion_function=tupak.gw.conversion.generate_all_bbh_parameters)
result.plot_corner()