35 #include <lal/CWMakeFakeData.h>
36 #include <lal/TimeSeries.h>
37 #include <lal/GeneratePulsarSignal.h>
38 #include <lal/TransientCW_utils.h>
39 #include <lal/LALString.h>
40 #include <lal/StringVector.h>
41 #include <lal/Units.h>
42 #include <lal/ConfigFile.h>
43 #include <lal/LFTandTSutils.h>
45 #include <lal/FFTWMutex.h>
46 #include <lal/ExtrapolatePulsarSpins.h>
47 #include <lal/ConfigFile.h>
52 #define SQ(x) ( (x) * (x) )
59 "This is a comma-separated list of file patterns for configuration files,\n"
60 "or else direct configuration strings in the following format:\n"
61 " * Enclose with curly braces ('{}').\n"
62 " * Give pulsar parameters as key=value pairs with a '=' separator.\n"
63 " * Separate each key=value pair with a semicolon (';').\n"
64 "Available parameters are:\n"
65 " * Required parameters: Alpha, Delta, Freq, refTime\n"
66 " * Optional parameters:\n"
67 " - Injection amplitudes: either (h0, cosi) or (aPlus, aCross), psi, phi0\n"
68 " - Higher-order spindowns: f1dot, f2dot, ... f6dot\n"
69 " - Binary sources: orbitTp, orbitArgp, orbitasini, orbitEcc, orbitPeriod\n"
70 " - Transient injections: transientWindowType, transientStartTime, transientTau\n"
72 " * '{Alpha=0; Delta=0; Freq=50; f1dot=1e-11; f2dot=0; refTime=1000000000; h0=1.00000000e-23; cosi=0; psi=0; phi0=0;}'\n"
73 " * 'file1.dat,someFiles*.txt,{Alpha=0;Delta=0;Freq=0;refTime=1000000000;},someOtherFiles[0-9].dat'\n\n";
127 if ( multiTseries != NULL ) {
136 if ( outMSFTs != NULL ) {
137 svp = &( outMSFTs->
data[X] );
139 if ( outMTS != NULL ) {
140 tsp = &( outMTS->
data[X] );
147 if ( multiTseries ) {
148 ( *multiTseries ) = outMTS;
152 ( *multiSFTs ) = outMSFTs;
206 if ( SFTvect != NULL ) {
207 UINT4 firstBinEff, numBinsEff;
210 REAL8 fBand_eff = ( numBinsEff - 1.0 ) /
Tsft;
213 REAL8 fMax_eff = fMin_eff + fBand_eff;
214 if ( ( fMin_eff !=
fMin ) || ( fBand_eff !=
fBand ) ) {
215 XLALPrintWarning(
"Caller asked for Band [%.16g, %.16g] Hz, effective SFT-Band produced is [%.16g, %.16g] Hz\n",
217 XLAL_CHECK(
dataParams->inputMultiTS == NULL,
XLAL_EINVAL,
"Cannot expand effective frequency band with input timeseries given. Timeseries seems inconsistent with SFTs\n" );
220 fSamp = 2.0 * fBand_eff;
235 if ( ( SFTvect != NULL ) && ( n1_fSamp != n0_fSamp ) ) {
237 XLALPrintWarning(
"GAPS: Initial SFT sampling frequency fSamp0= %d/%.0f = %g had to be increased to fSamp1 = %d/%.0f = %g\n",
238 n0_fSamp,
Tsft, fSamp, n1_fSamp,
Tsft, fSamp1 );
239 XLAL_CHECK(
dataParams->inputMultiTS == NULL,
XLAL_EINVAL,
"Cannot expand effective frequency band with input timeseries given. Timeseries seems inconsistent with SFT timestamps\n" );
251 firstGPS =
ts->epoch;
279 for (
UINT4 iInj = 0; iInj < numPulsars; iInj ++ ) {
289 if (
t0 <= firstGPS_REAL8 ) {
290 signalStartGPS = firstGPS;
291 }
else if (
t0 >= lastGPS_REAL8 ) {
292 signalStartGPS = lastGPS;
296 REAL8 offs0_aligned = round( (
t0 - firstGPS_REAL8 ) * fSamp ) / fSamp;
297 signalStartGPS = firstGPS;
303 if (
t1 >= lastGPS_REAL8 ) {
304 signalEndGPS = lastGPS;
305 }
else if (
t1 <= firstGPS_REAL8 ) {
306 signalEndGPS = firstGPS;
308 signalEndGPS.gpsSeconds =
t1;
311 REAL8 fCoverMin, fCoverMax;
315 memcpy( spinRange.fkdot, fkdot,
sizeof( spinRange.fkdot ) );
319 XLAL_CHECK( ( fCoverMin >=
fMin ) && ( fCoverMax <
fMin +
fBand ),
XLAL_EINVAL,
"Error: injection signal %d:'%s' needs frequency band [%f,%f]Hz, injecting into [%f,%f]Hz\n",
323 XLAL_CHECK( signalDuration >= 0,
XLAL_EFAILED,
"Something went wrong, got negative signal duration = %g\n", signalDuration );
324 if ( signalDuration > 0 ) {
331 REAL8 bin_mismatch = fabs( Delta_epoch / Tseries_sum->
deltaT );
332 REAL8 mismatch = fabs( bin_mismatch - round( bin_mismatch ) ) * Tseries_sum->
deltaT;
333 XLAL_CHECK( mismatch <= 1
e-9,
XLAL_EDATA,
"Incompatible start-times when adding signal time series %d of %d, bins misaligned by %g seconds (%g bins).\n", iInj, numPulsars, mismatch, bin_mismatch );
343 REAL8 noiseSigma = sqrtSn * sqrt( 0.5 * fSamp );
358 REAL8 bin_mismatch = fabs( Delta_epoch / outTS->
deltaT );
359 REAL8 mismatch = fabs( bin_mismatch - round( bin_mismatch ) ) * outTS->
deltaT;
360 XLAL_CHECK( mismatch <= 1
e-9,
XLAL_EDATA,
"Incompatible start-times when adding input noise time-series, bins misaligned by %g seconds (%g bins).\n", mismatch, bin_mismatch );
365 if ( SFTvect != NULL ) {
376 if ( Tseries != NULL ) {
377 ( *Tseries ) = outTS;
425 for (
UINT4 s = 0;
s < s_max;
s ++ ) {
438 params.pulsar.aPlus = aPlus;
439 params.pulsar.aCross = aCross;
443 params.pulsar.spindown = spindown;
451 params.fHeterodyne = fHet;
452 params.sourceDeltaT = sourceDeltaT;
457 params.samplingRate = fSamp;
482 const char *windowType,
495 UINT4 timestepsSFT = lround( timestepsSFT0 );
497 "Inconsistent sampling-step (dt=%g) and Tsft=%g: must be integer multiple Tsft/dt = %g >= %g\n",
502 if ( windowType != NULL ) {
509 UINT4 numSFTBins = timestepsSFT / 2 + 1;
510 fftw_complex *fftOut;
511 XLAL_CHECK_NULL( ( fftOut = fftw_malloc( numSFTBins *
sizeof( fftOut[0] ) ) ) != NULL,
XLAL_ENOMEM,
"fftw_malloc(%d*sizeof(complex)) failed\n", numSFTBins );
526 char buf1[256], buf2[256];
538 "XLALCreateSFTVector(numSFTs=%d, numBins=%d) failed.\n",
numSFTs, numSFTBins );
546 INT4 offsetBins = lround( offset /
dt );
549 memcpy( timeStretchCopy->
data, timeseries->
data->
data + offsetBins, timeStretchCopy->
length *
sizeof( timeStretchCopy->
data[0] ) );
552 REAL8 sigma_window = 1;
553 if ( window != NULL ) {
555 for (
UINT4 iBin = 0; iBin < timeStretchCopy->
length; iBin++ ) {
556 timeStretchCopy->
data[iBin] *= window->
data->
data[iBin];
561 fftw_execute( fftplan );
564 strcpy( thisSFT->
name, timeseries->
name );
566 thisSFT->
f0 = timeseries->
f0;
571 REAL8 norm =
dt / sigma_window;
572 for (
UINT4 k = 0;
k < numSFTBins ;
k ++ ) {
586 fftw_destroy_plan( fftplan );
629 XLAL_CHECK( TsftREAL == round( TsftREAL ),
XLAL_EDOM,
"Only exact integer-second Tsft allowed, got Tsft = %.16g s\n", TsftREAL );
646 INT4 nsdiff =
t1->gpsNanoSeconds -
t0->gpsNanoSeconds;
647 XLAL_CHECK( nsdiff == 0,
XLAL_EDOM,
"Only integer-second gaps allowed, found %d ns excess in gap between i=%d and i-1\n", nsdiff,
i );
649 INT4 gap_i0 =
t1->gpsSeconds -
t0->gpsSeconds;
650 XLAL_CHECK( gap_i0 > 0,
XLAL_EDOM,
"Timestamps must be sorted in increasing order, found negative gap %d s between i=%d and i-1\n", gap_i0,
i );
655 if ( ( (
INT8 )gap_i * nCur ) %
Tsft == 0 ) {
676 REAL8 Tg = TsftREAL / g;
677 UINT4 nNew = (
UINT4 ) ceil( nCur / Tg ) * Tg;
679 XLAL_CHECK( nNew > nCur,
XLAL_ETOL,
"This seems wrong: nNew = %d !> nCur = %d, but should be greater!\n", nNew, nCur );
680 XLALPrintInfo(
"Need to increase from fSamp = %d/Tsft = %g to fSamp = %d / Tsft = %g\n", nCur, nCur / TsftREAL, nNew, nNew / TsftREAL );
706 UINT4 next_numer, next_denom, rmdr;
710 while ( next_denom != 0 ) {
711 rmdr = next_numer % next_denom;
712 next_numer = next_denom;
728 if ( numPulsars > 0 ) {
742 if ( ppvect == NULL ) {
746 if ( ppvect->
data != NULL ) {
819 aPlus = 0.5 *
h0 * ( 1.0 +
SQ(
cosi ) );
829 pulsarParams->
Amp.
psi = psi;
843 if ( have_refTime ) {
845 }
else if ( refTimeDef != NULL ) {
917 REAL8 orbitasini = 0 ;
923 REAL8 orbitPeriod = 0;
927 if ( have_orbitasini || have_orbitEcc || have_orbitPeriod || have_orbitArgp || have_orbitTp ) {
929 XLAL_CHECK( ( orbitasini == 0 ) || ( have_orbitPeriod && ( orbitPeriod > 0 ) && have_orbitTp ),
XLAL_EINVAL,
"If orbitasini>0 then we also need 'orbitPeriod>0' and 'orbitTp'\n" );
944 char *transientWindowType = NULL;
945 BOOLEAN have_transientWindowType;
948 UINT4 transientStartTime = 0;
949 BOOLEAN have_transientStartTime;
950 XLAL_CHECK( XLALReadConfigUINT4Variable( &transientStartTime, cfgdata, secName,
"transientStartTime", &have_transientStartTime ) ==
XLAL_SUCCESS,
XLAL_EFUNC );
952 UINT4 transientTau = 0;
955 REAL8 transientTauDays = 0;
960 if ( have_transientWindowType ) {
967 XLAL_CHECK( have_transientStartTime && ( have_transientTau || have_transientTauDays ),
XLAL_EINVAL,
"For transientWindowType!=None, we also need transientStartTime and either transientTau (deprecated) or transientTau." );
968 XLAL_CHECK( !( have_transientTau && have_transientTauDays ),
XLAL_EINVAL,
"Cannot have both transientTau and transientTauDays; the latter is deprecated." );
971 if ( have_transientTauDays ) {
973 printf(
"Warning: Option transientTauDays is deprecated, please switch to using transientTau [UINT4, in seconds] instead.\n" );
980 XLAL_CHECK( !( have_transientStartTime || have_transientTau || have_transientTauDays ),
XLAL_EINVAL,
"Cannot use transientStartTime, transientTau or transientTauDays without transientWindowType!" );
1009 for (
UINT4 i = 0;
i < numPulsars;
i ++ ) {
1010 const char *sec_i = sections->
data[
i];
1012 if ( strcmp( sec_i,
"default" ) == 0 ) {
1043 if ( n_unread != NULL ) {
1044 XLALPrintError(
"ERROR: Pulsar params config file '%s' contained '%d' unknown entries:\n", fname, n_unread->
length );
1075 const char *thisInput = UserInput->
data[
l];
1077 if ( thisInput[0] !=
'{' ) {
1081 for (
UINT4 i = 0;
i < numFiles;
i ++ ) {
1093 UINT4 len = strlen( thisInput );
1094 XLAL_CHECK_NULL( ( thisInput[0] ==
'{' ) && ( thisInput[len - 1] ==
'}' ),
XLAL_EINVAL,
"Invalid file-content input:\n%s\n", thisInput );
1097 len = strlen( buf );
1108 strncpy( addSource->
data[0].
name,
"direct-string-input",
sizeof( addSource->
data[0].
name ) );
1134 if ( list == NULL ) {
1142 UINT4 newlen = oldlen + addlen;
1145 memcpy( ret->
data + oldlen, add->
data, addlen *
sizeof( ret->
data[0] ) );
1147 for (
UINT4 i = 0;
i < addlen;
i ++ ) {
1181 snprintf( col_name,
sizeof( col_name ),
"f%zudot [Hz/s^%zu]",
k,
k );
1194 for (
size_t i = 0;
i < list->
length; ++
i ) {
int XLALReadConfigSTRINGVariable(CHAR **varp, LALParsedDataFile *cfgdata, const CHAR *secName, const CHAR *varName, BOOLEAN *wasRead)
#define LAL_FFTW_WISDOM_UNLOCK
#define LAL_FFTW_WISDOM_LOCK
int XLALFITSTableWriteRow(FITSFile UNUSED *file, const void UNUSED *record)
int XLALFITSTableOpenWrite(FITSFile UNUSED *file, const CHAR UNUSED *name, const CHAR UNUSED *comment)
const WeaveSearchTimingDenominator denom
LIGOTimeGPSVector * timestamps
int XLALFindSmallestValidSamplingRate(UINT4 *n1, UINT4 n0, const LIGOTimeGPSVector *timestamps)
Find the smallest sampling rate of the form fsamp = n / Tsft, with n>=n0, such that all gap sizes Dt_...
int XLALCWMakeFakeData(SFTVector **SFTvect, REAL8TimeSeries **Tseries, const PulsarParamsVector *injectionSources, const CWMFDataParams *dataParams, UINT4 detectorIndex, const EphemerisData *edat)
Single-IFO version of XLALCWMakeFakeMultiData(), handling the actual work, but same input API.
PulsarParamsVector * XLALPulsarParamsVectorAppend(PulsarParamsVector *list, const PulsarParamsVector *add)
Append the given PulsarParamsVector 'add' to the vector 'list' ( which can be NULL),...
PulsarParamsVector * XLALCreatePulsarParamsVector(UINT4 numPulsars)
Create zero-initialized PulsarParamsVector for numPulsars.
int XLALReadPulsarParams(PulsarParams *pulsarParams, LALParsedDataFile *cfgdata, const CHAR *secName, const LIGOTimeGPS *refTimeDef)
Function to parse a config-file-type string (or section thereof) into a PulsarParams struct.
static UINT4 gcd(UINT4 numer, UINT4 denom)
int XLALCheckConfigFileWasFullyParsed(const char *fname, const LALParsedDataFile *cfgdata)
SFTVector * XLALMakeSFTsFromREAL8TimeSeries(const REAL8TimeSeries *timeseries, const LIGOTimeGPSVector *timestamps, const char *windowType, REAL8 windowParam)
Make SFTs from given REAL8TimeSeries at given timestamps, potentially applying a time-domain window o...
PulsarParamsVector * XLALPulsarParamsFromUserInput(const LALStringVector *UserInput, const LIGOTimeGPS *refTimeDef)
Function to determine the PulsarParamsVector input from a user-input defining CW sources.
const char *const InjectionSourcesHelpString
void XLALDestroyPulsarParamsVector(PulsarParamsVector *ppvect)
Destructor for PulsarParamsVector type.
int XLALFITSWritePulsarParamsVector(FITSFile *file, const CHAR *tableName, const PulsarParamsVector *list)
Write a PulsarParamsVector to a FITS file.
#define SQ(x)
Functions to generate 'fake' data containing CW signals and/or Gaussian noise. These basically presen...
void XLALDestroyCWMFDataParams(CWMFDataParams *params)
Destructor for a CWMFDataParams type.
REAL4TimeSeries * XLALGenerateCWSignalTS(const PulsarParams *pulsarParams, const LALDetector *site, LIGOTimeGPS startTime, REAL8 duration, REAL8 fSamp, REAL8 fHet, const EphemerisData *edat, REAL8 sourceDeltaT)
Generate a (heterodyned) REAL4 timeseries of a CW signal for given pulsarParams, site,...
PulsarParamsVector * XLALPulsarParamsFromFile(const char *fname, const LIGOTimeGPS *refTimeDef)
Parse a given 'CWsources' config file for PulsarParams, return vector of all pulsar definitions found...
int XLALCWMakeFakeMultiData(MultiSFTVector **multiSFTs, MultiREAL8TimeSeries **multiTseries, const PulsarParamsVector *injectionSources, const CWMFDataParams *dataParams, const EphemerisData *edat)
Generate fake 'CW' data, returned either as SFTVector or REAL4Timeseries or both, for given CW-signal...
UINT4Vector * XLALConfigFileGetUnreadEntries(const LALParsedDataFile *cfgdata)
int XLALParseDataFileContent(LALParsedDataFile **cfgdata, const CHAR *string)
int XLALParseDataFile(LALParsedDataFile **cfgdata, const CHAR *fname)
LALStringVector * XLALListConfigFileSections(const LALParsedDataFile *cfgdata)
void XLALDestroyParsedDataFile(LALParsedDataFile *cfgdata)
char * XLALGPSToStr(char *s, const LIGOTimeGPS *t)
#define XLAL_FITS_TABLE_COLUMN_BEGIN(record_type)
#define XLAL_FITS_TABLE_COLUMN_ADD_NAMED(file, type, field, col_name)
#define XLAL_FITS_TABLE_COLUMN_ADD_ARRAY(file, type, field)
struct tagFITSFile FITSFile
Representation of a FITS file.
REAL4TimeSeries * XLALGeneratePulsarSignal(const PulsarSignalParams *params)
Generate a time-series at the detector for a given pulsar.
void XLALDestroyMultiREAL8TimeSeries(MultiREAL8TimeSeries *multiTS)
Destroy a MultiREAL8TimeSeries, NULL-robust.
int XLALAddGaussianNoise(REAL4TimeSeries *inSeries, REAL4 sigma, INT4 seed)
Generate Gaussian noise with standard-deviation sigma, add it to inSeries.
int XLALcorrect_phase(SFTtype *sft, LIGOTimeGPS tHeterodyne)
Yousuke's phase-correction function, taken from makefakedata_v2.
#define XLAL_INIT_DECL(var,...)
void * XLALCalloc(size_t m, size_t n)
void * XLALRealloc(void *p, size_t n)
char char * XLALStringDuplicate(const char *s)
REAL8 PulsarSpins[PULSAR_MAX_SPINS]
Typedef for fixed-size array holding GW frequency and derivatives fk = d^k Freq/dt^k|(tau_ref)
#define PULSAR_MAX_SPINS
maximal number of spin-parameters (Freq + spindowns) we can handle
@ TRANSIENT_NONE
Note: in this case the window-parameters will be ignored, and treated as rect={data},...
void XLALDestroySFTVector(SFTVector *vect)
XLAL interface to destroy an 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,...
LALStringVector * XLALFindFiles(const CHAR *globstring)
Returns a list of filenames matching the input argument, which may be one of the following:
void XLALDestroyMultiTimestamps(MultiLIGOTimeGPSVector *multiTS)
Destroy a MultiLIGOTimeGPSVector timestamps vector.
SFTVector * XLALCreateSFTVector(UINT4 numSFTs, UINT4 numBins)
XLAL function to create an SFTVector of numSFT SFTs with SFTlen frequency-bins (which will be allocat...
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).
COORDINATESYSTEM_EQUATORIAL
void XLALDestroyStringVector(LALStringVector *vect)
REAL8TimeSeries * XLALConvertREAL4TimeSeriesToREAL8(const REAL4TimeSeries *series)
REAL8TimeSeries * XLALAddREAL8TimeSeries(REAL8TimeSeries *arg1, const REAL8TimeSeries *arg2)
void XLALDestroyREAL8TimeSeries(REAL8TimeSeries *series)
void XLALDestroyREAL4TimeSeries(REAL4TimeSeries *series)
REAL4TimeSeries * XLALAddREAL4TimeSeries(REAL4TimeSeries *arg1, const REAL4TimeSeries *arg2)
REAL4TimeSeries * XLALCreateREAL4TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
int XLALParseTransientWindowName(const char *windowName)
Parse a transient window name string into the corresponding transientWindowType.
int XLALGetTransientWindowTimespan(UINT4 *t0, UINT4 *t1, transientWindow_t transientWindow)
Helper-function to determine the total timespan of a transient CW window, ie.
int XLALApplyTransientWindow(REAL4TimeSeries *series, transientWindow_t transientWindow)
apply a "transient CW window" described by TransientWindowParams to the given timeseries
const LALUnit lalStrainUnit
REAL8Vector * XLALCreateREAL8Vector(UINT4 length)
void XLALDestroyREAL8Vector(REAL8Vector *vector)
REAL8Window * XLALCreateNamedREAL8Window(const char *windowName, REAL8 beta, UINT4 length)
void XLALDestroyREAL8Window(REAL8Window *window)
int int int XLALPrintInfo(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_CHECK(assertion,...)
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
int int XLALPrintWarning(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_CHECK_NULL(assertion,...)
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
REAL8 XLALGPSGetREAL8(const LIGOTimeGPS *epoch)
REAL8 XLALGPSDiff(const LIGOTimeGPS *t1, const LIGOTimeGPS *t0)
A vector of COMPLEX8FrequencySeries.
COMPLEX8FrequencySeries * data
Pointer to the data array.
Struct controlling all the aspects of the fake data (time-series + SFTs) to be produced by XLALCWMake...
This structure contains all information about the center-of-mass positions of the Earth and Sun,...
A vector of 'timestamps' of type LIGOTimeGPS.
REAL8 deltaT
'length' of each timestamp (e.g.
LIGOTimeGPS * data
array of timestamps
UINT4 length
number of timestamps
A collection of (multi-IFO) LIGOTimeGPSVector time-stamps vectors.
A collection of (multi-IFO) REAL8 time-series.
REAL8TimeSeries ** data
Pointer to the data array.
UINT4 length
Number of elements in array.
A collection of SFT vectors – one for each IFO in a multi-IFO search.
UINT4 length
number of ifos
SFTVector ** data
sftvector for each ifo
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, ...
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
CHAR name[LALNameLength]
'name' for this sources, can be an empty string
transientWindow_t Transient
Transient window-parameters (start-time, duration, window-type)
PulsarDopplerParams Doppler
'Phase-evolution parameters': {skypos, fkdot, orbital params }
Straightforward vector type of N PulsarParams structs.
PulsarParams * data
array of pulsar-signal parameters
UINT4 length
number of pulsar-signals
Input parameters to GeneratePulsarSignal(), defining the source and the time-series.
Contains a "spin-range", ie spins and corresponding bands at a given (SSB) reference GPS-time .
UINT4 t0
GPS start-time 't0'.
UINT4 tau
transient timescale tau in seconds
transientWindowType_t type
window-type: none, rectangular, exponential, ....