LALPulsar  6.1.0.1-89842e6
test_heterodyned_model.py
Go to the documentation of this file.
1 # Copyright (C) 2019 Matthew Pitkin
2 #
3 # This program is free software; you can redistribute it and/or modify it
4 # under the terms of the GNU General Public License as published by the
5 # Free Software Foundation; either version 2 of the License, or (at your
6 # option) any later version.
7 #
8 # This program is distributed in the hope that it will be useful, but
9 # WITHOUT ANY WARRANTY; without even the implied warranty of
10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11 # Public License for more details.
12 #
13 # You should have received a copy of the GNU General Public License along
14 # with this program; if not, write to the Free Software Foundation, Inc.,
15 # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
16 
17 """
18 Test script for the heterodyned pulsar model functions.
19 
20 These functions are tested against fiducial outputs from the reviewed
21 lalpulsar_parameter_estimation_nested code.
22 """
23 
24 from __future__ import print_function, division
25 
26 import os
27 import sys
28 
29 import numpy as np
30 from numpy.testing import assert_allclose, assert_equal
31 
32 import pytest
33 
34 import lal
35 
36 import lalpulsar
37 from lalpulsar.PulsarParametersWrapper import PulsarParametersPy
38 from lalpulsar.simulateHeterodynedCW import HeterodynedCWSimulator
39 
40 """
41 The first test output is to create a signal model for a source assuming that
42 the heterodyne precisely matches the signal, i.e., it only varies due to the
43 antenna pattern. lalpulsar_parameter_estimation_nested has been run with
44 the following pulsar parameter "inj.par" file:
45 
46 PSRJ TEST
47 RAJ 01:23:34.5
48 DECJ -45:01:23.4
49 F0 123.456789
50 F1 -9.87654321e-12
51 PEPOCH 58000
52 H0 5.6e-26
53 COSIOTA -0.2
54 PSI 0.4
55 PHI0 2.3
56 EPHEM DE405
57 UNITS TCB
58 
59 lalpulsar_parameter_estimation_nested --par-file inj.par --outfile test.hdf --harmonics 2 --inject-file inj.par --inject-only --inject-output inj.txt --fake-data H1,L1 --fake-starts 1000000000,1000000000 --fake-lengths 86400,86400 --fake-dt 3600,3600
60 
61 The outputs of "inj.txt_H1_2.0_signal_only" and "inj.txt_L1_2.0_signal_only"
62 are given below.
63 """
64 
65 t1output = {
66  "H1": np.array(
67  [
68  [1000000000.0, -4.365015078242e-27, -4.216860020522e-27],
69  [1000003600.0, -3.123840987225e-27, -7.303345528632e-27],
70  [1000007200.0, -1.392894194743e-27, -8.343779137631e-27],
71  [1000010800.0, 3.907780411954e-28, -7.371234643743e-27],
72  [1000014400.0, 1.803192456526e-27, -4.935204629063e-27],
73  [1000018000.0, 2.543802234234e-27, -1.933962480155e-27],
74  [1000021600.0, 2.510964710484e-27, 6.441088928009e-28],
75  [1000025200.0, 1.822605255722e-27, 1.999009861526e-27],
76  [1000028800.0, 7.770155400880e-28, 1.741944189379e-27],
77  [1000032400.0, -2.352775438513e-28, 1.624267692628e-30],
78  [1000036000.0, -8.441413982604e-28, -2.614474888096e-27],
79  [1000039600.0, -8.062469575161e-28, -5.193209436791e-27],
80  [1000043200.0, -7.603746349482e-29, -6.776100558983e-27],
81  [1000046800.0, 1.178187576640e-27, -6.635532319641e-27],
82  [1000050400.0, 2.617633050484e-27, -4.491492229213e-27],
83  [1000054000.0, 3.824370674045e-27, -6.087150052676e-28],
84  [1000057600.0, 4.415959728830e-27, 4.253199979849e-27],
85  [1000061200.0, 4.152430242107e-27, 9.024733948035e-27],
86  [1000064800.0, 3.006434304359e-27, 1.259817232053e-26],
87  [1000068400.0, 1.177337218278e-27, 1.411375249637e-26],
88  [1000072000.0, -9.549597694961e-28, 1.318431214347e-26],
89  [1000075600.0, -2.924662644700e-27, 9.998085875592e-27],
90  [1000079200.0, -4.298094465681e-27, 5.272226410224e-27],
91  [1000082800.0, -4.783899643376e-27, 6.937006559007e-29],
92  ]
93  ),
94  "L1": np.array(
95  [
96  [1000000000.0, 2.354774544244e-27, 2.989220110321e-27],
97  [1000003600.0, 1.102264220569e-27, 4.794624473630e-27],
98  [1000007200.0, -3.784783044226e-28, 4.901350258365e-27],
99  [1000010800.0, -1.726096919602e-27, 3.558770900936e-27],
100  [1000014400.0, -2.638664308831e-27, 1.375241458539e-27],
101  [1000018000.0, -2.950862387455e-27, -8.626400678749e-28],
102  [1000021600.0, -2.672975553462e-27, -2.415588636184e-27],
103  [1000025200.0, -1.981652053751e-27, -2.800065246359e-27],
104  [1000028800.0, -1.165367276485e-27, -1.923079648198e-27],
105  [1000032400.0, -5.396781233382e-28, -1.064506005032e-28],
106  [1000036000.0, -3.556885019889e-28, 2.005730381569e-27],
107  [1000039600.0, -7.267224431594e-28, 3.631191457942e-27],
108  [1000043200.0, -1.593447257663e-27, 4.074861141635e-27],
109  [1000046800.0, -2.737073733735e-27, 2.933672841813e-27],
110  [1000050400.0, -3.837271308977e-27, 2.246600367826e-28],
111  [1000054000.0, -4.559321270321e-27, -3.599805614295e-27],
112  [1000057600.0, -4.646886209647e-27, -7.754852546804e-27],
113  [1000061200.0, -3.995330812549e-27, -1.131770424887e-26],
114  [1000064800.0, -2.685441479913e-27, -1.346283635637e-26],
115  [1000068400.0, -9.681500316052e-28, -1.367509409753e-26],
116  [1000072000.0, 7.960559086451e-28, -1.188425812350e-26],
117  [1000075600.0, 2.227270382987e-27, -8.484795106371e-27],
118  [1000079200.0, 3.021952545474e-27, -4.235613315750e-27],
119  [1000082800.0, 3.029100989189e-27, -6.642344227215e-29],
120  ]
121  ),
122 }
123 
124 
125 """
126 The second test output is to create a signal model for a source assuming that
127 the heterodyne and injected signal do not precisely match.
128 lalpulsar_parameter_estimation_nested has been run with the following
129 pulsar parameter "inj.par" file:
130 
131 PSRJ TEST
132 RAJ 01:23:34.5
133 DECJ -45:01:23.4
134 F0 123.456789
135 F1 -9.87654321e-12
136 PEPOCH 58000
137 H0 5.6e-26
138 COSIOTA -0.2
139 PSI 0.4
140 PHI0 2.3
141 EPHEM DE405
142 UNITS TCB
143 
144 and "heterodyne" "het.par" file:
145 
146 PSRJ TEST
147 RAJ 01:23:34.6
148 DECJ -45:01:23.5
149 F0 123.4567
150 F1 -9.876e-12
151 PEPOCH 58000
152 H0 5.6e-26
153 COSIOTA -0.2
154 PSI 0.4
155 PHI0 2.3
156 EPHEM DE405
157 UNITS TCB
158 
159 lalpulsar_parameter_estimation_nested --par-file het.par --outfile test.hdf --harmonics 2 --inject-file inj.par --inject-only --inject-output inj.txt --fake-data H1 --fake-starts 1000000000 --fake-lengths 600 --fake-dt 60
160 
161 The output of "inj.txt_H1_2.0_signal_only" is given below.
162 """
163 
164 t2output = np.array(
165  [
166  [1000000000.0, -5.941104012924e-27, 1.561111119568e-27],
167  [1000000060.0, -6.078674025313e-27, 1.108663867367e-27],
168  [1000000120.0, -6.181603080801e-27, 6.451136978413e-28],
169  [1000000180.0, -6.249044892549e-27, 1.731192817987e-28],
170  [1000000240.0, -6.280357307660e-27, -3.045989420574e-28],
171  [1000000300.0, -6.275106653741e-27, -7.852744278380e-28],
172  [1000000360.0, -6.233070880130e-27, -1.266110481373e-27],
173  [1000000420.0, -6.154241482804e-27, -1.744296420398e-27],
174  [1000000480.0, -6.038824191671e-27, -2.217023814819e-27],
175  [1000000540.0, -5.887238436697e-27, -2.681502811182e-27],
176  ]
177 )
178 
179 
180 """
181 The third test output is to create a signal model for a source assuming that
182 the heterodyne and injected signal do not precisely matches, but with emission
183 at both once and twice the rotation frequency.
184 lalpulsar_parameter_estimation_nested has been run with
185 the following pulsar parameter "inj.par" file:
186 
187 PSRJ TEST
188 RAJ 01:23:34.5
189 DECJ -45:01:23.4
190 F0 123.456789
191 F1 -9.87654321e-12
192 PEPOCH 58000
193 C22 5.6e-26
194 C21 1.4e-25
195 COSIOTA -0.2
196 PSI 0.4
197 PHI21 2.3
198 PHI22 4.5
199 EPHEM DE405
200 UNITS TCB
201 
202 and "heterodyne" "het.par" file:
203 
204 PSRJ TEST
205 RAJ 01:23:34.6
206 DECJ -45:01:23.5
207 F0 123.4567
208 F1 -9.876e-12
209 PEPOCH 58000
210 C22 5.6e-26
211 C21 1.4e-25
212 COSIOTA -0.2
213 PSI 0.4
214 PHI21 2.3
215 PHI22 4.5
216 EPHEM DE405
217 UNITS TCB
218 
219 lalpulsar_parameter_estimation_nested --par-file het.par --outfile test.hdf --harmonics 1,2 --inject-file inj.par --inject-only --inject-output inj.txt --fake-data H1 --fake-starts 1000000000,1000000000 --fake-lengths 600,600 --fake-dt 60,60
220 
221 The outputs of "inj.txt_H1_1.0_signal_only" and "inj.txt_H1_2.0_signal_only"
222 are given below.
223 """
224 
225 t3output = {
226  "1": np.array(
227  [
228  [1000000000.0, -2.281620347081e-26, -7.267059337744e-27],
229  [1000000060.0, -2.242012648444e-26, -8.024920903732e-27],
230  [1000000120.0, -2.199802615415e-26, -8.763862546799e-27],
231  [1000000180.0, -2.155075521735e-26, -9.482964910868e-27],
232  [1000000240.0, -2.107919965732e-26, -1.018134596439e-26],
233  [1000000300.0, -2.058427708343e-26, -1.085816221335e-26],
234  [1000000360.0, -2.006693507621e-26, -1.151260986472e-26],
235  [1000000420.0, -1.952814948068e-26, -1.214392590347e-26],
236  [1000000480.0, -1.896892272393e-26, -1.275138910868e-26],
237  [1000000540.0, -1.839028196029e-26, -1.333432099023e-26],
238  ]
239  ),
240  "2": np.array(
241  [
242  [1000000000.0, 1.151114436476e-26, -4.292865557392e-27],
243  [1000000060.0, 1.187524854552e-26, -3.419959925105e-27],
244  [1000000120.0, 1.217263381782e-26, -2.518042744682e-27],
245  [1000000180.0, 1.240108521541e-26, -1.592235817764e-27],
246  [1000000240.0, 1.255878166730e-26, -6.478246234004e-28],
247  [1000000300.0, 1.264430777434e-26, 3.097719790379e-28],
248  [1000000360.0, 1.265666324682e-26, 1.275032881006e-27],
249  [1000000420.0, 1.259526996162e-26, 2.242366499355e-27],
250  [1000000480.0, 1.245997657263e-26, 3.206142957363e-27],
251  [1000000540.0, 1.225106070777e-26, 4.160726677163e-27],
252  ]
253  ),
254 }
255 
256 
257 """
258 The fourth test output is to create a signal model for a source assuming that
259 the heterodyne and injected signal do not precisely matches, adding in some
260 binary system parameters.
261 lalpulsar_parameter_estimation_nested has been run with
262 the following pulsar parameter "inj.par" file:
263 
264 PSRJ TEST
265 RAJ 01:23:34.5
266 DECJ -45:01:23.4
267 F0 123.456789
268 F1 -9.87654321e-12
269 PEPOCH 58000
270 H0 5.6e-26
271 COSIOTA -0.2
272 PSI 0.4
273 PHI0 2.3
274 BINARY BT
275 T0 58121.3
276 ECC 0.0001
277 OM 1.2
278 A1 8.9
279 PB 0.54
280 EPHEM DE405
281 UNITS TCB
282 
283 and "heterodyne" "het.par" file:
284 
285 PSRJ TEST
286 RAJ 01:23:34.6
287 DECJ -45:01:23.5
288 F0 123.4567
289 F1 -9.876e-12
290 PEPOCH 58000
291 H0 5.6e-26
292 COSIOTA -0.2
293 PSI 0.4
294 PHI0 2.3
295 BINARY BT
296 T0 58121.3
297 ECC 0.0001
298 OM 2.2
299 A1 8.9
300 PB 0.54
301 EPHEM DE405
302 UNITS TCB
303 
304 lalpulsar_parameter_estimation_nested --par-file het.par --outfile test.hdf --harmonics 2 --inject-file inj.par --inject-only --inject-output inj.txt --fake-data H1 --fake-starts 1000000000 --fake-lengths 600 --fake-dt 60
305 
306 The output of "inj.txt_H1_2.0_signal_only" is given below.
307 """
308 
309 t4output = np.array(
310  [
311  [1000000000.0, -2.005382682171e-27, -5.806222964895e-27],
312  [1000000060.0, 6.157890500690e-27, 5.097038927998e-28],
313  [1000000120.0, -2.949373857456e-27, 5.470793560413e-27],
314  [1000000180.0, -3.871712822495e-27, -4.908194390501e-27],
315  [1000000240.0, 6.071607014385e-27, -1.634397959568e-27],
316  [1000000300.0, -8.720674118593e-28, 6.263634557667e-27],
317  [1000000360.0, -5.468817484054e-27, -3.247498063718e-27],
318  [1000000420.0, 5.125788386822e-27, -3.826689384362e-27],
319  [1000000480.0, 1.601886955222e-27, 6.230292962298e-27],
320  [1000000540.0, -6.411302862101e-27, -8.632666722102e-28],
321  ]
322 )
323 
324 
325 """
326 The fifth test output is to create a signal model for a source assuming that
327 the heterodyne precisely matches the signal, i.e., it only varies due to the
328 antenna pattern, but including vector and scalar polarisation modes.
329 lalpulsar_parameter_estimation_nested has been run with
330 the following pulsar parameter "inj.par" file:
331 
332 PSRJ TEST
333 RAJ 01:23:34.5
334 DECJ -45:01:23.4
335 F0 123.456789
336 F1 -9.87654321e-12
337 PEPOCH 58000
338 HPLUS 5.6e-26
339 HCROSS 1.3e-26
340 HVECTORX 1.4e-26
341 HVECTORY 2.3e-26
342 HSCALARB 4.5e-26
343 HSCALARL 3.1e-26
344 PHI0TENSOR 0.4
345 PSITENSOR 1.2
346 PHI0SCALAR 3.1
347 PSISCALAR 0.2
348 PHI0VECTOR 4.5
349 PSIVECTOR 2.4
350 EPHEM DE405
351 UNITS TCB
352 
353 lalpulsar_parameter_estimation_nested --par-file inj.par --outfile test.hdf --harmonics 2 --inject-file inj.par --inject-only --inject-output inj.txt --fake-data H1 --fake-starts 1000000000 --fake-lengths 86400 --fake-dt 3600 --nonGR --inject-nonGR
354 
355 The output of "inj.txt_H1_2.0_signal_only" is given below.
356 """
357 
358 t5output = np.array(
359  [
360  [1000000000.0, 2.413832756525e-26, 1.185008925479e-26],
361  [1000003600.0, 2.261015869048e-26, 1.353795950522e-26],
362  [1000007200.0, 1.693819370778e-26, 1.205591484118e-26],
363  [1000010800.0, 8.826543665516e-27, 7.844031741297e-27],
364  [1000014400.0, 4.864166183254e-28, 2.043552377844e-27],
365  [1000018000.0, -5.961350254474e-27, -3.810915453418e-27],
366  [1000021600.0, -9.049037164685e-27, -8.203008361034e-27],
367  [1000025200.0, -8.338836580606e-27, -1.003879189105e-26],
368  [1000028800.0, -4.514037100223e-27, -8.935501871933e-27],
369  [1000032400.0, 8.389816591977e-28, -5.316862125252e-27],
370  [1000036000.0, 5.696076643445e-27, -2.902907800008e-28],
371  [1000039600.0, 8.180535124244e-27, 4.660428187861e-27],
372  [1000043200.0, 7.107590901248e-27, 8.083846851900e-27],
373  [1000046800.0, 2.338787027244e-27, 8.959804096635e-27],
374  [1000050400.0, -5.151227246717e-27, 6.981092071373e-27],
375  [1000054000.0, -1.351620286624e-26, 2.640990277497e-27],
376  [1000057600.0, -2.052378588485e-26, -2.896887794574e-27],
377  [1000061200.0, -2.415548913847e-26, -8.111602201827e-27],
378  [1000064800.0, -2.315941209325e-26, -1.153669175848e-26],
379  [1000068400.0, -1.740609207988e-26, -1.215933456752e-26],
380  [1000072000.0, -7.950868000824e-27, -9.699490143059e-27],
381  [1000075600.0, 3.216714991178e-27, -4.692952735182e-27],
382  [1000079200.0, 1.367058676109e-26, 1.644309606653e-27],
383  [1000082800.0, 2.116179045441e-26, 7.734541471774e-27],
384  ]
385 )
386 
387 
388 """
389 The sixth test output is to create a signal model phase evolution for a source
390 assuming that heterodyne and an updated signal do not precisely match, adding
391 in some binary system parameters, and glitch parameters for the updated signal.
392 The first glitch is at GPS 1000000600 (in TDB, which puts it into the data set
393 time, or MJD 55818.08161090822), second glitch is at GPS 1000000700 (or MJD
394 55818.08276831563)
395 lalpulsar_heterodyne has been run with the following pulsar parameter
396 "inj.par" file:
397 
398 PSRJ TEST
399 RAJ 01:23:34.5
400 DECJ -45:01:23.4
401 F0 123.456789
402 F1 -9.87654321e-12
403 PEPOCH 58000
404 BINARY BT
405 T0 58121.3
406 ECC 0.0001
407 OM 1.2
408 A1 8.9
409 PB 0.54
410 EPHEM DE405
411 UNITS TCB
412 GLPH_1 0.3
413 GLEP_1 55818.08161090822
414 GLF0_1 5.4e-6
415 GLF1_1 -3.2e-13
416 GLF0D_1 1.2e-5
417 GLTD_1 0.31
418 GLPH_2 0.7
419 GLEP_2 55818.08276831563
420 GLF0_2 3.4e-7
421 GLF1_2 -1.2e-14
422 GLF0D_2 -0.4e-6
423 GLTD_2 0.45
424 
425 and "heterodyne" "het.par" file:
426 
427 PSRJ TEST
428 RAJ 01:23:34.6
429 DECJ -45:01:23.5
430 F0 123.4567
431 F1 -9.876e-12
432 PEPOCH 58000
433 EPHEM DE405
434 UNITS TCB
435 
436 lalpulsar_heterodyne -i H1 --pulsar J0000+0000 -z 2 -f het.par -g inj.par -s 1/60 -r 1/60 -P -o testing.txt -l segments.txt -d inj.txt -L -e earth00-40-DE405.dat.gz -S sun00-40-DE405.dat.gz -t te405_2000-2040.dat.gz
437 
438 where inj.txt contains 10 rows with times between 1000000000 and 1000000540,
439 and segments.txt contains one row with:
440 1000000000 1000000600
441 
442 The output of "phase.txt" is given below.
443 """
444 
445 t6output = np.array(
446  [
447  0.962333260,
448  6.103230186,
449  4.076170449,
450  1.165825399,
451  3.656722577,
452  2.766150856,
453  6.024421923,
454  5.892164017,
455  4.884626261,
456  3.003294698,
457  ]
458 )
459 
460 
461 """
462 The seventh test output is to create a signal model phase evolution for a source
463 assuming that heterodyne and an updated signal do not precisely match, adding
464 in some binary system parameters, glitch parameters and timing noise (FITWAVES)
465 parameters for the updated signal.
466 lalpulsar_heterodyne has been run with the following pulsar parameter
467 "inj.par" file:
468 
469 PSRJ TEST
470 RAJ 04:23:34.5
471 DECJ -05:01:23.4
472 F0 153.456789
473 F1 -2.87654321e-11
474 PEPOCH 55810
475 BINARY BT
476 T0 58121.3
477 ECC 0.0002
478 OM 7.2
479 A1 14.9
480 PB 1.03
481 EPHEM DE405
482 UNITS TCB
483 GLPH_1 0.3
484 GLEP_1 55818.08161090822
485 GLF0_1 7.4e-6
486 GLF1_1 -3.2e-12
487 GLF0D_1 1.2e-5
488 GLTD_1 0.41
489 GLPH_2 0.91
490 GLEP_2 55818.08276831563
491 GLF0_2 3.4e-7
492 GLF1_2 -1.2e-14
493 GLF0D_2 -0.4e-6
494 GLTD_2 1.45
495 WAVEEPOCH 55818.0
496 WAVE_OM 0.005
497 WAVE1 0.098 0.056
498 WAVE2 0.078 -0.071
499 WAVE3 -0.03 -0.12
500 
501 
502 and "heterodyne" "het.par" file:
503 
504 PSRJ TEST
505 RAJ 04:23:34.6
506 DECJ -05:01:23.5
507 F0 153.4567
508 F1 -2.876e-11
509 PEPOCH 55810
510 EPHEM DE405
511 UNITS TCB
512 
513 lalpulsar_heterodyne -i H1 --pulsar J0000+0000 -z 2 -f het.par -g inj.par -s 1/60 -r 1/60 -P -o testing.txt -l segments.txt -d inj.txt -L -e earth00-40-DE405.dat.gz -S sun00-40-DE405.dat.gz -t te405_2000-2040.dat.gz
514 
515 where inj.txt contains 10 rows with times between 1000000000 and 1000000540,
516 and segments.txt contains one row with:
517 1000000000 1000000600
518 
519 The output of "phase.txt" is given below.
520 """
521 
522 t7output = np.array(
523  [
524  2.420723736,
525  4.682767534,
526  0.313020652,
527  1.879459964,
528  3.100512978,
529  3.977798875,
530  4.512941892,
531  4.707573756,
532  2.060503777,
533  0.462653250,
534  ]
535 )
536 
537 
538 # set ephemeris files
539 earthephem = os.path.join(os.environ["LAL_TEST_PKGDATADIR"], "earth00-40-DE405.dat.gz")
540 sunephem = os.path.join(os.environ["LAL_TEST_PKGDATADIR"], "sun00-40-DE405.dat.gz")
541 timefile = os.path.join(os.environ["LAL_TEST_PKGDATADIR"], "te405_2000-2040.dat.gz")
542 
543 # get ephemeris files
544 edat = lalpulsar.InitBarycenter(earthephem, sunephem)
545 tdat = lalpulsar.InitTimeCorrections(timefile)
546 
547 
548 @pytest.mark.parametrize("det", sorted(t1output.keys()))
549 def test_one(det):
550  par = PulsarParametersPy()
551  par["F"] = [123.456789, -9.87654321e-12] # set frequency
552  par["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
553  par["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
554  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
555  par["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
556  par["H0"] = 5.6e-26
557  par["COSIOTA"] = -0.2
558  par["PSI"] = 0.4
559  par["PHI0"] = 2.3
560 
561  freqfactor = 2.0 # set frequency factor
562 
563  # convert into GPS times
564  gpstimes = lalpulsar.CreateTimestampVector(len(t1output[det]))
565  for i, time in enumerate(t1output[det][:, 0]):
566  gpstimes.data[i] = lal.LIGOTimeGPS(time)
567 
568  detector = lalpulsar.GetSiteInfo(det)
569 
570  # set the response function look-up table
571  dt = t1output[det][1, 0] - t1output[det][0, 0] # time step
572  resp = lalpulsar.DetResponseLookupTable(
573  t1output[det][0, 0], detector, par["RAJ"], par["DECJ"], 2880, dt
574  )
575 
576  # get the heterodyned file SSB delay
577  hetSSBdelay = lalpulsar.HeterodynedPulsarGetSSBDelay(
578  par.PulsarParameters(),
579  gpstimes,
580  detector,
581  edat,
582  tdat,
583  lalpulsar.TIMECORRECTION_TCB,
584  )
585 
586  fullsignal = lalpulsar.HeterodynedPulsarGetModel(
587  par.PulsarParameters(),
588  par.PulsarParameters(),
589  freqfactor,
590  1,
591  0,
592  0,
593  gpstimes,
594  hetSSBdelay,
595  0,
596  None,
597  0,
598  None,
599  0,
600  None,
601  0,
602  resp,
603  edat,
604  tdat,
605  lalpulsar.TIMECORRECTION_TCB,
606  )
607 
608  # check output matches that from lalpulsar_parameter_estimation_nested
609  assert_allclose(fullsignal.data.data.real, t1output[det][:, 1])
610  assert_allclose(fullsignal.data.data.imag, t1output[det][:, 2])
611 
612 
613 def test_two():
614  parhet = PulsarParametersPy()
615  parhet["F"] = [123.4567, -9.876e-12] # set frequency
616  parhet["RAJ"] = lal.TranslateHMStoRAD("01:23:34.6") # set right ascension
617  parhet["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.5") # set declination
618  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
619  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
620  parhet["H0"] = 5.6e-26
621  parhet["COSIOTA"] = -0.2
622  parhet["PSI"] = 0.4
623  parhet["PHI0"] = 2.3
624 
625  parinj = PulsarParametersPy()
626  parinj["F"] = [123.456789, -9.87654321e-12] # set frequency
627  parinj["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
628  parinj["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
629  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
630  parinj["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
631  parinj["H0"] = 5.6e-26
632  parinj["COSIOTA"] = -0.2
633  parinj["PSI"] = 0.4
634  parinj["PHI0"] = 2.3
635 
636  freqfactor = 2.0 # set frequency factor
637  det = "H1" # the detector
638 
639  # convert into GPS times
640  gpstimes = lalpulsar.CreateTimestampVector(len(t2output))
641  for i, time in enumerate(t2output[:, 0]):
642  gpstimes.data[i] = lal.LIGOTimeGPS(time)
643 
644  detector = lalpulsar.GetSiteInfo(det)
645 
646  # set the response function look-up table
647  dt = t2output[1, 0] - t2output[0, 0] # time step
648  resp = lalpulsar.DetResponseLookupTable(
649  t2output[0, 0], detector, parhet["RAJ"], parhet["DECJ"], 2880, dt
650  )
651 
652  # get the heterodyned file SSB delay
653  hetSSBdelay = lalpulsar.HeterodynedPulsarGetSSBDelay(
654  parhet.PulsarParameters(),
655  gpstimes,
656  detector,
657  edat,
658  tdat,
659  lalpulsar.TIMECORRECTION_TCB,
660  )
661 
662  fullsignal = lalpulsar.HeterodynedPulsarGetModel(
663  parinj.PulsarParameters(),
664  parhet.PulsarParameters(),
665  freqfactor,
666  1,
667  0,
668  0,
669  gpstimes,
670  hetSSBdelay,
671  1,
672  None,
673  0,
674  None,
675  0,
676  None,
677  0,
678  resp,
679  edat,
680  tdat,
681  lalpulsar.TIMECORRECTION_TCB,
682  )
683 
684  # check output matches that from lalpulsar_parameter_estimation_nested
685  assert_allclose(fullsignal.data.data.real, t2output[:, 1])
686  assert_allclose(fullsignal.data.data.imag, t2output[:, 2])
687 
688 
689 @pytest.mark.parametrize("harmonic", sorted(t3output.keys()))
690 def test_three(harmonic):
691  parhet = PulsarParametersPy()
692  parhet["F"] = [123.4567, -9.876e-12] # set frequency
693  parhet["RAJ"] = lal.TranslateHMStoRAD("01:23:34.6") # set right ascension
694  parhet["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.5") # set declination
695  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
696  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
697  parhet["C22"] = 5.6e-26
698  parhet["C21"] = 1.4e-25
699  parhet["COSIOTA"] = -0.2
700  parhet["PSI"] = 0.4
701  parhet["PHI21"] = 2.3
702  parhet["PHI22"] = 4.5
703 
704  parinj = PulsarParametersPy()
705  parinj["F"] = [123.456789, -9.87654321e-12] # set frequency
706  parinj["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
707  parinj["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
708  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
709  parinj["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
710  parinj["C22"] = 5.6e-26
711  parinj["C21"] = 1.4e-25
712  parinj["COSIOTA"] = -0.2
713  parinj["PSI"] = 0.4
714  parinj["PHI21"] = 2.3
715  parinj["PHI22"] = 4.5
716 
717  det = "H1" # the detector
718  detector = lalpulsar.GetSiteInfo(det)
719 
720  freqfactor = float(harmonic) # set frequency factor
721 
722  # convert into GPS times
723  gpstimes = lalpulsar.CreateTimestampVector(len(t3output[harmonic]))
724  for i, time in enumerate(t3output[harmonic][:, 0]):
725  gpstimes.data[i] = lal.LIGOTimeGPS(time)
726 
727  # set the response function look-up table
728  dt = t3output[harmonic][1, 0] - t3output[harmonic][0, 0] # time step
729  resp = lalpulsar.DetResponseLookupTable(
730  t3output[harmonic][0, 0], detector, parhet["RAJ"], parhet["DECJ"], 2880, dt
731  )
732 
733  # get the heterodyned file SSB delay
734  hetSSBdelay = lalpulsar.HeterodynedPulsarGetSSBDelay(
735  parhet.PulsarParameters(),
736  gpstimes,
737  detector,
738  edat,
739  tdat,
740  lalpulsar.TIMECORRECTION_TCB,
741  )
742 
743  fullsignal = lalpulsar.HeterodynedPulsarGetModel(
744  parinj.PulsarParameters(),
745  parhet.PulsarParameters(),
746  freqfactor,
747  1,
748  0,
749  0,
750  gpstimes,
751  hetSSBdelay,
752  1,
753  None,
754  0,
755  None,
756  0,
757  None,
758  0,
759  resp,
760  edat,
761  tdat,
762  lalpulsar.TIMECORRECTION_TCB,
763  )
764 
765  # check output matches that from lalpulsar_parameter_estimation_nested
766  assert_allclose(fullsignal.data.data.real, t3output[harmonic][:, 1])
767  assert_allclose(fullsignal.data.data.imag, t3output[harmonic][:, 2])
768 
769 
770 def test_four():
771  parhet = PulsarParametersPy()
772  parhet["F"] = [123.4567, -9.876e-12] # set frequency
773  parhet["RAJ"] = lal.TranslateHMStoRAD("01:23:34.6") # set right ascension
774  parhet["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.5") # set declination
775  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
776  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
777  parhet["H0"] = 5.6e-26
778  parhet["COSIOTA"] = -0.2
779  parhet["PSI"] = 0.4
780  parhet["PHI0"] = 2.3
781  parhet["BINARY"] = "BT"
782  T0 = lal.TranslateStringMJDTTtoGPS("58121.3")
783  parhet["T0"] = T0.gpsSeconds + 1e-9 * T0.gpsNanoSeconds
784  parhet["OM"] = np.deg2rad(2.2)
785  parhet["A1"] = 8.9
786  parhet["PB"] = 0.54 * 86400.0
787  parhet["ECC"] = 0.0001
788 
789  parinj = PulsarParametersPy()
790  parinj["F"] = [123.456789, -9.87654321e-12] # set frequency
791  parinj["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
792  parinj["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
793  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
794  parinj["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
795  parinj["H0"] = 5.6e-26
796  parinj["COSIOTA"] = -0.2
797  parinj["PSI"] = 0.4
798  parinj["PHI0"] = 2.3
799  parinj["BINARY"] = "BT"
800  T0 = lal.TranslateStringMJDTTtoGPS("58121.3")
801  parinj["T0"] = T0.gpsSeconds + 1e-9 * T0.gpsNanoSeconds
802  parinj["OM"] = np.deg2rad(1.2)
803  parinj["A1"] = 8.9
804  parinj["PB"] = 0.54 * 86400.0
805  parinj["ECC"] = 0.0001
806 
807  freqfactor = 2.0 # set frequency factor
808  det = "H1" # the detector
809 
810  # convert into GPS times
811  gpstimes = lalpulsar.CreateTimestampVector(len(t4output))
812  for i, time in enumerate(t4output[:, 0]):
813  gpstimes.data[i] = lal.LIGOTimeGPS(time)
814 
815  detector = lalpulsar.GetSiteInfo(det)
816 
817  # set the response function look-up table
818  dt = t4output[1, 0] - t4output[0, 0] # time step
819  resp = lalpulsar.DetResponseLookupTable(
820  t4output[0, 0], detector, parhet["RAJ"], parhet["DECJ"], 2880, dt
821  )
822 
823  # get the heterodyned file SSB delay
824  hetSSBdelay = lalpulsar.HeterodynedPulsarGetSSBDelay(
825  parhet.PulsarParameters(),
826  gpstimes,
827  detector,
828  edat,
829  tdat,
830  lalpulsar.TIMECORRECTION_TCB,
831  )
832 
833  # get the heterodyned file BSB delay
834  hetBSBdelay = lalpulsar.HeterodynedPulsarGetBSBDelay(
835  parhet.PulsarParameters(), gpstimes, hetSSBdelay, edat
836  )
837 
838  fullsignal = lalpulsar.HeterodynedPulsarGetModel(
839  parinj.PulsarParameters(),
840  parhet.PulsarParameters(),
841  freqfactor,
842  1, # phase is varying between par files
843  0, # not using ROQ
844  0, # not using non-tensorial modes
845  gpstimes,
846  hetSSBdelay,
847  1, # the SSB delay should be updated compared to hetSSBdelay
848  hetBSBdelay,
849  1, # the BSB delay should be updated compared to hetBSBdelay
850  None,
851  0,
852  None,
853  0,
854  resp,
855  edat,
856  tdat,
857  lalpulsar.TIMECORRECTION_TCB,
858  )
859 
860  # check output matches that from lalpulsar_parameter_estimation_nested
861  assert_allclose(fullsignal.data.data.real, t4output[:, 1])
862  assert_allclose(fullsignal.data.data.imag, t4output[:, 2])
863 
864 
865 def test_five():
866  par = PulsarParametersPy()
867  par["F"] = [123.456789, -9.87654321e-12] # set frequency
868  par["DELTAF"] = [0.0, 0.0] # frequency difference
869  par["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
870  par["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
871  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
872  par["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
873  par["HPLUS"] = 5.6e-26
874  par["HCROSS"] = 1.3e-26
875  par["HVECTORX"] = 1.4e-26
876  par["HVECTORY"] = 2.3e-26
877  par["HSCALARB"] = 4.5e-26
878  par["HSCALARL"] = 3.1e-26
879  par["PHI0TENSOR"] = 0.4
880  par["PSITENSOR"] = 1.2
881  par["PHI0SCALAR"] = 3.1
882  par["PSISCALAR"] = 0.2
883  par["PHI0VECTOR"] = 4.5
884  par["PSIVECTOR"] = 2.4
885  par["PHI0"] = 2.3
886 
887  freqfactor = 2.0 # set frequency factor
888 
889  det = "H1" # detector
890  detector = lalpulsar.GetSiteInfo(det)
891 
892  gpstimes = lalpulsar.CreateTimestampVector(len(t5output))
893  for i, time in enumerate(t5output[:, 0]):
894  gpstimes.data[i] = lal.LIGOTimeGPS(time)
895 
896  # set the response function look-up table
897  dt = t5output[1, 0] - t5output[0, 0] # time step
898  resp = lalpulsar.DetResponseLookupTable(
899  t5output[0, 0], detector, par["RAJ"], par["DECJ"], 2880, dt
900  )
901 
902  # get the heterodyned file SSB delay
903  hetSSBdelay = lalpulsar.HeterodynedPulsarGetSSBDelay(
904  par.PulsarParameters(),
905  gpstimes,
906  detector,
907  edat,
908  tdat,
909  lalpulsar.TIMECORRECTION_TCB,
910  )
911 
912  fullsignal = lalpulsar.HeterodynedPulsarGetModel(
913  par.PulsarParameters(),
914  par.PulsarParameters(),
915  freqfactor,
916  1,
917  0,
918  1, # use non-GR modes
919  gpstimes,
920  hetSSBdelay,
921  0,
922  None,
923  0,
924  None,
925  0,
926  None,
927  0,
928  resp,
929  edat,
930  tdat,
931  lalpulsar.TIMECORRECTION_TCB,
932  )
933 
934  # check output matches that from lalpulsar_parameter_estimation_nested
935  assert_allclose(fullsignal.data.data.real, t5output[:, 1])
936  assert_allclose(fullsignal.data.data.imag, t5output[:, 2])
937 
938 
939 def test_six():
940  parhet = PulsarParametersPy()
941  parhet["F"] = [123.4567, -9.876e-12] # set frequency
942  parhet["RAJ"] = lal.TranslateHMStoRAD("01:23:34.6") # set right ascension
943  parhet["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.5") # set declination
944  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
945  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
946 
947  parinj = PulsarParametersPy()
948  parinj["F"] = [123.456789, -9.87654321e-12] # set frequency
949  parinj["RAJ"] = lal.TranslateHMStoRAD("01:23:34.5") # set right ascension
950  parinj["DECJ"] = lal.TranslateDMStoRAD("-45:01:23.4") # set declination
951  pepoch = lal.TranslateStringMJDTTtoGPS("58000")
952  parinj["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
953  parinj["BINARY"] = "BT"
954  T0 = lal.TranslateStringMJDTTtoGPS("58121.3")
955  parinj["T0"] = T0.gpsSeconds + 1e-9 * T0.gpsNanoSeconds
956  parinj["OM"] = np.deg2rad(1.2)
957  parinj["A1"] = 8.9
958  parinj["PB"] = 0.54 * 86400.0
959  parinj["ECC"] = 0.0001
960  parinj["GLF0"] = [5.4e-6, 3.4e-7]
961  parinj["GLF1"] = [-3.2e-13, -1.2e-14]
962  parinj["GLF0D"] = [1.2e-5, -0.4e-6]
963  parinj["GLTD"] = [0.31 * 86400, 0.45 * 86400]
964  parinj["GLPH"] = [0.3, 0.7]
965  glph1 = lal.TranslateStringMJDTTtoGPS("55818.08161090822")
966  glph2 = lal.TranslateStringMJDTTtoGPS("55818.08276831563")
967  parinj["GLEP"] = [
968  glph1.gpsSeconds + 1e-9 * glph1.gpsNanoSeconds,
969  glph2.gpsSeconds + 1e-9 * glph2.gpsNanoSeconds,
970  ]
971 
972  freqfactor = 2.0 # set frequency factor
973  det = "H1" # the detector
974 
975  # convert into GPS times
976  gpstimes = lalpulsar.CreateTimestampVector(len(t6output))
977  for i, time in enumerate(np.linspace(1000000000.0, 1000000540.0, 10)):
978  gpstimes.data[i] = lal.LIGOTimeGPS(time)
979 
980  detector = lalpulsar.GetSiteInfo(det)
981 
982  # replicate coarse heterodyne in which no SSB/BSB delay is applied
983  hetSSBdelay = lal.CreateREAL8Vector(len(t6output))
984  hetBSBdelay = lal.CreateREAL8Vector(len(t6output))
985  for i in range(len(t6output)):
986  hetSSBdelay.data[i] = 0.0
987  hetBSBdelay.data[i] = 0.0
988 
989  # get the heterodyne glitch phase (which should be zero)
990  glphase = lalpulsar.HeterodynedPulsarGetGlitchPhase(
991  parhet.PulsarParameters(), gpstimes, hetSSBdelay, hetBSBdelay
992  )
993 
994  fullphase = lalpulsar.HeterodynedPulsarPhaseDifference(
995  parinj.PulsarParameters(),
996  parhet.PulsarParameters(),
997  gpstimes,
998  freqfactor,
999  hetSSBdelay,
1000  1, # the SSB delay should be updated compared to hetSSBdelay
1001  hetBSBdelay,
1002  1, # the BSB delay should be updated compared to hetBSBdelay
1003  glphase,
1004  1, # the glitch phase should be updated compared to glphase
1005  None,
1006  0,
1007  detector,
1008  edat,
1009  tdat,
1010  lalpulsar.TIMECORRECTION_TCB,
1011  )
1012 
1013  # check output matches that from lalpulsar_heterodyne
1014  assert_allclose(2.0 * np.pi * np.fmod(fullphase.data, 1.0), t6output, rtol=1e-4)
1015 
1016 
1017 def test_seven():
1018  parhet = PulsarParametersPy()
1019  parhet["F"] = [153.4567, -2.876e-11] # set frequency
1020  parhet["RAJ"] = lal.TranslateHMStoRAD("04:23:34.6") # set right ascension
1021  parhet["DECJ"] = lal.TranslateDMStoRAD("-05:01:23.5") # set declination
1022  pepoch = lal.TranslateStringMJDTTtoGPS("55810")
1023  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
1024 
1025  parinj = PulsarParametersPy()
1026  parinj["F"] = [153.456789, -2.87654321e-11] # set frequency
1027  parinj["RAJ"] = lal.TranslateHMStoRAD("04:23:34.5") # set right ascension
1028  parinj["DECJ"] = lal.TranslateDMStoRAD("-05:01:23.4") # set declination
1029  pepoch = lal.TranslateStringMJDTTtoGPS("55810")
1030  parinj["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
1031  parinj["BINARY"] = "BT"
1032  T0 = lal.TranslateStringMJDTTtoGPS("58121.3")
1033  parinj["T0"] = T0.gpsSeconds + 1e-9 * T0.gpsNanoSeconds
1034  parinj["OM"] = np.deg2rad(7.2)
1035  parinj["A1"] = 14.9
1036  parinj["PB"] = 1.03 * 86400.0
1037  parinj["ECC"] = 0.0002
1038  parinj["GLF0"] = [7.4e-6, 3.4e-7]
1039  parinj["GLF1"] = [-3.2e-12, -1.2e-14]
1040  parinj["GLF0D"] = [1.2e-5, -0.4e-6]
1041  parinj["GLTD"] = [0.41 * 86400, 1.45 * 86400]
1042  parinj["GLPH"] = [0.3, 0.91]
1043  glep1 = lal.TranslateStringMJDTTtoGPS("55818.08161090822")
1044  glep2 = lal.TranslateStringMJDTTtoGPS("55818.08276831563")
1045  parinj["GLEP"] = [
1046  glep1.gpsSeconds + 1e-9 * glep1.gpsNanoSeconds,
1047  glep2.gpsSeconds + 1e-9 * glep2.gpsNanoSeconds,
1048  ]
1049  waveep = lal.TranslateStringMJDTTtoGPS("55818.0")
1050  parinj["WAVEEPOCH"] = waveep.gpsSeconds + 1e-9 * waveep.gpsNanoSeconds
1051  parinj["WAVE_OM"] = 0.005
1052  parinj["WAVESIN"] = [0.098, 0.078, -0.03]
1053  parinj["WAVECOS"] = [0.056, -0.071, -0.12]
1054 
1055  freqfactor = 2.0 # set frequency factor
1056  det = "H1" # the detector
1057 
1058  # convert into GPS times
1059  gpstimes = lalpulsar.CreateTimestampVector(len(t7output))
1060  for i, time in enumerate(np.linspace(1000000000.0, 1000000540.0, 10)):
1061  gpstimes.data[i] = lal.LIGOTimeGPS(time)
1062 
1063  detector = lalpulsar.GetSiteInfo(det)
1064 
1065  # replicate coarse heterodyne in which no SSB/BSB delay is applied
1066  hetSSBdelay = lal.CreateREAL8Vector(len(t6output))
1067  hetBSBdelay = lal.CreateREAL8Vector(len(t6output))
1068  for i in range(len(t6output)):
1069  hetSSBdelay.data[i] = 0.0
1070  hetBSBdelay.data[i] = 0.0
1071 
1072  # get the heterodyne glitch phase (which should be zero)
1073  glphase = lalpulsar.HeterodynedPulsarGetGlitchPhase(
1074  parhet.PulsarParameters(), gpstimes, hetSSBdelay, hetBSBdelay
1075  )
1076 
1077  assert_equal(glphase.data, np.zeros(len(t7output)))
1078 
1079  # get the FITWAVES phase (which should be zero)
1080  fwphase = lalpulsar.HeterodynedPulsarGetFITWAVESPhase(
1081  parhet.PulsarParameters(), gpstimes, hetSSBdelay, parhet["F0"]
1082  )
1083 
1084  assert_equal(fwphase.data, np.zeros(len(t7output)))
1085 
1086  fullphase = lalpulsar.HeterodynedPulsarPhaseDifference(
1087  parinj.PulsarParameters(),
1088  parhet.PulsarParameters(),
1089  gpstimes,
1090  freqfactor,
1091  hetSSBdelay,
1092  1, # the SSB delay should be updated compared to hetSSBdelay
1093  hetBSBdelay,
1094  1, # the BSB delay should be updated compared to hetBSBdelay
1095  glphase,
1096  1, # the glitch phase should be updated compared to glphase
1097  fwphase,
1098  1, # the FITWAVES phase should be updated compare to fwphase
1099  detector,
1100  edat,
1101  tdat,
1102  lalpulsar.TIMECORRECTION_TCB,
1103  )
1104 
1105  # check output matches that from lalpulsar_heterodyne
1106  assert_allclose(2.0 * np.pi * np.fmod(fullphase.data, 1.0), t7output, rtol=1e-3)
1107 
1108 
1109 @pytest.mark.parametrize("det", ["H1", "L1", "V1"])
1110 def test_transient(det):
1111  """
1112  Test the transient signal models with a rectangular and exponential decay.
1113  """
1114 
1115  parhet = PulsarParametersPy()
1116  parhet["F"] = [153.4567, -2.876e-11] # set frequency
1117  parhet["RAJ"] = lal.TranslateHMStoRAD("04:23:34.6") # set right ascension
1118  parhet["DECJ"] = lal.TranslateDMStoRAD("-05:01:23.5") # set declination
1119  pepoch = lal.TranslateStringMJDTTtoGPS("55810")
1120  parhet["PEPOCH"] = pepoch.gpsSeconds + 1e-9 * pepoch.gpsNanoSeconds
1121 
1122  parhet["H0"] = 1e-24
1123  parhet["COSIOTA"] = 0.4
1124  parhet["PSI"] = 1.3
1125  parhet["PHI0"] = 2.9
1126 
1127  times = np.arange(1000000000, 1000000000 + 8 * 86400, 600)
1128 
1129  parhet["TRANSIENTSTARTTIME"] = times[0] + 86400.0
1130  parhet["TRANSIENTTAU"] = 86400.0 * 2
1131 
1132  # model with no transient window
1133  hetfull = HeterodynedCWSimulator(
1134  parhet,
1135  det,
1136  times=times,
1137  earth_ephem=earthephem,
1138  sun_ephem=sunephem,
1139  time_corr=timefile,
1140  )
1141  fullmodel = hetfull.model(usephase=True)
1142 
1143  # model with rectangular window
1144  parhet["TRANSIENTWINDOWTYPE"] = "RECT"
1145  hetrect = HeterodynedCWSimulator(
1146  parhet,
1147  det,
1148  times=times,
1149  earth_ephem=earthephem,
1150  sun_ephem=sunephem,
1151  time_corr=timefile,
1152  )
1153  rectmodel = hetrect.model()
1154 
1155  # model with exponential decay window
1156  parhet["TRANSIENTWINDOWTYPE"] = "EXP"
1157  hetexp = HeterodynedCWSimulator(
1158  parhet,
1159  det,
1160  times=times,
1161  earth_ephem=earthephem,
1162  sun_ephem=sunephem,
1163  time_corr=timefile,
1164  )
1165  expmodel = hetexp.model()
1166 
1167  # check rectangular window
1168  wintimes = (times >= parhet["TRANSIENTSTARTTIME"]) & (
1169  times <= parhet["TRANSIENTSTARTTIME"] + parhet["TRANSIENTTAU"]
1170  )
1171  assert_equal(fullmodel[wintimes], rectmodel[wintimes])
1172  assert_equal(rectmodel[~wintimes], np.zeros_like(rectmodel[~wintimes]))
1173 
1174  # check exponential window
1175  for i, time in enumerate(
1176  np.arange(parhet["TRANSIENTSTARTTIME"], times[-1], parhet["TRANSIENTTAU"])
1177  ):
1178  idx = int((time - times[0]) / 600)
1179  assert_allclose([expmodel[idx]], [fullmodel[idx] / (np.e**i)], rtol=1e-3)
1180 
1181 
1182 if __name__ == "__main__":
1183  args = sys.argv[1:] or ["-v", "-rs", "--junit-xml=junit-heterodyned-model.xml"]
1184  sys.exit(pytest.main(args=[__file__] + args))
lalpulsar.PulsarParametersWrapper.PulsarParametersPy
A class to wrap the SWIG-wrapped lalpulsar.PulsarParameters structure.
Definition: PulsarParametersWrapper.py:256
lalpulsar_run_crosscorr_v2.int
int
Definition: lalpulsar_run_crosscorr_v2.py:41
test_flat_prior.range
range
Definition: test_flat_prior.py:105
test_heterodyned_model.test_six
def test_six()
Definition: test_heterodyned_model.py:939
test_heterodyned_model.test_five
def test_five()
Definition: test_heterodyned_model.py:865
test_heterodyned_model.test_three
def test_three(harmonic)
Definition: test_heterodyned_model.py:690
test_heterodyned_model.test_seven
def test_seven()
Definition: test_heterodyned_model.py:1017
test_heterodyned_model.test_two
def test_two()
Definition: test_heterodyned_model.py:613
test_heterodyned_model.test_one
def test_one(det)
Definition: test_heterodyned_model.py:549
test_heterodyned_model.test_transient
def test_transient(det)
Test the transient signal models with a rectangular and exponential decay.
Definition: test_heterodyned_model.py:1113
test_heterodyned_model.test_four
def test_four()
Definition: test_heterodyned_model.py:770