Note
Click here to download the full example code
Long transient MCMC search
MCMC search for a long transient CW signal.
7 import os
8
9 import numpy as np
10 import PyFstat_example_make_data_for_long_transient_search as data
11
12 import pyfstat
13 from pyfstat.helper_functions import get_predict_fstat_parameters_from_dict
14
15 if not os.path.isdir(data.outdir) or not np.any(
16 [f.endswith(".sft") for f in os.listdir(data.outdir)]
17 ):
18 raise RuntimeError(
19 "Please first run PyFstat_example_make_data_for_long_transient_search.py !"
20 )
21
22 tstart = data.tstart
23 duration = data.duration
24
25 inj = {
26 "tref": data.tstart,
27 "F0": data.F0,
28 "F1": data.F1,
29 "F2": data.F2,
30 "Alpha": data.Alpha,
31 "Delta": data.Delta,
32 "transient_tstart": data.transient_tstart,
33 "transient_duration": data.transient_duration,
34 }
35
36 DeltaF0 = 6e-7
37 DeltaF1 = 1e-13
38
39 # to make the search cheaper, we exactly target the transientStartTime
40 # to the injected value and only search over TransientTau
41 theta_prior = {
42 "F0": {
43 "type": "unif",
44 "lower": inj["F0"] - DeltaF0 / 2.0,
45 "upper": inj["F0"] + DeltaF0 / 2.0,
46 },
47 "F1": {
48 "type": "unif",
49 "lower": inj["F1"] - DeltaF1 / 2.0,
50 "upper": inj["F1"] + DeltaF1 / 2.0,
51 },
52 "F2": inj["F2"],
53 "Alpha": inj["Alpha"],
54 "Delta": inj["Delta"],
55 "transient_tstart": tstart + 0.25 * duration,
56 "transient_duration": {
57 "type": "halfnorm",
58 "loc": 0.001 * duration,
59 "scale": 0.5 * duration,
60 },
61 }
62
63 ntemps = 2
64 log10beta_min = -1
65 nwalkers = 100
66 nsteps = [100, 100]
67
68 transientWindowType = "rect"
69
70 mcmc = pyfstat.MCMCTransientSearch(
71 label="transient_search",
72 outdir=data.outdir,
73 sftfilepattern=os.path.join(data.outdir, "*simulated_transient_signal*sft"),
74 theta_prior=theta_prior,
75 tref=inj["tref"],
76 nsteps=nsteps,
77 nwalkers=nwalkers,
78 ntemps=ntemps,
79 log10beta_min=log10beta_min,
80 transientWindowType=transientWindowType,
81 )
82 mcmc.run(walker_plot_args={"plot_det_stat": True, "injection_parameters": inj})
83 mcmc.print_summary()
84 mcmc.plot_corner(add_prior=True, truths=inj)
85 mcmc.plot_prior_posterior(injection_parameters=inj)
86
87 # plot cumulative 2F, first building a dict as required for PredictFStat
88 d, maxtwoF = mcmc.get_max_twoF()
89 for key, val in mcmc.theta_prior.items():
90 if key not in d:
91 d[key] = val
92 d["h0"] = data.h0
93 d["cosi"] = data.cosi
94 d["psi"] = data.psi
95 PFS_input = get_predict_fstat_parameters_from_dict(
96 d, transientWindowType=transientWindowType
97 )
98 mcmc.plot_cumulative_max(PFS_input=PFS_input)
Total running time of the script: ( 0 minutes 0.000 seconds)