From 556e8d820c2c96d406dbe4b9cf5d6cf7f4b4779f Mon Sep 17 00:00:00 2001 From: Gregory Ashton Date: Thu, 14 Dec 2017 16:43:26 +0100 Subject: [PATCH 1/3] Adds new option 'SSBPREC_DMOFF' to SSB-precision in XLALGetSSBtimes() - Issue originally discussed on https://bugs.ligo.org/redmine/issues/5514 - only adds the 'DMOFF' options, other variants have been postponed. --- lalpulsar/src/SSBtimes.c | 48 ++++++++++++++++++++++++++++++++++++++++ lalpulsar/src/SSBtimes.h | 1 + 2 files changed, 49 insertions(+) diff --git a/lalpulsar/src/SSBtimes.c b/lalpulsar/src/SSBtimes.c index 1a254bcd66..9865ca327c 100644 --- a/lalpulsar/src/SSBtimes.c +++ b/lalpulsar/src/SSBtimes.c @@ -637,6 +637,54 @@ XLALGetSSBtimes ( const DetectorStateSeries *DetectorStates, /**< [in] detector- break; + case SSBPREC_DMOFF: /* switch off all demodulation terms */ + + baryinput.site = DetectorStates->detector; + baryinput.site.location[0] /= LAL_C_SI; + baryinput.site.location[1] /= LAL_C_SI; + baryinput.site.location[2] /= LAL_C_SI; + + baryinput.alpha = alpha; + baryinput.delta = delta; + baryinput.dInv = 0; + + for ( UINT4 i = 0; i < numSteps; i++ ) + { + EmissionTime emit; + DetectorState *state = &(DetectorStates->data[i]); + baryinput.tgps = state->tGPS; + + REAL8 tgps[2]; + tgps[0] = baryinput.tgps.gpsSeconds; + tgps[1] = baryinput.tgps.gpsNanoSeconds; + + if ( XLALBarycenterOpt ( &emit, &baryinput, &(state->earthState), &bBuffer ) != XLAL_SUCCESS ) + XLAL_ERROR_NULL ( XLAL_EFUNC, "XLALBarycenterOpt() failed with xlalErrno = %d\n", xlalErrno ); + + emit.deltaT = 0; + emit.tDot = 1; + INT4 deltaTint = floor ( emit.deltaT ); + + if ( ( 1e-9 * tgps[1] + emit.deltaT - deltaTint ) >= 1.e0 ) + { + emit.te.gpsSeconds = baryinput.tgps.gpsSeconds + deltaTint + 1; + emit.te.gpsNanoSeconds = floor ( 1e9 * ( tgps[1] * 1e-9 + emit.deltaT - deltaTint - 1.0 ) ); + } + else + { + emit.te.gpsSeconds = baryinput.tgps.gpsSeconds + deltaTint; + emit.te.gpsNanoSeconds = floor ( 1e9 * ( tgps[1] * 1e-9 + emit.deltaT - deltaTint ) ); + } + + ret->DeltaT->data[i] = XLALGPSGetREAL8 ( &emit.te ) - refTimeREAL8; + ret->Tdot->data[i] = emit.tDot; + + } /* for i < numSteps */ + // free buffer memory + XLALFree ( bBuffer ); + + break; + default: XLAL_ERROR_NULL (XLAL_EFAILED, "\n?? Something went wrong.. this should never be called!\n\n" ); break; diff --git a/lalpulsar/src/SSBtimes.h b/lalpulsar/src/SSBtimes.h index 05e51f703f..a3c68b32bc 100644 --- a/lalpulsar/src/SSBtimes.h +++ b/lalpulsar/src/SSBtimes.h @@ -46,6 +46,7 @@ typedef enum tagSSBprecision { SSBPREC_NEWTONIAN, /**< simple Newtonian: \f$\tau = t + \vec{r}\cdot\vec{n}/c\f$ */ SSBPREC_RELATIVISTIC, /**< detailed relativistic: \f$\tau=\tau(t; \vec{n}, \vec{r})\f$ */ SSBPREC_RELATIVISTICOPT, /**< optimized relativistic, numerically equivalent to #SSBPREC_RELATIVISTIC, but faster */ + SSBPREC_DMOFF, /**< switch off all demodulatoin terms */ SSBPREC_LAST /**< end marker */ } SSBprecision; -- GitLab From 1eca2f0ac0d046e2c710ec0ee026fb276a044724 Mon Sep 17 00:00:00 2001 From: Reinhard Prix Date: Sun, 28 Jan 2018 11:33:14 +0100 Subject: [PATCH 2/3] amend previous 'DMoff' patch: [to be squashed] - simplify 'SSBPREC_DMOFF' implementation - add 'DMoff' option to --SSBprecision code arguments by extending SSBprecisionChoices --- lalpulsar/src/SSBtimes.c | 43 +++------------------------------------- lalpulsar/src/SSBtimes.h | 2 +- 2 files changed, 4 insertions(+), 41 deletions(-) diff --git a/lalpulsar/src/SSBtimes.c b/lalpulsar/src/SSBtimes.c index 9865ca327c..2c39070be7 100644 --- a/lalpulsar/src/SSBtimes.c +++ b/lalpulsar/src/SSBtimes.c @@ -42,6 +42,7 @@ const UserChoices SSBprecisionChoices = { { SSBPREC_NEWTONIAN, "newtonian" }, { SSBPREC_RELATIVISTIC, "relativistic" }, { SSBPREC_RELATIVISTICOPT, "relativisticopt" }, + { SSBPREC_DMOFF, "DMoff"}, }; /*---------- internal prototypes ----------*/ @@ -639,50 +640,12 @@ XLALGetSSBtimes ( const DetectorStateSeries *DetectorStates, /**< [in] detector- case SSBPREC_DMOFF: /* switch off all demodulation terms */ - baryinput.site = DetectorStates->detector; - baryinput.site.location[0] /= LAL_C_SI; - baryinput.site.location[1] /= LAL_C_SI; - baryinput.site.location[2] /= LAL_C_SI; - - baryinput.alpha = alpha; - baryinput.delta = delta; - baryinput.dInv = 0; - for ( UINT4 i = 0; i < numSteps; i++ ) { - EmissionTime emit; DetectorState *state = &(DetectorStates->data[i]); - baryinput.tgps = state->tGPS; - - REAL8 tgps[2]; - tgps[0] = baryinput.tgps.gpsSeconds; - tgps[1] = baryinput.tgps.gpsNanoSeconds; - - if ( XLALBarycenterOpt ( &emit, &baryinput, &(state->earthState), &bBuffer ) != XLAL_SUCCESS ) - XLAL_ERROR_NULL ( XLAL_EFUNC, "XLALBarycenterOpt() failed with xlalErrno = %d\n", xlalErrno ); - - emit.deltaT = 0; - emit.tDot = 1; - INT4 deltaTint = floor ( emit.deltaT ); - - if ( ( 1e-9 * tgps[1] + emit.deltaT - deltaTint ) >= 1.e0 ) - { - emit.te.gpsSeconds = baryinput.tgps.gpsSeconds + deltaTint + 1; - emit.te.gpsNanoSeconds = floor ( 1e9 * ( tgps[1] * 1e-9 + emit.deltaT - deltaTint - 1.0 ) ); - } - else - { - emit.te.gpsSeconds = baryinput.tgps.gpsSeconds + deltaTint; - emit.te.gpsNanoSeconds = floor ( 1e9 * ( tgps[1] * 1e-9 + emit.deltaT - deltaTint ) ); - } - - ret->DeltaT->data[i] = XLALGPSGetREAL8 ( &emit.te ) - refTimeREAL8; - ret->Tdot->data[i] = emit.tDot; - + ret->DeltaT->data[i] = XLALGPSGetREAL8 ( &state->tGPS ) - refTimeREAL8; + ret->Tdot->data[i] = 1.0; } /* for i < numSteps */ - // free buffer memory - XLALFree ( bBuffer ); - break; default: diff --git a/lalpulsar/src/SSBtimes.h b/lalpulsar/src/SSBtimes.h index a3c68b32bc..0e49ebb14f 100644 --- a/lalpulsar/src/SSBtimes.h +++ b/lalpulsar/src/SSBtimes.h @@ -46,7 +46,7 @@ typedef enum tagSSBprecision { SSBPREC_NEWTONIAN, /**< simple Newtonian: \f$\tau = t + \vec{r}\cdot\vec{n}/c\f$ */ SSBPREC_RELATIVISTIC, /**< detailed relativistic: \f$\tau=\tau(t; \vec{n}, \vec{r})\f$ */ SSBPREC_RELATIVISTICOPT, /**< optimized relativistic, numerically equivalent to #SSBPREC_RELATIVISTIC, but faster */ - SSBPREC_DMOFF, /**< switch off all demodulatoin terms */ + SSBPREC_DMOFF, /**< switch off all demodulatoin terms */ SSBPREC_LAST /**< end marker */ } SSBprecision; -- GitLab From e47b4249178397027169f7c4cb3007ec94e1dd6d Mon Sep 17 00:00:00 2001 From: Heinz-Bernd Eggenstein Date: Mon, 29 Jan 2018 13:33:14 +0100 Subject: [PATCH 3/3] Add an option --SortToplist=7 to lalapps_HierarchSearchGCT - this option generates three toplists, sorted by 2F, BSGLtL, BtSGLtL respectively --- .../pulsar/EinsteinAtHome/hs_boinc_extras.c | 4 +- lalapps/src/pulsar/GCT/GCTtoplist.h | 1 + lalapps/src/pulsar/GCT/HierarchSearchGCT.c | 60 +++++++++++++---- lalapps/src/pulsar/GCT/testGCT.sh | 67 ++++++++++++++++++- 4 files changed, 113 insertions(+), 19 deletions(-) diff --git a/lalapps/src/pulsar/EinsteinAtHome/hs_boinc_extras.c b/lalapps/src/pulsar/EinsteinAtHome/hs_boinc_extras.c index 69a9d9ffc0..facda5bfa3 100644 --- a/lalapps/src/pulsar/EinsteinAtHome/hs_boinc_extras.c +++ b/lalapps/src/pulsar/EinsteinAtHome/hs_boinc_extras.c @@ -1143,8 +1143,8 @@ static int worker (void) { second_outfile = -1; } - /* record if there will be a second output file */ - else if (!strcmp("--SortToplist=6",argv[arg])) { + /* record if there will be triple output files */ + else if (!strcmp("--SortToplist=6",argv[arg]) || !strcmp("--SortToplist=7",argv[arg])) { rargv[rarg] = argv[arg]; LogPrintf(LOG_DEBUG,"BSGL output files\n"); bsgl_outfiles = -1; diff --git a/lalapps/src/pulsar/GCT/GCTtoplist.h b/lalapps/src/pulsar/GCT/GCTtoplist.h index 49b607a3f9..7fce454857 100644 --- a/lalapps/src/pulsar/GCT/GCTtoplist.h +++ b/lalapps/src/pulsar/GCT/GCTtoplist.h @@ -77,6 +77,7 @@ typedef enum SORTBY_BSGLtL = 4, //< sort by transient-line robust statistic B_S/GLtL SORTBY_BtSGLtL = 5, //< sort by transient-CW line robust statistic B_tS/GLtL SORTBY_TRIPLE_BStSGLtL = 6, //< triple toplists: one sorted by B_S/GL, one by B_S/GLtL, one by B_tS/GLtL + SORTBY_F_BSGLtL_BtSGLtL = 7,//< triple toplists: one sorted by 2F, one by B_S/GLtL, one by B_tS/GLtL SORTBY_LAST //< end-marker } SortBy_t; diff --git a/lalapps/src/pulsar/GCT/HierarchSearchGCT.c b/lalapps/src/pulsar/GCT/HierarchSearchGCT.c index 8f4aaa7bf1..c14a42c1d4 100644 --- a/lalapps/src/pulsar/GCT/HierarchSearchGCT.c +++ b/lalapps/src/pulsar/GCT/HierarchSearchGCT.c @@ -213,6 +213,7 @@ void PrintStackInfo( LALStatus *status, const SFTCatalogSequence *catalogSeq, FI void UpdateSemiCohToplists ( LALStatus *status, toplist_t *list1, toplist_t *list2, toplist_t *list3, FineGrid *in, REAL8 f1dot_fg, REAL8 f2dot_fg, REAL8 f3dot_fg, UsefulStageVariables *usefulparams, REAL4 NSegmentsInv, REAL4 *NSegmentsInvX, BOOLEAN have_f3dot ); void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, + SortBy_t toplist_sortby, toplist_t *list1, toplist_t *list2, toplist_t *list3, @@ -350,8 +351,8 @@ int MAIN( int argc, char *argv[]) { /* fstat candidate structure for candidate toplist*/ toplist_t *semiCohToplist=NULL; - toplist_t *semiCohToplist2=NULL; // only used for SORTBY_DUAL_F_BSGL or SORTBY_TRIPLE_BStSGLtL - toplist_t *semiCohToplist3=NULL; // only used for SORTBY_TRIPLE_BStSGLtL + toplist_t *semiCohToplist2=NULL; // only used for SORTBY_DUAL_F_BSGL, SORTBY_TRIPLE_BStSGLtL or SORTBY_F_BSGLtL_BtSGLtL + toplist_t *semiCohToplist3=NULL; // only used for SORTBY_TRIPLE_BStSGLtL and SORTBY_F_BSGLtL_BtSGLtL /* template and grid variables */ static DopplerSkyScanInit scanInit; /* init-structure for DopperScanner */ @@ -549,7 +550,7 @@ int MAIN( int argc, char *argv[]) { XLAL_CHECK_MAIN( XLALRegisterNamedUvar( &uvar_oLGX, "oLGX", STRINGVector, 0, OPTIONAL, "BSGL: prior per-detector line-vs-Gauss odds 'oLGX' (Defaults to oLGX=1/Ndet)") == XLAL_SUCCESS, XLAL_EFUNC); XLAL_CHECK_MAIN( XLALRegisterNamedUvar( &uvar_BSGLlogcorr, "BSGLlogcorr", BOOLEAN, 0, DEVELOPER, "BSGL: include log-correction terms (slower) or not (faster)") == XLAL_SUCCESS, XLAL_EFUNC); XLAL_CHECK_MAIN( XLALRegisterNamedUvar( &uvar_getMaxFperSeg, "getMaxFperSeg", BOOLEAN, 0, OPTIONAL, "Compute and output maximum F and FX over segments") == XLAL_SUCCESS, XLAL_EFUNC); - XLAL_CHECK_MAIN( XLALRegisterNamedUvar( &uvar_SortToplist, "SortToplist", INT4, 0, OPTIONAL, "Sort toplist by: 0=Fstat, 1=nc, 2=B_S/GL, 3='Fstat + B_S/GL', 4=B_S/GLtL, 5=B_tS/GLtL, 6='B_S/GL + B_S/GLtL + B_tS/GLtL'") == XLAL_SUCCESS, XLAL_EFUNC); + XLAL_CHECK_MAIN( XLALRegisterNamedUvar( &uvar_SortToplist, "SortToplist", INT4, 0, OPTIONAL, "Sort toplist by: 0=Fstat, 1=nc, 2=B_S/GL, 3='Fstat + B_S/GL', 4=B_S/GLtL, 5=B_tS/GLtL, 6='B_S/GL + B_S/GLtL + B_tS/GLtL', 7='Fstat + B_S/GLtL + B_tS/GLtL' ") == XLAL_SUCCESS, XLAL_EFUNC); // -------------------------------------------- XLAL_CHECK_MAIN( XLALRegisterNamedUvarAuxData( &uvar_FstatMethod, "FstatMethod", UserEnum, XLALFstatMethodChoices(), 0, OPTIONAL, "F-statistic method to use" ) == XLAL_SUCCESS, XLAL_EFUNC); @@ -637,12 +638,14 @@ int MAIN( int argc, char *argv[]) { } if ( (uvar_SortToplist == SORTBY_BSGL || uvar_SortToplist == SORTBY_DUAL_F_BSGL || uvar_SortToplist == SORTBY_BSGLtL || uvar_SortToplist == SORTBY_BtSGLtL || - uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL ) && !uvar_computeBSGL ) { - fprintf(stderr, "Toplist sorting by BSGL only possible if --computeBSGL given.\n"); + uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL || + uvar_SortToplist == SORTBY_F_BSGLtL_BtSGLtL) && !uvar_computeBSGL ) { + fprintf(stderr, "Toplist sorting by BSGL[tL] only possible if --computeBSGL given.\n"); return( HIERARCHICALSEARCH_EBAD ); } if ( ( uvar_SortToplist == SORTBY_BSGLtL || uvar_SortToplist == SORTBY_BtSGLtL || - uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL ) && !uvar_getMaxFperSeg ) { + uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL || + uvar_SortToplist == SORTBY_F_BSGLtL_BtSGLtL ) && !uvar_getMaxFperSeg ) { fprintf(stderr, "Toplist sorting by B[t]SGLtL only possible if --getMaxFperSeg given.\n"); return( HIERARCHICALSEARCH_EBAD ); } @@ -665,6 +668,15 @@ int MAIN( int argc, char *argv[]) { XLAL_CHECK ( 0 == create_gctFstat_toplist ( &semiCohToplist3, uvar_nCand1, SORTBY_BtSGLtL ), XLAL_EFUNC, "create_gctFstat_toplist() failed for nCand=%d and sortBy=%d\n", uvar_nCand1, SORTBY_BtSGLtL ); } + else if ( uvar_SortToplist == SORTBY_F_BSGLtL_BtSGLtL )// special treatement of 'triple' toplists: 1st one sorted by 'F', 2nd one by 'B_S/GLtL', 3rd by 'B_tS/GLtL' + { + XLAL_CHECK ( 0 == create_gctFstat_toplist ( &semiCohToplist, uvar_nCand1, SORTBY_F ), + XLAL_EFUNC, "create_gctFstat_toplist() failed for nCand=%d and sortBy=%d\n", uvar_nCand1, SORTBY_F ); + XLAL_CHECK ( 0 == create_gctFstat_toplist ( &semiCohToplist2, uvar_nCand1, SORTBY_BSGLtL ), + XLAL_EFUNC, "create_gctFstat_toplist() failed for nCand=%d and sortBy=%d\n", uvar_nCand1, SORTBY_BSGLtL ); + XLAL_CHECK ( 0 == create_gctFstat_toplist ( &semiCohToplist3, uvar_nCand1, SORTBY_BtSGLtL ), + XLAL_EFUNC, "create_gctFstat_toplist() failed for nCand=%d and sortBy=%d\n", uvar_nCand1, SORTBY_BtSGLtL ); + } else // 'normal' single-sorting toplist cases (sortby 'F', 'nc' or 'BSGL') { XLAL_CHECK ( 0 == create_gctFstat_toplist ( &semiCohToplist, uvar_nCand1, uvar_SortToplist), @@ -1761,8 +1773,9 @@ int MAIN( int argc, char *argv[]) { /* this is necessary here, because UpdateSemiCohToplists() might set a checkpoint that needs some information from here */ - if(uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL && finegrid.sumTwoFX && finegrid.maxTwoFXl ) { - LAL_CALL( UpdateSemiCohToplistsOptimTriple (&status, semiCohToplist, semiCohToplist2, semiCohToplist3, &finegrid, f1dot_fg, f2dot_fg, f3dot_fg, &usefulParams, NSegmentsInv, usefulParams.NSegmentsInvX, XLALUserVarWasSet(&uvar_f3dot) ), &status); + if(( uvar_SortToplist == SORTBY_TRIPLE_BStSGLtL || uvar_SortToplist == SORTBY_F_BSGLtL_BtSGLtL ) && + finegrid.sumTwoFX && finegrid.maxTwoFXl ) { + LAL_CALL( UpdateSemiCohToplistsOptimTriple (&status, uvar_SortToplist, semiCohToplist, semiCohToplist2, semiCohToplist3, &finegrid, f1dot_fg, f2dot_fg, f3dot_fg, &usefulParams, NSegmentsInv, usefulParams.NSegmentsInvX, XLALUserVarWasSet(&uvar_f3dot) ), &status); } else { LAL_CALL( UpdateSemiCohToplists (&status, semiCohToplist, semiCohToplist2, semiCohToplist3, &finegrid, f1dot_fg, f2dot_fg, f3dot_fg, &usefulParams, NSegmentsInv, usefulParams.NSegmentsInvX, XLALUserVarWasSet(&uvar_f3dot) ), &status); } @@ -1941,7 +1954,8 @@ int MAIN( int argc, char *argv[]) { strcpy(t3_suffix,""); break; }; - case SORTBY_TRIPLE_BStSGLtL : { + case SORTBY_TRIPLE_BStSGLtL: + case SORTBY_F_BSGLtL_BtSGLtL : { strcpy(t1_suffix,""); strcpy(t2_suffix,"-BSGLtL"); strcpy(t3_suffix,"-BtSGLtL"); @@ -2674,6 +2688,7 @@ void PrintStackInfo( LALStatus *status, * This function allows for inserting candidates into up to 3 toplists at once, which might be sorted differently! */ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, + SortBy_t toplists_sortby, toplist_t *list1, toplist_t *list2, //< optional (can be NULL): insert candidate into this 2nd toplist as well toplist_t *list3, //< optional (can be NULL): insert candidate into this 3rd toplist as well @@ -2691,6 +2706,7 @@ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, REAL8 freq_fg; UINT4 ifreq_fg; GCTtopOutputEntry line; + BOOLEAN delay_compute_BSGL = (toplists_sortby == SORTBY_F_BSGLtL_BtSGLtL); INITSTATUS(status); ATTATCHSTATUSPTR (status); @@ -2699,9 +2715,11 @@ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, ASSERT ( in != NULL, status, HIERARCHICALSEARCH_ENULL, HIERARCHICALSEARCH_MSGENULL ); ASSERT ( usefulparams != NULL, status, HIERARCHICALSEARCH_ENULL, HIERARCHICALSEARCH_MSGENULL ); - /* Optimized version for triple toplist case: + /* Optimized version for triple toplist cases: BSGL, BSGLtL,BtSGLtL or + 2F, BSGLtL,BtSGLtL - First compute all three detection metrics by which any of the toplists is sorted. + First compute all detection metrics by which any of the toplists is sorted (even as seconary + comparison criterion). Next test if the candidate is loud enough to make it into at least one of the toplists. If so, build the candidate structure. @@ -2730,7 +2748,15 @@ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, } xlalErrno = 0; - line.log10BSGL = XLALComputeBSGL ( sumTwoF, sumTwoFX, usefulparams->BSGLsetup ); + /* Note: sorting a toplist by 2F takes BSGL as a secondary sorting criterion, so even if the 1st toplist is + sorted by 2F, we still provide an overestimate of BSGL for a test if it can make it into the toplist in + the first place */ + + if(delay_compute_BSGL) { + line.log10BSGL = LAL_REAL4_MAX; + } else { + line.log10BSGL = XLALComputeBSGL ( sumTwoF, sumTwoFX, usefulparams->BSGLsetup ); + } if ( xlalErrno != 0 ) { XLALPrintError ("%s line %d : XLALComputeBSGL() failed with xlalErrno = %d.\n\n", __func__, __LINE__, xlalErrno ); ABORT ( status, HIERARCHICALSEARCH_EXLAL, HIERARCHICALSEARCH_MSGEXLAL ); @@ -2769,8 +2795,8 @@ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, /* take F-stat averages over segments */ line.avTwoF = sumTwoF*NSegmentsInv; /* average multi-2F by full number of segments */ -/* now test if this candidate makes it into any of the toplists, if not we don't need - to do any more copying of data from finegrid to toplist entry structure */ + /* now test if this candidate makes it into any of the toplists, if not we don't need + to do any more copying and computing of data from finegrid to toplist entry structure */ int isIncludedToplists = TEST_FSTAT_TOPLIST_INCLUSION( list1, &line) || @@ -2780,6 +2806,12 @@ void UpdateSemiCohToplistsOptimTriple ( LALStatus *status, if(likely(! isIncludedToplists)) continue; + + /* if we delayed the computation of BSGL for a 2F sorted first toplist, now is the time to actually compute it */ + if(delay_compute_BSGL) { + line.log10BSGL = XLALComputeBSGL ( sumTwoF, sumTwoFX, usefulparams->BSGLsetup ); + } + line.numDetectors = in->numDetectors; line.avTwoFrecalc = -1.0; /* initialise this to -1.0, so that it only gets written out by print_gctFstatline_to_str if later overwritten in recalcToplistStats step */ line.log10BSGLrecalc = -LAL_REAL4_MAX; /* for now, block field with minimal value, needed for output checking in print_gctFstatline_to_str() */ diff --git a/lalapps/src/pulsar/GCT/testGCT.sh b/lalapps/src/pulsar/GCT/testGCT.sh index b474e2bf83..a654d90490 100755 --- a/lalapps/src/pulsar/GCT/testGCT.sh +++ b/lalapps/src/pulsar/GCT/testGCT.sh @@ -389,8 +389,8 @@ fi echo echo "----------------------------------------------------------------------------------------------------" -echo " STEP 7: run HierarchSearchGCT using Resampling and triple toplist (BSGL, BSGLtL and BtSGLtL), " -echo " compared with separate runs for each of these 3 toplist rankings" +echo " STEP 7: run HierarchSearchGCT using Resampling and triple toplists (BSGL, BSGLtL and BtSGLtL), " +echo " (2F, BSGLtL and BtSGLtL), compared with separate runs for each of these 3 toplist rankings" echo "----------------------------------------------------------------------------------------------------" echo @@ -398,6 +398,8 @@ rm -f checkpoint.cpt # delete checkpoint to start correctly outfile_GCT_RS_triple="${testDir}/GCT_RS_triple.dat" timingsfile_RS_triple="${testDir}/timing_RS_triple.dat" +outfile_GCT_RS_triple2="${testDir}/GCT_RS_triple2.dat" +timingsfile_RS_triple2="${testDir}/timing_RS_triple2.dat" # set plan mode so that results should be deterministic, easier to compare results of different runs this way @@ -419,7 +421,42 @@ if ! eval "$cmdline"; then fi -## re-run, but now create the three toplists one by one, then compare results +# second variant of triple toplists, with 2F sorted first toplist + +cmdline="$gct_code $gct_CL_common --FstatMethod=ResampGeneric --fnameout='$outfile_GCT_RS_triple2' --outputTiming='$timingsfile_RS_triple2' ${BSGL_flags} --getMaxFperSeg --loudestSegOutput --SortToplist=7" +if [ -n "$DEBUG" ]; then + cmdline="$cmdline" +else + cmdline="$cmdline &> /dev/null" +fi + +echo $cmdline +if ! eval "$cmdline"; then + echo "Error.. something failed when running '$gct_code' ..." + exit 1 +fi + + + + +## re-run, but now create the in total four toplists one by one, then compare results + +outfile_GCT_RS_triple="${testDir}/GCT_RS_triple_0.dat" + +cmdline="$gct_code $gct_CL_common --FstatMethod=ResampGeneric --fnameout='$outfile_GCT_RS_triple' ${BSGL_flags} --getMaxFperSeg --loudestSegOutput --SortToplist=0" +if [ -n "$DEBUG" ]; then + cmdline="$cmdline" +else + cmdline="$cmdline &> /dev/null" +fi + +echo $cmdline +if ! eval "$cmdline"; then + echo "Error.. something failed when running '$gct_code' ..." + exit 1 +fi + + outfile_GCT_RS_triple="${testDir}/GCT_RS_triple_1.dat" @@ -469,6 +506,7 @@ fi # filter out comments +egrep -v "^%" ${testDir}/GCT_RS_triple_0.dat > ${testDir}/GCT_RS_triple_0.txt egrep -v "^%" ${testDir}/GCT_RS_triple_1.dat > ${testDir}/GCT_RS_triple_1.txt egrep -v "^%" ${testDir}/GCT_RS_triple_2.dat > ${testDir}/GCT_RS_triple_2.txt egrep -v "^%" ${testDir}/GCT_RS_triple_3.dat > ${testDir}/GCT_RS_triple_3.txt @@ -478,6 +516,11 @@ egrep -v "^%" ${testDir}/GCT_RS_triple.dat > ${testDir}/GCT_RS_triple.txt egrep -v "^%" ${testDir}/GCT_RS_triple.dat-BSGLtL > ${testDir}/GCT_RS_triple-BSGLtL.txt egrep -v "^%" ${testDir}/GCT_RS_triple.dat-BtSGLtL > ${testDir}/GCT_RS_triple-BtSGLtL.txt +egrep -v "^%" ${testDir}/GCT_RS_triple2.dat > ${testDir}/GCT_RS_triple2.txt +egrep -v "^%" ${testDir}/GCT_RS_triple2.dat-BSGLtL > ${testDir}/GCT_RS_triple2-BSGLtL.txt +egrep -v "^%" ${testDir}/GCT_RS_triple2.dat-BtSGLtL > ${testDir}/GCT_RS_triple2-BtSGLtL.txt + + if ! eval "diff ${testDir}/GCT_RS_triple_1.txt ${testDir}/GCT_RS_triple.txt"; then echo "Error: tripple toplists do not match separately generated toplists (1) " exit 1 @@ -493,6 +536,24 @@ if ! eval "diff ${testDir}/GCT_RS_triple_3.txt ${testDir}/GCT_RS_triple-BtSGLtL. exit 1 fi +if ! eval "diff ${testDir}/GCT_RS_triple_0.txt ${testDir}/GCT_RS_triple2.txt"; then + echo "Error: tripple toplists do not match separately generated toplists (4) " + exit 1 +fi + +if ! eval "diff ${testDir}/GCT_RS_triple_2.txt ${testDir}/GCT_RS_triple2-BSGLtL.txt"; then + echo "Error: tripple toplists do not match separately generated toplists (5) " + exit 1 +fi + +if ! eval "diff ${testDir}/GCT_RS_triple_3.txt ${testDir}/GCT_RS_triple2-BtSGLtL.txt"; then + echo "Error: tripple toplists do not match separately generated toplists (6) " + exit 1 +fi + + + + ## ---------- compute relative differences and check against tolerance -------------------- -- GitLab