36 #ifndef _READPULSARPARFILE_H
37 #define _READPULSARPARFILE_H
41 #include <lal/LALStdlib.h>
42 #include <lal/StringVector.h>
43 #include <lal/LALBarycenter.h>
45 #include <lal/LALHashTbl.h>
52 #define READPULSARPARFILEH_ENULLOUTPUT 1
54 #define READPULSARPARFILEH_MSGENULLOUTPUT "Output was Null"
57 #define PULSAR_HASHTABLE_SIZE 512
58 #define PULSAR_PARNAME_MAX 128
60 #define GPS0MJD 44244.0
63 #define TDT_TAI 32.184
65 #define GPS_TDT (TDT_TAI + XLAL_EPOCH_GPS_TAI_UTC)
69 #define LALPULSAR_TEMPO2_MTSUN_SI 4.925490947e-6
70 #define LALPULSAR_TEMPO2_MSUN_SI (LALPULSAR_TEMPO2_MTSUN_SI * LAL_MPL_SI / LAL_TPL_SI)
73 typedef enum tagPulsarParamType {
90 typedef struct tagPulsarParam {
96 struct tagPulsarParam *
next;
106 typedef struct tagPulsarParameters {
117 typedef struct tagBinaryPulsarParams {
442 SWIGLAL( OWNS_THIS_ARG(
const CHAR *, value ) );
458 SWIGLAL_CLEAR( OWNS_THIS_ARG(
const CHAR *, value ) );
462 SWIGLAL( OWNS_THIS_ARG(
const REAL8Vector *, value ) );
475 SWIGLAL_CLEAR( OWNS_THIS_ARG(
const REAL8Vector *, value ) );
595 CHAR *pulsarAndPath );
void PulsarSetParam(PulsarParameters *pars, const CHAR *name, const void *value)
Set the value of a parameter in the PulsarParameters structure.
void ParConvKpcToMetres(const CHAR *in, void *out)
Convert the input string from kiloparsecs to metres.
void * PulsarGetParamErr(const PulsarParameters *pars, const CHAR *name)
Get the required parameter error value from the PulsarParameters structure.
PulsarParameters * XLALReadTEMPOParFile(const CHAR *pulsarAndPath)
Read in the parameters from a TEMPO(2) parameter file into a PulsarParameters structure.
void ParConvToString(const CHAR *in, void *out)
Copy the input string into the output pointer.
void ParConvDegsToRads(const CHAR *in, void *out)
Convert the input string from degrees to radians.
LALStringVector * XLALReadTEMPOCorFile(REAL8Array *cormat, CHAR *corfile)
This function will read in a TEMPO-style parameter correlation matrix.
void ParConvInvArcsecsToInvRads(const CHAR *in, void *out)
Convert the input string from 1/acrseconds to 1/radians.
void ParConvDegPerYrToRadPerSec(const CHAR *in, void *out)
Convert the input string from degrees per year to radians per second.
void PulsarSetREAL8ParamErr(PulsarParameters *pars, const CHAR *name, REAL8 value, UINT4 fitFlag)
Set the error value for a REAL8 parameter.
void PulsarAddREAL8Param(PulsarParameters *pars, const CHAR *name, REAL8 value)
Add a REAL8 parameter to the PulsarParameters structure.
void PulsarCopyParams(PulsarParameters *origin, PulsarParameters *target)
Function to copy a PulsarParameters structure.
const UINT4 * PulsarGetParamFitFlag(const PulsarParameters *pars, const CHAR *name)
Get the fit flag array for a given parameter from the PulsarParameters structure.
void ParConvSecsToRads(const CHAR *in, void *out)
Convert the input string from seconds to radians.
REAL8 PulsarGetREAL8VectorParamErrIndividual(const PulsarParameters *pars, const CHAR *name)
Return an individual REAL8 value from a REAL8Vector parameter.
void ParConvBinaryUnits(const CHAR *in, void *out)
Convert the binary system parameter from a string to a double, but make the check (as performed by TE...
REAL8 PulsarGetREAL8ParamErr(const PulsarParameters *pars, const CHAR *name)
Return a REAL8 parameter error value.
const REAL8Vector * PulsarGetREAL8VectorParam(const PulsarParameters *pars, const CHAR *name)
Return a REAL8Vector parameter.
REAL8 PulsarGetREAL8ParamOrZero(const PulsarParameters *pars, const CHAR *name)
Return a REAL8 parameter if it exists, otherwise return zero.
int PulsarCheckParam(const PulsarParameters *pars, const CHAR *name)
Check for the existence of the parameter name in the PulsarParameters structure.
void PulsarSetREAL8VectorParamErr(PulsarParameters *pars, const CHAR *name, const REAL8Vector *value, const UINT4 *fitFlag)
Set the error values for a REAL8Vector parameter.
REAL8 XLALTCBMJDtoGPS(REAL8 MJD)
If you have an MJD arrival time on the Earth then this will convert it to the equivalent GPS time in ...
const UINT4Vector * PulsarGetParamFitFlagAsVector(const PulsarParameters *pars, const CHAR *name)
Get the fit flag array for a given parameter from the PulsarParameters structure.
void PulsarClearParams(PulsarParameters *pars)
Free all the parameters from a PulsarParameters structure.
void PulsarSetParamErr(PulsarParameters *pars, const CHAR *name, void *value, const UINT4 *fitFlag, UINT4 len)
Set the value of the error of a parameter in the PulsarParameters structure.
void ParConvDaysToSecs(const CHAR *in, void *out)
Convert the input string from days to seconds.
const CHAR * PulsarGetStringParam(const PulsarParameters *pars, const CHAR *name)
Return a string parameter.
void PrintPulsarParameters(BinaryPulsarParams params)
function to print out all the pulsar parameters read in from a par file
const REAL8Vector * PulsarGetREAL8VectorParamErr(const PulsarParameters *pars, const CHAR *name)
Return a REAL8Vector parameter error value.
void PulsarAddParam(PulsarParameters *pars, const CHAR *name, void *value, PulsarParamType type)
Add a parameter and value to the PulsarParameters structure.
void ParConvSolarMassToKg(const CHAR *in, void *out)
Convert the input string from solar masses to kilograms.
REAL8 PulsarGetREAL8Param(const PulsarParameters *pars, const CHAR *name)
Return a REAL8 parameter.
void PulsarRemoveParam(PulsarParameters *pars, const CHAR *name)
Remove a given parameter from a PulsarParameters structure.
REAL8 XLALTDBMJDtoGPS(REAL8 MJD)
If you have an MJD arrival time on the Earth then this will convert it to the equivalent GPS time in ...
UINT4 PulsarGetUINT4Param(const PulsarParameters *pars, const CHAR *name)
Return a UINT4 parameter.
REAL8 XLALTTMJDtoGPS(REAL8 MJD)
Functions for converting times given in Terrestrial time TT, TDB, or TCB in MJD to times in GPS - thi...
void PulsarAddStringParam(PulsarParameters *pars, const CHAR *name, const CHAR *value)
Add a string parameter to the PulsarParameters structure.
PulsarParamType
An enumerated type for denoting the type of a variable.
@ PULSARTYPE_REAL8Vector_t
void ParConvMicrosecToSec(const CHAR *in, void *out)
Convert an input string from microseconds into seconds.
void PulsarAddREAL8VectorParam(PulsarParameters *pars, const CHAR *name, const REAL8Vector *value)
Add a REAL8Vector parameter to the PulsarParameters structure.
REAL8 PulsarGetREAL8VectorParamIndividual(const PulsarParameters *pars, const CHAR *name)
Return an individual REAL8 value from a REAL8Vector parameter.
void ParConvMasPerYrToRadPerSec(const CHAR *in, void *out)
Convert the input string from milliarcsecs per year to radians per second.
void ParConvMasToRads(const CHAR *in, void *out)
Convert the input string from milliarcsecs to radians.
void PulsarFreeParams(PulsarParameters *par)
Function to free memory from pulsar parameters.
void * PulsarGetParam(const PulsarParameters *pars, const CHAR *name)
Get the required parameter value from the PulsarParameters structure.
void XLALReadTEMPOParFileOrig(BinaryPulsarParams *output, CHAR *pulsarAndPath)
function to read in a TEMPO parameter file
void ParConvToInt(const CHAR *in, void *out)
Convert the input string into a unsigned integer number.
void ParConvMJDToGPS(const CHAR *in, void *out)
Convert the input string from a TT MJD value into a GPS time.
void ParConvRAToRads(const CHAR *in, void *out)
Convert a right ascension input string in the format "hh:mm:ss.s" into radians.
void ParConvToFloat(const CHAR *in, void *out)
Conversion functions from units used in TEMPO parameter files into SI units.
PulsarParamType PulsarGetParamType(const PulsarParameters *pars, const char *name)
Get the required parameter's type.
void ParConvDecToRads(const CHAR *in, void *out)
Convert a declination input string in the format "dd:mm:ss.s" into radians.
#define PULSAR_PARNAME_MAX
void ParConvArcsecsToRads(const CHAR *in, void *out)
Convert the input string from arcseconds to radians.
UINT4 PulsarGetUINT4ParamOrZero(const PulsarParameters *pars, const CHAR *name)
Return a UINT4 parameter if it exists, otherwise return zero.
void PulsarAddUINT4Param(PulsarParameters *pars, const CHAR *name, UINT4 value)
Add a UINT4 parameter to the PulsarParameters structure.
A structure to contain all pulsar parameters and associated errors.
REAL8 phi0Vector
initial phase for vector modes
REAL8 iota
inclination angle
REAL8 w0
logitude of periastron (radians)
REAL8 hPlus
amplitude for the tensor plus polarisation
REAL8 hCross
amplitude for the tensor cross polarisation
REAL8 px
pulsar parallax (converted from milliarcsecs to rads)
REAL8 m2
companion mass (in kg)
REAL8 DM
dispersion measure
REAL8 f7
frequency seventh derivative (Hz/s^7)
REAL8 psiVector
phase polarisation angle for vector polarisation mode
REAL8 f6
frequency sixth derivative (Hz/s^6)
REAL8 pepoch
period/frequency epoch
REAL8 daop
parameter for the Kopeikin annual orbital parallax
REAL8 edot
rate of change of e
REAL8 f5
frequency fifth derivative (Hz/s^5)
REAL8 Q22
gravitational wave l=m=2 mass quadrupole moment
REAL8 x
projected semi-major axis/speed of light (light secs)
REAL8 f9
frequency ninth derivative (Hz/s^9)
REAL8 e
orbital eccentricity
REAL8 lambdapin
this is a longitude like angle between pinning axis and line of sight
REAL8 f2
frequency second derivative (Hz/s^2)
REAL8 hScalarL
amplitude for scalar longitudinal polarisation mode
CHAR * ephem
The JPL solar system ephemeris used e.g.
INT4 nfb
the number of fb coefficients
REAL8 Aplus
0.5*h0*(1+cos^2iota)
UINT4 nEll
set to zero if have eps time derivs (default) set to 1 if have wdot
REAL8 M
M = m1 + m2 (m1 = pulsar mass, m2 = companion mass) (in kg)
REAL8 ra
right ascension (rads)
REAL8 s
Shapiro 'shape' parameter sin i.
REAL8 I21
parameter for pinsf model.
REAL8 xdot
rate of change of x(=asini/c) - optional
REAL8 f3
frequency third derivative (Hz/s^3)
REAL8 f8
frequency eighth derivative (Hz/s^8)
REAL8 I31
parameter for pinsf model.
REAL8 hVectorY
amplitude for vector "y" polarisation mode
REAL8 hScalarB
amplitude for scalar breathing polarisation mode
REAL8 psiTensor
phase polarisation angle for tensor polarisation mode
REAL8 T0
time of orbital perisastron as measured in TDB (MJD)
REAL8 * fb
orbital frequency coefficients for BTX model
REAL8 wave_om
fundamental frequency timing noise terms
REAL8 cosiota
cosine of the pulsars inclination angle
INT4 daopset
set if daop is set from the par file
REAL8 h0
gravitational wave amplitude
REAL8 f4
frequency fourth derivative (Hz/s^4)
REAL8 hVectorX
amplitude for vector "x" polarisation mode
REAL8 Pb
orbital period (seconds)
REAL8 psiScalar
phase polarisation angle for scalar polarisation mode
REAL8 phi0Tensor
initial phase for tensor modes (to be used instead of phi0 or phi22)
CHAR * model
TEMPO binary model e.g.
CHAR * units
The time system used e.g.
REAL8 cgw
The speed of gravitational waves as a fraction of the speed of light LAL_C_SI
REAL8 pmdec
proper motion in dec (rads/s)
REAL8 pmra
proper motion in RA (rads/s)
REAL8 Tasc
time of the ascending node (used rather than T0)
REAL8 Pbdot
rate of change of Pb (dimensionless)
CHAR * jname
pulsar J name
REAL8 DM1
first derivative of dispersion measure
REAL8 startTime
start of parfile applicable time
REAL8 posepoch
position epoch
REAL8 f10
frequency tenth derivative (Hz/s^10)
REAL8 dec
declination (rads)
REAL8 wdot
precesion of longitude of periastron w = w0 + wdot(tb-T0) (degs/year)
REAL8 finishTime
finish of parfile applicable time
REAL8 f0
spin frequency (Hz)
REAL8 f1
frequency first derivative (Hz/s)
REAL8 phi0Scalar
initial phase for the scalar modes
REAL8 gamma
gravitational redshift and time dilation parameter (s)
REAL8 kin
binary inclination angle
REAL8 costheta
cosine of angle between rotation axis and pinning axis
REAL8 psi
polarisation angle
REAL8 dist
pulsar distance (in kiloparsecs)
CHAR * bname
pulsar B name
The PulsarParam list node structure.
void * value
Parameter value.
struct tagPulsarParam * next
void * err
Parameter error/uncertainty.
UINT4Vector * fitFlag
Set to 1 if the parameter has been fit in the par file.
PulsarParamType type
Parameter type e.g.
The PulsarParameters structure to contain a set of pulsar parameters.
PulsarParam * head
A linked list of PulsarParam structures.
LALHashTbl * hash_table
Hash table of parameters.
INT4 nparams
The total number of parameters in the structure.