Follow up example

Multi-stage MCMC follow up of a CW signal produced by an isolated source using a ladder of coherent times.

  8 import os
  9
 10 import matplotlib.pyplot as plt
 11 import numpy as np
 12
 13 import pyfstat
 14
 15 label = "PyFstat_example_semi_coherent_directed_follow_up"
 16 outdir = os.path.join("PyFstat_example_data", label)
 17
 18 # Properties of the GW data
 19 data_parameters = {
 20     "sqrtSX": 1e-23,
 21     "tstart": 1000000000,
 22     "duration": 100 * 86400,
 23     "detectors": "H1",
 24 }
 25 tend = data_parameters["tstart"] + data_parameters["duration"]
 26 mid_time = 0.5 * (data_parameters["tstart"] + tend)
 27
 28 # Properties of the signal
 29 depth = 40
 30 signal_parameters = {
 31     "F0": 30.0,
 32     "F1": -1e-10,
 33     "F2": 0,
 34     "Alpha": np.radians(83.6292),
 35     "Delta": np.radians(22.0144),
 36     "tref": mid_time,
 37     "h0": data_parameters["sqrtSX"] / depth,
 38     "cosi": 1.0,
 39 }
 40
 41 data = pyfstat.Writer(
 42     label=label, outdir=outdir, **data_parameters, **signal_parameters
 43 )
 44 data.make_data()
 45
 46 # The predicted twoF, given by lalapps_predictFstat can be accessed by
 47 twoF = data.predict_fstat()
 48 print("Predicted twoF value: {}\n".format(twoF))
 49
 50 # Search
 51 VF0 = VF1 = 1e5
 52 DeltaF0 = np.sqrt(VF0) * np.sqrt(3) / (np.pi * data_parameters["duration"])
 53 DeltaF1 = np.sqrt(VF1) * np.sqrt(180) / (np.pi * data_parameters["duration"] ** 2)
 54 theta_prior = {
 55     "F0": {
 56         "type": "unif",
 57         "lower": signal_parameters["F0"] - DeltaF0 / 2.0,
 58         "upper": signal_parameters["F0"] + DeltaF0 / 2,
 59     },
 60     "F1": {
 61         "type": "unif",
 62         "lower": signal_parameters["F1"] - DeltaF1 / 2.0,
 63         "upper": signal_parameters["F1"] + DeltaF1 / 2,
 64     },
 65 }
 66 for key in "F2", "Alpha", "Delta":
 67     theta_prior[key] = signal_parameters[key]
 68
 69
 70 ntemps = 3
 71 log10beta_min = -0.5
 72 nwalkers = 100
 73 nsteps = [100, 100]
 74
 75 mcmc = pyfstat.MCMCFollowUpSearch(
 76     label=label,
 77     outdir=outdir,
 78     sftfilepattern=os.path.join(outdir, "*{}*sft".format(label)),
 79     theta_prior=theta_prior,
 80     tref=mid_time,
 81     minStartTime=data_parameters["tstart"],
 82     maxStartTime=tend,
 83     nwalkers=nwalkers,
 84     nsteps=nsteps,
 85     ntemps=ntemps,
 86     log10beta_min=log10beta_min,
 87 )
 88
 89 NstarMax = 1000
 90 Nsegs0 = 100
 91 walkers_fig, walkers_axes = plt.subplots(nrows=2, figsize=(3.4, 3.5))
 92 mcmc.run(
 93     NstarMax=NstarMax,
 94     Nsegs0=Nsegs0,
 95     plot_walkers=True,
 96     walker_plot_args={
 97         "labelpad": 0.01,
 98         "plot_det_stat": False,
 99         "fig": walkers_fig,
100         "axes": walkers_axes,
101         "injection_parameters": signal_parameters,
102     },
103 )
104 walkers_fig.savefig(os.path.join(outdir, label + "_walkers.png"))
105 plt.close(walkers_fig)
106
107 mcmc.print_summary()
108 mcmc.plot_corner(add_prior=True, truths=signal_parameters)
109 mcmc.plot_prior_posterior(injection_parameters=signal_parameters)

Total running time of the script: ( 0 minutes 0.000 seconds)

Gallery generated by Sphinx-Gallery