Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALPulsar 7.1.1.1-3a66518
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
makefakedata_v4.c
Go to the documentation of this file.
1/*
2* Copyright (C) 2008, 2010, 2022 Karl Wette
3* Copyright (C) 2008 Chris Messenger
4* Copyright (C) 2007 Badri Krishnan, Reinhard Prix
5*
6* This program is free software; you can redistribute it and/or modify
7* it under the terms of the GNU General Public License as published by
8* the Free Software Foundation; either version 2 of the License, or
9* (at your option) any later version.
10*
11* This program is distributed in the hope that it will be useful,
12* but WITHOUT ANY WARRANTY; without even the implied warranty of
13* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14* GNU General Public License for more details.
15*
16* You should have received a copy of the GNU General Public License
17* along with with program; see the file COPYING. If not, write to the
18* Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
19* MA 02110-1301 USA
20*/
21
22/**
23 * \file
24 * \ingroup lalpulsar_bin_Tools
25 * \author R. Prix, M.A. Papa, X. Siemens, B. Allen, C. Messenger
26 */
27
28/*-----------------------------------------------------------------------
29 *
30 * File Name: makefakedata_v4.c
31 *
32 * Authors: R. Prix, M.A. Papa, X. Siemens, B. Allen, C. Messenger
33 *
34 * This code is a descendant of an earlier implementation 'makefakedata_v2.[ch]'
35 * by Badri Krishnan, Bruce Allen, Maria Alessandra Papa, Reinhard Prix, Xavier Siemens, Yousuke Itoh
36 *
37 *-----------------------------------------------------------------------
38 */
39
40/* ---------- includes ---------- */
41#include "config.h"
42
43#include <sys/stat.h>
44
45#include <lal/AVFactories.h>
46#include <lal/SeqFactories.h>
47#include <lal/FrequencySeries.h>
48#include <lal/LALInitBarycenter.h>
49#include <gsl/gsl_math.h>
50
51#include <lal/LALString.h>
52#include <lal/UserInput.h>
53#include <lal/SFTfileIO.h>
54#include <lal/GeneratePulsarSignal.h>
55#include <lal/SimulatePulsarSignal.h>
56#include <lal/TimeSeries.h>
57#include <lal/BinaryPulsarTiming.h>
58#include <lal/Window.h>
59#include <lal/TranslateAngles.h>
60#include <lal/TranslateMJD.h>
61#include <lal/TransientCW_utils.h>
62#include <lal/LALPulsarVCSInfo.h>
63
64#ifdef HAVE_LIBLALFRAME
65#include <lal/LALFrameIO.h>
66#endif
67
68/***************************************************/
69#define SQ(x) ( (x) * (x) )
70
71#define TRUE 1
72#define FALSE 0
73
74/*----------------------------------------------------------------------*/
75/** configuration-variables derived from user-variables */
76typedef struct {
77 PulsarParams pulsar; /**< pulsar signal-parameters (amplitude + doppler */
78 EphemerisData *edat; /**< ephemeris-data */
79 LALDetector site; /**< detector-site info */
80
81 LIGOTimeGPS startTimeGPS; /**< start-time of observation */
82 UINT4 duration; /**< total duration of observation in seconds */
83
84 LIGOTimeGPSVector *timestamps;/**< a vector of timestamps to generate time-series/SFTs for */
85
86 REAL8 fmin_eff; /**< 'effective' fmin: round down such that fmin*Tsft = integer! */
87 REAL8 fBand_eff; /**< 'effective' frequency-band such that samples/SFT is integer */
88 REAL8Vector *spindown; /**< vector of frequency-derivatives of GW signal */
89
90 SFTVector *noiseSFTs; /**< vector of noise-SFTs to be added to signal */
91 REAL8 noiseSigma; /**< sigma for Gaussian noise to be added */
92
93 REAL4Window *window; /**< window function for the time series */
94 CHAR *window_type; /**< name of window function */
95 REAL8 window_param; /**< parameter of window function */
96
97 COMPLEX8FrequencySeries *transfer; /**< detector's transfer function for use in hardware-injection */
98
99 transientWindow_t transientWindow; /**< properties of transient-signal window */
100 CHAR *VCSInfoString; /**< Git version string */
102
103typedef enum {
104 GENERATE_ALL_AT_ONCE = 0, /**< generate whole time-series at once before turning into SFTs */
105 GENERATE_PER_SFT, /**< generate timeseries individually for each SFT */
106 GENERATE_LAST /**< end-marker */
108
109// ----- User variables
110typedef struct {
111 /* output */
112 CHAR *outSFTbname; /**< Path and basefilename of output SFT files */
113 BOOLEAN outSingleSFT; /**< use to output a single concatenated SFT */
114
115 CHAR *TDDfile; /**< Filename for ASCII output time-series */
116 CHAR *TDDframedir; /**< directory for frame file output time-series */
117 CHAR *frameDesc; /**< description field entry in the frame filename */
118
119 BOOLEAN hardwareTDD; /**< Binary output timeseries in chunks of Tsft for hardware injections. */
120
121 CHAR *logfile; /**< name of logfile */
122
123 /* specify start + duration */
124 CHAR *timestampsFile; /**< Timestamps file */
125 LIGOTimeGPS startTime; /**< Start-time of requested signal in detector-frame (GPS seconds) */
126 INT4 duration; /**< Duration of requested signal in seconds */
127
128 /* generation mode of timer-series: all-at-once or per-sft */
129 INT4 generationMode; /**< How to generate the timeseries: all-at-once or per-sft */
130
131 /* time-series sampling + heterodyning frequencies */
132 REAL8 fmin; /**< Lowest frequency in output SFT (= heterodyning frequency) */
133 REAL8 Band; /**< bandwidth of output SFT in Hz (= 1/2 sampling frequency) */
134
135 /* SFT params */
136 REAL8 Tsft; /**< SFT time baseline Tsft */
137 REAL8 SFToverlap; /**< overlap SFTs by this many seconds */
138
139 /* noise to add [OPTIONAL] */
140 CHAR *noiseSFTs; /**< Glob-like pattern specifying noise-SFTs to be added to signal */
141 REAL8 noiseSqrtSh; /**< single-sided sqrt(Sh) for Gaussian noise */
142
143 /* Window function [OPTIONAL] */
144 CHAR *window; /**< Windowing function for the time series */
145 REAL8 windowParam; /**< Hann fraction of Tukey window (0.0=rect; 1,0=han; 0.5=default */
146
147 /* Detector and ephemeris */
148 CHAR *IFO; /**< Detector: H1, L1, H2, V1, ... */
149
150 CHAR *actuation; /**< filename containg detector actuation function */
151 REAL8 actuationScale; /**< Scale-factor to be applied to actuation-function */
152
153 CHAR *ephemEarth; /**< Earth ephemeris file to use */
154 CHAR *ephemSun; /**< Sun ephemeris file to use */
155
156 /* pulsar parameters [REQUIRED] */
157 LIGOTimeGPS refTime; /**< Pulsar reference epoch tRef in SSB ('0' means: use startTime converted to SSB) */
158
159 REAL8 h0; /**< overall signal amplitude h0 */
160 REAL8 cosi; /**< cos(iota) of inclination angle iota */
161 REAL8 aPlus; /**< ALTERNATIVE to {h0,cosi}: Plus polarization amplitude aPlus */
162 REAL8 aCross; /**< ALTERNATIVE to {h0,cosi}: Cross polarization amplitude aCross */
163 REAL8 psi; /**< Polarization angle psi */
164 REAL8 phi0; /**< Initial phase phi */
165
166 REAL8 Alpha; /**< Right ascension [radians] alpha of pulsar */
167 REAL8 Delta; /**< Declination [radians] delta of pulsar */
168
170 REAL8 f1dot; /**< First spindown parameter f' */
171 REAL8 f2dot; /**< Second spindown parameter f'' */
172 REAL8 f3dot; /**< Third spindown parameter f''' */
173
174 /* orbital parameters [OPTIONAL] */
175
176 REAL8 orbitasini; /**< Projected orbital semi-major axis in seconds (a/c) */
177 REAL8 orbitEcc; /**< Orbital eccentricity */
178 LIGOTimeGPS orbitTp; /**< 'true' epoch of periapsis passage */
179 REAL8 orbitPeriod; /**< Orbital period (seconds) */
180 REAL8 orbitArgp; /**< Argument of periapsis (radians) */
181
182 /* precision-level of signal-generation */
183 BOOLEAN exactSignal; /**< generate signal timeseries as exactly as possible (slow) */
184 BOOLEAN lineFeature; /**< generate a monochromatic line instead of a pulsar-signal */
185
186 INT4 randSeed; /**< allow user to specify random-number seed for reproducible noise-realizations */
187
188 CHAR *parfile; /** option .par file path */
189 CHAR *transientWindowType; /**< name of transient window ('rect', 'exp',...) */
190 REAL8 transientStartTime; /**< GPS start-time of transient window */
191 REAL8 transientTauDays; /**< time-scale in days of transient window */
192
193 REAL8 sourceDeltaT; /**< source-frame sampling period. '0' implies previous internal defaults */
194
196
197// ----- global variables ----------
198
199// ---------- local prototypes ----------
200int XLALInitUserVars( UserVariables_t *uvar, int argc, char *argv[] );
202int XLALWriteMFDlog( const char *logfile, const ConfigVars_t *cfg );
204int XLALFreeMem( ConfigVars_t *cfg );
205
206BOOLEAN is_directory( const CHAR *fname );
207int XLALIsValidDescriptionField( const char *desc );
208
209/*----------------------------------------------------------------------
210 * main function
211 *----------------------------------------------------------------------*/
212int
213main( int argc, char *argv[] )
214{
217 REAL4TimeSeries *Tseries = NULL;
218 UINT4 i_chunk, numchunks;
219 FILE *fpSingleSFT = NULL;
220 size_t len;
222
223
224 /* ------------------------------
225 * read user-input and set up shop
226 *------------------------------*/
227 XLAL_CHECK( XLALInitUserVars( &uvar, argc, argv ) == XLAL_SUCCESS, XLAL_EFUNC );
228
230
231 /*----------------------------------------
232 * fill the PulsarSignalParams struct
233 *----------------------------------------*/
234 /* pulsar params */
235 params.pulsar.refTime = GV.pulsar.Doppler.refTime;
236 params.pulsar.position.system = COORDINATESYSTEM_EQUATORIAL;
237 params.pulsar.position.longitude = GV.pulsar.Doppler.Alpha;
238 params.pulsar.position.latitude = GV.pulsar.Doppler.Delta;
239 params.pulsar.aPlus = GV.pulsar.Amp.aPlus;
240 params.pulsar.aCross = GV.pulsar.Amp.aCross;
241 params.pulsar.phi0 = GV.pulsar.Amp.phi0;
242 params.pulsar.psi = GV.pulsar.Amp.psi;
243
244 params.pulsar.f0 = GV.pulsar.Doppler.fkdot[0];
245 params.pulsar.spindown = GV.spindown;
246 params.orbit.tp = GV.pulsar.Doppler.tp;
247 params.orbit.argp = GV.pulsar.Doppler.argp;
248 params.orbit.asini = GV.pulsar.Doppler.asini;
249 params.orbit.ecc = GV.pulsar.Doppler.ecc;
250 params.orbit.period = GV.pulsar.Doppler.period;
251
252 params.sourceDeltaT = uvar.sourceDeltaT;
253
254 /* detector params */
255 params.transfer = GV.transfer; /* detector transfer function (NULL if not used) */
256 params.site = &( GV.site );
257 params.ephemerides = GV.edat;
258
259 /* characterize the output time-series */
260 if ( ! uvar.exactSignal ) { /* GeneratePulsarSignal() uses 'idealized heterodyning' */
261 params.samplingRate = 2.0 * GV.fBand_eff; /* sampling rate of time-series (=2*frequency-Band) */
262 params.fHeterodyne = GV.fmin_eff; /* heterodyning frequency for output time-series */
263 } else { /* in the exact-signal case: don't do heterodyning, sample at least twice highest frequency */
264 params.samplingRate = fmax( 2.0 * ( params.pulsar.f0 + 2 ), 2 * ( GV.fmin_eff + GV.fBand_eff ) );
265 params.fHeterodyne = 0;
266 }
267
268 /* set-up main-loop according to 'generation-mode' (all-at-once' or 'per-sft') */
269 switch ( uvar.generationMode ) {
271 params.duration = GV.duration;
272 numchunks = 1;
273 break;
274 case GENERATE_PER_SFT:
275 params.duration = ( UINT4 ) ceil( uvar.Tsft );
276 numchunks = GV.timestamps->length;
277 break;
278 default:
279 XLAL_ERROR( XLAL_EINVAL, "Illegal value for generationMode %d\n\n", uvar.generationMode );
280 break;
281 } /* switch generationMode */
282
283 /* if user requesting single concatenated SFT */
284 if ( uvar.outSingleSFT ) {
285 /* check that user isn't giving a directory */
286 if ( is_directory( uvar.outSFTbname ) ) {
287 XLAL_ERROR( XLAL_ETYPE, "'%s' is a directory, but --outSingleSFT expects a filename!\n", uvar.outSFTbname );
288 }
289
290 /* open concatenated SFT file for writing */
291 if ( ( fpSingleSFT = fopen( uvar.outSFTbname, "wb" ) ) == NULL ) {
292 XLAL_ERROR( XLAL_EIO, "Failed to open file '%s' for writing: %s\n\n", uvar.outSFTbname, strerror( errno ) );
293 }
294 } // if outSingleSFT
295
296 /* ----------
297 * Main loop: produce time-series and turn it into SFTs,
298 * either all-at-once or per-sft
299 * ----------*/
300 for ( i_chunk = 0; i_chunk < numchunks; i_chunk++ ) {
301 params.startTimeGPS = GV.timestamps->data[i_chunk];
302
303 /*----------------------------------------
304 * generate the signal time-series
305 *----------------------------------------*/
306 if ( uvar.exactSignal ) {
307 XLAL_CHECK( ( Tseries = XLALSimulateExactPulsarSignal( &params ) ) != NULL, XLAL_EFUNC );
308 } else if ( uvar.lineFeature ) {
309 XLAL_CHECK( ( Tseries = XLALGenerateLineFeature( &params ) ) != NULL, XLAL_EFUNC );
310 } else {
311 XLAL_CHECK( ( Tseries = XLALGeneratePulsarSignal( &params ) ) != NULL, XLAL_EFUNC );
312 }
313
314 XLAL_CHECK( XLALApplyTransientWindow( Tseries, GV.transientWindow ) == XLAL_SUCCESS, XLAL_EFUNC );
315
316 /* for HARDWARE-INJECTION:
317 * before the first chunk we send magic number and chunk-length to stdout
318 */
319 if ( uvar.hardwareTDD && ( i_chunk == 0 ) ) {
320 REAL4 magic = 1234.5;
321 UINT4 length = Tseries->data->length;
322 if ( ( 1 != fwrite( &magic, sizeof( magic ), 1, stdout ) ) || ( 1 != fwrite( &length, sizeof( INT4 ), 1, stdout ) ) ) {
323 perror( "Failed to write to stdout" );
325 }
326 } /* if hardware-injection and doing first chunk */
327
328
329 /* add Gaussian noise if requested */
330 if ( GV.noiseSigma > 0 ) {
331 // NOTE: seed=0 means randomize seed from /dev/urandom, otherwise we'll have to increment it for each chunk here
332 XLAL_CHECK( XLALAddGaussianNoise( Tseries, GV.noiseSigma, ( uvar.randSeed == 0 ) ? 0 : ( uvar.randSeed + i_chunk ) ) == XLAL_SUCCESS, XLAL_EFUNC );
333 }
334
335 /* output ASCII time-series if requested */
336 if ( uvar.TDDfile ) {
337 CHAR *fname = XLALCalloc( 1, len = strlen( uvar.TDDfile ) + 10 );
338 XLAL_CHECK( fname != NULL, XLAL_ENOMEM, "XLALCalloc(1,%zu) failed\n", len );
339 sprintf( fname, "%s.%02d", uvar.TDDfile, i_chunk );
341 XLALFree( fname );
342 } /* if outputting ASCII time-series */
343
344 /* output time-series to frames if requested */
345 if ( uvar.TDDframedir ) {
346#ifndef HAVE_LIBLALFRAME
347 XLAL_ERROR( XLAL_EINVAL, "--TDDframedir option not supported, code has to be compiled with lalframe\n" );
348#else
349 /* use standard frame output filename format */
351 len = strlen( uvar.TDDframedir ) + strlen( uvar.frameDesc ) + 100;
352 char *fname;
353 char IFO[2] = { Tseries->name[0], Tseries->name[1] };
354 XLAL_CHECK( ( fname = LALCalloc( 1, len ) ) != NULL, XLAL_ENOMEM );
355 size_t written = snprintf( fname, len, "%s/%c-%c%c_%s-%d-%d.gwf",
356 uvar.TDDframedir, IFO[0], IFO[0], IFO[1], uvar.frameDesc, params.startTimeGPS.gpsSeconds, ( int )params.duration );
357 XLAL_CHECK( written < len, XLAL_ESIZE, "Frame-filename exceeds expected maximal length (%zu): '%s'\n", len, fname );
358
359 /* define the output frame */
360 LALFrameH *outFrame;
361 XLAL_CHECK( ( outFrame = XLALFrameNew( &( params.startTimeGPS ), params.duration, uvar.frameDesc, 1, 0, 0 ) ) != NULL, XLAL_EFUNC );
362
363 /* add timeseries to the frame - make sure to change the timeseries name since this is used as the channel name */
364 char buffer[LALNameLength];
365 written = snprintf( buffer, LALNameLength, "%s:%s", Tseries->name, uvar.frameDesc );
366 XLAL_CHECK( written < LALNameLength, XLAL_ESIZE, "Updated frame name exceeds max length (%d): '%s'\n", LALNameLength, buffer );
367 strcpy( Tseries->name, buffer );
368
370
371 /* Here's where we add extra information into the frame - first we add the command line args used to generate it */
373 XLALFrameAddFrHistory( outFrame, __FILE__, hist );
374
375 /* then we add the version string */
376 XLALFrameAddFrHistory( outFrame, __FILE__, GV.VCSInfoString );
377
378 /* output the frame to file - compression level 1 (higher values make no difference) */
379 XLAL_CHECK( ( XLALFrameWrite( outFrame, fname ) == 0 ), XLAL_EFUNC );
380
381 /* free the frame, frame file name and history memory */
382 XLALFrameFree( outFrame );
383 LALFree( fname );
384 LALFree( hist );
385#endif
386 } /* if outputting time-series to frames */
387
388
389 /* if hardware injection: send timeseries in binary-format to stdout */
390 if ( uvar.hardwareTDD ) {
391 UINT4 length = Tseries->data->length;
392 REAL4 *datap = Tseries->data->data;
393
394 if ( length != fwrite( datap, sizeof( datap[0] ), length, stdout ) ) {
395 perror( "Fatal error in writing binary data to stdout\n" );
396 XLAL_ERROR( XLAL_EIO, "fwrite() failed\n" );
397 }
398 fflush( stdout );
399
400 } /* if hardware injections */
401
402 /*----------------------------------------
403 * last step: turn this timeseries into SFTs
404 * and output them to disk
405 *----------------------------------------*/
406 SFTVector *SFTs = NULL;
407 if ( uvar.outSFTbname ) {
408 SFTParams XLAL_INIT_DECL( sftParams );
410 SFTVector noise;
411
412 sftParams.Tsft = uvar.Tsft;
413
414 switch ( uvar.generationMode ) {
416 sftParams.timestamps = GV.timestamps;
417 sftParams.noiseSFTs = GV.noiseSFTs;
418 break;
419 case GENERATE_PER_SFT:
420 ts.length = 1;
421 ts.data = &( GV.timestamps->data[i_chunk] );
422 sftParams.timestamps = &( ts );
423 sftParams.noiseSFTs = NULL;
424
425 if ( GV.noiseSFTs ) {
426 noise.length = 1;
427 noise.data = &( GV.noiseSFTs->data[i_chunk] );
428 sftParams.noiseSFTs = &( noise );
429 }
430
431 break;
432
433 default:
434 XLAL_ERROR( XLAL_EINVAL, "Invalid Value --generationMode=%d\n", uvar.generationMode );
435 break;
436 }
437
438 /* Enter the window function into the SFTparams struct */
439 sftParams.window = GV.window;
440
441 /* get SFTs from timeseries */
442 XLAL_CHECK( ( SFTs = XLALSignalToSFTs( Tseries, &sftParams ) ) != NULL, XLAL_EFUNC );
443
444 /* extract requested band */
445 {
446 SFTVector *outSFTs;
447 XLAL_CHECK( ( outSFTs = XLALExtractStrictBandFromSFTVector( SFTs, uvar.fmin, uvar.Band ) ) != NULL, XLAL_EFUNC );
448 XLALDestroySFTVector( SFTs );
449 SFTs = outSFTs;
450 }
451
452 /* generate comment string */
453 CHAR *logstr;
455 char *comment = XLALCalloc( 1, len = strlen( logstr ) + strlen( GV.VCSInfoString ) + 512 );
456 XLAL_CHECK( comment != NULL, XLAL_ENOMEM, "XLALCalloc(1,%zu) failed.\n", len );
457 sprintf( comment, "Generated by:\n%s\n%s\n", logstr, GV.VCSInfoString );
458
459 /* if user requesting single concatenated SFT */
461 XLAL_CHECK( XLALFillSFTFilenameSpecStrings( &spec, uvar.outSFTbname, NULL, NULL, GV.window_type, "mfdv4", NULL, NULL ) == XLAL_SUCCESS, XLAL_EFUNC );
462 spec.window_param = GV.window_param;
463 if ( uvar.outSingleSFT ) {
464 /* write all SFTs to concatenated file */
465 for ( UINT4 k = 0; k < SFTs->length; k++ ) {
466 int ret = XLALWriteSFT2FilePointer( &( SFTs->data[k] ), fpSingleSFT, spec.window_type, spec.window_param, comment );
467 XLAL_CHECK( ret == XLAL_SUCCESS, XLAL_EFUNC, "XLALWriteSFT2FilePointer() failed to write SFT %d / %d to '%s'!\n", k + 1, SFTs->length, uvar.outSFTbname );
468 } // for k < numSFTs
469 } // if outSingleSFT
470 else {
471 XLAL_CHECK( is_directory( uvar.outSFTbname ), XLAL_EINVAL, "ERROR: the --outSFTbname '%s' is not a directory!\n", uvar.outSFTbname );
473 }
474 XLALFree( logstr );
475 XLALFree( comment );
476
477 } /* if outSFTbname */
478
479 /* free memory */
480 if ( Tseries ) {
482 Tseries = NULL;
483 }
484
485 if ( SFTs ) {
486 XLALDestroySFTVector( SFTs );
487 SFTs = NULL;
488 }
489
490 } /* for i_chunk < numchunks */
491
492 /* if user requesting single concatenated SFT */
493 if ( uvar.outSingleSFT ) {
494
495 /* close concatenated SFT */
496 fclose( fpSingleSFT );
497
498 }
499
500 XLALFreeMem( &GV ); /* free the config-struct */
501
503
504 return 0;
505} /* main */
506
507/**
508 * Handle user-input and set up shop accordingly, and do all
509 * consistency-checks on user-input.
510 */
511int
513{
514 XLAL_CHECK( cfg != NULL, XLAL_EINVAL, "Invalid NULL input 'cfg'\n" );
515 XLAL_CHECK( uvar != NULL, XLAL_EINVAL, "Invalid NULL input 'uvar'\n" );
516
518 XLAL_CHECK( cfg->VCSInfoString != NULL, XLAL_EFUNC, "XLALVCSInfoString failed.\n" );
519
520 BOOLEAN have_parfile = XLALUserVarWasSet( &uvar->parfile );
521 BinaryPulsarParams pulparams;
522
523 char *window_type_from_noiseSFTs = NULL;
524 REAL8 window_param_from_noiseSFTs = 0;
525
526 /* read in par file parameters if given */
527 if ( have_parfile ) {
528 XLALReadTEMPOParFileOrig( &pulparams, uvar->parfile );
529 XLAL_CHECK( xlalErrno == XLAL_SUCCESS, XLAL_EFUNC, "XLALReadTEMPOParFileOrig() failed for parfile = '%s', xlalErrno = %d\n", uvar->parfile, xlalErrno );
530 XLAL_CHECK( pulparams.f0 > 0, XLAL_EINVAL, "Invalid .par file values, need f0 > 0!\n" );
531 XLAL_CHECK( ( pulparams.pepoch > 0 ) || ( pulparams.posepoch > 0 ), XLAL_EINVAL, "Invalid .par file values, need PEPOCH or POSEPOCH!\n" );
532 }
533
534 /* if requested, log all user-input and code-versions */
535 if ( uvar->logfile ) {
536 XLAL_CHECK( XLALWriteMFDlog( uvar->logfile, cfg ) == XLAL_SUCCESS, XLAL_EFUNC, "XLALWriteMFDlog() failed with xlalErrno = %d\n", xlalErrno );
537 }
538
539 { /* ========== translate user-input into 'PulsarParams' struct ========== */
540 BOOLEAN have_h0 = XLALUserVarWasSet( &uvar->h0 );
541 BOOLEAN have_cosi = XLALUserVarWasSet( &uvar->cosi );
542 BOOLEAN have_aPlus = XLALUserVarWasSet( &uvar->aPlus );
543 BOOLEAN have_aCross = XLALUserVarWasSet( &uvar->aCross );
544 BOOLEAN have_Alpha = XLALUserVarWasSet( &uvar->Alpha );
545 BOOLEAN have_Delta = XLALUserVarWasSet( &uvar->Delta );
546
547 /*check .par file for gw parameters*/
548 if ( have_parfile ) {
549 if ( pulparams.h0 != 0 ) {
550 uvar->h0 = pulparams.h0;
551 uvar->cosi = pulparams.cosiota;
552 uvar->phi0 = pulparams.phi0;
553 uvar->psi = pulparams.psi;
554 have_h0 = 1; /*Set to TRUE as uvar->h0 not declared on command line -- same for rest*/
555 have_cosi = 1;
556 } else {
557 uvar->aPlus = pulparams.Aplus;
558 uvar->aCross = pulparams.Across;
559 uvar->phi0 = pulparams.phi0;
560 uvar->psi = pulparams.psi;
561 have_aPlus = 1;
562 have_aCross = 1;
563 }
564 uvar->Freq = 2.*pulparams.f0;
565 uvar->Alpha = pulparams.ra;
566 uvar->Delta = pulparams.dec;
567 have_Alpha = 1;
568 have_Delta = 1;
569 }
570
571 /* ----- {h0,cosi} or {aPlus,aCross} ----- */
572 if ( ( have_aPlus || have_aCross ) && ( have_h0 || have_cosi ) ) {
573 XLAL_ERROR( XLAL_EINVAL, "Need to specify EITHER {h0,cosi} OR {aPlus, aCross}!\n\n" );
574 }
575 if ( ( have_h0 ^ have_cosi ) ) {
576 XLAL_ERROR( XLAL_EINVAL, "Need BOTH --h0 and --cosi!\n\n" );
577 }
578 if ( ( have_aPlus ^ have_aCross ) ) {
579 XLAL_ERROR( XLAL_EINVAL, "Need BOTH --aPlus and --aCross !\n\n" );
580 }
581 if ( have_h0 && have_cosi ) {
582 /* translate {h_0, cosi} into A_{+,x} */
583 /* assume at 2f */
584 REAL8 h0 = uvar->h0;
585 REAL8 cosi = uvar->cosi;
586 cfg->pulsar.Amp.aPlus = 0.5 * h0 * ( 1.0 + SQ( cosi ) );
587 cfg->pulsar.Amp.aCross = h0 * cosi;
588 } else if ( have_aPlus && have_aCross ) {
589 cfg->pulsar.Amp.aPlus = uvar->aPlus;
590 cfg->pulsar.Amp.aCross = uvar->aCross;
591 } else {
592 cfg->pulsar.Amp.aPlus = 0.0;
593 cfg->pulsar.Amp.aCross = 0.0;
594 }
595 cfg->pulsar.Amp.phi0 = uvar->phi0;
596 cfg->pulsar.Amp.psi = uvar->psi;
597
598 /* ----- signal Frequency ----- */
599 cfg->pulsar.Doppler.fkdot[0] = uvar->Freq;
600
601 /* ----- skypos ----- */
602 if ( ( have_Alpha && !have_Delta ) || ( !have_Alpha && have_Delta ) ) {
603 XLAL_ERROR( XLAL_EINVAL, "\nSpecify skyposition: need BOTH --Alpha and --Delta!\n\n" );
604 }
605 cfg->pulsar.Doppler.Alpha = uvar->Alpha;
606 cfg->pulsar.Doppler.Delta = uvar->Delta;
607
608 } /* Pulsar signal parameters */
609
610 /* ---------- prepare vector of spindown parameters ---------- */
611 {
612 UINT4 msp = 0; /* number of spindown-parameters */
613 if ( have_parfile ) {
614 uvar->f1dot = 2.*pulparams.f1;
615 uvar->f2dot = 2.*pulparams.f2;
616 uvar->f3dot = 2.*pulparams.f3;
617 }
618 cfg->pulsar.Doppler.fkdot[1] = uvar->f1dot;
619 cfg->pulsar.Doppler.fkdot[2] = uvar->f2dot;
620 cfg->pulsar.Doppler.fkdot[3] = uvar->f3dot;
621
622 if ( uvar->f3dot != 0 ) {
623 msp = 3; /* counter number of spindowns */
624 } else if ( uvar->f2dot != 0 ) {
625 msp = 2;
626 } else if ( uvar->f1dot != 0 ) {
627 msp = 1;
628 } else {
629 msp = 0;
630 }
631 if ( msp ) {
632 /* memory not initialized, but ALL alloc'ed entries will be set below! */
633 cfg->spindown = XLALCreateREAL8Vector( msp );
634 XLAL_CHECK( cfg->spindown != NULL, XLAL_EFUNC, "XLALCreateREAL8Vector ( %d ) failed.\n", msp );
635 }
636 switch ( msp ) {
637 case 3:
638 cfg->spindown->data[2] = uvar->f3dot;
639#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
640 __attribute__( ( fallthrough ) );
641#endif
642 case 2:
643 cfg->spindown->data[1] = uvar->f2dot;
644#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
645 __attribute__( ( fallthrough ) );
646#endif
647 case 1:
648 cfg->spindown->data[0] = uvar->f1dot;
649 break;
650 case 0:
651 break;
652 default:
653 XLAL_ERROR( XLAL_EINVAL, "\nmsp = %d makes no sense to me...\n\n", msp );
654 } /* switch(msp) */
655
656 } /* END: prepare spindown parameters */
657
658 CHAR *channelName = NULL;
659 /* ---------- prepare detector ---------- */
660 {
661 const LALDetector *site;
662
663 site = XLALGetSiteInfo( uvar->IFO );
664 XLAL_CHECK( site != NULL, XLAL_EFUNC, "XLALGetSiteInfo('%s') failed\n", uvar->IFO );
665 channelName = XLALGetChannelPrefix( uvar->IFO );
666 XLAL_CHECK( channelName != NULL, XLAL_EFUNC, "XLALGetChannelPrefix('%s') failed.\n", uvar->IFO );
667
668 cfg->site = ( *site );
669 }
670
671 /* check for negative fmin and Band, which would break the fmin_eff, fBand_eff calculation below */
672 XLAL_CHECK( uvar->fmin >= 0, XLAL_EDOM, "Invalid negative frequency fmin=%f!\n\n", uvar->fmin );
673 XLAL_CHECK( uvar->Band >= 0, XLAL_EDOM, "Invalid negative frequency band Band=%f!\n\n", uvar->Band );
674
675 /* ---------- for SFT output: calculate effective fmin and Band ---------- */
676 // Note: this band is only used for internal data operations; ultimately SFTs covering
677 // the half-open interval uvar->[fmin,fmin+Band) are returned to the user using
678 // XLALExtractStrictBandFromSFTVector()
679 if ( XLALUserVarWasSet( &uvar->outSFTbname ) ) {
680 UINT4 firstBin, numBins;
681 /* calculate "effective" fmin from uvar->fmin:
682 * make sure that fmin_eff * Tsft = integer, such that freqBinIndex corresponds
683 * to a frequency-index of the non-heterodyned signal.
684 */
685
686 XLAL_CHECK( XLALFindCoveringSFTBins( &firstBin, &numBins, uvar->fmin, uvar->Band, uvar->Tsft ) == XLAL_SUCCESS, XLAL_EFUNC );
687
688 /* Adjust Band correspondingly */
689 REAL8 dFreq = 1.0 / uvar->Tsft;
690 cfg->fmin_eff = firstBin * dFreq;
691 cfg->fBand_eff = ( numBins - 1 ) * dFreq;
692
693 if ( lalDebugLevel ) {
694 if ( fabs( cfg->fmin_eff - uvar->fmin ) > LAL_REAL8_EPS
695 || fabs( cfg->fBand_eff - uvar->Band ) > LAL_REAL8_EPS )
696 printf( "\nWARNING: for SFT-creation we had to adjust (fmin,Band) to"
697 " fmin_eff=%.15g and Band_eff=%.15g\n\n", cfg->fmin_eff, cfg->fBand_eff );
698 }
699
700 } /* END: SFT-specific corrections to fmin and Band */
701 else {
702 /* producing pure time-series output ==> no adjustments necessary */
703 cfg->fmin_eff = uvar->fmin;
704 cfg->fBand_eff = uvar->Band;
705 }
706
707
708 /* ---------- determine timestamps to produce signal for ---------- */
709 {
710 /* check input consistency: *uvar->timestampsFile, uvar->startTime, uvar->duration */
711 BOOLEAN haveStart, haveDuration, haveTimestampsFile, haveOverlap;
712
713 haveStart = XLALUserVarWasSet( &uvar->startTime );
714 haveDuration = XLALUserVarWasSet( &uvar->duration );
715 haveTimestampsFile = ( uvar->timestampsFile != NULL );
716 haveOverlap = ( uvar->SFToverlap > 0 );
717
718 if ( ( haveDuration && !haveStart ) || ( !haveDuration && haveStart ) ) {
719 XLAL_ERROR( XLAL_EINVAL, "Need BOTH --startTime AND --duration if you give one of them !\n\n" );
720 }
721
722 /* don't allow using --SFToverlap with anything other than pure (--startTime,--duration) */
723 if ( haveOverlap && ( uvar->noiseSFTs || haveTimestampsFile ) ) {
724 XLAL_ERROR( XLAL_EINVAL, "I can't combine --SFToverlap with --noiseSFTs or --timestampsFile, only use with (--startTime, --duration)!\n\n" );
725 }
726
727 /* don't allow --SFToverlap with generationMode==1 (one timeseries per SFT), as this would
728 * result in independent noise realizations in the overlapping regions
729 */
730 if ( haveOverlap && ( uvar->generationMode != 0 ) ) {
731 XLAL_ERROR( XLAL_EINVAL, "--SFToverlap can only be used with --generationMode=0, otherwise we'll get overlapping independent noise!\n\n" );
732 }
733
734 /*-------------------- check special case: Hardware injection ---------- */
735 /* don't allow timestamps-file, noise-SFTs or SFT-output */
736 if ( uvar->hardwareTDD ) {
737 if ( haveTimestampsFile || uvar->noiseSFTs ) {
738 XLAL_ERROR( XLAL_EINVAL, "\nHardware injection: don't specify --timestampsFile or --noiseSFTs\n\n" );
739 }
740 if ( !haveStart || !haveDuration ) {
741 XLAL_ERROR( XLAL_EINVAL, "Hardware injection: need to specify --startTime and --duration!\n\n" );
742 }
743 if ( XLALUserVarWasSet( &uvar->outSFTbname ) ) {
744 XLAL_ERROR( XLAL_EINVAL, "Hardware injection mode is incompatible with producing SFTs\n\n" );
745 }
746 } /* if hardware-injection */
747
748 /* ----- load timestamps from file if given */
749 if ( haveTimestampsFile ) {
750 if ( haveStart || haveDuration || haveOverlap ) {
751 XLAL_ERROR( XLAL_EINVAL, "Using --timestampsFile is incompatible with either of --startTime, --duration or --SFToverlap\n\n" );
752 }
753 if ( ( cfg->timestamps = XLALReadTimestampsFile( uvar->timestampsFile ) ) == NULL ) {
754 XLAL_ERROR( XLAL_EFUNC, "XLALReadTimestampsFile() failed to read timestamps file '%s'\n", uvar->timestampsFile );
755 }
756 if ( ( cfg->timestamps->length == 0 ) || ( cfg->timestamps->data == NULL ) ) {
757 XLAL_ERROR( XLAL_EINVAL, "Got empty timestamps-list from file '%s'\n", uvar->timestampsFile );
758 }
759
760 } /* if haveTimestampsFile */
761
762 /* ----- if real noise-SFTs given: load them now using EITHER (start,start+duration) OR timestamps
763 * as constraints if given, otherwise load all of them. Also require window option to be given
764 */
765 if ( uvar->noiseSFTs ) {
766 REAL8 fMin, fMax;
767 SFTConstraints XLAL_INIT_DECL( constraints );
768 LIGOTimeGPS minStartTime, maxStartTime;
769
770 XLALPrintWarning( "\nWARNING: only SFTs corresponding to the noiseSFT-timestamps will be produced!\n" );
771
772 /* use all additional constraints user might have given */
773 if ( haveStart && haveDuration ) {
774 minStartTime = maxStartTime = uvar->startTime;
775 XLALGPSAdd( &maxStartTime, uvar->duration );
776 constraints.minStartTime = &minStartTime;
777 constraints.maxStartTime = &maxStartTime;
778 char bufGPS1[32], bufGPS2[32];
779 XLALPrintWarning( "\nWARNING: only noise-SFTs between GPS [%s, %s] will be used!\n", XLALGPSToStr( bufGPS1, &minStartTime ), XLALGPSToStr( bufGPS2, &maxStartTime ) );
780 } /* if start+duration given */
781 if ( cfg->timestamps ) {
782 constraints.timestamps = cfg->timestamps;
783 XLALPrintWarning( "\nWARNING: only noise-SFTs corresponding to given timestamps '%s' will be used!\n", uvar->timestampsFile );
784 } /* if we have timestamps already */
785
786 /* use detector-constraint */
787 constraints.detector = channelName ;
788
789 SFTCatalog *catalog = NULL;
790 XLAL_CHECK( ( catalog = XLALSFTdataFind( uvar->noiseSFTs, &constraints ) ) != NULL, XLAL_EFUNC );
791
792 /* check if anything matched */
793 if ( catalog->length == 0 ) {
794 XLAL_ERROR( XLAL_EFAILED, "No noise-SFTs matching the constraints (IFO, start+duration, timestamps) were found!\n" );
795 }
796
797 /* extract SFT window function */
798 window_type_from_noiseSFTs = XLALStringDuplicate( catalog->data[0].window_type );
799 XLAL_CHECK( window_type_from_noiseSFTs != NULL, XLAL_ENOMEM );
800 window_param_from_noiseSFTs = catalog->data[0].window_param;
801
802 /* load effective frequency-band from noise-SFTs */
803 fMin = cfg->fmin_eff;
804 fMax = fMin + cfg->fBand_eff;
805
806 cfg->noiseSFTs = XLALLoadSFTs( catalog, fMin, fMax );
807 XLAL_CHECK( cfg->noiseSFTs != NULL, XLAL_EFUNC, "XLALLoadSFTs() failed\n" );
808
809 XLALDestroySFTCatalog( catalog );
810
811 /* get timestamps from the loaded noise SFTs */
812 if ( cfg->timestamps ) {
814 }
816 XLAL_CHECK( cfg->timestamps != NULL, XLAL_EFUNC, "XLALExtractTimestampsFromSFTs() failed\n" );
817
818 } /* if uvar->noiseSFTs */
819
820 /* have we got our timestamps yet?: If not, we must get them from (start, duration) user-input */
821 if ( ! cfg->timestamps ) {
822 if ( !haveStart || !haveDuration ) {
823 XLAL_ERROR( XLAL_EINVAL, "Need to have either --timestampsFile OR (--startTime,--duration) OR --noiseSFTs\n\n" );
824 }
825 if ( uvar->SFToverlap > uvar->Tsft ) {
826 XLAL_ERROR( XLAL_EINVAL, "--SFToverlap cannot be larger than --Tsft!\n\n" );
827 }
828
829 /* internally always use timestamps, so we generate them */
830 XLAL_CHECK( ( cfg->timestamps = XLALMakeTimestamps( uvar->startTime, uvar->duration, uvar->Tsft, uvar->SFToverlap ) ) != NULL, XLAL_EFUNC );
831
832 } /* if !cfg->timestamps */
833
834 /* ----- figure out start-time and duration ----- */
835 {
838
839 t0 = cfg->timestamps->data[0];
840 t1 = cfg->timestamps->data[cfg->timestamps->length - 1 ];
841
842 duration = XLALGPSDiff( &t1, &t0 );
843 duration += uvar->Tsft;
844
845 cfg->startTimeGPS = cfg->timestamps->data[0];
846 cfg->duration = ( UINT4 )ceil( duration ); /* round up to seconds */
847 }
848
849 if ( cfg->duration < uvar->Tsft ) {
850 XLAL_ERROR( XLAL_EINVAL, "Requested duration of %d sec is less than minimal chunk-size of Tsft =%.0f sec.\n\n", uvar->duration, uvar->Tsft );
851 }
852
853 } /* END: setup signal start + duration */
854
855 /*----------------------------------------------------------------------*/
856 /* currently there are only two modes: [uvar->generationMode]
857 * 1) produce the whole timeseries at once [default],
858 * [GENERATE_ALL_AT_ONCE]
859 * 2) produce timeseries for each SFT,
860 * [GENERATE_PER_SFT]
861 *
862 * Mode 2 which is useful if 1) is too memory-intensive (e.g. for hardware-injections)
863 *
864 * intermediate modes would require a bit more work because the
865 * case with specified timestamps (allowing for gaps) would be
866 * a bit more tricky ==> this is left for future extensions if found useful
867 */
868 if ( ( uvar->generationMode < 0 ) || ( uvar->generationMode >= GENERATE_LAST ) ) {
869 XLAL_ERROR( XLAL_EINVAL, "Illegal input for 'generationMode': must lie within [0, %d]\n\n", GENERATE_LAST - 1 );
870 }
871
872 /* this is a violation of the UserInput-paradigm that no user-variables
873 * should by modified by the code, but this is the easiest way to do
874 * this here, and the log-output of user-variables will not be changed
875 * by this [it's already been written], so in this case it should be safe..
876 */
877 if ( uvar->hardwareTDD ) {
879 }
880
881 /*--------------------- Prepare windowing of time series ---------------------*/
882 cfg->window = NULL;
883 cfg->window_type = NULL;
884 cfg->window_param = 0;
885
886 {
887 const BOOLEAN have_window = XLALUserVarWasSet( &uvar->window );
888 const BOOLEAN have_param = XLALUserVarWasSet( &uvar->windowParam );
889 const CHAR *window_type_from_uvar = uvar->window;
890 const REAL8 window_param_from_uvar = have_param ? uvar->windowParam : 0;
891 if ( have_window ) {
892 XLAL_CHECK( XLALCheckNamedWindow( window_type_from_uvar, have_param ) == XLAL_SUCCESS, XLAL_EFUNC );
893 }
894 if ( uvar->noiseSFTs ) {
895 const BOOLEAN have_window_type_from_noiseSFTs = ( XLALStringCaseCompare( window_type_from_noiseSFTs, "unknown" ) != 0 );
896 /* need window info from noise SFTs or user input - if we have both, check for consistency */
897 if ( have_window_type_from_noiseSFTs && have_window ) {
898 XLAL_CHECK( XLALCompareSFTWindows( window_type_from_noiseSFTs, window_param_from_noiseSFTs, window_type_from_uvar, window_param_from_uvar ) == XLAL_SUCCESS, XLAL_EFUNC, "Inconsistent SFT window between --noiseSFTs and user input: [%s,%f] vs [%s,%f].", window_type_from_noiseSFTs, window_param_from_noiseSFTs, window_type_from_uvar, window_param_from_uvar );
899 cfg->window_type = XLALStringDuplicate( window_type_from_noiseSFTs );
900 cfg->window_param = window_param_from_noiseSFTs;
901 } else if ( have_window_type_from_noiseSFTs ) {
902 cfg->window_type = XLALStringDuplicate( window_type_from_noiseSFTs );
903 cfg->window_param = window_param_from_noiseSFTs;
904 } else if ( have_window ) {
905 cfg->window_type = XLALStringDuplicate( window_type_from_uvar );
906 cfg->window_param = window_param_from_uvar;
907 } else {
908 /* EITHER noise SFTs must have a known window OR user must specify the window function */
909 XLAL_ERROR( XLAL_EINVAL, "When --noiseSFTs is given and have unknown window, --window is required.\n" );
910 }
911 } else if ( have_window ) {
912 cfg->window_type = XLALStringDuplicate( window_type_from_uvar );
913 XLAL_CHECK( cfg->window_type != NULL, XLAL_ENOMEM );
914 cfg->window_param = window_param_from_uvar;
915 }
916 if ( cfg->window_type ) {
917 /* NOTE: a timeseries of length N*dT has no timestep at N*dT !! (convention) */
918 UINT4 lengthOfTimeSeries = ( UINT4 )round( uvar->Tsft * 2 * cfg->fBand_eff );
919
920 cfg->window = XLALCreateNamedREAL4Window( cfg->window_type, cfg->window_param, lengthOfTimeSeries );
921 XLAL_CHECK( cfg->window != NULL, XLAL_EFUNC, "XLALCreateNamedREAL4Window('%s', %g, %d) failed\n", cfg->window_type, cfg->window_param, lengthOfTimeSeries );
922 }
923 }
924
925 /* Init ephemerides */
926 XLAL_CHECK( ( cfg->edat = XLALInitBarycenter( uvar->ephemEarth, uvar->ephemSun ) ) != NULL, XLAL_EFUNC );
927
928 /* -------------------- handle binary orbital params if given -------------------- */
929
930 /* Consistency check: if any orbital parameters specified, we need all of them (except for nano-seconds)! */
931 {
932 if ( have_parfile ) {
933 if ( pulparams.model != NULL ) {
934 uvar->orbitasini = pulparams.x;
935 uvar->orbitPeriod = pulparams.Pb * 86400;
936 if ( strstr( pulparams.model, "ELL1" ) != NULL ) {
937 REAL8 w, e, eps1, eps2;
938 eps1 = pulparams.eps1;
939 eps2 = pulparams.eps2;
940 w = atan2( eps1, eps2 );
941 e = sqrt( eps1 * eps1 + eps2 * eps2 );
942 uvar->orbitArgp = w;
943 uvar->orbitEcc = e;
944 } else {
945 uvar->orbitArgp = pulparams.w0;
946 uvar->orbitEcc = pulparams.e;
947 }
948 if ( strstr( pulparams.model, "ELL1" ) != NULL ) {
949 REAL8 fe, uasc, Dt;
950 fe = sqrt( ( 1.0 - uvar->orbitEcc ) / ( 1.0 + uvar->orbitEcc ) );
951 uasc = 2.0 * atan( fe * tan( uvar->orbitArgp / 2.0 ) );
952 Dt = ( uvar->orbitPeriod / LAL_TWOPI ) * ( uasc - uvar->orbitEcc * sin( uasc ) );
953 pulparams.T0 = pulparams.Tasc + Dt;
954 }
955 uvar->orbitTp.gpsSeconds = ( UINT4 )floor( pulparams.T0 );
956 uvar->orbitTp.gpsNanoSeconds = ( UINT4 )floor( ( pulparams.T0 - uvar->orbitTp.gpsSeconds ) * 1e9 );
957 }
958 }
959 BOOLEAN set1 = XLALUserVarWasSet( &uvar->orbitasini );
960 BOOLEAN set2 = XLALUserVarWasSet( &uvar->orbitEcc );
961 BOOLEAN set3 = XLALUserVarWasSet( &uvar->orbitPeriod );
962 BOOLEAN set4 = XLALUserVarWasSet( &uvar->orbitArgp );
963 BOOLEAN set5 = XLALUserVarWasSet( &uvar->orbitTp );
964
965 if ( set1 || set2 || set3 || set4 || set5 ) {
966 if ( ( uvar->orbitasini > 0 ) && !( set1 && set2 && set3 && set4 && set5 ) ) {
967 XLAL_ERROR( XLAL_EINVAL, "\nPlease either specify ALL orbital parameters or NONE!\n\n" );
968 }
969 if ( ( uvar->orbitEcc < 0 ) || ( uvar->orbitEcc > 1 ) ) {
970 XLAL_ERROR( XLAL_EINVAL, "\nEccentricity = %g has to lie within [0, 1]\n\n", uvar->orbitEcc );
971 }
972
973 cfg->pulsar.Doppler.tp = uvar->orbitTp;
974
975 /* fill in orbital parameter structure */
976 cfg->pulsar.Doppler.period = uvar->orbitPeriod;
977 cfg->pulsar.Doppler.asini = uvar->orbitasini;
978 cfg->pulsar.Doppler.argp = uvar->orbitArgp;
979 cfg->pulsar.Doppler.ecc = uvar->orbitEcc;
980
981 } /* if one or more orbital parameters were set */
982 else {
983 cfg->pulsar.Doppler.asini = 0 /* isolated pulsar */;
984 }
985 } /* END: binary orbital params */
986
987
988 /* -------------------- handle NOISE params -------------------- */
989 cfg->noiseSigma = uvar->noiseSqrtSh * sqrt( cfg->fBand_eff ); /* convert Sh -> sigma */
990
991
992 /* ----- set "pulsar reference time", i.e. SSB-time at which pulsar params are defined ---------- */
993 if ( XLALUserVarWasSet( &uvar->parfile ) ) {
994 XLALGPSSetREAL8( &( uvar->refTime ), pulparams.pepoch ); /*XLALReadTEMPOParFileOrig converted pepoch to REAL8 */
995 XLALGPSSetREAL8( &( cfg->pulsar.Doppler.refTime ), pulparams.pepoch );
996 } else if ( XLALUserVarWasSet( &uvar->refTime ) ) {
997 cfg->pulsar.Doppler.refTime = uvar->refTime;
998 } else {
999 cfg->pulsar.Doppler.refTime = cfg->timestamps->data[0]; /* internal startTime always found in here*/
1000 }
1001
1002
1003 /* ---------- has the user specified an actuation-function file ? ---------- */
1004 if ( uvar->actuation ) {
1005 /* currently we only allow using a transfer-function for hardware-injections */
1006 if ( !uvar->hardwareTDD ) {
1007 XLAL_ERROR( XLAL_EINVAL, "Error: use of an actuation/transfer function restricted to hardare-injections\n\n" );
1008 } else {
1010 XLAL_CHECK( cfg->transfer != NULL, XLAL_EFUNC );
1011 }
1012
1013 } /* if uvar->actuation */
1014
1015 if ( !uvar->actuation && XLALUserVarWasSet( &uvar->actuationScale ) ) {
1016 XLAL_ERROR( XLAL_EINVAL, "Actuation-scale was specified without actuation-function file!\n\n" );
1017 }
1018
1019 XLALFree( channelName );
1020
1021 /* ----- handle transient-signal window if given ----- */
1022 int twtype;
1024 cfg->transientWindow.type = twtype;
1025
1027 cfg->transientWindow.tau = uvar->transientTauDays * 86400;
1028
1029 /* free memory */
1030 XLALFree( window_type_from_noiseSFTs );
1031
1032 return XLAL_SUCCESS;
1033
1034} /* XLALInitMakefakedata() */
1035
1036
1037/**
1038 * Register all "user-variables", and initialized them from command-line and config-files
1039 */
1040int
1041XLALInitUserVars( UserVariables_t *uvar, int argc, char *argv[] )
1042{
1043 int len;
1044
1045 XLAL_CHECK( uvar != NULL, XLAL_EINVAL, "Invalid NULL input 'uvar'\n" );
1046 XLAL_CHECK( argv != NULL, XLAL_EINVAL, "Invalid NULL input 'argv'\n" );
1047
1048 // ---------- set a few defaults ----------
1049 uvar->ephemEarth = XLALStringDuplicate( "earth00-40-DE405.dat.gz" );
1050 uvar->ephemSun = XLALStringDuplicate( "sun00-40-DE405.dat.gz" );
1051
1052 uvar->Tsft = 1800;
1053
1054 // per default we now generate a timeseries per SFT: slower, but avoids potential confusion about sft-"nudging"
1056
1057 uvar->actuationScale = + 1.0;
1058
1059#define DEFAULT_TRANSIENT "none"
1060 uvar->transientWindowType = XLALCalloc( 1, len = strlen( DEFAULT_TRANSIENT ) + 1 );
1061 XLAL_CHECK( uvar->transientWindowType != NULL, XLAL_ENOMEM, "XLALCalloc ( 1, %d ) failed.\n", len );
1062 strcpy( uvar->transientWindowType, DEFAULT_TRANSIENT );
1063
1064 // ---------- register all our user-variable ----------
1065 /* output options */
1066 XLALRegisterUvarMember( outSingleSFT, BOOLEAN, 's', OPTIONAL, "Write a single concatenated SFT (name given by --outSFTbname)" );
1067 XLALRegisterUvarMember( outSFTbname, STRING, 'n', OPTIONAL, "Output SFTs: target Directory (if --outSingleSFT=false) or filename (if --outSingleSFT=true)" );
1068
1069 XLALRegisterUvarMember( TDDfile, STRING, 't', OPTIONAL, "Basename for output of ASCII time-series" );
1070 XLALRegisterUvarMember( TDDframedir, STRING, 'F', OPTIONAL, "Directory to output frame time-series into" );
1071 XLALRegisterUvarMember( frameDesc, STRING, 0, OPTIONAL, "Description-field entry in frame filename" );
1072
1073 XLALRegisterUvarMember( logfile, STRING, 'l', OPTIONAL, "Filename for log-output" );
1074
1075 /* detector and ephemeris */
1076 XLALRegisterUvarMember( IFO, STRING, 'I', REQUIRED, "Detector: one of 'G1','L1','H1,'H2','V1', ..." );
1077
1078 XLALRegisterUvarMember( ephemEarth, STRING, 0, OPTIONAL, "Earth ephemeris file to use" );
1079 XLALRegisterUvarMember( ephemSun, STRING, 0, OPTIONAL, "Sun ephemeris file to use" );
1080
1081 /* start + duration of timeseries */
1082 XLALRegisterUvarMember( startTime, EPOCH, 'G', OPTIONAL, "Start-time of requested signal in detector-frame (format 'xx.yy[GPS|MJD]')" );
1083 XLALRegisterUvarMember( duration, INT4, 0, OPTIONAL, "Duration of requested signal in seconds" );
1084 XLALRegisterUvarMember( timestampsFile, STRING, 0, OPTIONAL, "ALTERNATIVE: File to read timestamps from (file-format: lines with <seconds> <nanoseconds>)" );
1085
1086 /* sampling and heterodyning frequencies */
1087 XLALRegisterUvarMember( fmin, REAL8, 0, REQUIRED, "Lowest frequency in output SFT (= heterodyning frequency)" );
1088 XLALRegisterUvarMember( Band, REAL8, 0, REQUIRED, "Bandwidth of output SFT in Hz (= 1/2 sampling frequency)" );
1089
1090 /* SFT properties */
1091 XLALRegisterUvarMember( Tsft, REAL8, 0, OPTIONAL, "Time baseline of one SFT in seconds" );
1092 XLALRegisterUvarMember( SFToverlap, REAL8, 0, OPTIONAL, "Overlap between successive SFTs in seconds (conflicts with --noiseSFTs or --timestampsFile)" );
1093 XLALRegisterUvarMember( window, STRING, 0, OPTIONAL, "Window function to apply to the SFTs ('rectangular', 'hann', 'tukey', etc.); when --noiseSFTs is given, required ONLY if noise SFTs have unknown window" );
1094 XLALRegisterUvarMember( windowParam, REAL8, 0, OPTIONAL, "Window parameter required for a few window-types (eg. 'tukey')" );
1095 XLALRegisterNamedUvar( NULL, "tukeyBeta", REAL8, 0, DEFUNCT, "Use " UVAR_STR( windowParam ) " instead" );
1096
1097 /* pulsar params */
1098 XLALRegisterUvarMember( refTime, EPOCH, 'S', OPTIONAL, "Pulsar SSB reference epoch: format 'xx.yy[GPS|MJD]' [default: startTime]" );
1099
1100 XLALRegisterUvarMember( Alpha, RAJ, 0, OPTIONAL, "Sky: equatorial J2000 right ascension (in radians or hours:minutes:seconds)" );
1101 XLALRegisterUvarMember( Delta, DECJ, 0, OPTIONAL, "Sky: equatorial J2000 declination (in radians or degrees:minutes:seconds)" );
1102
1103 XLALRegisterUvarMember( h0, REAL8, 0, OPTIONAL, "Overall signal-amplitude h0 (for emission at twice spin frequency only)" );
1104 XLALRegisterUvarMember( cosi, REAL8, 0, OPTIONAL, "cos(iota) of inclination-angle iota (for emission at twice spin frequency only)" );
1105 XLALRegisterUvarMember( aPlus, REAL8, 0, OPTIONAL, "ALTERNATIVE to {--h0,--cosi}: A_+ amplitude" );
1106 XLALRegisterUvarMember( aCross, REAL8, 0, OPTIONAL, "ALTERNATIVE to {--h0,--cosi}: A_x amplitude" );
1107
1108 XLALRegisterUvarMember( psi, REAL8, 0, OPTIONAL, "Polarization angle psi" );
1109 XLALRegisterUvarMember( phi0, REAL8, 0, OPTIONAL, "Initial GW phase phi" );
1110 XLALRegisterUvarMember( Freq, REAL8, 0, OPTIONAL, "Intrinsic GW-frequency at refTime" );
1111
1112 XLALRegisterUvarMember( f1dot, REAL8, 0, OPTIONAL, "First spindown parameter f' at refTime" );
1113 XLALRegisterUvarMember( f2dot, REAL8, 0, OPTIONAL, "Second spindown parameter f'' at refTime" );
1114 XLALRegisterUvarMember( f3dot, REAL8, 0, OPTIONAL, "Third spindown parameter f''' at refTime" );
1115
1116 /* binary-system orbital parameters */
1117 XLALRegisterUvarMember( orbitasini, REAL8, 0, OPTIONAL, "Projected orbital semi-major axis in seconds (a/c)" );
1118 XLALRegisterUvarMember( orbitEcc, REAL8, 0, OPTIONAL, "Orbital eccentricity" );
1119 XLALRegisterUvarMember( orbitTp, EPOCH, 0, OPTIONAL, "True epoch of periapsis passage: format 'xx.yy[GPS|MJD]'" );
1120 XLALRegisterUvarMember( orbitPeriod, REAL8, 0, OPTIONAL, "Orbital period (seconds)" );
1121 XLALRegisterUvarMember( orbitArgp, REAL8, 0, OPTIONAL, "Argument of periapsis (radians)" );
1122
1123 /* noise */
1124 XLALRegisterUvarMember( noiseSFTs, STRING, 'D', OPTIONAL, "Noise-SFTs to be added to signal. Possibilities are:\n"
1125 " - '<SFT file>;<SFT file>;...', where <SFT file> may contain wildcards\n - 'list:<file containing list of SFT files>'\n"
1126 "(Uses ONLY SFTs falling within (--startTime,--duration) or the set given in --timstampsFile, and ONLY within (--fmin,--Band).)" );
1127 XLALRegisterUvarMember( noiseSqrtSh, REAL8, 0, OPTIONAL, "Gaussian noise with single-sided PSD sqrt(Sh)" );
1128
1129 XLALRegisterUvarMember( lineFeature, BOOLEAN, 0, OPTIONAL, "Generate a monochromatic 'line' of amplitude h0 and frequency 'Freq'}" );
1130
1131 XLALRegisterUvarMember( parfile, STRING, 'p', OPTIONAL, "Directory path for optional .par files" ); /*registers .par file in mfd*/
1132
1133 /* transient signal window properties (name, start, duration) */
1134 XLALRegisterUvarMember( transientWindowType, STRING, 0, OPTIONAL, "Type of transient signal window to use. ('none', 'rect', 'exp')." );
1135 XLALRegisterUvarMember( transientStartTime, REAL8, 0, OPTIONAL, "GPS start-time 't0' of transient signal window." );
1136 XLALRegisterUvarMember( transientTauDays, REAL8, 0, OPTIONAL, "Timescale 'tau' of transient signal window in days." );
1137
1138 /* ----- 'expert-user/developer' options ----- */
1139 XLALRegisterUvarMember( sourceDeltaT, REAL8, 0, DEVELOPER, "Source-frame sampling period. '0' implies previous internal defaults" );
1140 XLALRegisterUvarMember( generationMode, INT4, 0, DEVELOPER, "How to generate timeseries: 0=all-at-once (faster), 1=per-sft (slower)" );
1141
1142 XLALRegisterUvarMember( hardwareTDD, BOOLEAN, 'b', DEVELOPER, "Hardware injection: output TDD in binary format (implies generationMode=1)" );
1143 XLALRegisterUvarMember( actuation, STRING, 0, DEVELOPER, "Filname containing actuation function of this detector" );
1144 XLALRegisterUvarMember( actuationScale, REAL8, 0, DEVELOPER, "(Signed) scale-factor to apply to the actuation-function." );
1145
1146 XLALRegisterUvarMember( exactSignal, BOOLEAN, 0, DEVELOPER, "Generate signal time-series as exactly as possible (slow)." );
1147 XLALRegisterUvarMember( randSeed, INT4, 0, DEVELOPER, "Specify random-number seed for reproducible noise (0 means use /dev/urandom for seeding)." );
1148
1149 /* read cmdline & cfgfile */
1150 BOOLEAN should_exit = 0;
1152 if ( should_exit ) {
1153 exit( 1 );
1154 }
1155
1156 return XLAL_SUCCESS;
1157
1158} /* XLALInitUserVars() */
1159
1160
1161/**
1162 * This routine frees up all the memory.
1163 */
1164int
1166{
1167 XLAL_CHECK( cfg != NULL, XLAL_EINVAL );
1168
1169 /* Free config-Variables and userInput stuff */
1171
1172 /* free timestamps if any */
1174
1175 /* free window if any */
1177 XLALFree( cfg->window_type );
1178
1179 /* free spindown-vector (REAL8) */
1181
1182 /* free noise-SFTs */
1184
1185 /* free transfer-function if we have one.. */
1187
1188 XLALFree( cfg->VCSInfoString );
1189
1190 /* Clean up earth/sun Ephemeris tables */
1192
1193 return XLAL_SUCCESS;
1194
1195} /* XLALFreeMem() */
1196
1197/**
1198 * Log the all relevant parameters of this run into a log-file.
1199 * The name of the log-file used is uvar_logfile
1200 * <em>NOTE:</em> Currently this function only logs the user-input and code-versions.
1201 */
1202int
1203XLALWriteMFDlog( const char *logfile, const ConfigVars_t *cfg )
1204{
1205 XLAL_CHECK( logfile != NULL, XLAL_EINVAL, "Invalid NULL input 'logfile'\n" );
1206 XLAL_CHECK( cfg != NULL, XLAL_EINVAL, "Invalid NULL input 'cfg'\n" );
1207
1208 FILE *fplog = fopen( logfile, "wb" );
1209 XLAL_CHECK( fplog != NULL, XLAL_EIO, "Failed to fopen log-file '%s' for writing ('wb').\n", logfile );
1210
1211 /* write out a log describing the complete user-input (in cfg-file format) */
1213 XLAL_CHECK( logstr != NULL, XLAL_EFUNC, "XLALUserVarGetLog(UVAR_LOGFMT_CFGFILE) failed.\n" );
1214
1215 fprintf( fplog, "## LOG-FILE of Makefakedata run\n\n" );
1216 fprintf( fplog, "# User-input: [formatted as config-file]\n" );
1217 fprintf( fplog, "# ----------------------------------------------------------------------\n\n" );
1218 fprintf( fplog, "%s", logstr );
1219 XLALFree( logstr );
1220
1221 /* write out a log describing the complete user-input (in commandline format) */
1223 XLAL_CHECK( logstr != NULL, XLAL_EFUNC, "XLALUserVarGetLog(UVAR_LOGFMT_CMDLINE) failed.\n" );
1224
1225 fprintf( fplog, "\n\n# User-input: [formatted as commandline]\n" );
1226 fprintf( fplog, "# ----------------------------------------------------------------------\n\n" );
1227 fprintf( fplog, "%s", logstr );
1228 XLALFree( logstr );
1229
1230 /* append an VCS-version string of the code used */
1231 fprintf( fplog, "\n\n# VCS-versions of executable:\n" );
1232 fprintf( fplog, "# ----------------------------------------------------------------------\n" );
1233 fprintf( fplog, "\n%s\n", cfg->VCSInfoString );
1234 fclose( fplog );
1235
1236 return XLAL_SUCCESS;
1237
1238} /* XLALWriteMFDLog() */
1239
1240/**
1241 * Reads an actuation-function in format (r,phi) from file 'fname',
1242 * and returns the associated transfer-function as a COMPLEX8FrequencySeries (Re,Im)
1243 * The transfer-function T is simply the inverse of the actuation A, so T=A^-1.
1244 */
1246XLALLoadTransferFunctionFromActuation( REAL8 actuationScale, /**< overall scale-factor to actuation */
1247 const CHAR *fname /**< file containing actuation-function */
1248 )
1249{
1250 int len;
1251
1252 XLAL_CHECK_NULL( fname != NULL, XLAL_EINVAL, "Invalid NULL input 'fname'\n" );
1253
1254 LALParsedDataFile *fileContents = NULL;
1255 XLAL_CHECK_NULL( XLALParseDataFile( &fileContents, fname ) == XLAL_SUCCESS, XLAL_EFUNC );
1256
1257
1258 /* skip first line if containing NaN's ... */
1259 UINT4 startline = 0;
1260 if ( strstr( fileContents->lines->tokens[startline], "NaN" ) != NULL ) {
1261 startline ++;
1262 }
1263
1264 UINT4 numLines = fileContents->lines->nTokens - startline;
1266 XLAL_CHECK_NULL( data != NULL, XLAL_EFUNC, "XLALCreateCOMPLEX8Vector(%d) failed\n", numLines );
1267
1268 COMPLEX8FrequencySeries *ret = XLALCalloc( 1, len = sizeof( *ret ) );
1269 XLAL_CHECK_NULL( ret != NULL, XLAL_ENOMEM, "Failed to XLALCalloc (1, %d)\n", len );
1270
1271 snprintf( ret->name, LALNameLength - 1, "Transfer-function from: %s", fname );
1272 ret->name[LALNameLength - 1] = 0;
1273
1274 /* initialize loop */
1275 REAL8 f0 = 0;
1276 REAL8 f1 = 0;
1277
1278 const CHAR *readfmt = "%" LAL_REAL8_FORMAT "%" LAL_REAL8_FORMAT "%" LAL_REAL8_FORMAT;
1279
1280 for ( UINT4 i = startline; i < fileContents->lines->nTokens; i++ ) {
1281 CHAR *thisline = fileContents->lines->tokens[i];
1282 REAL8 amp, phi;
1283
1284 f0 = f1;
1285 if ( 3 != sscanf( thisline, readfmt, &f1, &amp, &phi ) ) {
1286 XLAL_ERROR_NULL( XLAL_EIO, "Failed to read 3 floats from line %d of file %s\n\n", i, fname );
1287 }
1288
1289 if ( !gsl_finite( amp ) || !gsl_finite( phi ) ) {
1290 XLAL_ERROR_NULL( XLAL_EINVAL, "ERROR: non-finite number encountered in actuation-function at f!=0. Line=%d\n\n", i );
1291 }
1292
1293 /* first line: set f0 */
1294 if ( i == startline ) {
1295 ret->f0 = f1;
1296 }
1297
1298 /* second line: set deltaF */
1299 if ( i == startline + 1 ) {
1300 ret->deltaF = f1 - ret->f0;
1301 }
1302
1303 /* check constancy of frequency-step */
1304 if ( ( f1 - f0 != ret->deltaF ) && ( i > startline ) ) {
1305 XLAL_ERROR_NULL( XLAL_EINVAL, "Illegal frequency-step %f != %f in line %d of file %s\n\n", ( f1 - f0 ), ret->deltaF, i, fname );
1306 }
1307
1308 /* now convert into transfer-function and (Re,Im): T = A^-1 */
1309 data->data[i - startline] = crectf( cos( phi ) / ( amp * actuationScale ), -sin( phi ) / ( amp * actuationScale ) );
1310
1311 } /* for i < numlines */
1312
1313 XLALDestroyParsedDataFile( fileContents );
1314
1315 ret->data = data;
1316
1317 return ret;
1318
1319} /* XLALLoadTransferFunctionFromActuation() */
1320
1321/* determine if the given filepath is an existing directory or not */
1322BOOLEAN
1323is_directory( const CHAR *fname )
1324{
1325 struct stat stat_buf;
1326
1327 if ( stat( fname, &stat_buf ) ) {
1328 return 0;
1329 }
1330
1331 if ( ! S_ISDIR( stat_buf.st_mode ) ) {
1332 return 0;
1333 } else {
1334 return 1;
1335 }
1336
1337} /* is_directory() */
#define IFO
int k
void LALCheckMemoryLeaks(void)
#define LALCalloc(m, n)
#define LALFree(p)
const LALVCSInfoList lalPulsarVCSInfoList
NULL-terminated list of VCS and build information for LALPulsar and its dependencies
ConfigVariables GV
global container for various derived configuration settings
Definition: PredictFstat.c:70
#define STRING(a)
void XLALReadTEMPOParFileOrig(BinaryPulsarParams *output, CHAR *pulsarAndPath)
function to read in a TEMPO parameter file
const char * comment
Definition: SearchTiming.c:94
#define fprintf
double e
const double w
#define __attribute__(x)
int XLALParseDataFile(LALParsedDataFile **cfgdata, const CHAR *fname)
void XLALDestroyParsedDataFile(LALParsedDataFile *cfgdata)
char * XLALGPSToStr(char *s, const LIGOTimeGPS *t)
void XLALDestroyCOMPLEX8FrequencySeries(COMPLEX8FrequencySeries *series)
REAL4TimeSeries * XLALGeneratePulsarSignal(const PulsarSignalParams *params)
Generate a time-series at the detector for a given pulsar.
REAL4TimeSeries * XLALGenerateLineFeature(const PulsarSignalParams *params)
Generate a REAL4TimeSeries containing a sinusoid with amplitude 'h0', frequency 'Freq-fHeterodyne' an...
SFTVector * XLALSignalToSFTs(const REAL4TimeSeries *signalvec, const SFTParams *params)
Turn the given time-series into SFTs and add noise if given.
int XLALAddGaussianNoise(REAL4TimeSeries *inSeries, REAL4 sigma, INT4 seed)
Generate Gaussian noise with standard-deviation sigma, add it to inSeries.
EphemerisData * XLALInitBarycenter(const CHAR *earthEphemerisFile, const CHAR *sunEphemerisFile)
XLAL interface to reading ephemeris files 'earth' and 'sun', and return ephemeris-data in old backwar...
void XLALDestroyEphemerisData(EphemerisData *edat)
Destructor for EphemerisData struct, NULL robust.
#define LAL_REAL8_EPS
#define LAL_TWOPI
unsigned char BOOLEAN
double REAL8
#define XLAL_INIT_DECL(var,...)
#define crectf(re, im)
char CHAR
uint32_t UINT4
int32_t INT4
float REAL4
LALNameLength
#define lalDebugLevel
LALFrameUFrameH LALFrameH
int XLALFrameAddFrHistory(LALFrameH *frame, const char *name, const char *comment)
int XLALFrameWrite(LALFrameH *frame, const char *fname)
LALFrameH * XLALFrameNew(const LIGOTimeGPS *epoch, double duration, const char *project, int run, int frnum, INT8 detectorFlags)
int XLALFrameAddREAL4TimeSeriesProcData(LALFrameH *frame, const REAL4TimeSeries *series)
void XLALFrameFree(LALFrameH *frame)
void * XLALCalloc(size_t m, size_t n)
void XLALFree(void *p)
#define LAL_REAL8_FORMAT
char char * XLALStringDuplicate(const char *s)
int XLALStringCaseCompare(const char *s1, const char *s2)
char * XLALVCSInfoString(const LALVCSInfoList vcs_list, const int verbose, const char *prefix)
int XLALdumpREAL4TimeSeries(const char *fname, const REAL4TimeSeries *series)
const LALDetector * XLALGetSiteInfo(const CHAR *name)
Find the site geometry-information 'LALDetector' for given a detector name (or prefix).
Definition: SFTnaming.c:218
void XLALDestroySFTVector(SFTVector *vect)
XLAL interface to destroy an SFTVector.
Definition: SFTtypes.c:300
void XLALDestroySFTCatalog(SFTCatalog *catalog)
Free an 'SFT-catalogue'.
Definition: SFTcatalog.c:329
LIGOTimeGPSVector * XLALMakeTimestamps(LIGOTimeGPS tStart, REAL8 Tspan, REAL8 Tsft, REAL8 Toverlap)
Given a start-time, Tspan, Tsft and Toverlap, returns a list of timestamps covering this time-stretch...
LIGOTimeGPSVector * XLALExtractTimestampsFromSFTs(const SFTVector *sfts)
Extract timstamps-vector from the given SFTVector.
CHAR * XLALGetChannelPrefix(const CHAR *name)
Find the valid CW detector prefix.
Definition: SFTnaming.c:203
int XLALFindCoveringSFTBins(UINT4 *firstBin, UINT4 *numBins, REAL8 fMinIn, REAL8 BandIn, REAL8 Tsft)
Return the 'effective' frequency-band [fMinEff, fMaxEff] = [firstBin, lastBin] * 1/Tsft,...
Definition: SFTtypes.c:106
int XLALWriteSFT2FilePointer(const SFTtype *sft, FILE *fp, const CHAR *SFTwindowtype, const REAL8 SFTwindowparam, const CHAR *SFTcomment)
Write the given SFTtype to a FILE pointer.
Definition: SFTfileIO.c:488
int XLALCompareSFTWindows(const CHAR *type1, const REAL8 param1, const CHAR *type2, const REAL8 param2)
Check whether two SFT windows, each defined by a type name and parameter value, match.
Definition: SFTnaming.c:667
SFTVector * XLALLoadSFTs(const SFTCatalog *catalog, REAL8 fMin, REAL8 fMax)
Load the given frequency-band [fMin, fMax) (half-open) from the SFT-files listed in the SFT-'catalogu...
Definition: SFTfileIO.c:87
SFTCatalog * XLALSFTdataFind(const CHAR *file_pattern, const SFTConstraints *constraints)
Find the list of SFTs matching the file_pattern and satisfying the given constraints,...
Definition: SFTcatalog.c:71
int XLALCheckValidDescriptionField(const char *desc)
Check whether given string qualifies as a valid 'description' field of a FRAME filename (per ) or SFT...
Definition: SFTnaming.c:639
int XLALWriteSFTVector2StandardFile(const SFTVector *sftVect, SFTFilenameSpec *SFTfnspec, const CHAR *SFTcomment, const BOOLEAN merged)
Write the given SFTVector to SFT file(s) with a standard () filename(s).
Definition: SFTfileIO.c:755
LIGOTimeGPSVector * XLALReadTimestampsFile(const CHAR *fname)
backwards compatible wrapper to XLALReadTimestampsFileConstrained() without GPS-time constraints
int XLALFillSFTFilenameSpecStrings(SFTFilenameSpec *spec, const CHAR *path, const CHAR *extn, const CHAR *detector, const CHAR *window_type, const CHAR *privMisc, const CHAR *pubObsKind, const CHAR *pubChannel)
Convenience function for filling out the string fields in a SFTFilenameSpec.
Definition: SFTnaming.c:257
SFTVector * XLALExtractStrictBandFromSFTVector(const SFTVector *inSFTs, REAL8 fMin, REAL8 Band)
Return a copy of a vector of SFTs containing only the bins in [fMin, fMin+Band).
Definition: SFTtypes.c:658
void XLALDestroyTimestampVector(LIGOTimeGPSVector *vect)
De-allocate a LIGOTimeGPSVector.
Definition: SFTtimestamps.c:69
REAL4TimeSeries * XLALSimulateExactPulsarSignal(const PulsarSignalParams *params)
Simulate a pulsar signal to best accuracy possible.
COORDINATESYSTEM_EQUATORIAL
void XLALDestroyREAL4TimeSeries(REAL4TimeSeries *series)
int XLALParseTransientWindowName(const char *windowName)
Parse a transient window name string into the corresponding transientWindowType.
int XLALApplyTransientWindow(REAL4TimeSeries *series, transientWindow_t transientWindow)
apply a "transient CW window" described by TransientWindowParams to the given timeseries
int XLALUserVarReadAllInput(BOOLEAN *should_exit, int argc, char *argv[], const LALVCSInfoList vcs_list)
void XLALDestroyUserVars(void)
#define XLALRegisterUvarMember(name, type, option, category,...)
#define UVAR_STR(n)
void CHAR * XLALUserVarGetLog(UserVarLogFormat format)
#define XLALRegisterNamedUvar(cvar, name, type, option, category,...)
int XLALUserVarWasSet(const void *cvar)
UVAR_LOGFMT_CMDLINE
UVAR_LOGFMT_CFGFILE
COMPLEX8Vector * XLALCreateCOMPLEX8Vector(UINT4 length)
REAL8Vector * XLALCreateREAL8Vector(UINT4 length)
void XLALDestroyREAL8Vector(REAL8Vector *vector)
REAL4Window * XLALCreateNamedREAL4Window(const char *windowName, REAL8 beta, UINT4 length)
int XLALCheckNamedWindow(const char *windowName, const BOOLEAN haveBeta)
void XLALDestroyREAL4Window(REAL4Window *window)
#define XLAL_ERROR_NULL(...)
#define xlalErrno
#define XLAL_ERROR(...)
#define XLAL_CHECK(assertion,...)
int int XLALPrintWarning(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_CHECK_NULL(assertion,...)
XLAL_ENOMEM
XLAL_SUCCESS
XLAL_EFUNC
XLAL_EDOM
XLAL_ETYPE
XLAL_ESIZE
XLAL_EIO
XLAL_EINVAL
XLAL_EFAILED
LIGOTimeGPS * XLALGPSSetREAL8(LIGOTimeGPS *epoch, REAL8 t)
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
REAL8 XLALGPSDiff(const LIGOTimeGPS *t1, const LIGOTimeGPS *t0)
const char * actuation
Definition: hwinject.c:374
float data[BLOCKSIZE]
Definition: hwinject.c:360
int main(int argc, char *argv[])
int XLALIsValidDescriptionField(const char *desc)
GenerationMode
@ GENERATE_ALL_AT_ONCE
generate whole time-series at once before turning into SFTs
@ GENERATE_PER_SFT
generate timeseries individually for each SFT
@ GENERATE_LAST
end-marker
BOOLEAN is_directory(const CHAR *fname)
COMPLEX8FrequencySeries * XLALLoadTransferFunctionFromActuation(REAL8 actuationScale, const CHAR *fname)
Reads an actuation-function in format (r,phi) from file 'fname', and returns the associated transfer-...
int XLALInitUserVars(UserVariables_t *uvar, int argc, char *argv[])
Register all "user-variables", and initialized them from command-line and config-files.
int XLALInitMakefakedata(ConfigVars_t *cfg, UserVariables_t *uvar)
Handle user-input and set up shop accordingly, and do all consistency-checks on user-input.
int XLALWriteMFDlog(const char *logfile, const ConfigVars_t *cfg)
Log the all relevant parameters of this run into a log-file.
#define DEFAULT_TRANSIENT
#define SQ(x)
int XLALFreeMem(ConfigVars_t *cfg)
This routine frees up all the memory.
ts
double t0
A structure to contain all pulsar parameters and associated errors.
REAL8 w0
logitude of periastron (radians)
REAL8 pepoch
period/frequency epoch
REAL8 x
projected semi-major axis/speed of light (light secs)
REAL8 e
orbital eccentricity
REAL8 f2
frequency second derivative (Hz/s^2)
REAL8 Aplus
0.5*h0*(1+cos^2iota)
REAL8 Across
h0*cosiota
REAL8 ra
right ascension (rads)
REAL8 phi0
initial phase
REAL8 f3
frequency third derivative (Hz/s^3)
REAL8 T0
time of orbital perisastron as measured in TDB (MJD)
REAL8 cosiota
cosine of the pulsars inclination angle
REAL8 h0
gravitational wave amplitude
REAL8 Pb
orbital period (seconds)
CHAR * model
TEMPO binary model e.g.
REAL8 Tasc
time of the ascending node (used rather than T0)
REAL8 posepoch
position epoch
REAL8 dec
declination (rads)
REAL8 f0
spin frequency (Hz)
REAL8 f1
frequency first derivative (Hz/s)
REAL8 psi
polarisation angle
CHAR name[LALNameLength]
COMPLEX8Sequence * data
A vector of COMPLEX8FrequencySeries.
COMPLEX8FrequencySeries * data
Pointer to the data array.
UINT4 length
Number of elements in array.
CHAR * VCSInfoString
Git version string.
REAL8 refTime
reference time for pulsar phase definition
REAL8 Alpha
sky position alpha in radians
LIGOTimeGPSVector * timestamps
timestamps vector holding 1 element: timeGPS
REAL8 Delta
sky position delta in radians
EphemerisData * edat
ephemeris data
configuration-variables derived from user-variables
REAL8Vector * spindown
vector of frequency-derivatives of GW signal
UINT4 duration
total duration of observation in seconds
REAL8 noiseSigma
sigma for Gaussian noise to be added
REAL8 window_param
parameter of window function
PulsarParams pulsar
pulsar signal-parameters (amplitude + doppler
REAL8 fmin_eff
'effective' fmin: round down such that fmin*Tsft = integer!
CHAR * VCSInfoString
Git version string.
SFTVector * noiseSFTs
vector of noise-SFTs to be added to signal
LALDetector site
detector-site info
CHAR * window_type
name of window function
transientWindow_t transientWindow
properties of transient-signal window
LIGOTimeGPSVector * timestamps
a vector of timestamps to generate time-series/SFTs for
REAL4Window * window
window function for the time series
COMPLEX8FrequencySeries * transfer
detector's transfer function for use in hardware-injection
REAL8 fBand_eff
'effective' frequency-band such that samples/SFT is integer
LIGOTimeGPS startTimeGPS
start-time of observation
EphemerisData * edat
ephemeris-data
This structure contains all information about the center-of-mass positions of the Earth and Sun,...
TokenList * lines
INT4 gpsNanoSeconds
A vector of 'timestamps' of type LIGOTimeGPS.
Definition: SFTfileIO.h:188
LIGOTimeGPS * data
array of timestamps
Definition: SFTfileIO.h:193
UINT4 length
number of timestamps
Definition: SFTfileIO.h:192
REAL8 aCross
Signal amplitude (cross)
REAL8 psi
polarization angle psi
REAL8 aPlus
Signal amplitude (plus)
REAL8 phi0
initial signal-phase (at some reference time)
REAL8 period
Binary: orbital period (sec)
LIGOTimeGPS tp
Binary: time of observed periapsis passage (in SSB)
PulsarSpins fkdot
Intrinsic spins: [Freq, f1dot, f2dot, ... ] where fkdot = d^kFreq/dt^k.
REAL8 Delta
Sky position: DEC (latitude) in equatorial coords and radians.
LIGOTimeGPS refTime
Reference time of pulsar parameters (in SSB!)
REAL8 Alpha
Sky position: RA (longitude) in equatorial coords and radians.
REAL8 ecc
Binary: orbital eccentricity.
REAL8 asini
Binary: projected, normalized orbital semi-major axis (s).
REAL8 argp
Binary: argument of periapsis (radians)
Type defining the parameters of a pulsar-source of CW Gravitational waves.
PulsarAmplitudeParams Amp
'Amplitude-parameters': h0, cosi, phi0, psi
PulsarDopplerParams Doppler
'Phase-evolution parameters': {skypos, fkdot, orbital params }
Input parameters to GeneratePulsarSignal(), defining the source and the time-series.
CHAR name[LALNameLength]
REAL4Sequence * data
REAL4 * data
REAL8 * data
An "SFT-catalogue": a vector of SFTdescriptors, as returned by XLALSFTdataFind()
Definition: SFTfileIO.h:238
SFTDescriptor * data
array of data-entries describing matched SFTs
Definition: SFTfileIO.h:243
UINT4 length
number of SFTs in catalog
Definition: SFTfileIO.h:242
'Constraints' for SFT-matching: which detector, within which time-stretch and which timestamps exactl...
Definition: SFTfileIO.h:212
const char * window_type
window function applied to SFT
Definition: SFTfileIO.h:229
REAL8 window_param
parameter of window function, if required
Definition: SFTfileIO.h:230
Structure specifying an SFT file name, following the convention in .
Definition: SFTfileIO.h:266
Parameters defining the SFTs to be returned from LALSignalToSFTs().
CHAR ** tokens
UINT4 nTokens
user input variables
Definition: compareFstats.c:51
REAL8 Alpha
Right ascension [radians] alpha of pulsar.
REAL8 transientStartTime
GPS start-time of transient window.
CHAR * transientWindowType
option .par file path
INT4 duration
Duration of requested signal in seconds.
REAL8 orbitArgp
Argument of periapsis (radians)
REAL8 psi
Polarization angle psi.
CHAR * ephemSun
Sun ephemeris file to use.
REAL8 h0
overall signal amplitude h0
CHAR * frameDesc
description field entry in the frame filename
LIGOTimeGPS refTime
Pulsar reference epoch tRef in SSB ('0' means: use startTime converted to SSB)
INT4 randSeed
allow user to specify random-number seed for reproducible noise-realizations
CHAR * ephemEarth
Earth ephemeris file to use.
REAL8 transientTauDays
time-scale in days of transient window
REAL8 phi0
Initial phase phi.
REAL8 f2dot
Second spindown parameter f''.
CHAR * window
Windowing function for the time series.
LIGOTimeGPS orbitTp
'true' epoch of periapsis passage
REAL8 aCross
ALTERNATIVE to {h0,cosi}: Cross polarization amplitude aCross.
REAL8 Freq
physical start frequency to compute PSD for (excluding rngmed wings)
REAL8 fmin
Lowest frequency in output SFT (= heterodyning frequency)
REAL8 Band
bandwidth of output SFT in Hz (= 1/2 sampling frequency)
REAL8 noiseSqrtSh
single-sided sqrt(Sh) for Gaussian noise
REAL8 aPlus
ALTERNATIVE to {h0,cosi}: Plus polarization amplitude aPlus.
BOOLEAN exactSignal
generate signal timeseries as exactly as possible (slow)
BOOLEAN hardwareTDD
Binary output timeseries in chunks of Tsft for hardware injections.
LIGOTimeGPS startTime
Start-time of requested signal in detector-frame (GPS seconds)
REAL8 Delta
Declination [radians] delta of pulsar.
REAL8 SFToverlap
overlap SFTs by this many seconds
REAL8 orbitasini
Projected orbital semi-major axis in seconds (a/c)
CHAR * timestampsFile
Timestamps file.
CHAR * logfile
name of logfile
REAL8 sourceDeltaT
source-frame sampling period.
CHAR * outSFTbname
Path and basefilename of output SFT files.
REAL8 f1dot
First spindown parameter f'.
REAL8 actuationScale
Scale-factor to be applied to actuation-function.
CHAR * actuation
filename containg detector actuation function
REAL8 Tsft
SFT time baseline Tsft.
REAL8 orbitEcc
Orbital eccentricity.
BOOLEAN outSingleSFT
use to output a single concatenated SFT
BOOLEAN lineFeature
generate a monochromatic line instead of a pulsar-signal
CHAR * noiseSFTs
Glob-like pattern specifying noise-SFTs to be added to signal.
REAL8 f3dot
Third spindown parameter f'''.
INT4 generationMode
How to generate the timeseries: all-at-once or per-sft.
REAL8 orbitPeriod
Orbital period (seconds)
REAL8 cosi
cos(iota) of inclination angle iota
CHAR * TDDframedir
directory for frame file output time-series
REAL8 windowParam
Hann fraction of Tukey window (0.0=rect; 1,0=han; 0.5=default.
CHAR * TDDfile
Filename for ASCII output time-series.
CHAR * IFO
Detector: H1, L1, H2, V1, ...
double duration
double psi
Struct defining one transient window instance.
UINT4 t0
GPS start-time 't0'.
UINT4 tau
transient timescale tau in seconds
transientWindowType_t type
window-type: none, rectangular, exponential, ....
enum @4 site