Note
Click here to download the full example code
Software injection into pre-existing data filesΒΆ
Add a software injection into a set of SFTs.
In this case, the set of SFTs is generated using Makefakedata_v5, but the same procedure can be applied to any other set of SFTs (including real detector data).
12 import os
13
14 import numpy as np
15
16 import pyfstat
17
18 label = "PyFstat_example_injection_into_noise_sfts"
19 outdir = os.path.join("PyFstat_example_data", label)
20 logger = pyfstat.set_up_logger(label=label, outdir=outdir)
21
22 tstart = 1269734418
23 duration_Tsft = 100
24 Tsft = 1800
25 randSeed = 69420
26 IFO = "H1"
27 h0 = 1000
28 cosi = 0
29 F0 = 30
30 Alpha = 0
31 Delta = 0
32
33 Band = 2.0
34
35 # create sfts with a strong signal in them
36 # window options are optional here
37 noise_and_signal_writer = pyfstat.Writer(
38 label="test_noiseSFTs_noise_and_signal",
39 outdir=outdir,
40 h0=h0,
41 cosi=cosi,
42 F0=F0,
43 Alpha=Alpha,
44 Delta=Delta,
45 tstart=tstart,
46 duration=duration_Tsft * Tsft,
47 Tsft=Tsft,
48 Band=Band,
49 detectors=IFO,
50 randSeed=randSeed,
51 SFTWindowType="tukey",
52 SFTWindowBeta=0.001,
53 )
54 sftfilepattern = os.path.join(
55 noise_and_signal_writer.outdir,
56 "*{}*{}*sft".format(duration_Tsft, noise_and_signal_writer.label),
57 )
58
59 noise_and_signal_writer.make_data()
60
61 # compute Fstat
62 coherent_search = pyfstat.ComputeFstat(
63 tref=noise_and_signal_writer.tref,
64 sftfilepattern=sftfilepattern,
65 minCoverFreq=-0.5,
66 maxCoverFreq=-0.5,
67 )
68 FS_1 = coherent_search.get_fullycoherent_twoF(
69 noise_and_signal_writer.F0,
70 noise_and_signal_writer.F1,
71 noise_and_signal_writer.F2,
72 noise_and_signal_writer.Alpha,
73 noise_and_signal_writer.Delta,
74 )
75
76 # create noise sfts
77 # window options are again optional for this step
78 noise_writer = pyfstat.Writer(
79 label="test_noiseSFTs_only_noise",
80 outdir=outdir,
81 h0=0,
82 F0=F0,
83 tstart=tstart,
84 duration=duration_Tsft * Tsft,
85 Tsft=Tsft,
86 Band=Band,
87 detectors=IFO,
88 randSeed=randSeed,
89 SFTWindowType="tukey",
90 SFTWindowBeta=0.001,
91 )
92 noise_writer.make_data()
93
94 # then inject a strong signal
95 # window options *must* match those previously used for the noiseSFTs
96 add_signal_writer = pyfstat.Writer(
97 label="test_noiseSFTs_add_signal",
98 outdir=outdir,
99 F0=F0,
100 Alpha=Alpha,
101 Delta=Delta,
102 h0=h0,
103 cosi=cosi,
104 tstart=tstart,
105 duration=duration_Tsft * Tsft,
106 Tsft=Tsft,
107 Band=Band,
108 detectors=IFO,
109 sqrtSX=0,
110 noiseSFTs=os.path.join(
111 noise_writer.outdir, "*{}*{}*sft".format(duration_Tsft, noise_writer.label)
112 ),
113 SFTWindowType="tukey",
114 SFTWindowBeta=0.001,
115 )
116 sftfilepattern = os.path.join(
117 add_signal_writer.outdir,
118 "*{}*{}*sft".format(duration_Tsft, add_signal_writer.label),
119 )
120 add_signal_writer.make_data()
121
122 # compute Fstat
123 coherent_search = pyfstat.ComputeFstat(
124 tref=add_signal_writer.tref,
125 sftfilepattern=sftfilepattern,
126 minCoverFreq=-0.5,
127 maxCoverFreq=-0.5,
128 )
129 FS_2 = coherent_search.get_fullycoherent_twoF(
130 add_signal_writer.F0,
131 add_signal_writer.F1,
132 add_signal_writer.F2,
133 add_signal_writer.Alpha,
134 add_signal_writer.Delta,
135 )
136
137 logger.info("Base case Fstat: {}".format(FS_1))
138 logger.info("Noise + Signal Fstat: {}".format(FS_2))
139 logger.info("Relative Difference: {}".format(np.abs(FS_2 - FS_1) / FS_1))
Total running time of the script: ( 0 minutes 0.000 seconds)