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>
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>
64#ifdef HAVE_LIBLALFRAME
65#include <lal/LALFrameIO.h>
69#define SQ(x) ( (x) * (x) )
218 UINT4 i_chunk, numchunks;
219 FILE *fpSingleSFT = NULL;
239 params.pulsar.aPlus =
GV.pulsar.Amp.aPlus;
240 params.pulsar.aCross =
GV.pulsar.Amp.aCross;
241 params.pulsar.phi0 =
GV.pulsar.Amp.phi0;
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;
252 params.sourceDeltaT = uvar.sourceDeltaT;
260 if ( ! uvar.exactSignal ) {
261 params.samplingRate = 2.0 *
GV.fBand_eff;
264 params.samplingRate =
fmax( 2.0 * (
params.pulsar.f0 + 2 ), 2 * (
GV.fmin_eff +
GV.fBand_eff ) );
269 switch ( uvar.generationMode ) {
284 if ( uvar.outSingleSFT ) {
287 XLAL_ERROR(
XLAL_ETYPE,
"'%s' is a directory, but --outSingleSFT expects a filename!\n", uvar.outSFTbname );
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 ) );
300 for ( i_chunk = 0; i_chunk < numchunks; i_chunk++ ) {
306 if ( uvar.exactSignal ) {
308 }
else if ( uvar.lineFeature ) {
319 if ( uvar.hardwareTDD && ( i_chunk == 0 ) ) {
320 REAL4 magic = 1234.5;
322 if ( ( 1 != fwrite( &magic,
sizeof( magic ), 1, stdout ) ) || ( 1 != fwrite( &length,
sizeof(
INT4 ), 1, stdout ) ) ) {
323 perror(
"Failed to write to stdout" );
330 if (
GV.noiseSigma > 0 ) {
336 if ( uvar.TDDfile ) {
339 sprintf( fname,
"%s.%02d", uvar.TDDfile, i_chunk );
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" );
351 len = strlen( uvar.TDDframedir ) + strlen( uvar.frameDesc ) + 100;
353 char IFO[2] = { Tseries->
name[0], Tseries->
name[1] };
355 size_t written = snprintf( fname, len,
"%s/%c-%c%c_%s-%d-%d.gwf",
357 XLAL_CHECK( written < len,
XLAL_ESIZE,
"Frame-filename exceeds expected maximal length (%zu): '%s'\n", len, fname );
365 written = snprintf( buffer,
LALNameLength,
"%s:%s", Tseries->
name, uvar.frameDesc );
367 strcpy( Tseries->
name, buffer );
390 if ( uvar.hardwareTDD ) {
394 if ( length != fwrite( datap,
sizeof( datap[0] ), length, stdout ) ) {
395 perror(
"Fatal error in writing binary data to stdout\n" );
407 if ( uvar.outSFTbname ) {
412 sftParams.Tsft = uvar.Tsft;
414 switch ( uvar.generationMode ) {
417 sftParams.noiseSFTs =
GV.noiseSFTs;
422 sftParams.timestamps = &(
ts );
423 sftParams.noiseSFTs = NULL;
425 if (
GV.noiseSFTs ) {
427 noise.
data = &(
GV.noiseSFTs->data[i_chunk] );
428 sftParams.noiseSFTs = &( noise );
439 sftParams.window =
GV.window;
462 spec.window_param =
GV.window_param;
463 if ( uvar.outSingleSFT ) {
493 if ( uvar.outSingleSFT ) {
496 fclose( fpSingleSFT );
523 char *window_type_from_noiseSFTs = NULL;
524 REAL8 window_param_from_noiseSFTs = 0;
527 if ( have_parfile ) {
548 if ( have_parfile ) {
549 if ( pulparams.
h0 != 0 ) {
550 uvar->
h0 = pulparams.
h0;
553 uvar->
psi = pulparams.
psi;
560 uvar->
psi = pulparams.
psi;
564 uvar->
Freq = 2.*pulparams.
f0;
572 if ( ( have_aPlus || have_aCross ) && ( have_h0 || have_cosi ) ) {
575 if ( ( have_h0 ^ have_cosi ) ) {
578 if ( ( have_aPlus ^ have_aCross ) ) {
581 if ( have_h0 && have_cosi ) {
588 }
else if ( have_aPlus && have_aCross ) {
602 if ( ( have_Alpha && !have_Delta ) || ( !have_Alpha && have_Delta ) ) {
613 if ( have_parfile ) {
622 if ( uvar->
f3dot != 0 ) {
624 }
else if ( uvar->
f2dot != 0 ) {
626 }
else if ( uvar->
f1dot != 0 ) {
639#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
644#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
658 CHAR *channelName = NULL;
668 cfg->
site = ( *site );
696 printf(
"\nWARNING: for SFT-creation we had to adjust (fmin,Band) to"
711 BOOLEAN haveStart, haveDuration, haveTimestampsFile, haveOverlap;
718 if ( ( haveDuration && !haveStart ) || ( !haveDuration && haveStart ) ) {
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" );
731 XLAL_ERROR(
XLAL_EINVAL,
"--SFToverlap can only be used with --generationMode=0, otherwise we'll get overlapping independent noise!\n\n" );
737 if ( haveTimestampsFile || uvar->
noiseSFTs ) {
738 XLAL_ERROR(
XLAL_EINVAL,
"\nHardware injection: don't specify --timestampsFile or --noiseSFTs\n\n" );
740 if ( !haveStart || !haveDuration ) {
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" );
770 XLALPrintWarning(
"\nWARNING: only SFTs corresponding to the noiseSFT-timestamps will be produced!\n" );
773 if ( haveStart && haveDuration ) {
774 minStartTime = maxStartTime = uvar->
startTime;
776 constraints.minStartTime = &minStartTime;
777 constraints.maxStartTime = &maxStartTime;
778 char bufGPS1[32], bufGPS2[32];
787 constraints.detector = channelName ;
793 if ( catalog->
length == 0 ) {
794 XLAL_ERROR(
XLAL_EFAILED,
"No noise-SFTs matching the constraints (IFO, start+duration, timestamps) were found!\n" );
822 if ( !haveStart || !haveDuration ) {
823 XLAL_ERROR(
XLAL_EINVAL,
"Need to have either --timestampsFile OR (--startTime,--duration) OR --noiseSFTs\n\n" );
889 const CHAR *window_type_from_uvar = uvar->
window;
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 );
901 }
else if ( have_window_type_from_noiseSFTs ) {
904 }
else if ( have_window ) {
909 XLAL_ERROR(
XLAL_EINVAL,
"When --noiseSFTs is given and have unknown window, --window is required.\n" );
911 }
else if ( have_window ) {
932 if ( have_parfile ) {
933 if ( pulparams.
model != NULL ) {
936 if ( strstr( pulparams.
model,
"ELL1" ) != NULL ) {
938 eps1 = pulparams.
eps1;
939 eps2 = pulparams.
eps2;
940 w = atan2( eps1, eps2 );
941 e = sqrt( eps1 * eps1 + eps2 * eps2 );
948 if ( strstr( pulparams.
model,
"ELL1" ) != NULL ) {
951 uasc = 2.0 * atan( fe * tan( uvar->
orbitArgp / 2.0 ) );
953 pulparams.
T0 = pulparams.
Tasc + Dt;
965 if ( set1 || set2 || set3 || set4 || set5 ) {
966 if ( ( uvar->
orbitasini > 0 ) && !( set1 && set2 && set3 && set4 && set5 ) ) {
1007 XLAL_ERROR(
XLAL_EINVAL,
"Error: use of an actuation/transfer function restricted to hardare-injections\n\n" );
1030 XLALFree( window_type_from_noiseSFTs );
1059#define DEFAULT_TRANSIENT "none"
1067 XLALRegisterUvarMember( outSFTbname,
STRING,
'n', OPTIONAL,
"Output SFTs: target Directory (if --outSingleSFT=false) or filename (if --outSingleSFT=true)" );
1084 XLALRegisterUvarMember( timestampsFile,
STRING, 0, OPTIONAL,
"ALTERNATIVE: File to read timestamps from (file-format: lines with <seconds> <nanoseconds>)" );
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" );
1098 XLALRegisterUvarMember( refTime, EPOCH,
'S', OPTIONAL,
"Pulsar SSB reference epoch: format 'xx.yy[GPS|MJD]' [default: startTime]" );
1119 XLALRegisterUvarMember( orbitTp, EPOCH, 0, OPTIONAL,
"True epoch of periapsis passage: format 'xx.yy[GPS|MJD]'" );
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).)" );
1134 XLALRegisterUvarMember( transientWindowType,
STRING, 0, OPTIONAL,
"Type of transient signal window to use. ('none', 'rect', 'exp')." );
1140 XLALRegisterUvarMember( generationMode,
INT4, 0, DEVELOPER,
"How to generate timeseries: 0=all-at-once (faster), 1=per-sft (slower)" );
1142 XLALRegisterUvarMember( hardwareTDD,
BOOLEAN,
'b', DEVELOPER,
"Hardware injection: output TDD in binary format (implies generationMode=1)" );
1152 if ( should_exit ) {
1208 FILE *fplog = fopen( logfile,
"wb" );
1209 XLAL_CHECK( fplog != NULL,
XLAL_EIO,
"Failed to fopen log-file '%s' for writing ('wb').\n", logfile );
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 );
1225 fprintf( fplog,
"\n\n# User-input: [formatted as commandline]\n" );
1226 fprintf( fplog,
"# ----------------------------------------------------------------------\n\n" );
1227 fprintf( fplog,
"%s", logstr );
1231 fprintf( fplog,
"\n\n# VCS-versions of executable:\n" );
1232 fprintf( fplog,
"# ----------------------------------------------------------------------\n" );
1259 UINT4 startline = 0;
1260 if ( strstr( fileContents->
lines->
tokens[startline],
"NaN" ) != NULL ) {
1285 if ( 3 != sscanf( thisline, readfmt, &
f1, &, &phi ) ) {
1289 if ( !gsl_finite( amp ) || !gsl_finite( phi ) ) {
1294 if (
i == startline ) {
1299 if (
i == startline + 1 ) {
1304 if ( (
f1 -
f0 != ret->
deltaF ) && (
i > startline ) ) {
1309 data->data[
i - startline] =
crectf( cos( phi ) / ( amp * actuationScale ), -sin( phi ) / ( amp * actuationScale ) );
1325 struct stat stat_buf;
1327 if (
stat( fname, &stat_buf ) ) {
1331 if ( ! S_ISDIR( stat_buf.st_mode ) ) {
void LALCheckMemoryLeaks(void)
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
void XLALReadTEMPOParFileOrig(BinaryPulsarParams *output, CHAR *pulsarAndPath)
function to read in a TEMPO parameter file
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 XLAL_INIT_DECL(var,...)
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)
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).
void XLALDestroySFTVector(SFTVector *vect)
XLAL interface to destroy an SFTVector.
void XLALDestroySFTCatalog(SFTCatalog *catalog)
Free an 'SFT-catalogue'.
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.
int XLALFindCoveringSFTBins(UINT4 *firstBin, UINT4 *numBins, REAL8 fMinIn, REAL8 BandIn, REAL8 Tsft)
Return the 'effective' frequency-band [fMinEff, fMaxEff] = [firstBin, lastBin] * 1/Tsft,...
int XLALWriteSFT2FilePointer(const SFTtype *sft, FILE *fp, const CHAR *SFTwindowtype, const REAL8 SFTwindowparam, const CHAR *SFTcomment)
Write the given SFTtype to a FILE pointer.
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.
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...
SFTCatalog * XLALSFTdataFind(const CHAR *file_pattern, const SFTConstraints *constraints)
Find the list of SFTs matching the file_pattern and satisfying the given constraints,...
int XLALCheckValidDescriptionField(const char *desc)
Check whether given string qualifies as a valid 'description' field of a FRAME filename (per ) or SFT...
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).
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.
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).
void XLALDestroyTimestampVector(LIGOTimeGPSVector *vect)
De-allocate a LIGOTimeGPSVector.
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
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 XLAL_CHECK(assertion,...)
int int XLALPrintWarning(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_CHECK_NULL(assertion,...)
LIGOTimeGPS * XLALGPSSetREAL8(LIGOTimeGPS *epoch, REAL8 t)
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
REAL8 XLALGPSDiff(const LIGOTimeGPS *t1, const LIGOTimeGPS *t0)
int main(int argc, char *argv[])
int XLALIsValidDescriptionField(const char *desc)
@ 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
int XLALFreeMem(ConfigVars_t *cfg)
This routine frees up all the memory.
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 ra
right ascension (rads)
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
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,...
A vector of 'timestamps' of type LIGOTimeGPS.
LIGOTimeGPS * data
array of timestamps
UINT4 length
number of timestamps
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.
An "SFT-catalogue": a vector of SFTdescriptors, as returned by XLALSFTdataFind()
SFTDescriptor * data
array of data-entries describing matched SFTs
UINT4 length
number of SFTs in catalog
'Constraints' for SFT-matching: which detector, within which time-stretch and which timestamps exactl...
const char * window_type
window function applied to SFT
REAL8 window_param
parameter of window function, if required
Structure specifying an SFT file name, following the convention in .
Parameters defining the SFTs to be returned from LALSignalToSFTs().
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, ...
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, ....