Pulsar signal-generation routines for hardware- and software-injections.
This module also contains a few more general-purpose helper-functions:
params->pulsar.position
), the detector-site (params->site
) and the ephemeris-data (params->ephemerides
)are used from the PulsarSignalParams structure.LALGeneratePulsarSignal() is basically a wrapper for the two LAL-functions GenerateSpinOrbitCW() to produce the source-signal, and LALPulsarSimulateCoherentGW() to turn it into a time-series at the detector.
LALSignalToSFTs() uses LALForwardRealFFT() appropriately on the input-timeseries to produce the required output-SFTs.
LALSignalToSFTs() currently enforces the constraint of Tsft * Band
being an integer, so that the number of time-samples per SFT is even. This follows makefakedata_v2
.
Furthermore, if the timestamps for SFT-creation passed to LALSignalToSFTs() do not fit exactly on a time-step of the input time-series, it will be "nudged" to the closest one. If lalDebugLevel>0
a warning will be printed about this. The user could also see this effect in the actual timestamps of the SFTs returned.
The FFTW-"plan" is currently created using the ESTIMATE
flag, which is fast but only yields an approximate plan. Better results might be achieved by using MEASURE
and an appropriate buffering of the resulting plan (which doesnt change provided the SFT-length is the same). Measuring the plan takes longer but might lead to substantial speedups for many FFTs, which seems the most likely situation.
The functions LALComputeSkyAndZeroPsiAMResponse() and LALFastGeneratePulsarSFTs() use approximate analytic formulas to generate SFTs. This should be significantly faster than LALGeneratePulsarSignal() and LALSignalToSFTs(), which generate the time series data and then FFT it. Simple tests performed by the code in GeneratePulsarSignalTest.c indicate that the maximum modulus of the SFTs output by the approximate and exact codes differs by less that 10%. Since the tests are not exhaustive, the user should use caution and conduct their own test to compare LALComputeSkyAndZeroPsiAMResponse() and LALFastGeneratePulsarSFTs() with LALGeneratePulsarSignal() and LALSignalToSFTs().
The strain of a periodic signal at the detector is given by
\[ h(t) = F_+(t) A_+ {\rm cos}\Phi(t) + F_\times(t) A_\times {\rm sin}\Phi(t), \]
where \( F_+ \) and \( F_\times \) are the usual beam pattern response functions, \( A_+ \) and \( A_\times \) are the amplitudes of the gravitational wave for the plus and cross polarizations, and \( \Phi \) is the phase. The phase contains modulations from doppler shifts due to the relative motion between the source and the detector and the spin evolution of the source. (The functions discussed here support both isolated sources and those in binary systems. The binary case has not been tested.)
If we Taylor expand the phase out to first order about the time at the midpoint of an SFT and approximate \( F_+ \) and \( F_\times \) as constants, for one SFT we can write
\[ \Phi(t) \approx \Phi_{1/2} + 2\pi f_{1/2}(t - t_{1/2}). \]
The strain at discrete time \( t_j \) , measured from the start of the SFT, can thus be approximated as
\[ h_j \approx F_{+ 1/2} A_+ {\rm cos} [\Phi_{1/2} + 2\pi f_{1/2}(t_0 + t_j - t_{1/2})] + F_{\times 1/2} A_\times {\rm sin} [\Phi_{1/2} + 2\pi f_{1/2}(t_0 + t_j - t_{1/2})], \]
where \( t_0 \) is the time as the start of the SFT, and \( t_{1/2} - t_0 = T_{\rm sft}/2 \) , where \( T_{\rm sft} \) is the duration of one SFT. This simplifies to
\[ h_j \approx F_{+ 1/2} A_+ {\rm cos} (\Phi_0 + 2\pi f_{1/2}t_j) + F_{\times 1/2} A_\times {\rm sin} (\Phi_0 + 2\pi f_{1/2}t_j), \]
where \( \Phi_0 \) is the phase at the start of the SFT (not the initial phase at the start of the observation), i.e.,
\[ \Phi_0 = \Phi_{1/2} - 2 \pi f_{1/2} (T_{\rm sft} / 2). \]
Note that these are the same approximations used by LALDemod().
One can show that the Discrete Fourier Transform (DFT) of \( h_j \) above is:
\[ \tilde{h}_k = e^{i\Phi_0} \frac{1}{2} ( F_{+ 1/2} A_+ - i F_{\times 1/2} A_\times) \frac{ 1 - e^{2\pi i (\kappa - k)}}{1 - e^{2\pi i (\kappa - k)/N} } \\ + e^{-i\Phi_0} \frac{1}{2} ( F_{+ 1/2} A_+ + i F_{\times 1/2} A_\times) \frac{ 1 - e^{-2\pi i (\kappa + k)}}{ 1 - e^{-2\pi i (\kappa + k)/N} } \]
where \( N \) is the number of time samples used to find the DFT (i.e., the sample rate times \( T_{\rm sft} \) ), and
\[ \kappa \equiv f_{1/2} T_{\rm sft}, \]
is usually not an integer.
Note that the factor \( e^{\pm 2\pi i k} \) in the numerators of the equation for \( \tilde{h}_k \) equals 1. Furthermore, for \( 0 < \kappa < N/2 \) and \( |\kappa - k| << N \) the first term dominates and can be Taylor expanded to give:
\[ \tilde{h}_k = N e^{i\Phi_0} \frac{1}{2} ( F_{+ 1/2} A_+ - i F_{\times 1/2} A_\times) \left [ \, \frac{ {\rm sin} (2\pi\kappa)}{2 \pi (\kappa - k) } \, + \, i \frac{ 1 - {\rm cos} (2\pi\kappa)}{2 \pi (\kappa - k) } \, \right ] \]
Note that the last factor in square brackets is \( P_{\alpha k}^* \) and \( e^{i\Phi_0} = Q_{\alpha}^* \) used by LALDemod.
The structs used by LALComputeSkyAndZeroPsiAMResponse() and LALFastGeneratePulsarSFTs() are given in previous sections, and make use of those used by LALGeneratePulsarSignal() and LALSignalToSFTs() plus a small number of additional parameters. Thus it is fairly easy to change between the above approximate routines the exact routines. See GeneratePulsarSignalTest.c for an example implementation of the code.
Note that one needs to call LALComputeSkyAndZeroPsiAMResponse() once per sky position, and then call LALFastGeneratePulsarSFTs() for each set of signal parameters at that sky position. Thus, one could perform a Monte Carlo simulation, as shown by the pseudo code:
*outputSFTs
sent to LALFastGeneratePulsarSFTs() is NULL
then LALFastGeneratePulsarSFTs() allocates memory for the output SFTs; otherwise it assumes memory has already been allocated. Thus, the user does not have to deallocate memory for the SFTs until all calls to LALFastGeneratePulsarSFTs() are completed.fHeterodyne
and 0.5 * samplingRate
set in the PulsarSignalParams struct give the start frequency and frequency band of the SFTs output from LALFastGeneratePulsarSFTs().resTrig
is set to zero in the SFTandSignalParams struct, then the C math libary cos()
sin()
functions are called, else lookup tables (LUTs) are used for calls to trig functions. There may be a slight speedup in using LUTs.2*Dterms
centered on the signal frequency in each SFT. Dterms must be greater than zero and less than or equal to the number of frequency bins in the output SFTs. Note that Dterms is used the same way here as it is in LALDemod(). Nothing is done to the other bins, unless *outputSFTs
is NULL
; then, since memory is allocates for the output SFTs, the bins not in the 2*Dterms
band are initialized to zero. Prototypes | |
REAL4TimeSeries * | XLALGeneratePulsarSignal (const PulsarSignalParams *params) |
Generate a time-series at the detector for a given pulsar. More... | |
void | LALGeneratePulsarSignal (LALStatus *status, REAL4TimeSeries **signalvec, const PulsarSignalParams *params) |
SFTVector * | XLALSignalToSFTs (const REAL4TimeSeries *signalvec, const SFTParams *params) |
Turn the given time-series into SFTs and add noise if given. More... | |
void | LALSignalToSFTs (LALStatus *status, SFTVector **outputSFTs, const REAL4TimeSeries *signalvec, const SFTParams *params) |
void | LALComputeSkyAndZeroPsiAMResponse (LALStatus *status, SkyConstAndZeroPsiAMResponse *output, const SFTandSignalParams *params) |
/deprecated Move to attic? Wrapper for LALComputeSky() and LALComputeDetAMResponse() that finds the sky constants and \( F_+ \) and \( F_\times \) for use with LALFastGeneratePulsarSFTs(). More... | |
void | LALFastGeneratePulsarSFTs (LALStatus *status, SFTVector **outputSFTs, const SkyConstAndZeroPsiAMResponse *input, const SFTandSignalParams *params) |
Fast generation of Fake SFTs for a pure pulsar signal. More... | |
int | XLALConvertGPS2SSB (LIGOTimeGPS *SSBout, LIGOTimeGPS GPSin, const PulsarSignalParams *params) |
Convert earth-frame GPS time into barycentric-frame SSB time for given source. More... | |
int | XLALConvertSSB2GPS (LIGOTimeGPS *GPSout, LIGOTimeGPS SSBin, const PulsarSignalParams *params) |
Convert barycentric frame SSB time into earth-frame GPS time. More... | |
REAL4TimeSeries * | XLALGenerateLineFeature (const PulsarSignalParams *params) |
Generate a REAL4TimeSeries containing a sinusoid with amplitude 'h0', frequency 'Freq-fHeterodyne' and initial phase 'phi0'. More... | |
int | XLALAddGaussianNoise (REAL4TimeSeries *inSeries, REAL4 sigma, INT4 seed) |
Generate Gaussian noise with standard-deviation sigma, add it to inSeries. More... | |
void | XLALDestroyMultiREAL4TimeSeries (MultiREAL4TimeSeries *multiTS) |
Destroy a MultiREAL4TimeSeries, NULL-robust. More... | |
void | XLALDestroyMultiREAL8TimeSeries (MultiREAL8TimeSeries *multiTS) |
Destroy a MultiREAL8TimeSeries, NULL-robust. More... | |
static int | XLALcheck_timestamp_bounds (const LIGOTimeGPSVector *timestamps, LIGOTimeGPS t0, LIGOTimeGPS t1) |
Check that all timestamps given lie within the range [t0, t1]. More... | |
static int | XLALcheckNoiseSFTs (const SFTVector *sfts, REAL8 f0, REAL8 f1, REAL8 deltaF) |
Check if frequency-range and resolution of noiseSFTs is consistent with signal-band [f0, f1]. More... | |
int | XLALcorrect_phase (SFTtype *sft, LIGOTimeGPS tHeterodyne) |
Yousuke's phase-correction function, taken from makefakedata_v2. More... | |
Data Structures | |
struct | PulsarSignalParams |
Input parameters to GeneratePulsarSignal(), defining the source and the time-series. More... | |
struct | SFTParams |
Parameters defining the SFTs to be returned from LALSignalToSFTs(). More... | |
struct | SFTandSignalParams |
Parameters defining the pulsar signal and SFTs used by LALFastGeneratePulsarSFTs(). More... | |
struct | SkyConstAndZeroPsiAMResponse |
Sky Constants and beam pattern response functions used by LALFastGeneratePulsarSFTs(). More... | |
Files | |
file | GeneratePulsarSignalTest.c |
Error codes | |
#define | GENERATEPULSARSIGNALH_ENULL 1 |
Arguments contained an unexpected null pointer. More... | |
#define | GENERATEPULSARSIGNALH_ENONULL 2 |
Output pointer is not NULL. More... | |
#define | GENERATEPULSARSIGNALH_EMEM 3 |
Out of memory. More... | |
#define | GENERATEPULSARSIGNALH_ESAMPLING 4 |
Waveform sampling interval too large. More... | |
#define | GENERATEPULSARSIGNALH_ESSBCONVERT 5 |
SSB->GPS iterative conversion failed. More... | |
#define | GENERATEPULSARSIGNALH_ESYS 6 |
System error, probably while File I/O. More... | |
#define | GENERATEPULSARSIGNALH_ETIMEBOUND 7 |
Timestamp outside of allowed time-interval. More... | |
#define | GENERATEPULSARSIGNALH_ENUMSFTS 8 |
Inconsistent number of SFTs in timestamps and noise-SFTs. More... | |
#define | GENERATEPULSARSIGNALH_EINCONSBAND 9 |
Inconsistent values of sampling-rate and Tsft. More... | |
#define | GENERATEPULSARSIGNALH_ENOISEDELTAF 10 |
Frequency resolution of noise-SFTs inconsistent with signal. More... | |
#define | GENERATEPULSARSIGNALH_ENOISEBAND 11 |
Frequency band of noise-SFTs inconsistent with signal. More... | |
#define | GENERATEPULSARSIGNALH_ENOISEBINS 12 |
Frequency bins of noise-SFTs inconsistent with signal. More... | |
#define | GENERATEPULSARSIGNALH_EBADCOORDS 13 |
Current code requires sky position in equatorial coordinates. More... | |
#define | GENERATEPULSARSIGNALH_ELUTS 14 |
Lookup tables (LUTs) for trig functions must be defined on domain -2pi to 2pi inclusive. More... | |
#define | GENERATEPULSARSIGNALH_EDTERMS 15 |
Dterms must be greater than zero and less than or equal to half the number of SFT bins. More... | |
#define | GENERATEPULSARSIGNALH_EINPUT 16 |
Invalid input-arguments to function. More... | |
#define | GENERATEPULSARSIGNALH_EDETECTOR 17 |
Unknown detector-name. More... | |
REAL4TimeSeries * XLALGeneratePulsarSignal | ( | const PulsarSignalParams * | params | ) |
Generate a time-series at the detector for a given pulsar.
params | input params |
Definition at line 62 of file GeneratePulsarSignal.c.
void LALGeneratePulsarSignal | ( | LALStatus * | status, |
REAL4TimeSeries ** | signalvec, | ||
const PulsarSignalParams * | params | ||
) |
status | pointer to LALStatus structure |
signalvec | output time-series |
params | input params |
Definition at line 228 of file GeneratePulsarSignal.c.
SFTVector * XLALSignalToSFTs | ( | const REAL4TimeSeries * | signalvec, |
const SFTParams * | params | ||
) |
Turn the given time-series into SFTs and add noise if given.
signalvec | input time-series |
params | params for output-SFTs |
Definition at line 250 of file GeneratePulsarSignal.c.
void LALSignalToSFTs | ( | LALStatus * | status, |
SFTVector ** | outputSFTs, | ||
const REAL4TimeSeries * | signalvec, | ||
const SFTParams * | params | ||
) |
status | pointer to LALStatus structure | |
[out] | outputSFTs | SFT-vector |
signalvec | input time-series | |
params | params for output-SFTs |
Definition at line 436 of file GeneratePulsarSignal.c.
void LALComputeSkyAndZeroPsiAMResponse | ( | LALStatus * | status, |
SkyConstAndZeroPsiAMResponse * | output, | ||
const SFTandSignalParams * | params | ||
) |
/deprecated Move to attic? Wrapper for LALComputeSky() and LALComputeDetAMResponse() that finds the sky constants and \( F_+ \) and \( F_\times \) for use with LALFastGeneratePulsarSFTs().
Uses the output of LALComputeSkyAndZeroPsiAMResponse() and the same inputs as LALGeneratePulsarSignal() and LALSignalToSFTs(). This function used LALComputeSkyBinary() if params->pSigParams->orbit is not NULL, else it uses LALComputeSky() to find the skyConsts. NOTE THAT THIS FUNCTION COMPUTES \( F_+ \) and \( F_x \) for ZERO Psi!!! LALFastGeneratePulsarSFTs() used these to find \( F_+ \) and \( F_x \) for NONZERO Psi.
status | pointer to LALStatus structure |
output | output |
params | params |
Definition at line 475 of file GeneratePulsarSignal.c.
void LALFastGeneratePulsarSFTs | ( | LALStatus * | status, |
SFTVector ** | outputSFTs, | ||
const SkyConstAndZeroPsiAMResponse * | input, | ||
const SFTandSignalParams * | params | ||
) |
Fast generation of Fake SFTs for a pure pulsar signal.
Uses the output of LALComputeSkyAndZeroPsiAMResponse and the same inputs as LALGeneratePulsarSignal and LALSignalToSFTs. The fake signal is Taylor expanded to first order about the midpoint time of each SFT. Analytic expressions are used to find each SFT
Definition at line 573 of file GeneratePulsarSignal.c.
int XLALConvertGPS2SSB | ( | LIGOTimeGPS * | SSBout, |
LIGOTimeGPS | GPSin, | ||
const PulsarSignalParams * | params | ||
) |
Convert earth-frame GPS time into barycentric-frame SSB time for given source.
[out] | SSBout | arrival-time in SSB |
[in] | GPSin | GPS-arrival time at detector |
params | define source-location and detector |
Definition at line 824 of file GeneratePulsarSignal.c.
int XLALConvertSSB2GPS | ( | LIGOTimeGPS * | GPSout, |
LIGOTimeGPS | SSBin, | ||
const PulsarSignalParams * | params | ||
) |
Convert barycentric frame SSB time into earth-frame GPS time.
NOTE: this uses simply the inversion-routine used in the original makefakedata_v2
[out] | GPSout | GPS-arrival-time at detector |
[in] | SSBin | input: signal arrival time at SSB |
params | params defining source-location and detector |
Definition at line 870 of file GeneratePulsarSignal.c.
REAL4TimeSeries * XLALGenerateLineFeature | ( | const PulsarSignalParams * | params | ) |
Generate a REAL4TimeSeries containing a sinusoid with amplitude 'h0', frequency 'Freq-fHeterodyne' and initial phase 'phi0'.
Definition at line 941 of file GeneratePulsarSignal.c.
int XLALAddGaussianNoise | ( | REAL4TimeSeries * | inSeries, |
REAL4 | sigma, | ||
INT4 | seed | ||
) |
Generate Gaussian noise with standard-deviation sigma, add it to inSeries.
Definition at line 980 of file GeneratePulsarSignal.c.
void XLALDestroyMultiREAL4TimeSeries | ( | MultiREAL4TimeSeries * | multiTS | ) |
Destroy a MultiREAL4TimeSeries, NULL-robust.
Definition at line 1032 of file GeneratePulsarSignal.c.
void XLALDestroyMultiREAL8TimeSeries | ( | MultiREAL8TimeSeries * | multiTS | ) |
Destroy a MultiREAL8TimeSeries, NULL-robust.
Definition at line 1054 of file GeneratePulsarSignal.c.
|
static |
Check that all timestamps given lie within the range [t0, t1].
return: 0 if ok, ERROR if not
Definition at line 1084 of file GeneratePulsarSignal.c.
Check if frequency-range and resolution of noiseSFTs is consistent with signal-band [f0, f1].
return XLAL_SUCCESS if everything fine, error-code otherwise
Definition at line 1117 of file GeneratePulsarSignal.c.
int XLALcorrect_phase | ( | SFTtype * | sft, |
LIGOTimeGPS | tHeterodyne | ||
) |
Yousuke's phase-correction function, taken from makefakedata_v2.
Definition at line 1163 of file GeneratePulsarSignal.c.
#define GENERATEPULSARSIGNALH_ENULL 1 |
Arguments contained an unexpected null pointer.
Definition at line 236 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ENONULL 2 |
Output pointer is not NULL.
Definition at line 237 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EMEM 3 |
Out of memory.
Definition at line 238 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ESAMPLING 4 |
Waveform sampling interval too large.
Definition at line 239 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ESSBCONVERT 5 |
SSB->GPS iterative conversion failed.
Definition at line 240 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ESYS 6 |
System error, probably while File I/O.
Definition at line 241 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ETIMEBOUND 7 |
Timestamp outside of allowed time-interval.
Definition at line 242 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ENUMSFTS 8 |
Inconsistent number of SFTs in timestamps and noise-SFTs.
Definition at line 243 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EINCONSBAND 9 |
Inconsistent values of sampling-rate and Tsft.
Definition at line 244 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ENOISEDELTAF 10 |
Frequency resolution of noise-SFTs inconsistent with signal.
Definition at line 245 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ENOISEBAND 11 |
Frequency band of noise-SFTs inconsistent with signal.
Definition at line 246 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ENOISEBINS 12 |
Frequency bins of noise-SFTs inconsistent with signal.
Definition at line 247 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EBADCOORDS 13 |
Current code requires sky position in equatorial coordinates.
Definition at line 248 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_ELUTS 14 |
Lookup tables (LUTs) for trig functions must be defined on domain -2pi to 2pi inclusive.
Definition at line 249 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EDTERMS 15 |
Dterms must be greater than zero and less than or equal to half the number of SFT bins.
Definition at line 250 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EINPUT 16 |
Invalid input-arguments to function.
Definition at line 251 of file GeneratePulsarSignal.h.
#define GENERATEPULSARSIGNALH_EDETECTOR 17 |
Unknown detector-name.
Definition at line 252 of file GeneratePulsarSignal.h.