From fcd837b30670ea515226161dcd658243dddc5194 Mon Sep 17 00:00:00 2001 From: kevScheuer Date: Fri, 18 Sep 2026 08:19:49 -0400 Subject: [PATCH 1/3] Implement bootstrap sampling of signal data --- AmpTools/IUAmpTools/AmpToolsInterface.cc | 16 +++++++-- AmpTools/IUAmpTools/AmpToolsInterface.h | 13 +++++-- AmpTools/IUAmpTools/DataReader.h | 18 ++++++++++ AmpTools/IUAmpTools/LikelihoodCalculator.cc | 16 ++++++++- AmpTools/IUAmpTools/LikelihoodCalculator.h | 2 ++ AmpTools/IUAmpToolsMPI/DataReaderMPI.h | 34 +++++++++++++++++-- .../IUAmpToolsMPI/LikelihoodCalculatorMPI.cc | 16 +++++++++ .../IUAmpToolsMPI/LikelihoodCalculatorMPI.h | 7 ++++ .../IUAmpToolsMPI/LikelihoodManagerMPI.cc | 5 +++ AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h | 1 + 10 files changed, 121 insertions(+), 7 deletions(-) diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.cc b/AmpTools/IUAmpTools/AmpToolsInterface.cc index d4cb3c49..35be34f7 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.cc +++ b/AmpTools/IUAmpTools/AmpToolsInterface.cc @@ -170,7 +170,7 @@ AmpToolsInterface::resetConfigurationInfo(ConfigurationInfo* configurationInfo){ for (unsigned int i = 0; i < m_userDataReaders.size(); i++){ if (reaction->data().first == m_userDataReaders[i]->name()) m_dataReaderMap[reactionName] - = m_userDataReaders[i]->newDataReader(reaction->data().second); + = m_userDataReaders[i]->newDataReader(reaction->data().second); if (reaction->bkgnd().first == m_userDataReaders[i]->name()) m_bkgndReaderMap[reactionName] = m_userDataReaders[i]->newDataReader(reaction->bkgnd().second); @@ -199,7 +199,7 @@ AmpToolsInterface::resetConfigurationInfo(ConfigurationInfo* configurationInfo){ report( WARNING, kModule ) << "not creating a DataReader for accepted MC associated with reaction " << reactionName << endl; if( m_functionality == kFull ){ - + // ************************ // create a NormIntInterface // ************************ @@ -400,6 +400,18 @@ AmpToolsInterface::randomizeParameter( const string& parName, float min, float m invalidateAmps(); } +void +AmpToolsInterface::bootstrapSignalData(const string& reactionName){ + + LikelihoodCalculator* likCalc = likelihoodCalculator(reactionName); + if( likCalc == NULL ){ + report( ERROR, kModule ) << "no LikelihoodCalculator for reaction: " << reactionName << endl; + return; + } + likCalc->bootstrapSignalData(); + +} + void AmpToolsInterface::finalizeFit( const string& tag ){ diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.h b/AmpTools/IUAmpTools/AmpToolsInterface.h index 24c7b8af..00ac21fc 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.h +++ b/AmpTools/IUAmpTools/AmpToolsInterface.h @@ -259,10 +259,19 @@ class AmpToolsInterface{ */ void randomizeParameter( const string& parName, float min = 0, float max = 1 ); + + /** This function will call the dataReader.resample() method intended to randomly + * sample, with replacement, the signal events when getEvent() is called. This is + * particularly useful for fits that contain free parameters in the amplitudes + * themselves, as the bootstrap distributions will provide better uncertainty + * estimates. + */ + + void bootstrapSignalData(const string& reactionName); /** Print final fit results to a file. The tag can be used to * generate a unique name in the case that multiple results are - * written for a singele fit job. + * written for a single fit job. */ virtual void finalizeFit( const string& tag = "" ); @@ -270,7 +279,7 @@ class AmpToolsInterface{ /** For manual calculations: clear all events and calculations. * Call this before loading events from a new reaction or to start - * new calcuations. + * new calculations. * * \param[in] iDataSet used to index simultaneous manual calculations * diff --git a/AmpTools/IUAmpTools/DataReader.h b/AmpTools/IUAmpTools/DataReader.h index 19000fce..d951692a 100644 --- a/AmpTools/IUAmpTools/DataReader.h +++ b/AmpTools/IUAmpTools/DataReader.h @@ -105,6 +105,24 @@ class DataReader * \see DataReaderMPI */ virtual void resetSource() = 0; + + /** + * Reserved for DataReaders that support redrawing samples in place (e.g. + * bootstrap resampling with replacement) without reconstructing the object. + * The default implementation reports that no resampling is supported. + * Overridden by a reader that does. + * + * This method should be virtual in the user's class if the DataReaderMPI + * template is to be used. + * + * + * \see DataReaderMPI + */ + virtual void resample(){ + report( ERROR, kModule ) << name() << " does not implement resample(). " + << "Cannot bootstrap this data source." << endl; + assert( false ); + } /** * The user should override this function with one that returns the number diff --git a/AmpTools/IUAmpTools/LikelihoodCalculator.cc b/AmpTools/IUAmpTools/LikelihoodCalculator.cc index af659398..4b2132d2 100644 --- a/AmpTools/IUAmpTools/LikelihoodCalculator.cc +++ b/AmpTools/IUAmpTools/LikelihoodCalculator.cc @@ -230,7 +230,7 @@ SCOREP_USER_REGION_BEGIN( dataTerm, "dataTerm", SCOREP_USER_REGION_TYPE_COMMON ) m_numDataEvents = m_ampVecsSignal.m_iNTrueEvents; m_sumDataWeights = m_ampVecsSignal.m_dSumWeights; - + if( m_ampVecsSignal.m_hasNonUnityWeights && m_hasBackground ){ report( WARNING, kModule ) << "\n" @@ -305,6 +305,20 @@ SCOREP_USER_REGION_END( dataTerm ) return sumLnI; } +void +LikelihoodCalculator::bootstrapSignalData(){ +#ifdef SCOREP +SCOREP_USER_REGION_DEFINE( bootstrapSignalData ) +SCOREP_USER_REGION_BEGIN( bootstrapSignalData, "bootstrapSignalData", SCOREP_USER_REGION_TYPE_COMMON ) +#endif + m_ampVecsSignal.deallocAmpVecs(); + m_dataReaderSignal->resample(); + m_firstDataCalc = true; + #ifdef SCOREP + SCOREP_USER_REGION_END( bootstrapSignalData ) + #endif +} + void LikelihoodCalculator::invalidateTerms(){ diff --git a/AmpTools/IUAmpTools/LikelihoodCalculator.h b/AmpTools/IUAmpTools/LikelihoodCalculator.h index 62df8c0d..e5c9351a 100644 --- a/AmpTools/IUAmpTools/LikelihoodCalculator.h +++ b/AmpTools/IUAmpTools/LikelihoodCalculator.h @@ -87,6 +87,8 @@ class LikelihoodCalculator : public MIFunctionContribution virtual double numSignalEvents(); void invalidateTerms(); + + virtual void bootstrapSignalData(); protected: diff --git a/AmpTools/IUAmpToolsMPI/DataReaderMPI.h b/AmpTools/IUAmpToolsMPI/DataReaderMPI.h index a98b88e5..9cf1c1a6 100644 --- a/AmpTools/IUAmpToolsMPI/DataReaderMPI.h +++ b/AmpTools/IUAmpToolsMPI/DataReaderMPI.h @@ -31,8 +31,8 @@ struct KinStruct { * of the instances of this class on the follower nodes will behave as if * they have only a subset of the data. The instance of it on the leader * node will behave as if it has all of the data. Be sure that the getEvent, - * resetSource, and numEvents methods in the user-defined class are declared - * virtual. + * resetSource, numEvents, and (if used for bootstrapping) resample methods in + * the user-defined class are declared virtual. * * \ingroup IUAmpToolsMPI */ @@ -66,6 +66,7 @@ class DataReaderMPI : public T Kinematics* getEvent(); void resetSource(); + void resample(); unsigned int numEvents() const; @@ -229,6 +230,35 @@ void DataReaderMPI::resetSource() } } +template< class T > +void DataReaderMPI::resample() +{ + if( m_isLeader ){ + + report( DEBUG, kDRModule ) << "Resampling leader data source " << m_rank << endl; + + // redraw the bootstrap sample using the user's definition in the reader, then + // push the selection out to the follower's exactly as done in startup + T::resample(); + distributeData(); + } + else{ + + report( DEBUG, kDRModule ) << "Receiving resampled data on process with rank " << m_rank << endl; + + //discard the cached partition + for( vector::iterator ptrItr = m_ptrCache.begin(); + ptrItr != m_ptrCache.end(); + ++ptrItr ){ + + delete *ptrItr; + } + m_ptrCache.clear(); + + receiveData(); + } +} + template< class T > void DataReaderMPI::provideData() { diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc index cfcacf8d..e28b7ec8 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc +++ b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc @@ -215,6 +215,22 @@ LikelihoodCalculatorMPI::operator()() return -2 * lnL; } +void +LikelihoodCalculatorMPI::bootstrapSignalData(){ + + if( m_isLeader ){ + + int cmnd[2]; + cmnd[0] = m_thisId; + cmnd[1] = LikelihoodManagerMPI::kBootstrapSignalData; + MPI_Bcast( cmnd, MPI_INT, 0, MPI_COMM_WORLD); + } + + // calls resample() on DataReaderMPI for the leader, and redistributes to the + // followers, discarding the old cached partition + LikelihoodCalculator::bootstrapSignalData(); +} + double LikelihoodCalculatorMPI::numSignalEvents(){ diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h index eb88cd6b..1a22f50f 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h +++ b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h @@ -123,6 +123,13 @@ class LikelihoodCalculatorMPI : public LikelihoodCalculator */ double operator()(); + /** + * Triggers a bootstrap resampling of the signal data. On a leader it + * broadcasts the command to all followers and resamples the leader's own + * DataReaderMPI, pushing newly-drawn events to the followers. + */ + void bootstrapSignalData(); + /** * This sends the finalize fit flag to all follower jobs which breaks them * out of the deliverLikelihood method of the LikelihoodManagerMPI. diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc index 3582a22e..c8fb5c4b 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc +++ b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc @@ -106,6 +106,11 @@ LikelihoodManagerMPI::deliverLikelihood() likCalc->computeLikelihood(); break; + + case kBootstrapSignalData: + + likCalc->bootstrapSignalData(); + break; default: diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h index 5d1c79ed..eba6e7f7 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h +++ b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h @@ -60,6 +60,7 @@ class LikelihoodManagerMPI kComputeIntegrals, kUpdateParameters, kUpdateAmpParameter, + kBootstrapSignalData, kFinalizeFit, kExit }; From bf9444773c151cc00eab4b0e8c47ae667243cc74 Mon Sep 17 00:00:00 2001 From: kevScheuer Date: Mon, 21 Sep 2026 07:23:20 -0400 Subject: [PATCH 2/3] Add resample ability for bkgnd, and optional seed Background events can also be resampled with replacement via similar methods. An optional seed integer can be given to the dataReader, so that bootstrap samples can be reproduced. Other fixes: - Data reader report fixed for resample method - LikelihoodCalculator::dataTerm now checks for signal and background recalculations separately --- AmpTools/IUAmpTools/AmpToolsInterface.cc | 15 +++- AmpTools/IUAmpTools/AmpToolsInterface.h | 19 ++-- AmpTools/IUAmpTools/DataReader.cc | 1 + AmpTools/IUAmpTools/DataReader.h | 7 +- AmpTools/IUAmpTools/LikelihoodCalculator.cc | 90 +++++++++++-------- AmpTools/IUAmpTools/LikelihoodCalculator.h | 4 +- AmpTools/IUAmpToolsMPI/DataReaderMPI.h | 6 +- .../IUAmpToolsMPI/LikelihoodCalculatorMPI.cc | 30 ++++++- .../IUAmpToolsMPI/LikelihoodCalculatorMPI.h | 9 +- .../IUAmpToolsMPI/LikelihoodManagerMPI.cc | 12 ++- AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h | 1 + 11 files changed, 138 insertions(+), 56 deletions(-) diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.cc b/AmpTools/IUAmpTools/AmpToolsInterface.cc index 35be34f7..1348eacf 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.cc +++ b/AmpTools/IUAmpTools/AmpToolsInterface.cc @@ -401,17 +401,28 @@ AmpToolsInterface::randomizeParameter( const string& parName, float min, float m } void -AmpToolsInterface::bootstrapSignalData(const string& reactionName){ +AmpToolsInterface::bootstrapSignalData(const string& reactionName, unsigned int seed){ LikelihoodCalculator* likCalc = likelihoodCalculator(reactionName); if( likCalc == NULL ){ report( ERROR, kModule ) << "no LikelihoodCalculator for reaction: " << reactionName << endl; return; } - likCalc->bootstrapSignalData(); + likCalc->bootstrapSignalData( seed ); } +void +AmpToolsInterface::bootstrapBackgroundData(const string& reactionName, unsigned int seed){ + + LikelihoodCalculator* likCalc = likelihoodCalculator(reactionName); + if( likCalc == NULL ){ + report( ERROR, kModule ) << "no LikelihoodCalculator for reaction: " << reactionName << endl; + return; + } + likCalc->bootstrapBackgroundData( seed ); +} + void AmpToolsInterface::finalizeFit( const string& tag ){ diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.h b/AmpTools/IUAmpTools/AmpToolsInterface.h index 00ac21fc..b2352aa9 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.h +++ b/AmpTools/IUAmpTools/AmpToolsInterface.h @@ -260,14 +260,21 @@ class AmpToolsInterface{ void randomizeParameter( const string& parName, float min = 0, float max = 1 ); - /** This function will call the dataReader.resample() method intended to randomly - * sample, with replacement, the signal events when getEvent() is called. This is - * particularly useful for fits that contain free parameters in the amplitudes - * themselves, as the bootstrap distributions will provide better uncertainty - * estimates. + /** This function will call the dataReader.resample() method intended to + * randomly sample, with replacement, the signal events before getEvent() + * is called. This is particularly useful for fits that contain free + * parameters in the amplitudes themselves, as the bootstrap distributions + * will provide better uncertainty estimates. If no seed value is given, a + * random one is chosen. */ - void bootstrapSignalData(const string& reactionName); + void bootstrapSignalData(const string& reactionName, unsigned int seed = static_cast( time(NULL) ) ); + + /** This function will call the dataReader.resample() method intended to + * randomly sample, with replacement, the background events before getEvent() + * is called. If no seed value is given, a random one is chosen. + */ + void bootstrapBackgroundData(const string& reactionName, unsigned int seed = static_cast( time(NULL) ) ); /** Print final fit results to a file. The tag can be used to * generate a unique name in the case that multiple results are diff --git a/AmpTools/IUAmpTools/DataReader.cc b/AmpTools/IUAmpTools/DataReader.cc index 0e0a2405..8b3a44ab 100644 --- a/AmpTools/IUAmpTools/DataReader.cc +++ b/AmpTools/IUAmpTools/DataReader.cc @@ -39,6 +39,7 @@ #include #include "IUAmpTools/DataReader.h" +const char* DataReader::kModule = "DataReader"; using namespace std; diff --git a/AmpTools/IUAmpTools/DataReader.h b/AmpTools/IUAmpTools/DataReader.h index d951692a..72be429f 100644 --- a/AmpTools/IUAmpTools/DataReader.h +++ b/AmpTools/IUAmpTools/DataReader.h @@ -39,6 +39,9 @@ #include #include +#include + +#include "IUAmpTools/report.h" using namespace std; @@ -118,7 +121,7 @@ class DataReader * * \see DataReaderMPI */ - virtual void resample(){ + virtual void resample( unsigned int seed ){ report( ERROR, kModule ) << name() << " does not implement resample(). " << "Cannot bootstrap this data source." << endl; assert( false ); @@ -184,7 +187,7 @@ class DataReader vector m_args; - + static const char* kModule; }; #endif diff --git a/AmpTools/IUAmpTools/LikelihoodCalculator.cc b/AmpTools/IUAmpTools/LikelihoodCalculator.cc index 4b2132d2..2df5c4d5 100644 --- a/AmpTools/IUAmpTools/LikelihoodCalculator.cc +++ b/AmpTools/IUAmpTools/LikelihoodCalculator.cc @@ -67,6 +67,7 @@ m_normInt( normInt ), m_dataReaderSignal( dataReaderSignal ), m_dataReaderBkgnd( dataReaderBkgnd ), m_firstDataCalc( true ), +m_firstBkgndCalc( true ), m_firstNormIntCalc( true ), m_sumBkgWeights( 0 ), m_numBkgEvents( 0 ), @@ -106,7 +107,7 @@ LikelihoodCalculator::operator()(){ double LikelihoodCalculator::numSignalEvents(){ - if( m_firstDataCalc ) dataTerm(); + if( m_firstDataCalc || m_firstBkgndCalc ) dataTerm(); return numDataEvents() - sumBkgWeights(); } @@ -244,44 +245,45 @@ SCOREP_USER_REGION_BEGIN( dataTerm, "dataTerm", SCOREP_USER_REGION_TYPE_COMMON ) << "* contribution to the likelihood. *\n" << "****************************************************************\n" << endl; } + } - if( m_hasBackground ){ + if( m_firstBkgndCalc && m_hasBackground ){ - m_ampVecsBkgnd.loadData( m_dataReaderBkgnd, m_intenManager.needsUserVarsOnly() ); - m_ampVecsBkgnd.allocateTerms( m_intenManager, true ); + m_ampVecsBkgnd.loadData( m_dataReaderBkgnd, m_intenManager.needsUserVarsOnly() ); + m_ampVecsBkgnd.allocateTerms( m_intenManager, true ); - if( m_ampVecsBkgnd.m_hasMixedSignWeights ){ - report( NOTICE, kModule ) << "***************************************************************" << endl; - report( NOTICE, kModule ) << "* NOTICE: Weights with both positive and negative signs were *" << endl; - report( NOTICE, kModule ) << "* detected in the background file. This may be desirable for *" << endl; - report( NOTICE, kModule ) << "* some applications. Older versions of AmpTools (v0.10.x and *" << endl; - report( NOTICE, kModule ) << "* prior) will not properly handle this case and will also not *" << endl; - report( NOTICE, kModule ) << "* print this notification to the screen. *" << endl; - report( NOTICE, kModule ) << "***************************************************************" << endl; - } - - m_sumBkgWeights = m_ampVecsBkgnd.m_dSumWeights; - - // the extra boolean allows MPI jobs to suppress this check which - // may fail on one of the follower nodes if the background sample - // is sparse - if( m_sumBkgWeights < 0 && !suppressError ){ - report( ERROR, kModule ) << "****************************************************************" << endl; - report( ERROR, kModule ) << "* ERROR: The sum of all background weights is negative. This *" << endl; - report( ERROR, kModule ) << "* implies a negative background in the signal region, which *" << endl; - report( ERROR, kModule ) << "* unphysical. The weighted sum of the background events *" << endl; - report( ERROR, kModule ) << "* should represent the background contribution to the signal *" << endl; - report( ERROR, kModule ) << "* region. *" << endl; - report( ERROR, kModule ) << "****************************************************************" << endl; - - assert( false ); - } + if( m_ampVecsBkgnd.m_hasMixedSignWeights ){ + report( NOTICE, kModule ) << "***************************************************************" << endl; + report( NOTICE, kModule ) << "* NOTICE: Weights with both positive and negative signs were *" << endl; + report( NOTICE, kModule ) << "* detected in the background file. This may be desirable for *" << endl; + report( NOTICE, kModule ) << "* some applications. Older versions of AmpTools (v0.10.x and *" << endl; + report( NOTICE, kModule ) << "* prior) will not properly handle this case and will also not *" << endl; + report( NOTICE, kModule ) << "* print this notification to the screen. *" << endl; + report( NOTICE, kModule ) << "***************************************************************" << endl; + } + + m_sumBkgWeights = m_ampVecsBkgnd.m_dSumWeights; + + // the extra boolean allows MPI jobs to suppress this check which + // may fail on one of the follower nodes if the background sample + // is sparse + if( m_sumBkgWeights < 0 && !suppressError ){ + report( ERROR, kModule ) << "****************************************************************" << endl; + report( ERROR, kModule ) << "* ERROR: The sum of all background weights is negative. This *" << endl; + report( ERROR, kModule ) << "* implies a negative background in the signal region, which *" << endl; + report( ERROR, kModule ) << "* unphysical. The weighted sum of the background events *" << endl; + report( ERROR, kModule ) << "* should represent the background contribution to the signal *" << endl; + report( ERROR, kModule ) << "* region. *" << endl; + report( ERROR, kModule ) << "****************************************************************" << endl; - m_numBkgEvents = m_ampVecsBkgnd.m_iNTrueEvents; + assert( false ); } - report( DEBUG, kModule ) << "\tDone." << endl; + m_numBkgEvents = m_ampVecsBkgnd.m_iNTrueEvents; + m_firstBkgndCalc = false; } + + report( DEBUG, kModule ) << "\tDone." << endl; double sumLnI = m_intenManager.calcSumLogIntensity( m_ampVecsSignal ); @@ -306,17 +308,31 @@ SCOREP_USER_REGION_END( dataTerm ) } void -LikelihoodCalculator::bootstrapSignalData(){ +LikelihoodCalculator::bootstrapSignalData( unsigned int seed ){ #ifdef SCOREP SCOREP_USER_REGION_DEFINE( bootstrapSignalData ) SCOREP_USER_REGION_BEGIN( bootstrapSignalData, "bootstrapSignalData", SCOREP_USER_REGION_TYPE_COMMON ) #endif m_ampVecsSignal.deallocAmpVecs(); - m_dataReaderSignal->resample(); + m_dataReaderSignal->resample( seed ); m_firstDataCalc = true; - #ifdef SCOREP - SCOREP_USER_REGION_END( bootstrapSignalData ) - #endif +#ifdef SCOREP +SCOREP_USER_REGION_END( bootstrapSignalData ) +#endif +} + +void +LikelihoodCalculator::bootstrapBackgroundData( unsigned int seed ){ +#ifdef SCOREP +SCOREP_USER_REGION_DEFINE( bootstrapBackgroundData ) +SCOREP_USER_REGION_BEGIN( bootstrapBackgroundData, "bootstrapBackgroundData", SCOREP_USER_REGION_TYPE_COMMON ) +#endif + m_ampVecsBkgnd.deallocAmpVecs(); + m_dataReaderBkgnd->resample( seed ); + m_firstBkgndCalc = true; +#ifdef SCOREP +SCOREP_USER_REGION_END( bootstrapBackgroundData ) +#endif } void diff --git a/AmpTools/IUAmpTools/LikelihoodCalculator.h b/AmpTools/IUAmpTools/LikelihoodCalculator.h index e5c9351a..bc9cd0f5 100644 --- a/AmpTools/IUAmpTools/LikelihoodCalculator.h +++ b/AmpTools/IUAmpTools/LikelihoodCalculator.h @@ -88,7 +88,8 @@ class LikelihoodCalculator : public MIFunctionContribution void invalidateTerms(); - virtual void bootstrapSignalData(); + virtual void bootstrapSignalData( unsigned int seed ); + virtual void bootstrapBackgroundData( unsigned int seed); protected: @@ -123,6 +124,7 @@ class LikelihoodCalculator : public MIFunctionContribution bool m_firstNormIntCalc; bool m_firstDataCalc; + bool m_firstBkgndCalc; double* m_prodFactorArray; const double* m_normIntArray; diff --git a/AmpTools/IUAmpToolsMPI/DataReaderMPI.h b/AmpTools/IUAmpToolsMPI/DataReaderMPI.h index 9cf1c1a6..23e01d63 100644 --- a/AmpTools/IUAmpToolsMPI/DataReaderMPI.h +++ b/AmpTools/IUAmpToolsMPI/DataReaderMPI.h @@ -66,7 +66,7 @@ class DataReaderMPI : public T Kinematics* getEvent(); void resetSource(); - void resample(); + void resample( unsigned int seed ); unsigned int numEvents() const; @@ -231,7 +231,7 @@ void DataReaderMPI::resetSource() } template< class T > -void DataReaderMPI::resample() +void DataReaderMPI::resample( unsigned int seed ) { if( m_isLeader ){ @@ -239,7 +239,7 @@ void DataReaderMPI::resample() // redraw the bootstrap sample using the user's definition in the reader, then // push the selection out to the follower's exactly as done in startup - T::resample(); + T::resample( seed ); distributeData(); } else{ diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc index e28b7ec8..f71ad543 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc +++ b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.cc @@ -216,19 +216,43 @@ LikelihoodCalculatorMPI::operator()() } void -LikelihoodCalculatorMPI::bootstrapSignalData(){ +LikelihoodCalculatorMPI::bootstrapSignalData( unsigned int seed ){ if( m_isLeader ){ int cmnd[2]; cmnd[0] = m_thisId; cmnd[1] = LikelihoodManagerMPI::kBootstrapSignalData; - MPI_Bcast( cmnd, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast( cmnd, 2, MPI_INT, 0, MPI_COMM_WORLD); + + // second, command-specific broadcast carrying the randomized seed for bootstrapping + unsigned int seedBuf = seed; + MPI_Bcast( &seedBuf, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD ); + } + + // calls resample() on DataReaderMPI for the leader, and redistributes to the + // followers, discarding the old cached partition + LikelihoodCalculator::bootstrapSignalData( seed ); +} + +void +LikelihoodCalculatorMPI::bootstrapBackgroundData( unsigned int seed ){ + + if( m_isLeader ){ + + int cmnd[2]; + cmnd[0] = m_thisId; + cmnd[1] = LikelihoodManagerMPI::kBootstrapBackgroundData; + MPI_Bcast( cmnd, 2, MPI_INT, 0, MPI_COMM_WORLD); + + // second, command-specific broadcast carrying the randomized seed for bootstrapping + unsigned int seedBuf = seed; + MPI_Bcast( &seedBuf, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD ); } // calls resample() on DataReaderMPI for the leader, and redistributes to the // followers, discarding the old cached partition - LikelihoodCalculator::bootstrapSignalData(); + LikelihoodCalculator::bootstrapBackgroundData( seed ); } double diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h index 1a22f50f..3a4dab55 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h +++ b/AmpTools/IUAmpToolsMPI/LikelihoodCalculatorMPI.h @@ -128,7 +128,14 @@ class LikelihoodCalculatorMPI : public LikelihoodCalculator * broadcasts the command to all followers and resamples the leader's own * DataReaderMPI, pushing newly-drawn events to the followers. */ - void bootstrapSignalData(); + void bootstrapSignalData( unsigned int seed ); + + /** + * Triggers a bootstrap resampling of the background data. On a leader it + * broadcasts the command to all followers and resamples the leader's own + * DataReaderMPI, pushing newly-drawn events to the followers. + */ + void bootstrapBackgroundData( unsigned int seed ); /** * This sends the finalize fit flag to all follower jobs which breaks them diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc index c8fb5c4b..dd3aab23 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc +++ b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.cc @@ -72,6 +72,7 @@ LikelihoodManagerMPI::deliverLikelihood() int* fitFlag = &(cmnd[1]); LikelihoodCalculatorMPI* likCalc; + unsigned int seed; map< int, LikelihoodCalculatorMPI* >::iterator mapItr; MPI_Bcast( cmnd, 2, MPI_INT, 0, MPI_COMM_WORLD ); @@ -108,8 +109,17 @@ LikelihoodManagerMPI::deliverLikelihood() break; case kBootstrapSignalData: + + MPI_Bcast( &seed, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD ); + + likCalc->bootstrapSignalData( seed ); + break; - likCalc->bootstrapSignalData(); + case kBootstrapBackgroundData: + + MPI_Bcast( &seed, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD ); + + likCalc->bootstrapBackgroundData( seed ); break; default: diff --git a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h index eba6e7f7..045f9dfd 100644 --- a/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h +++ b/AmpTools/IUAmpToolsMPI/LikelihoodManagerMPI.h @@ -61,6 +61,7 @@ class LikelihoodManagerMPI kUpdateParameters, kUpdateAmpParameter, kBootstrapSignalData, + kBootstrapBackgroundData, kFinalizeFit, kExit }; From c0d527ad7dce6d886e292e40753e2b4108a8924d Mon Sep 17 00:00:00 2001 From: kevScheuer Date: Tue, 22 Sep 2026 04:49:32 -0400 Subject: [PATCH 3/3] Ensure bkgnd reader exists before bootstrapping Also added a missing include --- AmpTools/IUAmpTools/AmpToolsInterface.cc | 7 ++++++- AmpTools/IUAmpTools/AmpToolsInterface.h | 1 + AmpTools/IUAmpTools/LikelihoodCalculator.cc | 4 ++++ 3 files changed, 11 insertions(+), 1 deletion(-) diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.cc b/AmpTools/IUAmpTools/AmpToolsInterface.cc index 1348eacf..5a3a21ba 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.cc +++ b/AmpTools/IUAmpTools/AmpToolsInterface.cc @@ -409,17 +409,22 @@ AmpToolsInterface::bootstrapSignalData(const string& reactionName, unsigned int return; } likCalc->bootstrapSignalData( seed ); - + } void AmpToolsInterface::bootstrapBackgroundData(const string& reactionName, unsigned int seed){ LikelihoodCalculator* likCalc = likelihoodCalculator(reactionName); + DataReader* backgroundReader = bkgndReader(reactionName); if( likCalc == NULL ){ report( ERROR, kModule ) << "no LikelihoodCalculator for reaction: " << reactionName << endl; return; } + if( backgroundReader == NULL ){ + report( ERROR, kModule ) << "no background data for reaction: " << reactionName << endl; + return; + } likCalc->bootstrapBackgroundData( seed ); } diff --git a/AmpTools/IUAmpTools/AmpToolsInterface.h b/AmpTools/IUAmpTools/AmpToolsInterface.h index b2352aa9..914131a2 100644 --- a/AmpTools/IUAmpTools/AmpToolsInterface.h +++ b/AmpTools/IUAmpTools/AmpToolsInterface.h @@ -44,6 +44,7 @@ #include #include #include +#include #include "MinuitInterface/MinuitMinimizationManager.h" #include "IUAmpTools/IntensityManager.h" diff --git a/AmpTools/IUAmpTools/LikelihoodCalculator.cc b/AmpTools/IUAmpTools/LikelihoodCalculator.cc index 2df5c4d5..af6f3288 100644 --- a/AmpTools/IUAmpTools/LikelihoodCalculator.cc +++ b/AmpTools/IUAmpTools/LikelihoodCalculator.cc @@ -327,6 +327,10 @@ LikelihoodCalculator::bootstrapBackgroundData( unsigned int seed ){ SCOREP_USER_REGION_DEFINE( bootstrapBackgroundData ) SCOREP_USER_REGION_BEGIN( bootstrapBackgroundData, "bootstrapBackgroundData", SCOREP_USER_REGION_TYPE_COMMON ) #endif + if( !m_hasBackground ){ + report( ERROR, kModule ) << "Requested to bootstrap sample non-existant background dataset" << endl; + assert( false ); + } m_ampVecsBkgnd.deallocAmpVecs(); m_dataReaderBkgnd->resample( seed ); m_firstBkgndCalc = true;