Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALPulsar 7.1.1.1-3a66518
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
computeSignalDetector.c
Go to the documentation of this file.
1#include <lal/DopplerScan.h>
2#include <lal/LALInitBarycenter.h>
3#include <lal/UserInput.h>
4#include <lal/LALString.h>
5#include <gsl/gsl_sf_trig.h>
6#include <gsl/gsl_rng.h>
7#include "../TwoSpectSpecFunc.h"
8
9typedef struct {
10 REAL8 Tsft;
11 REAL8 SFToverlap;
13 REAL8 Tobs;
14 REAL8 cosi;
15 REAL8 psi;
21 CHAR *ephemEarth;
22 CHAR *ephemSun;
26
27INT4 InitUserVars( UserVariables_t *uvar, int argc, char *argv[] );
28
29int main( int argc, char *argv[] )
30{
32 XLAL_CHECK( InitUserVars( &uvar, argc, argv ) == XLAL_SUCCESS, XLAL_EFUNC );
33
35 XLAL_CHECK( ( detectors = XLALMalloc( sizeof( MultiLALDetector ) ) ) != NULL, XLAL_ENOMEM );
36 detectors->length = uvar.IFO->length;
37 for ( UINT4 ii = 0; ii < detectors->length; ii++ ) {
38 if ( strcmp( "H1", uvar.IFO->data[ii] ) == 0 ) {
40 } else if ( strcmp( "L1", uvar.IFO->data[ii] ) == 0 ) {
42 } else if ( strcmp( "V1", uvar.IFO->data[ii] ) == 0 ) {
44 } else if ( strcmp( "H2", uvar.IFO->data[ii] ) == 0 ) {
46 } else if ( strcmp( "H2r", uvar.IFO->data[ii] ) == 0 ) {
50 memset( &( H2.frDetector.name ), 0, sizeof( CHAR )*LALNameLength );
51 snprintf( H2.frDetector.name, LALNameLength, "%s", "LHO_2k_rotatedPiOver4" );
52 XLAL_CHECK( ( XLALCreateDetector( &( detectors->sites[ii] ), &( H2.frDetector ), LALDETECTORTYPE_IFODIFF ) ) != NULL, XLAL_EFUNC );
53 } else {
54 XLAL_ERROR( XLAL_EINVAL, "Not using valid interferometer! Expected 'H1', 'H2', 'H2r' (rotated H2), 'L1', or 'V1' not %s.\n", uvar.IFO->data[ii] );
55 }
56 }
57
58 EphemerisData *edat = NULL;
59 XLAL_CHECK( ( edat = XLALInitBarycenter( uvar.ephemEarth, uvar.ephemSun ) ) != NULL, XLAL_EFUNC );
60
61 LIGOTimeGPS tStart;
62 XLALGPSSetREAL8( &tStart, uvar.t0 );
63 XLAL_CHECK( xlalErrno == 0, XLAL_EFUNC, "XLALGPSSetREAL8 failed\n" );
64
66 XLAL_CHECK( ( multiTimestamps = XLALMakeMultiTimestamps( tStart, uvar.Tobs, uvar.Tsft, uvar.SFToverlap, detectors->length ) ) != NULL, XLAL_EFUNC );
67
68 LIGOTimeGPS refTime = multiTimestamps->data[0]->data[0];
69
70 MultiDetectorStateSeries *multiStateSeries = NULL;
71 XLAL_CHECK( ( multiStateSeries = XLALGetMultiDetectorStates( multiTimestamps, detectors, edat, uvar.SFToverlap ) ) != NULL, XLAL_EFUNC );
72
73 gsl_rng *rng = NULL;
74 XLAL_CHECK( ( rng = gsl_rng_alloc( gsl_rng_mt19937 ) ) != NULL, XLAL_EFUNC );
75 gsl_rng_set( rng, 0 );
76
77 FILE *OUTPUT;
78 XLAL_CHECK( ( OUTPUT = fopen( uvar.outfilename, "w" ) ) != NULL, XLAL_EIO, "Output file %s could not be opened\n", uvar.outfilename );
79
80 for ( INT4 n = 0; n < uvar.skylocations; n++ ) {
81 SkyPosition skypos;
82 if ( XLALUserVarWasSet( &( uvar.alpha ) ) && XLALUserVarWasSet( &( uvar.delta ) ) && n == 0 ) {
83 skypos.longitude = uvar.alpha;
84 skypos.latitude = uvar.delta;
86 } else {
87 skypos.longitude = LAL_TWOPI * gsl_rng_uniform( rng );
88 skypos.latitude = LAL_PI * gsl_rng_uniform( rng ) - LAL_PI_2;
90 }
91
92 REAL8 cosi0, psi0;
93 if ( XLALUserVarWasSet( &( uvar.cosi ) ) && n == 0 ) {
94 cosi0 = uvar.cosi;
95 } else {
96 cosi0 = 2.0 * gsl_rng_uniform( rng ) - 1.0;
97 }
98 if ( XLALUserVarWasSet( &( uvar.psi ) ) && n == 0 ) {
99 psi0 = uvar.psi;
100 } else {
101 psi0 = LAL_PI * gsl_rng_uniform( rng );
102 }
103
104 MultiAMCoeffs *multiAMcoefficients = NULL;
105 XLAL_CHECK( ( multiAMcoefficients = XLALComputeMultiAMCoeffs( multiStateSeries, NULL, skypos ) ) != NULL, XLAL_EFUNC );
106
107 MultiSSBtimes *multissb = NULL;
108 XLAL_CHECK( ( multissb = XLALGetMultiSSBtimes( multiStateSeries, skypos, refTime, SSBPREC_RELATIVISTICOPT ) ) != NULL, XLAL_EFUNC );
109
110 REAL8 frequency = 1000.0;
111 REAL8 frequency0 = frequency + ( gsl_rng_uniform( rng ) - 0.5 ) / uvar.Tsft;
112
113 for ( UINT4 ii = 0; ii < multiAMcoefficients->data[0]->a->length; ii++ ) {
114 REAL4 Fplus0 = multiAMcoefficients->data[0]->a->data[ii] * cos( 2.0 * psi0 ) + multiAMcoefficients->data[0]->b->data[ii] * sin( 2.0 * psi0 );
115 REAL4 Fcross0 = multiAMcoefficients->data[0]->b->data[ii] * cos( 2.0 * psi0 ) - multiAMcoefficients->data[0]->a->data[ii] * sin( 2.0 * psi0 );
116 REAL4 Fplus1 = multiAMcoefficients->data[1]->a->data[ii] * cos( 2.0 * psi0 ) + multiAMcoefficients->data[1]->b->data[ii] * sin( 2.0 * psi0 );
117 REAL4 Fcross1 = multiAMcoefficients->data[1]->b->data[ii] * cos( 2.0 * psi0 ) - multiAMcoefficients->data[1]->a->data[ii] * sin( 2.0 * psi0 );
118 COMPLEX16 RatioTerm0 = crect( 0.5 * Fplus1 * ( 1.0 + cosi0 * cosi0 ), Fcross1 * cosi0 ) / crect( 0.5 * Fplus0 * ( 1.0 + cosi0 * cosi0 ), Fcross0 * cosi0 ); //real det-sig ratio term
119
120 REAL4 detPhaseArg = 0.0, detPhaseMag = 0.0;
121 BOOLEAN loopbroken = 0;
122 for ( INT4 jj = 0; jj < 16 && !loopbroken; jj++ ) {
123 REAL4 psi = 0.0625 * jj * LAL_PI;
124 Fplus0 = multiAMcoefficients->data[0]->a->data[ii] * cos( 2.0 * psi ) + multiAMcoefficients->data[0]->b->data[ii] * sin( 2.0 * psi );
125 Fcross0 = multiAMcoefficients->data[0]->b->data[ii] * cos( 2.0 * psi ) - multiAMcoefficients->data[0]->a->data[ii] * sin( 2.0 * psi );
126 Fplus1 = multiAMcoefficients->data[1]->a->data[ii] * cos( 2.0 * psi ) + multiAMcoefficients->data[1]->b->data[ii] * sin( 2.0 * psi );
127 Fcross1 = multiAMcoefficients->data[1]->b->data[ii] * cos( 2.0 * psi ) - multiAMcoefficients->data[1]->a->data[ii] * sin( 2.0 * psi );
128 for ( INT4 kk = 0; kk < 21 && !loopbroken; kk++ ) {
129 REAL4 cosi = 1.0 - 2.0 * 0.05 * kk;
130 if ( !uvar.unrestrictedCosi ) {
131 if ( cosi0 < 0.0 ) {
132 cosi = -0.05 * kk;
133 } else {
134 cosi = 0.05 * kk;
135 }
136 }
137 COMPLEX16 complexnumerator = crect( 0.5 * Fplus1 * ( 1.0 + cosi * cosi ), Fcross1 * cosi );
138 COMPLEX16 complexdenominator = crect( 0.5 * Fplus0 * ( 1.0 + cosi * cosi ), Fcross0 * cosi );
139 if ( cabs( complexdenominator ) > 1.0e-6 ) {
140 COMPLEX16 complexval = complexnumerator / complexdenominator;
141 detPhaseMag += fmin( cabs( complexval ), 10.0 );
142 detPhaseArg += gsl_sf_angle_restrict_pos( carg( complexval ) );
143 } else {
144 loopbroken = 1;
145 detPhaseMag = 0.0;
146 detPhaseArg = 0.0;
147 }
148 }
149 }
150 detPhaseMag /= 336.0;
151 detPhaseArg /= 336.0;
152 COMPLEX16 RatioTerm = cpolar( detPhaseMag, detPhaseArg );
153
154 //Bin of interest
155 REAL8 signalFrequencyBin = round( multissb->data[0]->Tdot->data[ii] * frequency0 * uvar.Tsft ) - frequency * uvar.Tsft; //estimated nearest freq in ref IFO
156
157 REAL8 timediff0 = multissb->data[0]->DeltaT->data[ii] - 0.5 * uvar.Tsft * multissb->data[0]->Tdot->data[ii];
158 REAL8 timediff1 = multissb->data[1]->DeltaT->data[ii] - 0.5 * uvar.Tsft * multissb->data[1]->Tdot->data[ii];
159 REAL8 tau = timediff1 - timediff0;
160 REAL8 freqshift0 = -LAL_TWOPI * tau * frequency0; //real freq shift
161 REAL8 freqshift = -LAL_TWOPI * tau * ( round( multissb->data[0]->Tdot->data[ii] * frequency0 * uvar.Tsft ) / uvar.Tsft ); //estimated freq shift
162 COMPLEX16 phaseshift0 = cpolar( 1.0, freqshift0 );
163 COMPLEX16 phaseshift = cpolar( 1.0, freqshift );
164
165 REAL8 delta0_0 = ( multissb->data[0]->Tdot->data[ii] * frequency0 - frequency ) * uvar.Tsft - signalFrequencyBin;
166 REAL8 delta1_0 = ( multissb->data[1]->Tdot->data[ii] * frequency0 - frequency ) * uvar.Tsft - signalFrequencyBin;
167 REAL8 realSigBinDiff = round( delta0_0 ) - round( delta1_0 );
168 delta1_0 += realSigBinDiff;
169 REAL8 delta0 = round( multissb->data[0]->Tdot->data[ii] * frequency0 * uvar.Tsft ) * ( multissb->data[0]->Tdot->data[ii] - 1.0 );
170 REAL8 delta1 = round( multissb->data[0]->Tdot->data[ii] * frequency0 * uvar.Tsft ) * ( multissb->data[1]->Tdot->data[ii] - 1.0 );
171 REAL8 estSigBinDiff = round( delta0 ) - round( delta1 );
172 delta1 += estSigBinDiff;
173 COMPLEX16 dirichlet0;
174 if ( !uvar.rectWindow ) {
175 if ( fabsf( ( REAL4 )delta1_0 ) < ( REAL4 )1.0e-6 ) {
176 if ( fabsf( ( REAL4 )delta0_0 ) < ( REAL4 )1.0e-6 ) {
177 dirichlet0 = crect( 1.0, 0.0 );
178 } else if ( fabsf( ( REAL4 )( delta0_0 * delta0_0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
179 dirichlet0 = crect( -2.0, 0.0 );
180 } else if ( fabsf( ( REAL4 )( delta0_0 - roundf( delta0_0 ) ) ) < ( REAL4 )1.0e-6 ) {
181 dirichlet0 = crect( 0.0, 0.0 );
182 continue;
183 } else {
184 dirichlet0 = -LAL_PI * crectf( cos( LAL_PI * delta0_0 ), -sin( LAL_PI * delta0_0 ) ) * delta0_0 * ( delta0_0 * delta0_0 - 1.0 ) / sin( LAL_PI * delta0_0 );
185 }
186 } else if ( fabsf( ( REAL4 )( delta1_0 * delta1_0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
187 if ( fabsf( ( REAL4 )delta0_0 ) < ( REAL4 )1.0e-6 ) {
188 dirichlet0 = crect( -0.5, 0.0 );
189 } else if ( fabsf( ( REAL4 )( delta0_0 * delta0_0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
190 dirichlet0 = crect( 1.0, 0.0 );
191 } else if ( fabsf( ( REAL4 )( delta0_0 - roundf( delta0_0 ) ) ) < ( REAL4 )1.0e-6 ) {
192 dirichlet0 = crect( 0.0, 0.0 );
193 continue;
194 } else {
195 dirichlet0 = -LAL_PI_2 * crectf( -cos( LAL_PI * delta0_0 ), sin( LAL_PI * delta0_0 ) ) * delta0_0 * ( delta0_0 * delta0_0 - 1.0 ) / sin( LAL_PI * delta0_0 );
196 }
197 } else if ( fabsf( ( REAL4 )delta0_0 ) < ( REAL4 )1.0e-6 ) {
198 dirichlet0 = -LAL_1_PI * crectf( cos( LAL_PI * delta1_0 ), sin( LAL_PI * delta1_0 ) ) * sin( LAL_PI * delta1_0 ) / ( delta1_0 * ( delta1_0 * delta1_0 - 1.0 ) );
199 } else if ( fabsf( ( REAL4 )( delta0_0 * delta0_0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
200 dirichlet0 = LAL_2_PI * crectf( cos( LAL_PI * delta1_0 ), sin( LAL_PI * delta1_0 ) ) * sin( LAL_PI * delta1_0 ) / ( delta1_0 * ( delta1_0 * delta1_0 - 1.0 ) );
201 } else if ( fabsf( ( REAL4 )( delta0_0 - roundf( delta0_0 ) ) ) < ( REAL4 )1.0e-6 ) {
202 dirichlet0 = crect( 0.0, 0.0 );
203 continue;
204 } else {
205 dirichlet0 = sin( LAL_PI * delta1_0 ) / sin( LAL_PI * delta0_0 ) * ( delta0_0 * ( delta0_0 * delta0_0 - 1.0 ) ) / ( delta1_0 * ( delta1_0 * delta1_0 - 1.0 ) ) * cpolarf( 1.0, LAL_PI * ( delta1_0 - delta0_0 ) );
206 }
207 } else {
208 if ( fabsf( ( REAL4 )delta1_0 ) < ( REAL4 )1.0e-6 ) {
209 if ( fabsf( ( REAL4 )delta0_0 ) < ( REAL4 )1.0e-6 ) {
210 dirichlet0 = crect( 1.0, 0.0 );
211 } else if ( fabsf( ( REAL4 )( delta0_0 - roundf( delta0_0 ) ) ) < ( REAL4 )1.0e-6 ) {
212 dirichlet0 = crect( 0.0, 0.0 );
213 continue;
214 } else {
215 dirichlet0 = conj( 1.0 / ( ( cpolar( 1.0, LAL_TWOPI * delta0_0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta0_0 ) ) );
216 }
217 } else if ( fabsf( ( REAL4 )delta0_0 ) < ( REAL4 )1.0e-6 ) {
218 dirichlet0 = conj( ( cpolar( 1.0, LAL_TWOPI * delta1_0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta1_0 ) );
219 } else if ( fabsf( ( REAL4 )( delta0_0 - roundf( delta0_0 ) ) ) < ( REAL4 )1.0e-6 ) {
220 dirichlet0 = crect( 0.0, 0.0 );
221 continue;
222 } else {
223 dirichlet0 = conj( ( ( cpolar( 1.0, LAL_TWOPI * delta1_0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta1_0 ) ) / ( ( cpolar( 1.0, LAL_TWOPI * delta0_0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta0_0 ) ) );
224 }
225 }
226
227 COMPLEX16 dirichlet;
228 if ( !uvar.rectWindow ) {
229 if ( fabsf( ( REAL4 )delta1 ) < ( REAL4 )1.0e-6 ) {
230 if ( fabsf( ( REAL4 )delta0 ) < ( REAL4 )1.0e-6 ) {
231 dirichlet = crect( 1.0, 0.0 );
232 } else if ( fabsf( ( REAL4 )( delta0 * delta0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
233 dirichlet = crect( -2.0, 0.0 );
234 } else if ( fabsf( ( REAL4 )( delta0 - roundf( delta0 ) ) ) < ( REAL4 )1.0e-6 ) {
235 dirichlet = crect( 0.0, 0.0 );
236 } else {
237 dirichlet = -LAL_PI * crectf( cos( LAL_PI * delta0 ), -sin( LAL_PI * delta0 ) ) * delta0 * ( delta0 * delta0 - 1.0 ) / sin( LAL_PI * delta0 );
238 }
239 } else if ( fabsf( ( REAL4 )( delta1 * delta1 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
240 if ( fabsf( ( REAL4 )delta0 ) < ( REAL4 )1.0e-6 ) {
241 dirichlet = crect( -0.5, 0.0 );
242 } else if ( fabsf( ( REAL4 )( delta0 * delta0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
243 dirichlet = crect( 1.0, 0.0 );
244 } else if ( fabsf( ( REAL4 )( delta0 - roundf( delta0 ) ) ) < ( REAL4 )1.0e-6 ) {
245 dirichlet = crect( 0.0, 0.0 );
246 } else {
247 dirichlet = -LAL_PI_2 * crectf( -cos( LAL_PI * delta0 ), sin( LAL_PI * delta0 ) ) * delta0 * ( delta0 * delta0 - 1.0 ) / sin( LAL_PI * delta0 );
248 }
249 } else if ( fabsf( ( REAL4 )delta0 ) < ( REAL4 )1.0e-6 ) {
250 dirichlet = -LAL_1_PI * crectf( cos( LAL_PI * delta1 ), sin( LAL_PI * delta1 ) ) * sin( LAL_PI * delta1 ) / ( delta1 * ( delta1 * delta1 - 1.0 ) );
251 } else if ( fabsf( ( REAL4 )( delta0 * delta0 - 1.0 ) ) < ( REAL4 )1.0e-6 ) {
252 dirichlet = LAL_2_PI * crectf( cos( LAL_PI * delta1 ), sin( LAL_PI * delta1 ) ) * sin( LAL_PI * delta1 ) / ( delta1 * ( delta1 * delta1 - 1.0 ) );
253 } else if ( fabsf( ( REAL4 )( delta0 - roundf( delta0 ) ) ) < ( REAL4 )1.0e-6 ) {
254 dirichlet = crect( 0.0, 0.0 );
255 } else {
256 dirichlet = sin( LAL_PI * delta1 ) / sin( LAL_PI * delta0 ) * ( delta0 * ( delta0 * delta0 - 1.0 ) ) / ( delta1 * ( delta1 * delta1 - 1.0 ) ) * cpolarf( 1.0, LAL_PI * ( delta1 - delta0 ) );
257 }
258 } else {
259 if ( fabsf( ( REAL4 )delta1 ) < ( REAL4 )1.0e-6 ) {
260 if ( fabsf( ( REAL4 )delta0 ) < ( REAL4 )1.0e-6 ) {
261 dirichlet = crect( 1.0, 0.0 );
262 } else if ( fabsf( ( REAL4 )( delta0 - roundf( delta0 ) ) ) < ( REAL4 )1.0e-6 ) {
263 dirichlet = crect( 0.0, 0.0 );
264 } else {
265 dirichlet = conj( 1.0 / ( ( cpolar( 1.0, LAL_TWOPI * delta0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta0 ) ) );
266 }
267 } else if ( fabsf( ( REAL4 )delta0 ) < ( REAL4 )1.0e-6 ) {
268 dirichlet = conj( ( cpolar( 1.0, LAL_TWOPI * delta1 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta1 ) );
269 } else if ( fabsf( ( REAL4 )( delta0 - roundf( delta0 ) ) ) < ( REAL4 )1.0e-6 ) {
270 dirichlet = crect( 0.0, 0.0 );
271 } else {
272 dirichlet = conj( ( ( cpolar( 1.0, LAL_TWOPI * delta1 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta1 ) ) / ( ( cpolar( 1.0, LAL_TWOPI * delta0 ) - 1.0 ) / crect( 0.0, LAL_TWOPI * delta0 ) ) );
273 }
274 }
275 dirichlet = cpolar( 1.0, carg( dirichlet ) );
276 if ( fabs( floor( delta0 ) - floor( delta1 ) ) >= 1.0 ) {
277 dirichlet *= -1.0;
278 }
279
280 COMPLEX16 realRatio = RatioTerm0 * phaseshift0 * conj( dirichlet0 );
281 COMPLEX16 estRatio = RatioTerm * phaseshift * conj( dirichlet );
282
283 fprintf( OUTPUT, "%g %g %g %g\n", cabs( realRatio ), gsl_sf_angle_restrict_pos( carg( realRatio ) ), cabs( estRatio ), gsl_sf_angle_restrict_pos( carg( estRatio ) ) );
284 }
285
286 XLALDestroyMultiAMCoeffs( multiAMcoefficients );
287 XLALDestroyMultiSSBtimes( multissb );
288 }
289
290 fclose( OUTPUT );
291 gsl_rng_free( rng );
292 XLALDestroyMultiDetectorStateSeries( multiStateSeries );
297}
298
299INT4 InitUserVars( UserVariables_t *uvar, int argc, char *argv[] )
300{
301 XLAL_CHECK( uvar != NULL, XLAL_EINVAL, "Invalid NULL input 'uvar'\n" );
302 XLAL_CHECK( argv != NULL, XLAL_EINVAL, "Invalid NULL input 'argv'\n" );
303
304 uvar->ephemEarth = XLALStringDuplicate( "earth00-40-DE405.dat.gz" );
305 uvar->ephemSun = XLALStringDuplicate( "sun00-40-DE405.dat.gz" );
306 uvar->outfilename = XLALStringDuplicate( "output.dat" );
307 uvar->Tsft = 1800;
308 uvar->SFToverlap = 900;
309 uvar->skylocations = 1;
310 uvar->unrestrictedCosi = 0;
311 uvar->rectWindow = 0;
312
313 XLALRegisterUvarMember( Tsft, REAL8, 0, OPTIONAL, "SFT coherence time" );
314 XLALRegisterUvarMember( SFToverlap, REAL8, 0, OPTIONAL, "SFT overlap in seconds, usually Tsft/2" );
315 XLALRegisterUvarMember( t0, REAL8, 0, OPTIONAL, "GPS start time of the search" );
316 XLALRegisterUvarMember( Tobs, REAL8, 0, OPTIONAL, "Duration of the search (in seconds)" );
317 XLALRegisterUvarMember( cosi, REAL8, 0, OPTIONAL, "Cosine of NS inclinaiont angle" );
318 XLALRegisterUvarMember( psi, REAL8, 0, OPTIONAL, "Polarization angle of GW" );
319 XLALRegisterUvarMember( alpha, REAL8, 0, OPTIONAL, "Right ascension of source (in radians)" );
320 XLALRegisterUvarMember( delta, REAL8, 0, OPTIONAL, "Declination of source (in radians)" );
321 XLALRegisterUvarMember( skylocations, INT4, 0, OPTIONAL, "Number of sky locations" );
322 XLALRegisterUvarMember( IFO, STRINGVector, 0, REQUIRED, "CSV list of detectors, eg. \"H1,H2,L1,G1, ...\" " );
323 XLALRegisterUvarMember( outfilename, STRING, 0, OPTIONAL, "Output filename" );
324 XLALRegisterUvarMember( ephemEarth, STRING, 0, OPTIONAL, "Earth ephemeris file" );
325 XLALRegisterUvarMember( ephemSun, STRING, 0, OPTIONAL, "Sun ephemeris file" );
326 XLALRegisterUvarMember( unrestrictedCosi, BOOLEAN, 0, OPTIONAL, "Marginalize over cos(iota) from -1 to 1" );
327 XLALRegisterUvarMember( rectWindow, BOOLEAN, 0, OPTIONAL, "Use rectangular window function instead of Hann windowing" );
328
329 BOOLEAN should_exit = 0;
331 if ( should_exit ) {
332 exit( 1 );
333 }
334
335 return XLAL_SUCCESS;
336}
#define IFO
const LALVCSInfoList lalPulsarVCSInfoList
NULL-terminated list of VCS and build information for LALPulsar and its dependencies
static double double delta
#define STRING(a)
#define fprintf
int main(int argc, char *argv[])
INT4 InitUserVars(UserVariables_t *uvar, int argc, char *argv[])
LALDetector * XLALCreateDetector(LALDetector *detector, const LALFrDetector *frDetector, LALDetectorType type)
const LALDetector lalCachedDetectors[LAL_NUM_DETECTORS]
void XLALDestroyMultiDetectorStateSeries(MultiDetectorStateSeries *mdetStates)
Helper function to get rid of a multi-IFO DetectorStateSeries Note, this is "NULL-robust" in the sens...
MultiDetectorStateSeries * XLALGetMultiDetectorStates(const MultiLIGOTimeGPSVector *multiTS, const MultiLALDetector *multiIFO, const EphemerisData *edat, REAL8 tOffset)
Get the detector-time series for the given MultiLIGOTimeGPSVector.
EphemerisData * XLALInitBarycenter(const CHAR *earthEphemerisFile, const CHAR *sunEphemerisFile)
XLAL interface to reading ephemeris files 'earth' and 'sun', and return ephemeris-data in old backwar...
void XLALDestroyEphemerisData(EphemerisData *edat)
Destructor for EphemerisData struct, NULL robust.
void XLALDestroyMultiAMCoeffs(MultiAMCoeffs *multiAMcoef)
Destroy a MultiAMCoeffs structure.
Definition: LALComputeAM.c:469
MultiAMCoeffs * XLALComputeMultiAMCoeffs(const MultiDetectorStateSeries *multiDetStates, const MultiNoiseWeights *multiWeights, SkyPosition skypos)
Multi-IFO version of XLALComputeAMCoeffs().
Definition: LALComputeAM.c:379
#define LAL_PI_2
#define LAL_2_PI
#define LAL_PI
#define LAL_1_PI
#define LAL_TWOPI
#define crect(re, im)
unsigned char BOOLEAN
#define cpolar(r, th)
double complex COMPLEX16
#define cpolarf(r, th)
double REAL8
#define XLAL_INIT_DECL(var,...)
#define crectf(re, im)
char CHAR
uint32_t UINT4
int32_t INT4
float REAL4
LALNameLength
LALDETECTORTYPE_IFODIFF
LAL_LLO_4K_DETECTOR
LAL_VIRGO_DETECTOR
LAL_LHO_2K_DETECTOR
LAL_LHO_4K_DETECTOR
void * XLALMalloc(size_t n)
void XLALFree(void *p)
char char * XLALStringDuplicate(const char *s)
MultiLIGOTimeGPSVector * XLALMakeMultiTimestamps(LIGOTimeGPS tStart, REAL8 Tspan, REAL8 Tsft, REAL8 Toverlap, UINT4 numDet)
Same as XLALMakeTimestamps() just for several detectors, additionally specify the number of detectors...
void XLALDestroyMultiTimestamps(MultiLIGOTimeGPSVector *multiTS)
Destroy a MultiLIGOTimeGPSVector timestamps vector.
MultiSSBtimes * XLALGetMultiSSBtimes(const MultiDetectorStateSeries *multiDetStates, SkyPosition skypos, LIGOTimeGPS refTime, SSBprecision precision)
Multi-IFO version of XLALGetSSBtimes().
Definition: SSBtimes.c:657
void XLALDestroyMultiSSBtimes(MultiSSBtimes *multiSSB)
Destroy a MultiSSBtimes structure.
Definition: SSBtimes.c:845
@ SSBPREC_RELATIVISTICOPT
optimized relativistic, numerically equivalent to SSBPREC_RELATIVISTIC, but faster
Definition: SSBtimes.h:48
COORDINATESYSTEM_EQUATORIAL
int XLALUserVarReadAllInput(BOOLEAN *should_exit, int argc, char *argv[], const LALVCSInfoList vcs_list)
void XLALDestroyUserVars(void)
#define XLALRegisterUvarMember(name, type, option, category,...)
int XLALUserVarWasSet(const void *cvar)
#define xlalErrno
#define XLAL_ERROR(...)
#define XLAL_CHECK(assertion,...)
XLAL_ENOMEM
XLAL_SUCCESS
XLAL_EFUNC
XLAL_EIO
XLAL_EINVAL
LIGOTimeGPS * XLALGPSSetREAL8(LIGOTimeGPS *epoch, REAL8 t)
n
LALDetector detectors[LAL_NUM_DETECTORS]
double alpha
double t0
REAL4Vector * b
(weighted) per-SFT antenna-pattern function
Definition: LALComputeAM.h:65
REAL4Vector * a
(weighted) per-SFT antenna-pattern function
Definition: LALComputeAM.h:64
This structure contains all information about the center-of-mass positions of the Earth and Sun,...
LALFrDetector frDetector
REAL4 xArmAzimuthRadians
REAL4 yArmAzimuthRadians
CHAR name[LALNameLength]
Multi-IFO container for antenna-pattern coefficients and atenna-pattern matrix .
Definition: LALComputeAM.h:137
AMCoeffs ** data
noise-weighted AM-coeffs , and
Definition: LALComputeAM.h:142
Multi-IFO time-series of DetectorStates.
array of detectors definitions 'LALDetector'
A collection of (multi-IFO) LIGOTimeGPSVector time-stamps vectors.
Definition: SFTfileIO.h:198
Multi-IFO container for SSB timings.
Definition: SSBtimes.h:67
SSBtimes ** data
array of SSBtimes (pointers)
Definition: SSBtimes.h:72
REAL4 * data
REAL8 * data
REAL8Vector * Tdot
dT/dt : time-derivative of SSB-time wrt local time for SFT-alpha
Definition: SSBtimes.h:63
REAL8Vector * DeltaT
Time-difference of SFT-alpha - tau0 in SSB-frame.
Definition: SSBtimes.h:62
REAL8 longitude
REAL8 latitude
CoordinateSystem system
user input variables
Definition: compareFstats.c:51
CHAR * ephemSun
Sun ephemeris file to use.
CHAR * ephemEarth
Earth ephemeris file to use.
REAL8 SFToverlap
overlap SFTs by this many seconds
REAL8 Tsft
SFT time baseline Tsft.
LALStringVector * IFO