Commit 60305296 authored by Gregory Ashton's avatar Gregory Ashton
Browse files

Merge branch 'master' into add-process-method-for-mcmc

parents a01e581b 6e733890
......@@ -213,7 +213,7 @@ class InjectionCreator(Input):
uncertainty=self.deltaT / 2.0,
)
geocent_times.append(geocent_time)
inj_df["geocenter_times"] = geocent_times
inj_df["geocent_time"] = geocent_times
return inj_df
@staticmethod
......
......@@ -168,6 +168,7 @@ class DataGenerationInput(Input):
args.spline_calibration_phase_uncertainty_dict
)
self.spline_calibration_nodes = args.spline_calibration_nodes
self.calibration_prior_boundary = args.calibration_prior_boundary
# Marginalization
self.distance_marginalization = args.distance_marginalization
......@@ -520,6 +521,7 @@ class DataGenerationInput(Input):
waveform_generator=waveform_generator,
parameters=parameters,
det=ifo.name,
power_spectral_density=ifo.power_spectral_density,
outdir=outdir,
label=label,
)
......
......@@ -980,6 +980,7 @@ class Input(object):
maximum_frequency=self.maximum_frequency_dict[det],
n_nodes=self.spline_calibration_nodes,
label=det,
boundary=self.calibration_prior_boundary,
)
)
elif (
......
......@@ -18,9 +18,11 @@ from .overview import create_overview
def get_trigger_time_list(inputs):
""" Returns a list of GPS trigger times for each data segment """
if inputs.gaussian_noise and inputs.trigger_time is None:
if (inputs.gaussian_noise or inputs.zero_noise) and inputs.trigger_time is None:
trigger_times = [0] * inputs.n_simulation
elif inputs.gaussian_noise and isinstance(inputs.trigger_time, float):
elif (inputs.gaussian_noise or inputs.zero_noise) and isinstance(
inputs.trigger_time, float
):
trigger_times = [inputs.trigger_time] * inputs.n_simulation
elif inputs.trigger_time is not None:
trigger_times = [inputs.trigger_time]
......
#!/usr/bin/env python
"""
This a command line tool to prepare bilby results for use with ligo skymap
"""
import argparse
import os
import numpy as np
from h5py import File
from bilby.core.result import read_in_result
from .utils import logger
def get_args():
parser = argparse.ArgumentParser(
prog="bilby_pipe_to_ligo_skymap_samples", description=__doc__
)
parser.add_arg("input_file", type=str, help="The bilby file to convert")
parser.add_arg("-o", "--out", type=str, help="The output hdf5 filename")
parser.add_arg(
"-n", "--nsamples", type=int, help="The maximum number of samples", default=None
)
args = parser.parse_args()
return args
def main():
args = get_args()
all_keys = [
"ra",
"dec",
"luminosity_distance",
"time",
"mass_1",
"mass_2",
"spin1z",
"spin2z",
]
logger.info(f"Converting bilby result file {args.input_file}")
result = read_in_result(args.input_file)
posterior = result.posterior
if "geocent_time" in posterior:
posterior["time"] = posterior.pop("geocent_time")
keys = [key for key in all_keys if key in posterior]
samples = posterior[keys].to_records(index=False)
if os.path.exists(result.outdir) is False:
logger.info(
f"The directory {result.outdir} is not accessible, falling back to"
" the current working directory"
)
result.outdir = "."
if args.out:
output_file = args.out
else:
output_file = f"{result.outdir}/{result.label}_posterior_samples.hdf5"
n = len(samples)
if args.nsamples is not None and args.nsamples < n:
samples = np.random.choice(samples, args.nsamples)
n = len(samples)
logger.info(f"Writing {n} ligo-skymap ready samples to {output_file}")
with File(output_file, "w") as ff:
ff.create_dataset(
"posterior_samples",
shape=samples.shape,
dtype=samples.dtype,
data=samples,
)
......@@ -106,6 +106,12 @@ def create_parser(top_level=True):
"uncertainty model"
),
)
calibration_parser.add(
"--calibration-prior-boundary",
type=nonestr,
default="reflective",
help="Boundary methods for the calibration prior boundary",
)
if top_level is False:
parser.add("--idx", type=int, help="The level A job index", default=0)
......@@ -389,6 +395,7 @@ def create_parser(top_level=True):
injection_parser.add(
"--injection-numbers",
action="append",
type=nonestr,
default=None,
help=(
"Specific injections rows to use from the injection_file, e.g. "
......
......@@ -6,6 +6,7 @@ This a command line tool to convert XML injection files
import argparse
import json
import os
from math import pi
import lalsimulation as lalsim
import pandas as pd
......@@ -19,7 +20,7 @@ except ImportError:
raise ImportError("You do not have ligo.lw install: $ pip install python-liw-lw")
def xml_to_dataframe(prior_file, reference_frequency):
def xml_to_dataframe(prior_file, reference_frequency, convert_negative_ra=False):
table = Table.read(prior_file, format="ligolw", tablename="sim_inspiral")
injection_values = {
"mass_1": [],
......@@ -49,7 +50,10 @@ def xml_to_dataframe(prior_file, reference_frequency):
float(row["geocent_end_time"])
+ float(row["geocent_end_time_ns"]) * (10 ** -9)
)
injection_values["ra"].append(float(row["longitude"]))
if convert_negative_ra and float(row["longitude"]) < 0:
injection_values["ra"].append(float(row["longitude"]) + 2 * pi)
else:
injection_values["ra"].append(float(row["longitude"]))
injection_values["dec"].append(float(row["latitude"]))
args_list = [
......@@ -108,9 +112,18 @@ def main():
help=("The reference frequency to use for converting from xml"),
required=True,
)
parser.add_arg(
"--convert-negative-ra",
default=False,
help=("Convert (-pi,pi) RA range from (0,2pi)"),
action="store_true",
required=False,
)
args = parser.parse_args()
injection_values = xml_to_dataframe(args.xml_file, args.reference_frequency)
injection_values = xml_to_dataframe(
args.xml_file, args.reference_frequency, args.convert_negative_ra
)
basename = os.path.splitext(args.xml_file)[0]
path = basename + os.path.extsep + args.format
if args.format == "json":
......
......@@ -122,6 +122,7 @@ setup(
"bilby_pipe_gracedb=bilby_pipe.gracedb:main",
"bilby_pipe_write_default_ini=bilby_pipe.parser:main",
"bilby_pipe_process_mcmc=bilby_pipe.process_bilby_mcmc:main",
"bilby_pipe_to_ligo_skymap_samples=bilby_pipe.ligo_skymap:main",
]
},
classifiers=[
......
......@@ -130,7 +130,7 @@ class TestCreateInjections(unittest.TestCase):
len(df.columns.values) > 1, f"Column names: {df.columns.values}"
)
self.assertAlmostEqual(
df["geocenter_times"].iloc[0] / 100, gps_vals[0] / 100, places=1
df["geocent_time"].iloc[0] / 100, gps_vals[0] / 100, places=1
)
def test_create_injection_file_json(self):
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment