28#ifndef _HETERODYNE_PULSAR_H
29#define _HETERODYNE_PULSAR_H
41#ifndef HAVE_LOCALTIME_R
42#define localtime_r(timep, result) memcpy((result), localtime(timep), sizeof(struct tm))
45#include <lal/BinaryPulsarTiming.h>
46#include <lal/LALgetopt.h>
47#include <lal/LogPrintf.h>
48#include <lal/FileIO.h>
49#include <lal/LALStdlib.h>
50#include <lal/LALAtomicDatatypes.h>
51#include <lal/LALDatatypes.h>
52#include <lal/AVFactories.h>
53#include <lal/LALCache.h>
54#include <lal/LALFrStream.h>
55#include <lal/IIRFilter.h>
56#include <lal/ZPGFilter.h>
57#include <lal/LALBarycenter.h>
58#include <lal/LALInitBarycenter.h>
59#include <lal/SkyCoordinates.h>
60#include <lal/DetectorSite.h>
61#include <lal/DetResponse.h>
62#include <lal/BandPassTimeSeries.h>
63#include <lal/FrequencySeries.h>
64#include <lal/RealFFT.h>
65#include <lal/ComplexFFT.h>
66#include <lal/SFTfileIO.h>
67#include <lal/LALString.h>
69#include <lal/TimeSeries.h>
70#include <lal/XLALError.h>
71#include <lal/LALPulsarVCSInfo.h>
79"Usage: %s [options]\n\n"\
80" --help (-h) display this message\n"\
81" --verbose (-v) display all error messages\n"\
82" --ifo (-i) name of ifo e.g. L1, H1, H2, G1\n"\
83" --pulsar (-p) name of pulsar e.g. J0534+2200\n"\
84" --heterodyne-flag (-z) (int) coarse heterodyne 0, fine heterodyne 1,\n\
85 parameter update 2, full heterodyne in one go 3,\n\
86 re-heterodyne an already fine heterodyned file 4\n"\
87" --param-file (-f) name of file containing initial pulsar parameters\n\
89" --param-file-update (-g) name of file containing updated pulsar parameters\n"\
90" --manual-epoch (-M) a hardwired epoch for the pulsar frequency and\n\
91 position (for use when dealing with hardware\n\
92 injections when this should be set to 751680013.0)\n"\
93" --ephem-earth-file (-e) name of file containing the earth ephemeris\n"\
94" --ephem-sun-file (-S) name of file containing the sun ephemeris\n"\
95" --ephem-time-file (-t) name of file containing the Einstein delay\n\
97" --filter-knee (-k) knee frequency of low-pass filter (don't filter\n\
99" --sample-rate (-s) sample rate of input data\n"\
100" --resample-rate (-r) sample rate for output data (i.e. 1, 1/60,\n\
102" --data-file (-d) file containing list of frame files (in frame\n\
103 cache format) or previously heterodyned data\n\
105" --data-chunk-length (-D) maximum length of data (in seconds) to be read in\n\
106 at one time from a frame file, overriding\n\
108" --channel (-c) frame data channel (i.e. LSC-DARM_ERR)\n"\
109" --output-file (-o) full path and filename for the output data\n"\
110" --seg-file (-l) name of file containing science segment list\n"\
111" --calibrate (-A) if specified calibrate data (no argument)\n"\
112" --response-file (-R) name of file containing the response function\n"\
113" --coefficient-file (-C) name of file containing the calibration\n\
115" --sensing-function (-F) name of sensing function file (CAV_GAIN)\n"\
116" --open-loop-gain (-O) name of open loop gain file(OLOOP_GAIN)\n"\
117" --stddev-thresh (-T) standard deviation threshold veto for outliers\n"\
118" --freq-factor (-m) factor which mulitplies the pulsar spin frequency\n\
119 (default 2.0 i.e. a signal from a triaxial pulsar)\n"\
120" --scale-factor (-G) factor to scale the calibrated h(t) data by\n"\
121" --high-pass-freq (-H) high-pass frequency for calibrated h(t) data\n"\
122" --binary-input (-B) read in input data from binary file (for fine and\n\
123 update heterodynes only)\n"\
124" --binary-output (-b) output data to a binary file\n"\
125" --legacy-input (-L) if input heterodyned file is a legacy files produced\n\
126 before the introduction of headers, then set this\n\
128" --gzip-output (-Z) if not outputting in binary format then gzip the\n\
129 output file. This will be done by default if the\n\
130 --output-file specified has the \".gz\" suffix, but\n\
131 if not this suffix will be appended\n"\
132" --output-phase (-P) if set, output the phase evolution to a text file\n\
133 (for debugging purposes)\n"\
136#define MAXDATALENGTH 256
137#define MAXSTRLENGTH 1024
138#define MAXLISTLENGTH 25000
143#define FILTERFFTTIME 200
145#define HEADERSIZE 2048
149typedef struct tagCalibrationFiles {
157typedef struct tagFrameCache {
164typedef struct tagInputParams {
207typedef struct tagHeterodyneParams {
226typedef struct tagFilters {
235typedef struct tagFilterResponse {
TimeCorrectionType
Enumerated type denoting the time system type to be produced in the solar system barycentring routine...
void destroy_filter_response(FilterResponse *filtresp)
void set_filters(Filters *iirFilters, REAL8 filterKnee, REAL8 samplerate)
INT4 remove_outliers(COMPLEX16TimeSeries *data, REAL8Vector *times, REAL8 stddevthresh)
void calibrate(COMPLEX16TimeSeries *series, REAL8Vector *datatimes, CalibrationFiles calfiles, REAL8 frequency, CHAR *channel)
void heterodyne_data(COMPLEX16TimeSeries *data, REAL8Vector *times, HeterodyneParams hetParams, REAL8 freqfactor, FilterResponse *filtResp)
void get_calibration_values(REAL8 *magnitude, REAL8 *phase, CHAR *calibfilename, REAL8 frequency)
LALCache * set_frame_files(INT4 *starts, INT4 *stops, LALCache *cache, INT4 *position, INT4 maxchunklength)
FilterResponse * create_filter_response(REAL8 filterKnee)
INT4 get_segment_list(INT4Vector *starts, INT4Vector *stops, CHAR *seglistfile, INT4 heterodyneflag)
REAL8TimeSeries * get_frame_data(LALCache *framecache, CHAR *channel, REAL8 gpstime, REAL8 length, INT4 duration, REAL8 samplerate, REAL8 scalefac, REAL8 highpass)
void get_input_args(InputParams *inputParams, int argc, char *argv[])
void get_frame_times(CHAR *framefile, REAL8 *gpstime, INT4 *duration)
COMPLEX16TimeSeries * resample_data(COMPLEX16TimeSeries *data, REAL8Vector *times, INT4Vector *starts, INT4Vector *stops, REAL8 sampleRate, REAL8 resampleRate, INT4 heterodyneflag)
void filter_data(COMPLEX16TimeSeries *data, Filters *iirFilters)
def phase(point, coeffs, params, ignoreintcheck=False)
CHAR * calibcoefficientfile
CHAR * sensingfunctionfile
CHAR * responsefunctionfile
REAL8IIRFilter * filter2Im
REAL8IIRFilter * filter3Im
REAL8IIRFilter * filter3Re
REAL8IIRFilter * filter1Re
REAL8IIRFilter * filter1Im
REAL8IIRFilter * filter2Re
PulsarParameters * hetUpdate
The PulsarParameters structure to contain a set of pulsar parameters.