Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALPulsar 7.1.1.1-b246709
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
heterodyne_pulsar.h
Go to the documentation of this file.
1/*
2* Copyright (C) 2007 Matt Pitkin
3*
4* This program is free software; you can redistribute it and/or modify
5* it under the terms of the GNU General Public License as published by
6* the Free Software Foundation; either version 2 of the License, or
7* (at your option) any later version.
8*
9* This program is distributed in the hope that it will be useful,
10* but WITHOUT ANY WARRANTY; without even the implied warranty of
11* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12* GNU General Public License for more details.
13*
14* You should have received a copy of the GNU General Public License
15* along with with program; see the file COPYING. If not, write to the
16* Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17* MA 02110-1301 USA
18*/
19
20/*
21 Author: Pitkin, M. D.
22*/
23
24/* Matt Pitkin 09/02/06 -------------- heterodyne_pulsar.h */
25
26/* header file for heterodyne_pulsar.c */
27
28#ifndef _HETERODYNE_PULSAR_H
29#define _HETERODYNE_PULSAR_H
30
31#include "config.h"
32
33#include <stdio.h>
34#include <stdlib.h>
35#include <math.h>
36#include <string.h>
37#include <time.h>
38#include <sys/types.h>
39
40// definition taken from lal/lib/support/LogPrintf.c
41#ifndef HAVE_LOCALTIME_R
42#define localtime_r(timep, result) memcpy((result), localtime(timep), sizeof(struct tm))
43#endif
44
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>
68#include <lal/Units.h>
69#include <lal/TimeSeries.h>
70#include <lal/XLALError.h>
71#include <lal/LALPulsarVCSInfo.h>
72
73#ifdef __cplusplus
74extern "C" {
75#endif
76
77/* Usage format string. */
78#define USAGE \
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\
88 (.par file)\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\
96 correction values\n"\
97" --filter-knee (-k) knee frequency of low-pass filter (don't filter\n\
98 if = 0)\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\
101 8192 etc) \n"\
102" --data-file (-d) file containing list of frame files (in frame\n\
103 cache format) or previously heterodyned data\n\
104 file\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\
107 MAXDATALENGTH\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\
114 coefficients\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\
127 flag\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"\
134"\n"
135
136#define MAXDATALENGTH 256 /* maximum length of data to be read from frames */
137#define MAXSTRLENGTH 1024 /* maximum number of characters in a frame filename */
138#define MAXLISTLENGTH 25000 /* maximum length of a list of frames files */
139
140#define ALPHAMIN 0.0 /* minimum acceptable value of alpha calib coefficient */
141#define ALPHAMAX 2.0 /* maximum acceptable value of alpha calib coefficient */
142
143#define FILTERFFTTIME 200
144
145#define HEADERSIZE 2048 /* number of bytes in header for output files */
146
147/* define structures */
148
149typedef struct tagCalibrationFiles {
155
156/* structure to store data from a lal frame cache as output from LSCdataFind */
157typedef struct tagFrameCache {
158 CHAR **framelist; /* list of file names in frame cache file */
159 INT4 *duration; /* duration of each frame file */
160 INT4 *starttime; /* start time of each frame file */
161 UINT4 length; /* length of cache */
162} FrameCache;
163
164typedef struct tagInputParams {
167
169 CHAR paramfile[256];
170 CHAR paramfileupdate[256];
172
174
175 CHAR earthfile[256];
176 CHAR sunfile[256];
178
182
186
187 CHAR outputfile[256];
188 CHAR segfile[256];
189
192
194
197
199
206
207typedef struct tagHeterodyneParams {
210
212
214
218
219 CHAR earthfile[256];
220 CHAR sunfile[256];
225
226typedef struct tagFilters {
227 REAL8IIRFilter *filter1Re; /* filters for real and imaginary parts of heterodyed data */
233} Filters;
234
235typedef struct tagFilterResponse {
238
239 REAL8 srate; /* sample rate */
240 REAL8 deltaf; /* frequency step between successive points */
242
243/* define functions */
244void get_input_args( InputParams *inputParams, int argc, char *argv[] );
245
247 REAL8 freqfactor, FilterResponse *filtResp );
248
249void set_filters( Filters *iirFilters, REAL8 filterKnee, REAL8 samplerate );
250
251void filter_data( COMPLEX16TimeSeries *data, Filters *iirFilters );
252
254 *starts, INT4Vector *stops, REAL8 sampleRate, REAL8 resampleRate, INT4 heterodyneflag );
255
256/* function to extract the frame time and duration from the file name */
257void get_frame_times( CHAR *framefile, REAL8 *gpstime, INT4 *duration );
258
259/* reads in a time series from frames */
261 REAL8 length, INT4 duration, REAL8 samplerate, REAL8 scalefac,
262 REAL8 highpass );
263
264/* read in science segment list file - returns the number of segments */
265INT4 get_segment_list( INT4Vector *starts, INT4Vector *stops, CHAR *seglistfile, INT4 heterodyneflag );
266
267/* get frame data for partcular science segment */
268LALCache *set_frame_files( INT4 *starts, INT4 *stops, LALCache *cache, INT4 *position, INT4 maxchunklength );
269
270/* calibrate data */
271void calibrate( COMPLEX16TimeSeries *series, REAL8Vector *datatimes, CalibrationFiles calfiles,
273
274/* function to extract the correct calibration values from the files */
275void get_calibration_values( REAL8 *magnitude, REAL8 *phase, CHAR *calibfilename, REAL8
276 frequency );
277
278/* function to remove outliers above a certain standard deviation threshold - returns the number
279of outliers removed */
281
283
284/* free memory for filter response structure */
286
287#ifdef __cplusplus
288}
289#endif
290
291#endif /* _HETERODYNE_PULSAR_H */
TimeCorrectionType
Enumerated type denoting the time system type to be produced in the solar system barycentring routine...
Definition: LALBarycenter.h:72
double REAL8
char CHAR
uint32_t UINT4
int32_t INT4
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)
float data[BLOCKSIZE]
Definition: hwinject.c:360
int gpstime
Definition: hwinject.c:365
def phase(point, coeffs, params, ignoreintcheck=False)
REAL8Vector * freqResp
REAL8Vector * phaseResp
REAL8IIRFilter * filter2Im
REAL8IIRFilter * filter3Im
REAL8IIRFilter * filter3Re
REAL8IIRFilter * filter1Re
REAL8IIRFilter * filter1Im
REAL8IIRFilter * filter2Re
PulsarParameters * hetUpdate
PulsarParameters * het
TimeCorrectionType ttype
CalibrationFiles calibfiles
The PulsarParameters structure to contain a set of pulsar parameters.