Note
Click here to download the full example code
Randomly sampling parameter space points
Application of dedicated classes to sample software injection parameters according to the specified parameter space priors.
8 import os
9 import numpy as np
10 import matplotlib.pyplot as plt
11 from pyfstat import (
12 InjectionParametersGenerator,
13 AllSkyInjectionParametersGenerator,
14 Writer,
15 )
16
17 label = "PyFstat_example_InjectionParametersGenerator"
18 outdir = os.path.join("PyFstat_example_data", label)
19
20 # Properties of the GW data
21 gw_data = {
22 "sqrtSX": 1e-23,
23 "tstart": 1000000000,
24 "duration": 86400,
25 "detectors": "H1,L1",
26 "Band": 1,
27 "Tsft": 1800,
28 }
29
30 print("Drawing random signal parameters...")
31
32 # Draw random signal phase parameters.
33 # The AllSkyInjectionParametersGenerator covers [Alpha,Delta] priors automatically.
34 # The rest can be a mix of nontrivial prior distributions and fixed values.
35 phase_params_generator = AllSkyInjectionParametersGenerator(
36 priors={
37 "F0": {"uniform": {"low": 29.0, "high": 31.0}},
38 "F1": -1e-10,
39 "F2": 0,
40 },
41 seed=23,
42 )
43 phase_parameters = phase_params_generator.draw()
44 phase_parameters["tref"] = gw_data["tstart"]
45
46 # Draw random signal amplitude parameters.
47 # Here we use the plain InjectionParametersGenerator class.
48 amplitude_params_generator = InjectionParametersGenerator(
49 priors={
50 "h0": {"normal": {"loc": 1e-24, "scale": 1e-24}},
51 "cosi": {"uniform": {"low": 0.0, "high": 1.0}},
52 "phi": {"uniform": {"low": 0.0, "high": 2 * np.pi}},
53 "psi": {"uniform": {"low": 0.0, "high": np.pi}},
54 },
55 seed=42,
56 )
57 amplitude_parameters = amplitude_params_generator.draw()
58
59 # Now we can pass the parameter dictionaries to the Writer class and make SFTs.
60 data = Writer(
61 label=label,
62 outdir=outdir,
63 **gw_data,
64 **phase_parameters,
65 **amplitude_parameters,
66 )
67 data.make_data()
68
69 # Now we draw many phase parameters and check the sky distribution
70 Ndraws = 10000
71 phase_parameters = [phase_params_generator.draw() for n in range(Ndraws)]
72 Alphas = np.array([p["Alpha"] for p in phase_parameters])
73 Deltas = np.array([p["Delta"] for p in phase_parameters])
74 plotfile = os.path.join(outdir, label + "_allsky.png")
75 print(f"Plotting sky distribution of {Ndraws} points to file: {plotfile}")
76 plt.subplot(111, projection="aitoff")
77 plt.plot(Alphas - np.pi, Deltas, ".", markersize=1)
78 plt.savefig(plotfile, dpi=300)
79 plt.close()
80 plotfile = os.path.join(outdir, label + "_alpha_hist.png")
81 print(f"Plotting Alpha distribution of {Ndraws} points to file: {plotfile}")
82 plt.hist(Alphas, 50)
83 plt.xlabel("Alpha")
84 plt.ylabel("draws")
85 plt.savefig(plotfile, dpi=100)
86 plt.close()
87 plotfile = os.path.join(outdir, label + "_delta_hist.png")
88 print(f"Plotting Delta distribution of {Ndraws} points to file: {plotfile}")
89 plt.hist(Deltas, 50)
90 plt.xlabel("Delta")
91 plt.ylabel("draws")
92 plt.savefig(plotfile, dpi=100)
93 plt.close()
94 plotfile = os.path.join(outdir, label + "_sindelta_hist.png")
95 print(f"Plotting sin(Delta) distribution of {Ndraws} points to file: {plotfile}")
96 plt.hist(np.sin(Deltas), 50)
97 plt.xlabel("sin(Delta)")
98 plt.ylabel("draws")
99 plt.savefig(plotfile, dpi=100)
100 plt.close()
Total running time of the script: ( 0 minutes 0.000 seconds)