26 #ifndef _GENERATEPULSARSIGNAL_H
27 #define _GENERATEPULSARSIGNAL_H
29 #include <lal/LALDatatypes.h>
30 #include <lal/DetResponse.h>
31 #include <lal/DetectorSite.h>
32 #include <lal/GenerateSpinOrbitCW.h>
34 #include <lal/LALBarycenter.h>
35 #include <lal/PulsarDataTypes.h>
36 #include <lal/ComputeSky.h>
37 #include <lal/Window.h>
38 #include <lal/SFTfileIO.h>
236 #define GENERATEPULSARSIGNALH_ENULL 1
237 #define GENERATEPULSARSIGNALH_ENONULL 2
238 #define GENERATEPULSARSIGNALH_EMEM 3
239 #define GENERATEPULSARSIGNALH_ESAMPLING 4
240 #define GENERATEPULSARSIGNALH_ESSBCONVERT 5
241 #define GENERATEPULSARSIGNALH_ESYS 6
242 #define GENERATEPULSARSIGNALH_ETIMEBOUND 7
243 #define GENERATEPULSARSIGNALH_ENUMSFTS 8
244 #define GENERATEPULSARSIGNALH_EINCONSBAND 9
245 #define GENERATEPULSARSIGNALH_ENOISEDELTAF 10
246 #define GENERATEPULSARSIGNALH_ENOISEBAND 11
247 #define GENERATEPULSARSIGNALH_ENOISEBINS 12
248 #define GENERATEPULSARSIGNALH_EBADCOORDS 13
249 #define GENERATEPULSARSIGNALH_ELUTS 14
250 #define GENERATEPULSARSIGNALH_EDTERMS 15
251 #define GENERATEPULSARSIGNALH_EINPUT 16
252 #define GENERATEPULSARSIGNALH_EDETECTOR 17
256 #define GENERATEPULSARSIGNALH_MSGENULL "Arguments contained an unexpected null pointer"
257 #define GENERATEPULSARSIGNALH_MSGENONULL "Output pointer is not NULL"
258 #define GENERATEPULSARSIGNALH_MSGEMEM "Out of memory"
259 #define GENERATEPULSARSIGNALH_MSGESAMPLING "Waveform sampling interval too large."
260 #define GENERATEPULSARSIGNALH_MSGESSBCONVERT "SSB->GPS iterative conversion failed"
261 #define GENERATEPULSARSIGNALH_MSGESYS "System error, probably while File I/O"
262 #define GENERATEPULSARSIGNALH_MSGETIMEBOUND "Timestamp outside of allowed time-interval"
263 #define GENERATEPULSARSIGNALH_MSGENUMSFTS "Inconsistent number of SFTs in timestamps and noise-SFTs"
264 #define GENERATEPULSARSIGNALH_MSGEINCONSBAND "Inconsistent values of sampling-rate and Tsft"
265 #define GENERATEPULSARSIGNALH_MSGENOISEDELTAF "Frequency resolution of noise-SFTs inconsistent with signal"
266 #define GENERATEPULSARSIGNALH_MSGENOISEBAND "Frequency band of noise-SFTs inconsistent with signal"
267 #define GENERATEPULSARSIGNALH_MSGENOISEBINS "Frequency bins of noise-SFTs inconsistent with signal"
268 #define GENERATEPULSARSIGNALH_MSGEBADCOORDS "Current code requires sky position in equatorial coordinates"
269 #define GENERATEPULSARSIGNALH_MSGELUTS "Lookup tables (LUTs) for trig functions must be defined on domain -2pi to 2pi inclusive"
270 #define GENERATEPULSARSIGNALH_MSGEDTERMS "Dterms must be greater than zero and less than or equal to half the number of SFT bins"
271 #define GENERATEPULSARSIGNALH_MSGEINPUT "Invalid input-arguments to function"
272 #define GENERATEPULSARSIGNALH_MSGEDETECTOR "Unknown detector-name"
279 typedef struct tagPulsarSignalParams {
319 SWIGLAL( IMMUTABLE_MEMBERS( tagSFTParams,
timestamps, noiseSFTs, window ) );
321 typedef struct tagSFTParams {
333 typedef struct tagSFTandSignalParams {
349 typedef struct tagSkyConstAndZeroPsiAMResponse {
LIGOTimeGPSVector * timestamps
void XLALDestroyMultiREAL4TimeSeries(MultiREAL4TimeSeries *multiTS)
Destroy a MultiREAL4TimeSeries, NULL-robust.
int XLALConvertGPS2SSB(LIGOTimeGPS *SSBout, LIGOTimeGPS GPSin, const PulsarSignalParams *params)
Convert earth-frame GPS time into barycentric-frame SSB time for given source.
REAL4TimeSeries * XLALGeneratePulsarSignal(const PulsarSignalParams *params)
Generate a time-series at the detector for a given pulsar.
int XLALConvertSSB2GPS(LIGOTimeGPS *GPSout, LIGOTimeGPS SSBin, const PulsarSignalParams *params)
Convert barycentric frame SSB time into earth-frame GPS time.
void LALGeneratePulsarSignal(LALStatus *status, REAL4TimeSeries **signalvec, const PulsarSignalParams *params)
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.
void XLALDestroyMultiREAL8TimeSeries(MultiREAL8TimeSeries *multiTS)
Destroy a MultiREAL8TimeSeries, NULL-robust.
void LALFastGeneratePulsarSFTs(LALStatus *status, SFTVector **outputSFTs, const SkyConstAndZeroPsiAMResponse *input, const SFTandSignalParams *params)
Fast generation of Fake SFTs for a pure pulsar signal.
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 s...
int XLALAddGaussianNoise(REAL4TimeSeries *inSeries, REAL4 sigma, INT4 seed)
Generate Gaussian noise with standard-deviation sigma, add it to inSeries.
A vector of COMPLEX8FrequencySeries.
This structure contains all information about the center-of-mass positions of the Earth and Sun,...
A vector of 'timestamps' of type LIGOTimeGPS.
A collection of (multi-IFO) REAL4 time-series.
A collection of (multi-IFO) REAL8 time-series.
Input parameters to GeneratePulsarSignal(), defining the source and the time-series.
LIGOTimeGPS tp
time of observed periapsis passage (in SSB)
REAL8 asini
projected, normalized orbital semi-major axis (s)
SkyPosition position
source location (in radians)
const COMPLEX8FrequencySeries * transfer
detector transfer function (NULL if not used)
REAL8 period
orbital period (sec)
REAL8 sourceDeltaT
source-frame sampling period ('0' means use previous internal defaults)
const EphemerisData * ephemerides
Earth and Sun ephemerides.
UINT4 dtPolBy2
half-interval for the polarisation response look-up table for LALPulsarSimulateCoherentGW()
REAL8 fHeterodyne
heterodyning frequency for output time-series
REAL4 aPlus
plus-polarization amplitude at tRef
REAL4 psi
polarization angle (radians) at tRef
REAL8 argp
argument of periapsis (radians)
REAL8 phi0
initial phase (radians) at tRef
LIGOTimeGPS startTimeGPS
start time of output time series
UINT4 dtDelayBy2
half-interval for the Doppler delay look-up table for LALPulsarSimulateCoherentGW()
const LALDetector * site
detector location and orientation
REAL8Vector * spindown
wave-frequency spindowns at tRef (NOT f0-normalized!)
REAL8 ecc
orbital eccentricity
LIGOTimeGPS refTime
reference time of pulsar parameters (in SSB!)
REAL8 f0
WAVE-frequency(!) at tRef (in Hz)
REAL4 aCross
cross-polarization amplitude at tRef
REAL8 samplingRate
sampling rate of time-series (= 2 * frequency-Band)
UINT4 duration
length of time series in seconds
Parameters defining the SFTs to be returned from LALSignalToSFTs().
REAL8 Tsft
length of each SFT in seconds
const REAL4Window * window
window function for the time series (can be NULL, which is equivalent to a rectangular window)
const LIGOTimeGPSVector * timestamps
timestamps to produce SFTs for (can be NULL)
const SFTVector * noiseSFTs
noise SFTs to be added (can be NULL)
Parameters defining the pulsar signal and SFTs used by LALFastGeneratePulsarSFTs().
REAL8 * sinVal
sinVal holds lookup table (LUT) values for doing trig sin calls
INT4 nSamples
nsample from noise SFT header; 2x this equals effective number of time samples
PulsarSignalParams * pSigParams
REAL8 * cosVal
cosVal holds lookup table (LUT) values for doing trig cos calls
REAL8 * trigArg
array of arguments to hold lookup table (LUT) values for doing trig calls
INT4 Dterms
use this to fill in SFT bins with fake data as per LALDemod else fill in bin with zero
INT4 resTrig
length sinVal, cosVal; domain: -2pi to 2pi; resolution = 4pi/resTrig
Sky Constants and beam pattern response functions used by LALFastGeneratePulsarSFTs().
REAL8 * skyConst
vector of A and B sky constants
REAL4 * fPlusZeroPsi
vector of Fplus values for psi = 0 at midpoint of each SFT
REAL4 * fCrossZeroPsi
vector of Fcross values for psi = 0 at midpoint of each SFT