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
OutputResults.c
Go to the documentation of this file.
1//
2// Copyright (C) 2016, 2017 Karl Wette
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/// \file
22/// \ingroup lalpulsar_bin_Weave
23///
24
25#include "OutputResults.h"
26#include "ResultsToplist.h"
27#include "Statistics.h"
28
29#include <lal/UserInputPrint.h>
30
31///
32/// Histogram bin width of mean multi-F-statistics
33///
35
36///
37/// Histogram bin of mean multi-F-statistics
38///
39typedef struct {
40 // Lower bin boundary
42 // Upper bin boundary
44 // Bin count
47
48///
49/// Output results from a search
50///
52 /// Struct holding all parameters for which statistics to output and compute, when, and how
53 /// NOTE: this is the *owner* of WeaveStatisticsParams, which is where it will be freed at the end
54 /// while toplists will simply hold a reference-pointer
55 WeaveStatisticsParams *statistics_params;
56 /// Reference time at which search is conducted
58 /// Number of spindown parameters to output
59 size_t nspins;
60 // Maximum size of toplists
62 /// Number of output results toplists
63 size_t ntoplists;
64 /// Output result toplists
65 WeaveResultsToplist *toplists[8];
66 // Vector to store histogram of mean multi-F-statistics
68 // Number of mean multi-F-statistics below range of histogram
70 // Number of mean multi-F-statistics above range of histogram
72 // Temporary REAL4 vector for generating histogram of mean multi-F-statistics
74 // Temporary INT4 vector for generating histogram of mean multi-F-statistics
76};
77
78///
79/// \name Internal functions
80///
81/// @{
82
83/// @}
84
85///
86/// \name Functions for results toplist ranked by summed multi-detector F-statistic
87///
88/// @{
89
90static const REAL4 *toplist_results_sum2F( const WeaveSemiResults *semi_res )
91{
92 return semi_res->sum2F->data;
93}
94static REAL4 toplist_item_get_sum2F( const WeaveResultsToplistItem *item )
95{
96 return item->stage[0].sum2F;
97}
98static void toplist_item_set_sum2F( WeaveResultsToplistItem *item, const REAL4 value )
99{
100 item->stage[0].sum2F = value;
101}
102
103/// @}
104///
105/// \name Functions for results toplist ranked by mean multi-detector F-statistic
106///
107/// @{
108
109static const REAL4 *toplist_results_mean2F( const WeaveSemiResults *semi_res )
110{
111 return semi_res->mean2F->data;
112}
113static REAL4 toplist_item_get_mean2F( const WeaveResultsToplistItem *item )
114{
115 return item->stage[0].mean2F;
116}
117static void toplist_item_set_mean2F( WeaveResultsToplistItem *item, const REAL4 value )
118{
119 item->stage[0].mean2F = value;
120}
121
122/// @}
123///
124/// \name Functions for results toplist ranked by line-robust log10(B_S/GL) statistic
125///
126/// @{
127
128static const REAL4 *toplist_results_log10BSGL( const WeaveSemiResults *semi_res )
129{
130 return semi_res->log10BSGL->data;
131}
132static REAL4 toplist_item_get_log10BSGL( const WeaveResultsToplistItem *item )
133{
134 return item->stage[0].log10BSGL;
135}
136static void toplist_item_set_log10BSGL( WeaveResultsToplistItem *item, const REAL4 value )
137{
138 item->stage[0].log10BSGL = value;
139}
140
141/// @}
142///
143/// \name Functions for results toplist ranked by transient-line-robust log10(B_S/GLtL) statistic
144///
145/// @{
146
147static const REAL4 *toplist_results_log10BSGLtL( const WeaveSemiResults *semi_res )
148{
149 return semi_res->log10BSGLtL->data;
150}
151static REAL4 toplist_item_get_log10BSGLtL( const WeaveResultsToplistItem *item )
152{
153 return item->stage[0].log10BSGLtL;
154}
155static void toplist_item_set_log10BSGLtL( WeaveResultsToplistItem *item, const REAL4 value )
156{
157 item->stage[0].log10BSGLtL = value;
158}
159
160/// @}
161///
162/// \name Functions for results toplist ranked by transient-signal line-robust log10(B_tS/GLtL) statistic
163///
164/// @{
165
166static const REAL4 *toplist_results_log10BtSGLtL( const WeaveSemiResults *semi_res )
167{
168 return semi_res->log10BtSGLtL->data;
169}
170static REAL4 toplist_item_get_log10BtSGLtL( const WeaveResultsToplistItem *item )
171{
172 return item->stage[0].log10BtSGLtL;
173}
174static void toplist_item_set_log10BtSGLtL( WeaveResultsToplistItem *item, const REAL4 value )
175{
176 item->stage[0].log10BtSGLtL = value;
177}
178
179/// @}
180
181///
182/// Create output results
183///
185 const LIGOTimeGPS *ref_time,
186 const size_t nspins,
187 WeaveStatisticsParams *statistics_params,
188 const UINT4 toplist_limit,
189 const BOOLEAN mean2F_hgrm
190)
191{
192 // Check input
193 XLAL_CHECK_NULL( ref_time != NULL, XLAL_EFAULT );
194 XLAL_CHECK_NULL( statistics_params != NULL, XLAL_EFAULT );
195
196 // Allocate memory
197 WeaveOutputResults *out = XLALCalloc( 1, sizeof( *out ) );
198 XLAL_CHECK_NULL( out != NULL, XLAL_ENOMEM );
199
200 // Set fields
201 out->ref_time = *ref_time;
202 out->nspins = nspins;
203 out->toplist_limit = toplist_limit;
204 out->statistics_params = statistics_params;
205
206 WeaveStatisticType toplist_statistics = statistics_params->toplist_statistics;
207
208 // Initialise number of toplists
209 out->ntoplists = 0;
210
211 // Create a toplist which ranks results by mean multi-detector F-statistic
212 if ( toplist_statistics & WEAVE_STATISTIC_MEAN2F ) {
213 out->toplists[out->ntoplists] = XLALWeaveResultsToplistCreate( nspins, statistics_params, WEAVE_STATISTIC_NAME( WEAVE_STATISTIC_MEAN2F ), "average multi-detector F-statistic", toplist_limit, toplist_results_mean2F, toplist_item_get_mean2F, toplist_item_set_mean2F );
214 XLAL_CHECK_NULL( out->toplists[out->ntoplists] != NULL, XLAL_EFUNC );
215 XLAL_CHECK_NULL( out->ntoplists < XLAL_NUM_ELEM( out->toplists ), XLAL_EFAILED );
216 out->ntoplists++;
217 }
218
219 // Create a toplist which ranks results by summed multi-detector F-statistic
220 if ( toplist_statistics & WEAVE_STATISTIC_SUM2F ) {
221 out->toplists[out->ntoplists] = XLALWeaveResultsToplistCreate( nspins, statistics_params, WEAVE_STATISTIC_NAME( WEAVE_STATISTIC_SUM2F ), "summed multi-detector F-statistic", toplist_limit, toplist_results_sum2F, toplist_item_get_sum2F, toplist_item_set_sum2F );
222 XLAL_CHECK_NULL( out->toplists[out->ntoplists] != NULL, XLAL_EFUNC );
223 XLAL_CHECK_NULL( out->ntoplists < XLAL_NUM_ELEM( out->toplists ), XLAL_EFAILED );
224 out->ntoplists++;
225 }
226
227 // Create a toplist which ranks results by line-robust log10(B_S/GL) statistic
228 if ( toplist_statistics & WEAVE_STATISTIC_BSGL ) {
229 out->toplists[out->ntoplists] = XLALWeaveResultsToplistCreate( nspins, statistics_params, WEAVE_STATISTIC_NAME( WEAVE_STATISTIC_BSGL ), "line-robust log10BSGL statistic", toplist_limit, toplist_results_log10BSGL, toplist_item_get_log10BSGL, toplist_item_set_log10BSGL );
230 XLAL_CHECK_NULL( out->toplists[out->ntoplists] != NULL, XLAL_EFUNC );
231 XLAL_CHECK_NULL( out->ntoplists < XLAL_NUM_ELEM( out->toplists ), XLAL_EFAILED );
232 out->ntoplists++;
233 }
234
235 // Create a toplist which ranks results by transient-line-robust log10(B_S/GLtL) statistic
236 if ( toplist_statistics & WEAVE_STATISTIC_BSGLtL ) {
237 out->toplists[out->ntoplists] = XLALWeaveResultsToplistCreate( nspins, statistics_params, WEAVE_STATISTIC_NAME( WEAVE_STATISTIC_BSGLtL ), "transient line-robust log10BSGLtL statistic", toplist_limit, toplist_results_log10BSGLtL, toplist_item_get_log10BSGLtL, toplist_item_set_log10BSGLtL );
238 XLAL_CHECK_NULL( out->toplists[out->ntoplists] != NULL, XLAL_EFUNC );
239 XLAL_CHECK_NULL( out->ntoplists < XLAL_NUM_ELEM( out->toplists ), XLAL_EFAILED );
240 out->ntoplists++;
241 }
242
243 // Create a toplist which ranks results by transient-signal line-robust log10(B_tS/GLtL) statistic
244 if ( toplist_statistics & WEAVE_STATISTIC_BtSGLtL ) {
245 out->toplists[out->ntoplists] = XLALWeaveResultsToplistCreate( nspins, statistics_params, WEAVE_STATISTIC_NAME( WEAVE_STATISTIC_BtSGLtL ), "transient signal line-robust log10BtSGLtL statistic", toplist_limit, toplist_results_log10BtSGLtL, toplist_item_get_log10BtSGLtL, toplist_item_set_log10BtSGLtL );
246 XLAL_CHECK_NULL( out->toplists[out->ntoplists] != NULL, XLAL_EFUNC );
247 XLAL_CHECK_NULL( out->ntoplists < XLAL_NUM_ELEM( out->toplists ), XLAL_EFAILED );
248 out->ntoplists++;
249 }
250
251 // Comnsistency check on number of toplists
252 XLAL_CHECK_NULL( out->ntoplists == statistics_params->ntoplists, XLAL_EFAILED );
253
254 // Create histogram of mean multi-F-statistic
255 if ( mean2F_hgrm ) {
256 out->mean2F_hgrm_bins = XLALCreateUINT8Vector( 10000 );
257 XLAL_CHECK_NULL( out->mean2F_hgrm_bins != NULL, XLAL_EFUNC );
258 memset( out->mean2F_hgrm_bins->data, 0, sizeof( out->mean2F_hgrm_bins->data[0] ) * out->mean2F_hgrm_bins->length );
259 }
260
261 return out;
262
263}
264
265///
266/// Free output results
267///
269 WeaveOutputResults *out
270)
271{
272 if ( out != NULL ) {
273 XLALWeaveStatisticsParamsDestroy( out->statistics_params );
274 for ( size_t i = 0; i < out->ntoplists; ++i ) {
276 }
277 XLALDestroyUINT8Vector( out->mean2F_hgrm_bins );
278 XLALDestroyREAL4Vector( out->mean2F_hgrm_tmp_REAL4 );
279 XLALDestroyINT4Vector( out->mean2F_hgrm_tmp_INT4 );
280 XLALFree( out );
281 }
282}
283
284///
285/// Add semicoherent results to output
286///
288 WeaveOutputResults *out,
289 const WeaveSemiResults *semi_res,
290 const UINT4 semi_nfreqs
291)
292{
293
294 // Check input
295 XLAL_CHECK( out != NULL, XLAL_EFAULT );
296 XLAL_CHECK( semi_res != NULL, XLAL_EFAULT );
297
298 // Add results to toplists
299 for ( size_t i = 0; i < out->ntoplists; ++i ) {
300 XLAL_CHECK( XLALWeaveResultsToplistAdd( out->toplists[i], semi_res, semi_nfreqs ) == XLAL_SUCCESS, XLAL_EFUNC );
301 }
302
303 // Add to histogram of mean multi-F-statistics
304 if ( out->mean2F_hgrm_bins != NULL ) {
305
306 // Check input
307 XLAL_CHECK( semi_res->mean2F != NULL, XLAL_EINVAL );
308
309 // Allocate memory
310 if ( out->mean2F_hgrm_tmp_REAL4 == NULL || out->mean2F_hgrm_tmp_REAL4->length < semi_res->mean2F->length ) {
311 out->mean2F_hgrm_tmp_REAL4 = XLALResizeREAL4Vector( out->mean2F_hgrm_tmp_REAL4, semi_res->mean2F->length );
312 XLAL_CHECK( out->mean2F_hgrm_tmp_REAL4 != NULL, XLAL_EFUNC );
313 out->mean2F_hgrm_tmp_INT4 = XLALResizeINT4Vector( out->mean2F_hgrm_tmp_INT4, semi_res->mean2F->length );
314 XLAL_CHECK( out->mean2F_hgrm_tmp_INT4 != NULL, XLAL_EFUNC );
315 }
316
317 // Put mean multi-F-statistics into bins
318 XLAL_CHECK( XLALVectorScaleREAL4( out->mean2F_hgrm_tmp_REAL4->data, 1.0 / mean2F_hgrm_bin_width, semi_res->mean2F->data, semi_res->mean2F->length ) == XLAL_SUCCESS, XLAL_EFUNC );
319 XLAL_CHECK( XLALVectorINT4FromREAL4( out->mean2F_hgrm_tmp_INT4->data, out->mean2F_hgrm_tmp_REAL4->data, semi_res->mean2F->length ) == XLAL_SUCCESS, XLAL_EFUNC );
320
321 // Add to histogram bins
322 for ( size_t j = 0; j < semi_res->mean2F->length; ++j ) {
323 INT4 bin = out->mean2F_hgrm_tmp_INT4->data[j];
324 const INT4 bin_max = out->mean2F_hgrm_bins->length;
325 if ( bin < 0 ) {
326 ++out->mean2F_hgrm_underflow;
327 } else if ( bin >= bin_max ) {
328 ++out->mean2F_hgrm_overflow;
329 } else {
330 ++out->mean2F_hgrm_bins->data[bin];
331 }
332 }
333
334 }
335
336 return XLAL_SUCCESS;
337
338}
339
340///
341/// Compute all the missing 'completion-loop' statistics for all toplist entries
342///
344 WeaveOutputResults *out
345)
346{
347 // Check input
348 XLAL_CHECK( out != NULL, XLAL_EFAULT );
349
350 // Iterate over all toplists
351 for ( size_t i = 0; i < out->ntoplists; ++i ) {
353 }
354
355 return XLAL_SUCCESS;
356
357}
358
359///
360/// Write output results to a FITS file
361///
363 FITSFile *file,
364 const WeaveOutputResults *out
365)
366{
367
368 // Check input
369 XLAL_CHECK( file != NULL, XLAL_EFAULT );
370 XLAL_CHECK( out != NULL, XLAL_EFAULT );
371
372 // Write reference time
373 XLAL_CHECK( XLALFITSHeaderWriteGPSTime( file, "date-obs", &out->ref_time, "reference time" ) == XLAL_SUCCESS, XLAL_EFUNC );
374
375 // Write number of spindowns
376 XLAL_CHECK( XLALFITSHeaderWriteUINT4( file, "nspins", out->nspins, "number of spindowns" ) == XLAL_SUCCESS, XLAL_EFUNC );
377
378 // Write list of detectors (if outputting per-detector quantities
379 XLAL_CHECK( XLALFITSHeaderWriteStringVector( file, "detect", out->statistics_params->detectors, "list of detectors" ) == XLAL_SUCCESS, XLAL_EFUNC );
380
381 // Write number of segments
382 XLAL_CHECK( XLALFITSHeaderWriteUINT4( file, "segments", out->statistics_params->nsegments, "number of segments" ) == XLAL_SUCCESS, XLAL_EFUNC );
383
384 // Write names of selected toplist types (ie ranking statistics)
385 {
386 char *toplist_statistics = XLALPrintStringValueOfUserFlag( ( const int * ) & ( out->statistics_params->toplist_statistics ), &WeaveToplistChoices );
387 XLAL_CHECK( toplist_statistics != NULL, XLAL_EFUNC );
388 XLAL_CHECK( XLALFITSHeaderWriteString( file, "toplists", toplist_statistics, "names of selected toplist statistics" ) == XLAL_SUCCESS, XLAL_EFUNC );
389 XLALFree( toplist_statistics );
390 }
391
392 // Write names of all selected 'extra' output statistics
393 {
394 WeaveStatisticType extras = ( out->statistics_params->statistics_to_output[0] & ~out->statistics_params->toplist_statistics );
395 if ( extras == 0 ) {
396 extras = WEAVE_STATISTIC_NONE;
397 }
398 char *extras_names = XLALPrintStringValueOfUserFlag( ( const int * )&extras, &WeaveStatisticChoices );
399 XLAL_CHECK( extras_names != NULL, XLAL_EFUNC );
400 XLAL_CHECK( XLALFITSHeaderWriteString( file, "extras", extras_names, "names of additional selected output statistics" ) == XLAL_SUCCESS, XLAL_EFUNC );
401 XLALFree( extras_names );
402 }
403 // Write names of all requested 'recalc' statistics
404 {
405 char *recalc_names = XLALPrintStringValueOfUserFlag( ( const int * )&out->statistics_params->statistics_to_output[1], &WeaveStatisticChoices );
406 XLAL_CHECK( recalc_names != NULL, XLAL_EFUNC );
407 XLAL_CHECK( XLALFITSHeaderWriteString( file, "recalc", recalc_names, "names of selected recalc statistics" ) == XLAL_SUCCESS, XLAL_EFUNC );
408 XLALFree( recalc_names );
409 }
410
411 // Write maximum size of toplists
412 XLAL_CHECK( XLALFITSHeaderWriteUINT4( file, "toplimit", out->toplist_limit, "maximum size of toplists" ) == XLAL_SUCCESS, XLAL_EFUNC );
413
414 // Write whether a histogram of mean multi-F-statistics will be written
415 XLAL_CHECK( XLALFITSHeaderWriteBOOLEAN( file, "m2Fhgrm", out->mean2F_hgrm_bins != NULL, "mean 2F histogram?" ) == XLAL_SUCCESS, XLAL_EFUNC );
416
417 // Write toplists
418 for ( size_t i = 0; i < out->ntoplists; ++i ) {
420 }
421
422 // Write histogram of mean multi-F-statistics
423 if ( out->mean2F_hgrm_bins != NULL ) {
424
425 // Open and describe FITS table for writing histogram bins
426 XLAL_CHECK( XLALFITSTableOpenWrite( file, "mean2F_hgrm", "histogram of mean multi-F-statistics" ) == XLAL_SUCCESS, XLAL_EFUNC );
431
432 // Write histogram bins
434 if ( out->mean2F_hgrm_underflow > 0 ) {
435 bin.lower = GSL_NEGINF;
436 bin.upper = 0;
437 bin.count = out->mean2F_hgrm_underflow;
439 }
440 for ( size_t j = 0; j < out->mean2F_hgrm_bins->length; ++j ) {
441 const UINT4 bin_count = out->mean2F_hgrm_bins->data[j];
442 if ( bin_count > 0 ) {
444 bin.upper = mean2F_hgrm_bin_width * ( j + 1 );
445 bin.count = bin_count;
447 }
448 }
449 if ( out->mean2F_hgrm_underflow > 0 ) {
450 bin.lower = mean2F_hgrm_bin_width * out->mean2F_hgrm_bins->length;
451 bin.upper = GSL_POSINF;
452 bin.count = out->mean2F_hgrm_underflow;
454 }
455
456 }
457
458 return XLAL_SUCCESS;
459
460}
461
462///
463/// Read results from a FITS file and append to new/existing output results
464///
466 FITSFile *file,
467 WeaveOutputResults **out,
468 UINT4 toplist_limit
469)
470{
471
472 // Check input
473 XLAL_CHECK( file != NULL, XLAL_EFAULT );
474 XLAL_CHECK( out != NULL, XLAL_EFAULT );
475
476 // Read reference time
477 LIGOTimeGPS ref_time;
478 XLAL_CHECK( XLALFITSHeaderReadGPSTime( file, "date-obs", &ref_time ) == XLAL_SUCCESS, XLAL_EFUNC );
479
480 // Read number of spindowns
481 UINT4 nspins = 0;
483
484 // ----- read elements for 'statistics_params' struct ----------
485 WeaveStatisticsParams *statistics_params = XLALCalloc( 1, sizeof( *statistics_params ) );
486 XLAL_CHECK( statistics_params != NULL, XLAL_ENOMEM );
487
488 // Read list of detectors
489 XLAL_CHECK( XLALFITSHeaderReadStringVector( file, "detect", &( statistics_params->detectors ) ) == XLAL_SUCCESS, XLAL_EFUNC );
490
491 /// Number of segments
492 XLAL_CHECK( XLALFITSHeaderReadUINT4( file, "segments", &( statistics_params->nsegments ) ) == XLAL_SUCCESS, XLAL_EFUNC );
493
494 // Read names of selected toplist statistics
495 char *toplist_stats_names = NULL;
496 int toplist_stats = 0;
497 XLAL_CHECK( XLALFITSHeaderReadString( file, "toplists", &toplist_stats_names ) == XLAL_SUCCESS, XLAL_EFUNC );
498 XLAL_CHECK( XLALParseStringValueAsUserFlag( &toplist_stats, &WeaveToplistChoices, toplist_stats_names ) == XLAL_SUCCESS, XLAL_EFUNC );
499 XLALFree( toplist_stats_names );
500
501 // Read names of selected extra output stats
502 char *extras_names = NULL;
503 int extra_stats = 0;
504 XLAL_CHECK( XLALFITSHeaderReadString( file, "extras", &extras_names ) == XLAL_SUCCESS, XLAL_EFUNC );
506 XLALFree( extras_names );
507
508 // Read names of selected recalc stats
509 int recalc_stats = 0;
510 BOOLEAN exists = 0;
512 if ( exists ) {
513 char *recalc_names = NULL;
514 XLAL_CHECK( XLALFITSHeaderReadString( file, "recalc", &recalc_names ) == XLAL_SUCCESS, XLAL_EFUNC );
516 XLALFree( recalc_names );
517 }
518
519 // Compute and fill the full stats-dependency map
520 XLAL_CHECK( XLALWeaveStatisticsParamsSetDependencyMap( statistics_params, toplist_stats, extra_stats, recalc_stats ) == XLAL_SUCCESS, XLAL_EFUNC );
521
522 // Read whether a histogram of mean multi-F-statistics will be written
523 BOOLEAN mean2F_hgrm = 0;
524 {
525 exists = 0;
527 if ( exists ) {
528 XLAL_CHECK( XLALFITSHeaderReadBOOLEAN( file, "m2Fhgrm", &mean2F_hgrm ) == XLAL_SUCCESS, XLAL_EFUNC );
529 }
530 }
531
532 if ( *out == NULL ) {
533
534 // Create new output results
535 *out = XLALWeaveOutputResultsCreate( &ref_time, nspins, statistics_params, toplist_limit, mean2F_hgrm );
536 XLAL_CHECK( *out != NULL, XLAL_EFUNC );
537
538 } else {
539
540 // Check reference time
541 XLAL_CHECK( XLALGPSCmp( &ref_time, &( *out )->ref_time ) == 0, XLAL_EIO, "Inconsistent reference time: %" LAL_GPS_FORMAT " != %" LAL_GPS_FORMAT, LAL_GPS_PRINT( ref_time ), LAL_GPS_PRINT( ( *out )->ref_time ) );
542
543 // Check number of spindowns
544 XLAL_CHECK( ( size_t ) nspins == ( *out )->nspins, XLAL_EIO, "Inconsistent number of spindowns: %i != %zu", nspins, ( *out )->nspins );
545
546 // Check if list of detectors agrees
547 if ( statistics_params->detectors != NULL ) {
548 XLAL_CHECK( statistics_params->detectors->length == ( *out )->statistics_params->detectors->length, XLAL_EIO, "Inconsistent number of detectors: %u != %u", statistics_params->detectors->length, ( *out )->statistics_params->detectors->length );
549 for ( size_t i = 0; i < statistics_params->detectors->length; ++i ) {
550 XLAL_CHECK( strcmp( statistics_params->detectors->data[i], ( *out )->statistics_params->detectors->data[i] ) == 0, XLAL_EIO, "Inconsistent detectors: %s != %s", statistics_params->detectors->data[i], ( *out )->statistics_params->detectors->data[i] );
551 }
552 }
553
554 // Check if number of segments agrees
555 XLAL_CHECK( statistics_params->nsegments == ( *out )->statistics_params->nsegments, XLAL_EIO, "Inconsistent output per segment?: %i != %u", statistics_params->nsegments, ( *out )->statistics_params->nsegments );
556
557 // Check list of selected toplist statistics
558 if ( statistics_params->toplist_statistics != ( *out )->statistics_params->toplist_statistics ) {
559 char *toplists1, *toplists2;
560 toplists1 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( statistics_params->toplist_statistics ), &WeaveToplistChoices );
561 XLAL_CHECK( toplists1 != NULL, XLAL_EFUNC );
562 toplists2 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( ( *out )->statistics_params->toplist_statistics ), &WeaveToplistChoices );
563 XLAL_CHECK( toplists2 != NULL, XLAL_EFUNC );
564 XLALPrintError( "Inconsistent set of toplist statistics: %s != %s\n", toplists1, toplists2 );
565 XLALFree( toplists1 );
566 XLALFree( toplists2 );
568 }
569 // Check list of selected output statistics
570 for ( UINT4 istage = 0; istage < 2; ++ istage ) {
571 if ( statistics_params->statistics_to_output[istage] != ( *out )->statistics_params->statistics_to_output[istage] ) {
572 char *output1, *output2;
573 output1 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( statistics_params->statistics_to_output[istage] ), &WeaveStatisticChoices );
574 XLAL_CHECK( output1 != NULL, XLAL_EFUNC );
575 output2 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( ( *out )->statistics_params->statistics_to_output[istage] ), &WeaveStatisticChoices );
576 XLAL_CHECK( output2 != NULL, XLAL_EFUNC );
577 XLALPrintError( "Inconsistent set of stage-%d output statistics: {%s} != {%s}\n", istage, output1, output2 );
578 XLALFree( output1 );
579 XLALFree( output2 );
581 }
582 }
583
584 XLALWeaveStatisticsParamsDestroy( statistics_params ); // Not creating a new output, so we need to free this
585
586 // Check whether a histogram of mean multi-F-statistics will be written
587 XLAL_CHECK( !mean2F_hgrm == !( ( *out )->mean2F_hgrm_bins != NULL ), XLAL_EIO, "Inconsistent mean 2F histogram? %i != %i", mean2F_hgrm, ( *out )->mean2F_hgrm_bins != NULL );
588
589 }
590
591 // Read and append to toplists
592 for ( size_t i = 0; i < ( *out )->ntoplists; ++i ) {
594 }
595
596 // Read and append histogram of mean multi-F-statistics
597 if ( ( *out )->mean2F_hgrm_bins != NULL ) {
598
599 // Open and describe FITS table for writing histogram bins
600 UINT8 nrows = 0;
601 XLAL_CHECK( XLALFITSTableOpenRead( file, "mean2F_hgrm", &nrows ) == XLAL_SUCCESS, XLAL_EFUNC );
606
607 // Read histogram bins
609 const REAL4 bin_upper_max = mean2F_hgrm_bin_width * ( *out )->mean2F_hgrm_bins->length;
610 while ( nrows > 0 ) {
612 if ( bin.lower < 0 ) {
613 ( *out )->mean2F_hgrm_underflow += bin.count;
614 } else if ( bin.upper >= bin_upper_max ) {
615 ( *out )->mean2F_hgrm_overflow += bin.count;
616 } else {
617 const size_t j = bin.lower / mean2F_hgrm_bin_width;
618 ( *out )->mean2F_hgrm_bins->data[j] += bin.count;
619 }
620 }
621
622 }
623
624 return XLAL_SUCCESS;
625
626}
627
628///
629/// Compare two output results and return whether they are equal
630///
632 BOOLEAN *equal,
633 const WeaveSetupData *setup,
634 const BOOLEAN sort_by_semi_phys,
635 const UINT4 round_param_to_dp,
636 const UINT4 round_param_to_sf,
637 const REAL8 unmatched_item_tol,
638 const REAL8 param_tol_mism,
639 const VectorComparison *result_tol,
640 const UINT4 toplist_compare_limit,
641 const WeaveOutputResults *out_1,
642 const WeaveOutputResults *out_2
643)
644{
645
646 // Check input
647 XLAL_CHECK( equal != NULL, XLAL_EFAULT );
648 XLAL_CHECK( setup != NULL, XLAL_EFAULT );
649 XLAL_CHECK( unmatched_item_tol >= 0, XLAL_EINVAL );
650 XLAL_CHECK( param_tol_mism >= 0, XLAL_EINVAL );
651 XLAL_CHECK( result_tol != NULL, XLAL_EFAULT );
652 XLAL_CHECK( out_1 != NULL, XLAL_EFAULT );
653 XLAL_CHECK( out_2 != NULL, XLAL_EFAULT );
654
655 // Output results are assumed equal until we find otherwise
656 *equal = 1;
657
658 // Compare reference times
659 if ( XLALGPSCmp( &out_1->ref_time, &out_2->ref_time ) != 0 ) {
660 *equal = 0;
661 XLALPrintInfo( "%s: unequal reference times: %" LAL_GPS_FORMAT " != %" LAL_GPS_FORMAT "\n", __func__, LAL_GPS_PRINT( out_1->ref_time ), LAL_GPS_PRINT( out_2->ref_time ) );
662 return XLAL_SUCCESS;
663 }
664
665 // Compare number of spindowns
666 if ( out_1->nspins != out_2->nspins ) {
667 *equal = 0;
668 XLALPrintInfo( "%s: unequal number of spindowns: %zu != %zu\n", __func__, out_1->nspins, out_2->nspins );
669 return XLAL_SUCCESS;
670 }
671
672 // Compare list of detectors
673 {
674 const BOOLEAN have_det_1 = ( out_1->statistics_params->detectors != NULL ), have_det_2 = ( out_2->statistics_params->detectors != NULL );
675 if ( have_det_1 != have_det_2 ) {
676 *equal = 0;
677 XLALPrintInfo( "%s: unequal presence of detector lists: %i != %i\n", __func__, have_det_1, have_det_2 );
678 return XLAL_SUCCESS;
679 }
680 if ( have_det_1 ) {
681 if ( out_1->statistics_params->detectors->length != out_2->statistics_params->detectors->length ) {
682 *equal = 0;
683 XLALPrintInfo( "%s: unequal number of detectors: %u != %u\n", __func__, out_1->statistics_params->detectors->length, out_2->statistics_params->detectors->length );
684 return XLAL_SUCCESS;
685 }
686 for ( size_t i = 0; i < out_1->statistics_params->detectors->length; ++i ) {
687 if ( strcmp( out_1->statistics_params->detectors->data[i], out_2->statistics_params->detectors->data[i] ) != 0 ) {
688 *equal = 0;
689 XLALPrintInfo( "%s: unequal detectors: %s != %s\n", __func__, out_1->statistics_params->detectors->data[i], out_2->statistics_params->detectors->data[i] );
690 return XLAL_SUCCESS;
691 }
692 }
693 }
694 }
695
696 // Compare number of per-segment items
697 if ( out_1->statistics_params->nsegments != out_2->statistics_params->nsegments ) {
698 *equal = 0;
699 XLALPrintInfo( "%s: unequal number of segments: %u != %u\n", __func__, out_1->statistics_params->nsegments, out_2->statistics_params->nsegments );
700 return XLAL_SUCCESS;
701 }
702
703 // Compare toplist statistics
704 if ( out_1->statistics_params->toplist_statistics != out_2->statistics_params->toplist_statistics ) {
705 *equal = 0;
706 char *toplists1, *toplists2;
707 toplists1 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( out_1->statistics_params->toplist_statistics ), &WeaveToplistChoices );
708 XLAL_CHECK( toplists1 != NULL, XLAL_EFUNC );
709 toplists2 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( out_2->statistics_params->toplist_statistics ), &WeaveToplistChoices );
710 XLAL_CHECK( toplists2 != NULL, XLAL_EFUNC );
711 XLALPrintError( "%s: Inconsistent set of toplist statistics: {%s} != {%s}\n", __func__, toplists1, toplists2 );
712 XLALFree( toplists1 );
713 XLALFree( toplists2 );
714 return XLAL_SUCCESS;
715 }
716
717 // Compare statistics_to_output
718 for ( UINT4 istage = 0; istage < 2; ++ istage ) {
719 if ( out_1->statistics_params->statistics_to_output[istage] != out_2->statistics_params->statistics_to_output[istage] ) {
720 *equal = 0;
721 char *outputs1, *outputs2;
722 outputs1 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( out_1->statistics_params->statistics_to_output[istage] ), &WeaveStatisticChoices );
723 XLAL_CHECK( outputs1 != NULL, XLAL_EFUNC );
724 outputs2 = XLALPrintStringValueOfUserFlag( ( const int * ) & ( out_2->statistics_params->statistics_to_output[istage] ), &WeaveStatisticChoices );
725 XLAL_CHECK( outputs2 != NULL, XLAL_EFUNC );
726 XLALPrintError( "%s: Inconsistent set of stage-%d ouput statistics: {%s} != {%s}\n", __func__, istage, outputs1, outputs2 );
727 XLALFree( outputs1 );
728 XLALFree( outputs2 );
729 return XLAL_SUCCESS;
730 }
731 }
732
733 // Compare toplists
734 for ( size_t i = 0; i < out_1->ntoplists; ++i ) {
736 setup, sort_by_semi_phys,
737 round_param_to_dp, round_param_to_sf, unmatched_item_tol, param_tol_mism, result_tol, toplist_compare_limit,
738 out_1->toplists[i], out_2->toplists[i]
739 ) == XLAL_SUCCESS, XLAL_EFUNC );
740 if ( !*equal ) {
741 return XLAL_SUCCESS;
742 }
743 }
744
745
746 return XLAL_SUCCESS;
747
748}
749
750// Local Variables:
751// c-file-style: "linux"
752// c-basic-offset: 2
753// End:
#define __func__
log an I/O error, i.e.
int XLALFITSHeaderWriteUINT4(FITSFile UNUSED *file, const CHAR UNUSED *key, const UINT4 UNUSED value, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:762
int XLALFITSTableWriteRow(FITSFile UNUSED *file, const void UNUSED *record)
Definition: FITSFileIO.c:2550
int XLALFITSHeaderWriteStringVector(FITSFile UNUSED *file, const CHAR UNUSED *key, const LALStringVector UNUSED *values, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:1412
int XLALFITSHeaderWriteGPSTime(FITSFile UNUSED *file, const CHAR UNUSED *key, const LIGOTimeGPS UNUSED *value, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:1494
int XLALFITSHeaderWriteBOOLEAN(FITSFile UNUSED *file, const CHAR UNUSED *key, const BOOLEAN UNUSED value, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:635
int XLALFITSTableReadRow(FITSFile UNUSED *file, void UNUSED *record, UINT8 UNUSED *rem_nrows)
Definition: FITSFileIO.c:2621
int XLALFITSHeaderQueryKeyExists(FITSFile UNUSED *file, const CHAR UNUSED *key, BOOLEAN UNUSED *exists)
Definition: FITSFileIO.c:563
int XLALFITSHeaderReadBOOLEAN(FITSFile UNUSED *file, const CHAR UNUSED *key, BOOLEAN UNUSED *value)
Definition: FITSFileIO.c:668
int XLALFITSHeaderReadGPSTime(FITSFile UNUSED *file, const CHAR UNUSED *key, LIGOTimeGPS UNUSED *value)
Definition: FITSFileIO.c:1547
int XLALFITSTableOpenWrite(FITSFile UNUSED *file, const CHAR UNUSED *name, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:2162
int XLALFITSHeaderReadStringVector(FITSFile UNUSED *file, const CHAR UNUSED *key, LALStringVector UNUSED **values)
Definition: FITSFileIO.c:1456
int XLALFITSHeaderReadUINT4(FITSFile UNUSED *file, const CHAR UNUSED *key, UINT4 UNUSED *value)
Definition: FITSFileIO.c:797
int XLALFITSTableOpenRead(FITSFile UNUSED *file, const CHAR UNUSED *name, UINT8 UNUSED *nrows)
Definition: FITSFileIO.c:2198
int XLALFITSHeaderReadString(FITSFile UNUSED *file, const CHAR UNUSED *key, CHAR UNUSED **value)
Definition: FITSFileIO.c:1377
int XLALFITSHeaderWriteString(FITSFile UNUSED *file, const CHAR UNUSED *key, const CHAR UNUSED *value, const CHAR UNUSED *comment)
Definition: FITSFileIO.c:1339
int j
int XLALWeaveOutputResultsAdd(WeaveOutputResults *out, const WeaveSemiResults *semi_res, const UINT4 semi_nfreqs)
Add semicoherent results to output.
static const REAL4 * toplist_results_log10BtSGLtL(const WeaveSemiResults *semi_res)
static REAL4 toplist_item_get_log10BSGLtL(const WeaveResultsToplistItem *item)
static void toplist_item_set_mean2F(WeaveResultsToplistItem *item, const REAL4 value)
static void toplist_item_set_sum2F(WeaveResultsToplistItem *item, const REAL4 value)
Definition: OutputResults.c:98
int XLALWeaveOutputResultsCompare(BOOLEAN *equal, const WeaveSetupData *setup, const BOOLEAN sort_by_semi_phys, const UINT4 round_param_to_dp, const UINT4 round_param_to_sf, const REAL8 unmatched_item_tol, const REAL8 param_tol_mism, const VectorComparison *result_tol, const UINT4 toplist_compare_limit, const WeaveOutputResults *out_1, const WeaveOutputResults *out_2)
Compare two output results and return whether they are equal.
static void toplist_item_set_log10BtSGLtL(WeaveResultsToplistItem *item, const REAL4 value)
int XLALWeaveOutputResultsWrite(FITSFile *file, const WeaveOutputResults *out)
Write output results to a FITS file.
static const REAL4 * toplist_results_sum2F(const WeaveSemiResults *semi_res)
Definition: OutputResults.c:90
static const REAL4 * toplist_results_log10BSGLtL(const WeaveSemiResults *semi_res)
static REAL4 toplist_item_get_sum2F(const WeaveResultsToplistItem *item)
Definition: OutputResults.c:94
static REAL4 toplist_item_get_log10BSGL(const WeaveResultsToplistItem *item)
int XLALWeaveOutputResultsReadAppend(FITSFile *file, WeaveOutputResults **out, UINT4 toplist_limit)
Read results from a FITS file and append to new/existing output results.
static void toplist_item_set_log10BSGL(WeaveResultsToplistItem *item, const REAL4 value)
static const REAL4 * toplist_results_mean2F(const WeaveSemiResults *semi_res)
void XLALWeaveOutputResultsDestroy(WeaveOutputResults *out)
Free output results.
static void toplist_item_set_log10BSGLtL(WeaveResultsToplistItem *item, const REAL4 value)
int XLALWeaveOutputResultsCompletionLoop(WeaveOutputResults *out)
Compute all the missing 'completion-loop' statistics for all toplist entries.
static const REAL4 * toplist_results_log10BSGL(const WeaveSemiResults *semi_res)
static REAL4 toplist_item_get_log10BtSGLtL(const WeaveResultsToplistItem *item)
WeaveOutputResults * XLALWeaveOutputResultsCreate(const LIGOTimeGPS *ref_time, const size_t nspins, WeaveStatisticsParams *statistics_params, const UINT4 toplist_limit, const BOOLEAN mean2F_hgrm)
Create output results.
static REAL4 toplist_item_get_mean2F(const WeaveResultsToplistItem *item)
const REAL4 mean2F_hgrm_bin_width
Histogram bin width of mean multi-F-statistics.
Definition: OutputResults.c:34
Module which handles the output results.
int XLALWeaveResultsToplistAdd(WeaveResultsToplist *toplist, const WeaveSemiResults *semi_res, const UINT4 semi_nfreqs)
Add semicoherent results to toplist.
void XLALWeaveResultsToplistDestroy(WeaveResultsToplist *toplist)
Free results toplist.
WeaveResultsToplist * XLALWeaveResultsToplistCreate(const size_t nspins, WeaveStatisticsParams *statistics_params, const char *stat_name, const char *stat_desc, const UINT4 toplist_limit, WeaveResultsToplistRankingStats toplist_rank_stats_fcn, WeaveResultsToplistItemGetRankStat toplist_item_get_rank_stat_fcn, WeaveResultsToplistItemSetRankStat toplist_item_set_rank_stat_fcn)
Create results toplist.
int XLALWeaveResultsToplistReadAppend(FITSFile *file, WeaveResultsToplist *toplist)
Read results from a FITS file and append to existing results toplist.
int XLALWeaveResultsToplistCompare(BOOLEAN *equal, const WeaveSetupData *setup, const BOOLEAN sort_by_semi_phys, const UINT4 round_param_to_dp, const UINT4 round_param_to_sf, const REAL8 unmatched_item_tol, const REAL8 param_tol_mism, const VectorComparison *result_tol, const UINT4 toplist_compare_limit, const WeaveResultsToplist *toplist_1, const WeaveResultsToplist *toplist_2)
Compare two results toplists and return whether they are equal.
int XLALWeaveResultsToplistWrite(FITSFile *file, const WeaveResultsToplist *toplist)
Write results toplist to a FITS file.
int XLALWeaveResultsToplistCompletionLoop(WeaveResultsToplist *toplist)
Compute all missing 'extra' (non-toplist-ranking) statistics for all toplist entries.
Module which handles toplists of results.
enum tagWeaveStatisticType WeaveStatisticType
Definition: Weave.h:61
int XLALWeaveStatisticsParamsSetDependencyMap(WeaveStatisticsParams *statistics_params, const WeaveStatisticType toplist_stats, const WeaveStatisticType extra_output_stats, const WeaveStatisticType recalc_stats)
Fill StatisticsParams logic for given toplist and extra-output stats.
const UserChoices WeaveToplistChoices
User input choices for toplist ranking statistics.
const UserChoices WeaveStatisticChoices
User input choices for all supported statistics.
void XLALWeaveStatisticsParamsDestroy(WeaveStatisticsParams *statistics_params)
Destroy a StatisticsParams struct.
@ WEAVE_STATISTIC_BtSGLtL
(transient-)line robust log10(B_tS/GLtL) statistic
@ WEAVE_STATISTIC_MEAN2F
Multi-detector average (over segments) F-statistic.
@ WEAVE_STATISTIC_BSGLtL
(transient-)line robust log10(B_S/GLtL) statistic
@ WEAVE_STATISTIC_BSGL
Line-robust log10(B_S/GL) statistic.
@ WEAVE_STATISTIC_NONE
@ WEAVE_STATISTIC_SUM2F
Multi-detector sum (over segments) F-statistic.
#define WEAVE_STATISTIC_NAME(ws)
Names of all possible statistics.
#define LAL_GPS_FORMAT
#define LAL_GPS_PRINT(gps)
#define XLAL_FITS_TABLE_COLUMN_BEGIN(record_type)
Definition: FITSFileIO.h:239
#define XLAL_FITS_TABLE_COLUMN_ADD(file, type, field)
Definition: FITSFileIO.h:243
struct tagFITSFile FITSFile
Representation of a FITS file.
Definition: FITSFileIO.h:54
unsigned char BOOLEAN
#define XLAL_NUM_ELEM(x)
uint64_t UINT8
double REAL8
uint32_t UINT4
int32_t INT4
float REAL4
void * XLALCalloc(size_t m, size_t n)
void XLALFree(void *p)
int XLALParseStringValueAsUserFlag(int *valFlag, const UserChoices *flagData, const char *valString)
char * XLALPrintStringValueOfUserFlag(const int *valFlag, const UserChoices *flagData)
INT4Vector * XLALResizeINT4Vector(INT4Vector *vector, UINT4 length)
void XLALDestroyREAL4Vector(REAL4Vector *vector)
void XLALDestroyINT4Vector(INT4Vector *vector)
void XLALDestroyUINT8Vector(UINT8Vector *vector)
REAL4Vector * XLALResizeREAL4Vector(REAL4Vector *vector, UINT4 length)
UINT8Vector * XLALCreateUINT8Vector(UINT4 length)
int XLALVectorScaleREAL4(REAL4 *out, REAL4 scalar, const REAL4 *in, const UINT4 len)
int XLALVectorINT4FromREAL4(INT4 *out, const REAL4 *in, const UINT4 len)
int int int XLALPrintInfo(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_ERROR(...)
#define XLAL_CHECK(assertion,...)
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_CHECK_NULL(assertion,...)
XLAL_ENOMEM
XLAL_SUCCESS
XLAL_EFAULT
XLAL_EFUNC
XLAL_EIO
XLAL_EINVAL
XLAL_EFAILED
int XLALGPSCmp(const LIGOTimeGPS *t0, const LIGOTimeGPS *t1)
long long count
Definition: hwinject.c:371
out
Struct holding the results of comparing two floating-point vectors (real-valued or complex),...
Definition: LFTandTSutils.h:64
Histogram bin of mean multi-F-statistics.
Definition: OutputResults.c:39
Output results from a search.
Definition: OutputResults.c:51
REAL4Vector * mean2F_hgrm_tmp_REAL4
Definition: OutputResults.c:73
WeaveResultsToplist * toplists[8]
Output result toplists.
Definition: OutputResults.c:65
LIGOTimeGPS ref_time
Reference time at which search is conducted.
Definition: OutputResults.c:57
size_t nspins
Number of spindown parameters to output.
Definition: OutputResults.c:59
INT4Vector * mean2F_hgrm_tmp_INT4
Definition: OutputResults.c:75
size_t ntoplists
Number of output results toplists.
Definition: OutputResults.c:63
UINT8Vector * mean2F_hgrm_bins
Definition: OutputResults.c:67
WeaveStatisticsParams * statistics_params
Struct holding all parameters for which statistics to output and compute, when, and how NOTE: this is...
Definition: OutputResults.c:55