25#include <gsl/gsl_blas.h>
27#include <lal/SuperskyMetrics.h>
28#include <lal/LALStdlib.h>
29#include <lal/LALInitBarycenter.h>
30#include <lal/MetricUtils.h>
32#include <lal/GSLHelpers.h>
37#define REF_TIME { 900100100, 0 }
38#define FIDUCIAL_FREQ 100.0
42#define MAX_SKY_OFFSETS PULSAR_MAX_SPINS
57 { .
refTime =
REF_TIME, .Alpha = 0.00000000000000, .Delta = 0.000000000000000, .fkdot = {100.0000000000000, 0.00000000000000e-00} },
58 { .refTime =
REF_TIME, .Alpha = 4.88014010120016, .Delta = -0.954446475246007, .fkdot = { 99.9999983978492, 1.48957780038094e-09} },
59 { .refTime =
REF_TIME, .Alpha = 0.52587274931672, .Delta = 0.685297319257976, .fkdot = { 99.9999923150006, -1.41365319702693e-09} },
60 { .refTime =
REF_TIME, .Alpha = 3.53542175437611, .Delta = -1.502778038590950, .fkdot = {100.0000064863180, -1.28748375084384e-09} },
61 { .refTime =
REF_TIME, .Alpha = 1.36054903961191, .Delta = 0.241343663657163, .fkdot = { 99.9999901679571, 3.37107171004537e-10} },
62 { .refTime =
REF_TIME, .Alpha = 2.85470536965808, .Delta = 1.1575340928032900, .fkdot = {100.0000074463050, 2.46412240438217e-09} },
63 { .refTime =
REF_TIME, .Alpha = 1.82755817952460, .Delta = 0.667995269285982, .fkdot = { 99.9999897239871, 1.79900370692270e-10} },
64 { .refTime =
REF_TIME, .Alpha = 1.70734223243163, .Delta = -1.213787405673430, .fkdot = {100.0000026535270, -1.07122135891104e-09} },
65 { .refTime =
REF_TIME, .Alpha = 2.30597131157246, .Delta = 0.348657791621429, .fkdot = {100.0000133749770, -5.43309003215614e-10} },
66 { .refTime =
REF_TIME, .Alpha = 3.31129323970275, .Delta = -1.225892709583030, .fkdot = {100.0000062524320, 8.07713885739405e-10} }
70 { 5.069230179517075e+02, 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00},
71 { 0.000000000000000e+00, 6.547662357441416e+01, 0.000000000000000e+00, 0.000000000000000e+00},
72 { 0.000000000000000e+00, 0.000000000000000e+00, 7.064620283536166e+22, -8.077428080473405e+02},
73 { 0.000000000000000e+00, 0.000000000000000e+00, -8.077428080473405e+02, 2.210286062098707e+11},
78 { 6.568635205161581e+01, 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00},
79 { 0.000000000000000e+00, 6.237311686068988e+01, 0.000000000000000e+00, 0.000000000000000e+00},
80 { 0.000000000000000e+00, 0.000000000000000e+00, 1.058455565334740e+23, -1.527749726118672e+17},
81 { 0.000000000000000e+00, 0.000000000000000e+00, -1.527749726118672e+17, 2.210286062098809e+11},
83 { 6.568666765075774e+01, 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00},
84 { 0.000000000000000e+00, 6.236546718525190e+01, 0.000000000000000e+00, 0.000000000000000e+00},
85 { 0.000000000000000e+00, 0.000000000000000e+00, 2.474954556329869e+20, -1.208418861582654e+02},
86 { 0.000000000000000e+00, 0.000000000000000e+00, -1.208418861582654e+02, 2.210286062098680e+11},
88 { 6.568680817974267e+01, 0.000000000000000e+00, 0.000000000000000e+00, 0.000000000000000e+00},
89 { 0.000000000000000e+00, 6.235717845927233e+01, 0.000000000000000e+00, 0.000000000000000e+00},
90 { 0.000000000000000e+00, 0.000000000000000e+00, 1.058455565227733e+23, 1.527749726120594e+17},
91 { 0.000000000000000e+00, 0.000000000000000e+00, 1.527749726120594e+17, 2.210286062097871e+11},
101 { 9.363383835591375e-01, 3.276518345636219e-01, 1.261535048302427e-01},
102 {-3.309081001489509e-01, 9.436487575620593e-01, 5.181853663774516e-03},
103 {-1.173467542357819e-01, -4.659718509388588e-02, 9.919971983890140e-01},
107 {-2.134760488419512e-11, 1.830302603631687e-09, 7.303907936812212e-10},
108 { 9.866360627623324e-03, 4.620189300680953e-05, 1.748754269165573e-04},
119 { 7.014188729589919e-01, -7.127492992955380e-01, 3.179083552258733e-05},
120 { 7.127491812384910e-01, 7.014187811320953e-01, 5.460027715689418e-04},
121 {-4.114617819526091e-04, -3.603177566767635e-04, 9.999998504351469e-01},
125 {-1.518511649108543e-09, 9.119239045225907e-10, 7.652734977238319e-10},
126 { 4.238436305149103e-03, 8.987129848963495e-03, 1.432613602463570e-03},
134 { 7.928332170641094e-01, -6.094386653165459e-01, 5.600859050599027e-05},
135 { 6.094384042541772e-01, 7.928329562280915e-01, 8.572856860540078e-04},
136 {-5.668685006887667e-04, -6.455507823947125e-04, 9.999996309620771e-01},
140 {-1.538506600219701e-09, 9.036045301253465e-10, 7.329842343015285e-10},
141 { 5.337970249571395e-03, 8.252674722491856e-03, 1.420014590206618e-03},
149 { 8.689897105711624e-01, -4.948301493223842e-01, 7.901287635748977e-05},
150 { 4.948297433622081e-01, 8.689891816122729e-01, 1.152095883124570e-03},
151 {-6.387531126429617e-04, -9.620615466963718e-04, 9.999993332157983e-01},
155 {-1.560552363459718e-09, 9.006124427330534e-10, 6.868353365772801e-10},
156 { 6.434079710711734e-03, 7.519415923170119e-03, 1.434855908485180e-03},
162 {0.000000e+00, 2.842696e+07, 1.714926e+05, 2.821377e+07, 3.877800e+06, 2.839587e+07, 1.360228e+07, 2.139293e+07, 3.333989e+07, 4.206510e+07},
163 {2.842696e+07, 0.000000e+00, 2.419283e+07, 2.504598e+05, 1.251935e+07, 1.002015e+06, 3.282244e+06, 5.689312e+05, 3.277731e+05, 1.336220e+06},
164 {1.714926e+05, 2.419283e+07, 0.000000e+00, 2.398770e+07, 2.552860e+06, 2.429737e+07, 1.080643e+07, 1.773417e+07, 2.877145e+07, 3.688169e+07},
165 {2.821377e+07, 2.504598e+05, 2.398770e+07, 0.000000e+00, 1.312504e+07, 2.237378e+06, 4.014595e+06, 5.120383e+05, 9.922287e+05, 1.759271e+06},
166 {3.877800e+06, 1.251935e+07, 2.552860e+06, 1.312504e+07, 0.000000e+00, 1.153347e+07, 3.162766e+06, 8.455535e+06, 1.533802e+07, 2.178058e+07},
167 {2.839587e+07, 1.002015e+06, 2.429737e+07, 2.237378e+06, 1.153347e+07, 0.000000e+00, 2.696814e+06, 1.922870e+06, 6.601354e+05, 2.417185e+06},
168 {1.360228e+07, 3.282244e+06, 1.080643e+07, 4.014595e+06, 3.162766e+06, 2.696814e+06, 0.000000e+00, 1.740316e+06, 4.592892e+06, 8.444680e+06},
169 {2.139293e+07, 5.689312e+05, 1.773417e+07, 5.120383e+05, 8.455535e+06, 1.922870e+06, 1.740316e+06, 0.000000e+00, 1.695582e+06, 3.585868e+06},
170 {3.333989e+07, 3.277731e+05, 2.877145e+07, 9.922287e+05, 1.533802e+07, 6.601354e+05, 4.592892e+06, 1.695582e+06, 0.000000e+00, 6.169347e+05},
171 {4.206510e+07, 1.336220e+06, 3.688169e+07, 1.759271e+06, 2.178058e+07, 2.417185e+06, 8.444680e+06, 3.585868e+06, 6.169347e+05, 0.000000e+00},
176 {0.000000e+00, 3.024189e+07, 1.543644e+05, 2.361445e+07, 7.297336e+06, 4.341212e+07, 2.026653e+07, 1.998370e+07, 4.088115e+07, 4.591125e+07},
177 {3.024189e+07, 0.000000e+00, 2.607532e+07, 4.103494e+05, 7.831289e+06, 1.189821e+06, 9.965225e+05, 1.059082e+06, 8.007433e+05, 1.629541e+06},
178 {1.543644e+05, 2.607532e+07, 0.000000e+00, 1.995049e+07, 5.329463e+06, 3.838988e+07, 1.688386e+07, 1.662551e+07, 3.601172e+07, 4.074165e+07},
179 {2.361445e+07, 4.103494e+05, 1.995049e+07, 0.000000e+00, 4.662728e+06, 2.997510e+06, 1.322075e+05, 1.516108e+05, 2.356988e+06, 3.673699e+06},
180 {7.297336e+06, 7.831289e+06, 5.329463e+06, 4.662728e+06, 0.000000e+00, 1.511270e+07, 3.242171e+06, 3.132814e+06, 1.363639e+07, 1.660431e+07},
181 {4.341212e+07, 1.189821e+06, 3.838988e+07, 2.997510e+06, 1.511270e+07, 0.000000e+00, 4.355434e+06, 4.492228e+06, 3.929618e+04, 3.792766e+04},
182 {2.026653e+07, 9.965225e+05, 1.688386e+07, 1.322075e+05, 3.242171e+06, 4.355434e+06, 0.000000e+00, 3.558143e+03, 3.580300e+06, 5.172549e+06},
183 {1.998370e+07, 1.059082e+06, 1.662551e+07, 1.516108e+05, 3.132814e+06, 4.492228e+06, 3.558143e+03, 0.000000e+00, 3.701188e+06, 5.315678e+06},
184 {4.088115e+07, 8.007433e+05, 3.601172e+07, 2.356988e+06, 1.363639e+07, 3.929618e+04, 3.580300e+06, 3.701188e+06, 0.000000e+00, 1.462400e+05},
185 {4.591125e+07, 1.629541e+06, 4.074165e+07, 3.673699e+06, 1.660431e+07, 3.792766e+04, 5.172549e+06, 5.315678e+06, 1.462400e+05, 0.000000e+00},
187 {0.000000e+00, 2.876453e+07, 1.736967e+05, 2.845102e+07, 3.414811e+06, 2.743194e+07, 1.322922e+07, 2.165178e+07, 3.346665e+07, 4.253244e+07},
188 {2.876453e+07, 0.000000e+00, 2.447187e+07, 1.779353e+03, 1.236248e+07, 1.946008e+04, 2.981507e+06, 5.045797e+05, 1.784609e+05, 1.342001e+06},
189 {1.736967e+05, 2.447187e+07, 0.000000e+00, 2.418261e+07, 2.049651e+06, 2.324404e+07, 1.037386e+07, 1.795027e+07, 2.882252e+07, 3.727495e+07},
190 {2.845102e+07, 1.779353e+03, 2.418261e+07, 0.000000e+00, 1.215974e+07, 1.722632e+04, 2.883870e+06, 4.636396e+05, 2.061967e+05, 1.412060e+06},
191 {3.414811e+06, 1.236248e+07, 2.049651e+06, 1.215974e+07, 0.000000e+00, 1.149079e+07, 3.202415e+06, 7.874615e+06, 1.550436e+07, 2.184963e+07},
192 {2.743194e+07, 1.946008e+04, 2.324404e+07, 1.722632e+04, 1.149079e+07, 0.000000e+00, 2.561189e+06, 3.466390e+05, 3.014610e+05, 1.652932e+06},
193 {1.322922e+07, 2.981507e+06, 1.037386e+07, 2.883870e+06, 3.202415e+06, 2.561189e+06, 0.000000e+00, 1.035231e+06, 4.614106e+06, 8.322602e+06},
194 {2.165178e+07, 5.045797e+05, 1.795027e+07, 4.636396e+05, 7.874615e+06, 3.466390e+05, 1.035231e+06, 0.000000e+00, 1.282484e+06, 3.491914e+06},
195 {3.346665e+07, 1.784609e+05, 2.882252e+07, 2.061967e+05, 1.550436e+07, 3.014610e+05, 4.614106e+06, 1.282484e+06, 0.000000e+00, 5.430866e+05},
196 {4.253244e+07, 1.342001e+06, 3.727495e+07, 1.412060e+06, 2.184963e+07, 1.652932e+06, 8.322602e+06, 3.491914e+06, 5.430866e+05, 0.000000e+00},
198 {0.000000e+00, 2.627447e+07, 1.867214e+05, 3.257584e+07, 9.213898e+05, 1.434368e+07, 7.311169e+06, 2.254331e+07, 2.567187e+07, 3.775160e+07},
199 {2.627447e+07, 0.000000e+00, 2.203155e+07, 3.392506e+05, 1.736439e+07, 1.796882e+06, 5.868773e+06, 1.431319e+05, 4.129221e+03, 1.037116e+06},
200 {1.867214e+05, 2.203155e+07, 0.000000e+00, 2.783016e+07, 2.794658e+05, 1.125817e+07, 5.161574e+06, 1.862685e+07, 2.148010e+07, 3.262861e+07},
201 {3.257584e+07, 3.392506e+05, 2.783016e+07, 0.000000e+00, 2.255274e+07, 3.697478e+06, 9.027756e+06, 9.208643e+05, 4.135111e+05, 1.920540e+05},
202 {9.213898e+05, 1.736439e+07, 2.794658e+05, 2.255274e+07, 0.000000e+00, 7.996908e+06, 3.043713e+06, 1.435923e+07, 1.687331e+07, 2.688786e+07},
203 {1.434368e+07, 1.796882e+06, 1.125817e+07, 3.697478e+06, 7.996908e+06, 0.000000e+00, 1.173820e+06, 9.298057e+05, 1.639649e+06, 5.560758e+06},
204 {7.311169e+06, 5.868773e+06, 5.161574e+06, 9.027756e+06, 3.043713e+06, 1.173820e+06, 0.000000e+00, 4.182195e+06, 5.584268e+06, 1.183893e+07},
205 {2.254331e+07, 1.431319e+05, 1.862685e+07, 9.208643e+05, 1.435923e+07, 9.298057e+05, 4.182195e+06, 0.000000e+00, 1.030798e+05, 1.950011e+06},
206 {2.567187e+07, 4.129221e+03, 2.148010e+07, 4.135111e+05, 1.687331e+07, 1.639649e+06, 5.584268e+06, 1.030798e+05, 0.000000e+00, 1.161485e+06},
207 {3.775160e+07, 1.037116e+06, 3.262861e+07, 1.920540e+05, 2.688786e+07, 5.560758e+06, 1.183893e+07, 1.950011e+06, 1.161485e+06, 0.000000e+00},
211#define CHECK_RELERR(A, B, TOL) do { \
212 const double lhs = fabs( (A) - (B) ); \
213 const double tol = (TOL); \
214 const double rhs = GSL_MAX( 1.0, fabs( (A) + (B) ) ); \
215 XLALPrintInfo( #A"=%0.5e "#B"=%0.5e |"#A" - "#B"|=%0.5e tol=%0.5e |"#A" + "#B"|=%0.5e\n", A, B, lhs, tol, rhs ); \
216 XLAL_CHECK( lhs <= tol * rhs, XLAL_ETOL, "|"#A" - "#B"| = %0.5e > %0.5e = %0.5e * |"#A" + "#B"|", lhs, tol * rhs, tol ); \
231 const gsl_matrix *rssky_metric,
232 const double rssky_metric_ref[4][4],
233 const SuperskyTransformData *rssky_transf,
234 const SuperskyTransformData *rssky_transf_ref,
236 const double phys_mismatch_tol
242 gsl_matrix_const_view rssky_metric_ref_view = gsl_matrix_const_view_array( (
const double * )rssky_metric_ref, 4, 4 );
244 XLAL_CHECK( err <= err_tol, XLAL_ETOL, "'rssky_metric' check failed: err = %0.3e > %0.3
e = err_tol
", err, err_tol );
247 XLAL_CHECK( rssky_transf->ndim == 4, XLAL_ESIZE );
248 XLAL_CHECK( rssky_transf->nsky_offsets == 2, XLAL_ESIZE );
249 const double err_tol = 1e-5;
250 for ( size_t i = 0; i < 3; ++i ) {
251 for ( size_t j = 0; j < 3; ++j ) {
252 const double align_sky_ij = rssky_transf->align_sky[i][j];
253 const double align_sky_ref_ij = rssky_transf_ref->align_sky[i][j];
254 CHECK_RELERR( align_sky_ij, align_sky_ref_ij, err_tol );
257 for ( size_t i = 0; i < rssky_transf->nsky_offsets; ++i ) {
258 for ( size_t j = 0; j < 3; ++j ) {
259 const double sky_offsets_ij = rssky_transf->sky_offsets[i][j];
260 const double sky_offsets_ref_ij = rssky_transf_ref->sky_offsets[i][j];
261 CHECK_RELERR( sky_offsets_ij, sky_offsets_ref_ij, err_tol );
266 // Check round-trip conversions of each test point
268 gsl_matrix *GAMAT( rssky_points, 4, NUM_POINTS );
269 for ( size_t j = 0; j < NUM_POINTS; ++j ) {
270 gsl_vector_view rssky_point = gsl_matrix_column( rssky_points, j );
271 PulsarDopplerParams XLAL_INIT_DECL( new_phys_point );
272 XLAL_CHECK( XLALConvertPhysicalToSuperskyPoint( &rssky_point.vector, &phys_points[j], rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
273 XLAL_CHECK( XLALConvertSuperskyToPhysicalPoint( &new_phys_point, &rssky_point.vector, NULL, rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
274 XLAL_CHECK( CompareDoppler( &phys_points[j], &new_phys_point ) == EXIT_SUCCESS, XLAL_EFUNC );
276 gsl_matrix *intm_phys_points = NULL;
277 gsl_matrix *new_rssky_points = NULL;
278 XLAL_CHECK( XLALConvertSuperskyToPhysicalPoints( &intm_phys_points, rssky_points, rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
279 XLAL_CHECK( XLALConvertPhysicalToSuperskyPoints( &new_rssky_points, intm_phys_points, rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
280 const double err_tol = 1e-6;
281 for ( size_t i = 0; i < 4; ++i ) {
282 for ( size_t j = 0; j < NUM_POINTS; ++j ) {
283 const double rssky_points_ij = gsl_matrix_get( rssky_points, i, j );
284 const double new_rssky_points_ij = gsl_matrix_get( new_rssky_points, i, j );
285 CHECK_RELERR( rssky_points_ij, new_rssky_points_ij, err_tol );
288 GFMAT( rssky_points, intm_phys_points, new_rssky_points );
291 // Check mismatches between pairs of points
293 gsl_vector *GAVEC( rssky_point_i, 4 );
294 gsl_vector *GAVEC( rssky_point_j, 4 );
295 gsl_vector *GAVEC( temp, 4 );
296 for ( size_t i = 0; i < NUM_POINTS; ++i ) {
297 XLAL_CHECK( XLALConvertPhysicalToSuperskyPoint( rssky_point_i, &phys_points[i], rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
298 for ( size_t j = 0; j < NUM_POINTS; ++j ) {
299 XLAL_CHECK( XLALConvertPhysicalToSuperskyPoint( rssky_point_j, &phys_points[j], rssky_transf ) == XLAL_SUCCESS, XLAL_EFUNC );
300 gsl_vector_sub( rssky_point_j, rssky_point_i );
301 gsl_blas_dgemv( CblasNoTrans, 1.0, rssky_metric, rssky_point_j, 0.0, temp );
302 double mismatch = 0.0;
303 gsl_blas_ddot( rssky_point_j, temp, &mismatch );
304 CHECK_RELERR( mismatch, phys_mismatch[i][j], phys_mismatch_tol );
307 GFVEC( rssky_point_i, rssky_point_j, temp );
317 // Load ephemeris data
318 EphemerisData *edat = XLALInitBarycenter( TEST_PKG_DATA_DIR "earth00-40-DE405.dat.gz
",
319 TEST_PKG_DATA_DIR "sun00-40-DE405.dat.gz
" );
320 XLAL_CHECK_MAIN( edat != NULL, XLAL_EFUNC );
322 // Create segment list
324 XLAL_CHECK_MAIN( XLALSegListInit( &segments ) == XLAL_SUCCESS, XLAL_EFUNC );
325 for ( size_t n = 0; n < NUM_SEGS; ++n ) {
326 const double Tspan = 3 * 86400;
327 const double deltat[NUM_SEGS] = { -8 * 86400, 0, 8 * 86400 };
329 LIGOTimeGPS start_time = REF_TIME, end_time = REF_TIME;
330 XLALGPSAdd( &start_time, deltat[n] - 0.5 * Tspan );
331 XLALGPSAdd( &end_time, deltat[n] + 0.5 * Tspan );
332 XLAL_CHECK_MAIN( XLALSegSet( &segment, &start_time, &end_time, 0 ) == XLAL_SUCCESS, XLAL_EFUNC );
333 XLAL_CHECK_MAIN( XLALSegListAppend( &segments, &segment ) == XLAL_SUCCESS, XLAL_EFUNC );
336 // Compute supersky metrics
337 SuperskyMetrics *metrics = NULL;
339 const LIGOTimeGPS ref_time = REF_TIME;
340 const MultiLALDetector detectors = { .length = 1, .sites = { lalCachedDetectors[LAL_LLO_4K_DETECTOR] } };
341 metrics = XLALComputeSuperskyMetrics( SUPERSKY_METRIC_TYPE, 1, &ref_time, &segments, FIDUCIAL_FREQ, &detectors, NULL, DETMOTION_SPIN | DETMOTION_PTOLEORBIT, edat );
343 XLAL_CHECK_MAIN( metrics != NULL, XLAL_EFUNC );
345 // Check coherent and semicoherent metrics
346 for ( size_t n = 0; n < NUM_SEGS; ++n ) {
347 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
348 metrics->coh_rssky_metric[n], coh_rssky_metric_refs[n],
349 metrics->coh_rssky_transf[n], &coh_rssky_transf_refs[n],
350 coh_phys_mismatches[n], 1e-2
351 ) == XLAL_SUCCESS, XLAL_EFUNC );
353 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
354 metrics->semi_rssky_metric, semi_rssky_metric_ref,
355 metrics->semi_rssky_transf, &semi_rssky_transf_ref,
356 semi_phys_mismatch, 3e-2
357 ) == XLAL_SUCCESS, XLAL_EFUNC );
359 // Check semicoherent metric after round-trip frequency rescaling
360 XLAL_CHECK_MAIN( XLALScaleSuperskyMetricsFiducialFreq( metrics, 257.52 ) == XLAL_SUCCESS, XLAL_EFUNC );
362 double semi_rssky_metric_rescale[4][4];
363 memcpy( semi_rssky_metric_rescale, semi_rssky_metric_ref, sizeof( semi_rssky_metric_ref ) );
364 gsl_matrix_view semi_rssky_metric_rescale_view = gsl_matrix_view_array( ( double * )semi_rssky_metric_rescale, 4, 4 );
365 gsl_matrix_view sky_sky = gsl_matrix_submatrix( &semi_rssky_metric_rescale_view.matrix, 0, 0, 2, 2 );
366 gsl_matrix_scale( &sky_sky.matrix, ( 257.52 / FIDUCIAL_FREQ ) * ( 257.52 / FIDUCIAL_FREQ ) );
367 const double err = XLALCompareMetrics( metrics->semi_rssky_metric, &semi_rssky_metric_rescale_view.matrix ), err_tol = 1e-6;
368 XLAL_CHECK_MAIN( err <= err_tol, XLAL_ETOL, "'rssky_metric' check failed:
err = %0.3e > %0.3e = err_tol
", err, err_tol );
370 XLAL_CHECK_MAIN( XLALScaleSuperskyMetricsFiducialFreq( metrics, FIDUCIAL_FREQ ) == XLAL_SUCCESS, XLAL_EFUNC );
371 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
372 metrics->semi_rssky_metric, semi_rssky_metric_ref,
373 metrics->semi_rssky_transf, &semi_rssky_transf_ref,
374 semi_phys_mismatch, 3e-2
375 ) == XLAL_SUCCESS, XLAL_EFUNC );
377 // Check that XLALEqualizeReducedSuperskyMetricsFreqSpacing() with equal mismatches does nothing
378 XLAL_CHECK_MAIN( XLALEqualizeReducedSuperskyMetricsFreqSpacing( metrics, 0.1, 0.1 ) == XLAL_SUCCESS, XLAL_EFUNC );
379 for ( size_t n = 0; n < NUM_SEGS; ++n ) {
380 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
381 metrics->coh_rssky_metric[n], coh_rssky_metric_refs[n],
382 metrics->coh_rssky_transf[n], &coh_rssky_transf_refs[n],
383 coh_phys_mismatches[n], 1e-2
384 ) == XLAL_SUCCESS, XLAL_EFUNC );
386 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
387 metrics->semi_rssky_metric, semi_rssky_metric_ref,
388 metrics->semi_rssky_transf, &semi_rssky_transf_ref,
389 semi_phys_mismatch, 3e-2
390 ) == XLAL_SUCCESS, XLAL_EFUNC );
392 // Check coherent and semicoherent metrics
393 SuperskyMetrics *copy_metrics = XLALCopySuperskyMetrics( metrics );
394 XLAL_CHECK_MAIN( copy_metrics != NULL, XLAL_EFUNC );
395 for ( size_t n = 0; n < NUM_SEGS; ++n ) {
396 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
397 copy_metrics->coh_rssky_metric[n], coh_rssky_metric_refs[n],
398 copy_metrics->coh_rssky_transf[n], &coh_rssky_transf_refs[n],
399 coh_phys_mismatches[n], 1e-2
400 ) == XLAL_SUCCESS, XLAL_EFUNC );
402 XLAL_CHECK_MAIN( CheckSuperskyMetrics(
403 copy_metrics->semi_rssky_metric, semi_rssky_metric_ref,
404 copy_metrics->semi_rssky_transf, &semi_rssky_transf_ref,
405 semi_phys_mismatch, 3e-2
406 ) == XLAL_SUCCESS, XLAL_EFUNC );
409 // Check writing/reading SuperskyMetrics to a FITS file
410#if !defined(HAVE_LIBCFITSIO)
411 fprintf( stderr, "CFITSIO library is not available; skipping FITS writing/reading test\n
" );
414 FITSFile *file = XLALFITSFileOpenWrite( "SuperskyMetricsTest.fits
" );
415 XLAL_CHECK_MAIN( file != NULL, XLAL_EFUNC );
416 XLAL_CHECK_MAIN( XLALFITSWriteSuperskyMetrics( file, metrics ) == XLAL_SUCCESS, XLAL_EFUNC );
417 XLALFITSFileClose( file );
420 FITSFile *file = XLALFITSFileOpenRead( "SuperskyMetricsTest.fits
" );
421 XLAL_CHECK_MAIN( file != NULL, XLAL_EFUNC );
422 SuperskyMetrics *metrics_from_fits = NULL;
423 XLAL_CHECK_MAIN( XLALFITSReadSuperskyMetrics( file, &metrics_from_fits ) == XLAL_SUCCESS, XLAL_EFUNC );
424 XLALFITSFileClose( file );
425 for ( size_t n = 0; n < NUM_SEGS; ++n ) {
426 XLAL_CHECK_MAIN( XLALCompareMetrics( metrics->coh_rssky_metric[n], metrics_from_fits->coh_rssky_metric[n] ) == 0,
427 XLAL_EIO, "Failed to read coh_rssky_metric[%zu] from FITS
file", n );
428 XLAL_CHECK_MAIN( memcmp( metrics->coh_rssky_transf[n], metrics_from_fits->coh_rssky_transf[n], sizeof( *metrics->coh_rssky_transf[n] ) ) == 0,
429 XLAL_EIO, "Failed to read coh_rssky_transf[%zu] from FITS
file", n );
431 XLAL_CHECK_MAIN( XLALCompareMetrics( metrics->semi_rssky_metric, metrics_from_fits->semi_rssky_metric ) == 0,
432 XLAL_EIO, "Failed to read semi_rssky_metric from FITS
file" );
433 XLAL_CHECK_MAIN( memcmp( metrics->semi_rssky_transf, metrics_from_fits->semi_rssky_transf, sizeof( *metrics->semi_rssky_transf ) ) == 0,
434 XLAL_EIO, "Failed to read semi_rssky_transf from FITS
file" );
435 XLALDestroySuperskyMetrics( metrics_from_fits );
440 XLALDestroyEphemerisData( edat );
441 XLALSegListClear( &segments );
442 XLALDestroySuperskyMetrics( metrics );
443 XLALDestroySuperskyMetrics( copy_metrics );
444 LALCheckMemoryLeaks();
451// c-file-style: "linux
"
const double semi_phys_mismatch[NUM_POINTS][NUM_POINTS]
const double coh_phys_mismatches[NUM_SEGS][NUM_POINTS][NUM_POINTS]
const SuperskyTransformData semi_rssky_transf_ref
const double coh_rssky_metric_refs[NUM_SEGS][4][4]
#define CHECK_RELERR(A, B, TOL)
const double semi_rssky_metric_ref[4][4]
const SuperskyTransformData coh_rssky_transf_refs[NUM_SEGS]
const PulsarDopplerParams phys_points[NUM_POINTS]
static int CompareDoppler(const PulsarDopplerParams *a, const PulsarDopplerParams *b)
static int CheckSuperskyMetrics(const gsl_matrix *rssky_metric, const double rssky_metric_ref[4][4], const SuperskyTransformData *rssky_transf, const SuperskyTransformData *rssky_transf_ref, const double phys_mismatch[NUM_POINTS][NUM_POINTS], const double phys_mismatch_tol)
REAL8 XLALCompareMetrics(const gsl_matrix *g1_ij, const gsl_matrix *g2_ij)
Flexible comparison function between two metrics and .
#define XLAL_CHECK(assertion,...)
int XLALGPSCmp(const LIGOTimeGPS *t0, const LIGOTimeGPS *t1)
Type containing the 'Doppler-parameters' affecting the time-evolution of the phase.
PulsarSpins fkdot
Intrinsic spins: [Freq, f1dot, f2dot, ... ] where fkdot = d^kFreq/dt^k.
REAL8 Delta
Sky position: DEC (latitude) in equatorial coords and radians.
LIGOTimeGPS refTime
Reference time of pulsar parameters (in SSB!)
REAL8 Alpha
Sky position: RA (longitude) in equatorial coords and radians.