33 gPARMS->SetDefaultParameter(
"ANALYSIS:DEBUG_LEVEL",
dDebugLevel);
39 vector<const DMCThrown*> locMCThrowns;
40 locEventLoop->Get(locMCThrowns);
47 for(
size_t loc_i = 0; loc_i < locReactions.size(); ++loc_i)
49 const DReaction* locReaction = locReactions[loc_i];
53 for(
size_t loc_j = 0; loc_j < locNumActions; ++loc_j)
57 cout <<
"Initialize Action # " << loc_j + 1 <<
": " << locAnalysisAction->
Get_ActionName() <<
" of reaction: " << locReaction->
Get_ReactionName() << endl;
61 if(locMCThrowns.empty())
65 bool locExactMatchFlag =
true;
67 locExactMatchFlag =
false;
68 else if(!locReactions[loc_i]->Get_MissingPIDs().empty())
69 locExactMatchFlag =
false;
92 set<string> locReactionNames;
93 set<string> locDuplicateReactionNames;
94 for(
auto& locReaction : locReactions)
96 string locReactionName = locReaction->Get_ReactionName();
97 if(locReactionNames.find(locReactionName) == locReactionNames.end())
98 locReactionNames.insert(locReactionName);
100 locDuplicateReactionNames.insert(locReactionName);
103 if(locDuplicateReactionNames.empty())
106 cout <<
"ERROR: MULTIPLE DREACTIONS WITH THE SAME NAME(S): " << endl;
107 for(
auto& locReactionName : locDuplicateReactionNames)
108 cout << locReactionName <<
", ";
110 cout <<
"ABORTING" << endl;
116 string locHistName, locHistTitle;
122 TDirectory* locCurrentDir = gDirectory;
124 for(
size_t loc_i = 0; loc_i < locReactions.size(); ++loc_i)
126 const DReaction* locReaction = locReactions[loc_i];
132 vector<string> locActionNames;
133 bool locPreKinFitFlag =
true;
134 for(
auto& locAction : locActions)
136 if(locPreKinFitFlag && locAction->Get_UseKinFitResultsFlag() && (locKinFitType !=
d_NoFit))
138 locPreKinFitFlag =
false;
139 locActionNames.push_back(
"KinFit Convergence");
141 locActionNames.push_back(locAction->Get_ActionName());
143 if(locPreKinFitFlag && (locKinFitType !=
d_NoFit))
144 locActionNames.push_back(
"KinFit Convergence");
146 string locDirName = locReactionName;
147 string locDirTitle = locReactionName;
150 TDirectoryFile* locDirectoryFile =
static_cast<TDirectoryFile*
>(locCurrentDir->GetDirectory(locDirName.c_str()));
151 if(locDirectoryFile == NULL)
152 locDirectoryFile =
new TDirectoryFile(locDirName.c_str(), locDirTitle.c_str());
153 locDirectoryFile->cd();
155 locHistName =
"NumParticleCombos";
156 loc1DHist =
static_cast<TH1D*
>(locDirectoryFile->Get(locHistName.c_str()));
157 if(loc1DHist == NULL)
159 double* locBinArray =
new double[55];
160 for(
unsigned int loc_j = 0; loc_j < 6; ++loc_j)
162 for(
unsigned int loc_k = 1; loc_k <= 9; ++loc_k)
163 locBinArray[loc_j*9 + loc_k - 1] =
double(loc_k)*pow(10.0,
double(loc_j));
165 locBinArray[54] = 1.0E6;
166 locHistTitle = locReactionName +
string(
";# Particle Combinations;# Events");
167 loc1DHist =
new TH1D(locHistName.c_str(), locHistTitle.c_str(), 54, locBinArray);
168 delete[] locBinArray;
172 locHistName =
"NumEventsSurvivedAction";
173 loc1DHist =
static_cast<TH1D*
>(locDirectoryFile->Get(locHistName.c_str()));
174 if(loc1DHist == NULL)
176 locHistTitle = locReactionName +
string(
";;# Events Survived Action");
177 loc1DHist =
new TH1D(locHistName.c_str(), locHistTitle.c_str(), locActionNames.size() + 2, -0.5, locActionNames.size() + 2.0 - 0.5);
178 loc1DHist->GetXaxis()->SetBinLabel(1,
"Input");
179 loc1DHist->GetXaxis()->SetBinLabel(2,
"Has Particle Combos");
180 for(
size_t loc_j = 0; loc_j < locActionNames.size(); ++loc_j)
181 loc1DHist->GetXaxis()->SetBinLabel(3 + loc_j, locActionNames[loc_j].c_str());
187 locHistName =
"NumEventsWhereTrueComboSurvivedAction";
188 loc1DHist =
static_cast<TH1D*
>(locDirectoryFile->Get(locHistName.c_str()));
189 if(loc1DHist == NULL)
191 locHistTitle = locReactionName +
string(
";;# Events Where True Combo Survived Action");
192 loc1DHist =
new TH1D(locHistName.c_str(), locHistTitle.c_str(), locActionNames.size() + 1, -0.5, locActionNames.size() + 1.0 - 0.5);
193 loc1DHist->GetXaxis()->SetBinLabel(1,
"Has Particle Combos");
194 for(
size_t loc_j = 0; loc_j < locActionNames.size(); ++loc_j)
195 loc1DHist->GetXaxis()->SetBinLabel(2 + loc_j, locActionNames[loc_j].c_str());
200 locHistName =
"NumCombosSurvivedAction";
201 loc2DHist =
static_cast<TH2D*
>(locDirectoryFile->Get(locHistName.c_str()));
202 if(loc2DHist == NULL)
204 double* locBinArray =
new double[55];
205 for(
unsigned int loc_j = 0; loc_j < 6; ++loc_j)
207 for(
unsigned int loc_k = 1; loc_k <= 9; ++loc_k)
208 locBinArray[loc_j*9 + loc_k - 1] =
double(loc_k)*pow(10.0,
double(loc_j));
210 locBinArray[54] = 1.0E6;
212 locHistTitle = locReactionName +
string(
";;# Particle Combos Survived Action");
213 loc2DHist =
new TH2D(locHistName.c_str(), locHistTitle.c_str(), locActionNames.size() + 1, -0.5, locActionNames.size() + 1 - 0.5, 54, locBinArray);
214 delete[] locBinArray;
215 loc2DHist->GetXaxis()->SetBinLabel(1,
"Has Particle Combos");
216 for(
size_t loc_j = 0; loc_j < locActionNames.size(); ++loc_j)
217 loc2DHist->GetXaxis()->SetBinLabel(2 + loc_j, locActionNames[loc_j].c_str());
221 locHistName =
"NumCombosSurvivedAction1D";
222 loc1DHist =
static_cast<TH1D*
>(locDirectoryFile->Get(locHistName.c_str()));
223 if(loc1DHist == NULL)
225 locHistTitle = locReactionName +
string(
";;# Particle Combos Survived Action");
226 loc1DHist =
new TH1D(locHistName.c_str(), locHistTitle.c_str(), locActionNames.size() + 1, -0.5, locActionNames.size() + 1 - 0.5);
227 loc1DHist->GetXaxis()->SetBinLabel(1,
"Combos Constructed");
228 for(
size_t loc_j = 0; loc_j < locActionNames.size(); ++loc_j)
229 loc1DHist->GetXaxis()->SetBinLabel(2 + loc_j, locActionNames[loc_j].c_str());
233 locDirectoryFile->cd(
"..");
246 VT_TRACER(
"DAnalysisResults_factory::evnt()");
250 cout <<
"Analyze event: " << eventnumber << endl;
254 locEventLoop->GetSingle(locTrigger);
271 cout <<
"# DReactions: " << locReactions.size() << endl;
274 vector<const DReactionVertexInfo*> locReactionVertexInfos;
275 locEventLoop->Get(locReactionVertexInfos);
277 for(
auto& locReactionVertexInfo : locReactionVertexInfos)
281 cout <<
"Build combos for reaction: " << locReactionVertexInfo->Get_Reaction()->Get_ReactionName() << endl;
285 for(
auto& locReactionComboPair : locReactionComboMap)
287 auto& locReaction = locReactionComboPair.first;
288 auto& locCombos = locReactionComboPair.second;
293 auto locTrueParticleCombo =
Find_TrueCombo(locEventLoop, locReaction, locCombos);
294 int locLastActionTrueComboSurvives = (locTrueParticleCombo !=
nullptr) ? -1 : -2;
298 locAnalysisResults->Set_Reaction(locReaction);
301 auto locActions = locReaction->Get_AnalysisActions();
302 for(
auto& locAction : locActions)
303 locAction->Reset_NewEvent();
306 cout <<
"Combos built, execute actions for reaction: " << locReaction->Get_ReactionName() << endl;
309 auto locIsKinFit = (locReaction->Get_KinFitType() !=
d_NoFit);
310 auto locNumActionsForHist = locIsKinFit ? locActions.size() + 2 : locActions.size() + 1;
311 vector<size_t> locNumCombosSurvived(locNumActionsForHist, 0);
312 locNumCombosSurvived[0] = locCombos.size();
313 for(
auto& locCombo : locCombos)
316 size_t locActionIndex = 0;
317 if(!
Execute_Actions(locEventLoop, locIsKinFit, locCombo, locTrueParticleCombo,
true, locActions, locActionIndex, locNumCombosSurvived, locLastActionTrueComboSurvives))
321 auto locPostKinFitCombo =
Handle_ComboFit(locReactionVertexInfo, locCombo, locReaction);
322 if(locPostKinFitCombo ==
nullptr)
325 ++(locNumCombosSurvived[locActionIndex + 1]);
328 if(!
Execute_Actions(locEventLoop, locIsKinFit, locPostKinFitCombo, locTrueParticleCombo,
false, locActions, locActionIndex, locNumCombosSurvived, locLastActionTrueComboSurvives))
332 locAnalysisResults->Add_PassedParticleCombo(locPostKinFitCombo);
336 cout <<
"Action loop completed, # surviving combos: " << locAnalysisResults->Get_NumPassedParticleCombos() << endl;
339 japp->WriteLock(
"DAnalysisResults");
342 if(locNumCombosSurvived[0] > 0)
344 for(
size_t loc_j = 0; loc_j < locNumCombosSurvived.size(); ++loc_j)
346 if(locNumCombosSurvived[loc_j] == 0)
350 if(locNumCombosSurvived[loc_j] > 0)
358 for(
int loc_j = -1; loc_j <= locLastActionTrueComboSurvives; ++loc_j)
362 japp->Unlock(
"DAnalysisResults");
365 _data.push_back(locAnalysisResults);
374 for(; locActionIndex < locActions.size(); ++locActionIndex)
376 auto locAction = locActions[locActionIndex];
377 if(locPreKinFitFlag && locAction->Get_UseKinFitResultsFlag())
380 cout <<
"Execute action " << locActionIndex <<
": " << locAction->Get_ActionName() << endl;
381 if(!(*locAction)(locEventLoop, locCombo))
384 auto locComboSurvivedIndex = (locIsKinFit && !locPreKinFitFlag) ? locActionIndex + 2 : locActionIndex + 1;
385 ++(locNumCombosSurvived[locComboSurvivedIndex]);
386 if(locCombo == locTrueCombo)
387 locLastActionTrueComboSurvives = locActionIndex;
398 for(
auto& locCombo : locCombos)
400 if((*locAction)(locEventLoop, locCombo))
410 return locParticleCombo;
415 auto locComboKinFitTuple = std::make_tuple(locParticleCombo, locKinFitType, locUpdateCovMatricesFlag, locNoConstrainMassSteps);
420 return locComboIterator->second;
424 cout <<
"Do kinfit" << endl;
425 auto locKinFitResultsPair =
Fit_Kinematics(locReactionVertexInfo, locReaction, locParticleCombo, locKinFitType, locUpdateCovMatricesFlag);
426 if(locKinFitResultsPair.second ==
nullptr)
431 cout <<
"Create new combo" << endl;
434 return locNewParticleCombo;
443 vector<shared_ptr<DKinFitConstraint_Vertex>> locSortedVertexConstraints;
444 auto locConstraints =
dKinFitUtils->
Create_Constraints(locReactionVertexInfo, locReaction, locParticleCombo, locKinFitChain, locKinFitType, locSortedVertexConstraints);
445 if(locConstraints.empty())
446 return pair<shared_ptr<const DKinFitChain>,
const DKinFitResults*>(
nullptr,
nullptr);
449 auto locResultPair = std::make_pair(locConstraints, locUpdateCovMatricesFlag);
455 if(locKinFitResults !=
nullptr)
460 return std::make_pair(locOutputKinFitChain, locKinFitResults);
464 return pair<shared_ptr<const DKinFitChain>,
const DKinFitResults*>(
nullptr,
nullptr);
479 locKinFitResults =
Build_KinFitResults(locParticleCombo, locKinFitType, locOutputKinFitChain);
481 return std::make_pair(locOutputKinFitChain, locKinFitResults);
487 return pair<shared_ptr<const DKinFitChain>,
const DKinFitResults*>(
nullptr,
nullptr);
493 locKinFitResults->Set_KinFitType(locKinFitType);
508 locKinFitResults->Add_OutputKinFitParticles(locOutputKinFitParticles);
514 map<const JObject*, map<DKinFitPullType, double> > locPulls_JObject;
517 map<shared_ptr<DKinFitParticle>, map<DKinFitPullType, double> > locPulls_KinFitParticle;
521 auto locMapIterator = locPulls_KinFitParticle.begin();
522 for(; locMapIterator != locPulls_KinFitParticle.end(); ++locMapIterator)
524 auto locOutputKinFitParticle = locMapIterator->first;
528 locPulls_JObject[locSourceJObject] = locMapIterator->second;
531 locKinFitResults->Set_Pulls(locPulls_JObject);
534 for(
auto& locKinFitParticle : locOutputKinFitParticles)
536 if(locKinFitParticle ==
nullptr)
540 if(locSourceJObject != NULL)
542 locKinFitResults->Add_ParticleMapping_SourceToOutput(locSourceJObject, locKinFitParticle);
547 if(locInputKinFitParticle != NULL)
550 if(locSourceJObject != NULL)
551 locKinFitResults->Add_ParticleMapping_SourceToOutput(locSourceJObject, locKinFitParticle);
555 return locKinFitResults;
void Reset_NewEvent(void)
bool Get_IsFirstStepBeam(const DReaction *locReaction)
void Set_UpdateCovarianceMatricesFlag(bool locUpdateCovarianceMatricesFlag)
set< size_t > Get_NoConstrainMassSteps(const DReaction *locReaction)
double dMinThrownMatchFOM
jerror_t brun(JEventLoop *locEventLoop, int32_t runnumber)
Called everytime a new run number is detected.
bool dRequireKinFitConvergence
void Check_ReactionNames(vector< const DReaction * > &locReactions) const
void Reset_NewEvent(void)
map< tuple< const DParticleCombo *, DKinFitType, bool, set< size_t > >, const DParticleCombo * > dPreToPostKinFitComboMap
DKinFitType Get_KinFitType(void) const
const DParticleCombo * Handle_ComboFit(const DReactionVertexInfo *locReactionVertexInfo, const DParticleCombo *locParticleCombo, const DReaction *locReaction)
DKinFitResults * Get_KinFitResultsResource(void)
vector< const DReaction * > Get_Reactions(JEventLoop *locEventLoop)
DCombosByReaction Build_ParticleCombos(const DReactionVertexInfo *locReactionVertexInfo)
unordered_map< const DReaction *, TH1 * > dHistMap_NumEventsWhereTrueComboSurvivedAction
const DParticleCombo * Create_KinFitCombo_NewCombo(const DParticleCombo *locOrigCombo, const DReaction *locReaction, const DKinFitResults *locKinFitResults, const shared_ptr< const DKinFitChain > &locKinFitChain)
void Make_ControlHistograms(vector< const DReaction * > &locReactions)
unordered_map< const DReaction *, DCutAction_TrueCombo * > dTrueComboCuts
DApplication * dApplication
unordered_map< const DReaction *, bool > dMCReactionExactMatchFlags
unordered_map< const DReaction *, TH1 * > dHistMap_NumEventsSurvivedAction_All
bool Execute_Actions(JEventLoop *locEventLoop, bool locIsKinFit, const DParticleCombo *locCombo, const DParticleCombo *locTrueCombo, bool locPreKinFitFlag, const vector< DAnalysisAction * > &locActions, size_t &locActionIndex, vector< size_t > &locNumCombosSurvived, int &locLastActionTrueComboSurvives)
void Recycle_LastFitMemory(void)
unordered_map< const DReaction *, TH1 * > dHistMap_NumCombosSurvivedAction1D
const DParticleCombo * Find_TrueCombo(JEventLoop *locEventLoop, const DReaction *locReaction, const vector< const DParticleCombo * > &locCombos)
bool Get_IsPhysicsEvent(void) const
DParticleComboCreator * Get_ParticleComboCreator(void) const
vector< DKinFitResults * > dCreatedKinFitResults
size_t Get_NumAnalysisActions(void) const
map< pair< set< shared_ptr< DKinFitConstraint > >, bool >, DKinFitResults * > dConstraintResultsMap
DParticleComboCreator * dParticleComboCreator
unsigned int Get_NumConstraintEquations(void) const
jerror_t init(void)
Called once at program start.
unsigned int dKinFitDebugLevel
set< shared_ptr< DKinFitParticle > > Get_OutputKinFitParticles(void) const
shared_ptr< const DKinFitChain > Build_OutputKinFitChain(const shared_ptr< const DKinFitChain > &locInputKinFitChain, set< shared_ptr< DKinFitParticle >> &locKinFitOutputParticles)
string Get_ReactionName(void) const
DResourcePool< DKinFitResults > dResourcePool_KinFitResults
set< shared_ptr< DKinFitConstraint > > Get_KinFitConstraints(void) const
unsigned int Get_NDF(void) const
h_means GetYaxis() -> SetTitleOffset(1.8)
void Set_DebugLevel(int locDebugLevel)
double Get_ChiSq(void) const
unordered_map< const DReaction *, TH2 * > dHistMap_NumCombosSurvivedAction
DKinFitResults * Build_KinFitResults(const DParticleCombo *locParticleCombo, DKinFitType locKinFitType, const shared_ptr< const DKinFitChain > &locKinFitChain)
void Add_Constraints(const set< shared_ptr< DKinFitConstraint >> &locKinFitConstraints)
set< shared_ptr< DKinFitConstraint > > Create_Constraints(const DReactionVertexInfo *locReactionVertexInfo, const DReaction *locReaction, const DParticleCombo *locParticleCombo, const shared_ptr< const DKinFitChain > &locKinFitChain, DKinFitType locKinFitType, vector< shared_ptr< DKinFitConstraint_Vertex >> &locSortedVertexConstraints)
void Get_Pulls(map< shared_ptr< DKinFitParticle >, map< DKinFitPullType, double > > &locPulls) const
const JObject * Get_SourceJObject(const shared_ptr< DKinFitParticle > &locInputKinFitParticle) const
void Reset_NewEvent(JEventLoop *locEventLoop)
const TMatrixDSym & Get_VXi(void)
set< shared_ptr< DKinFitParticle > > Get_KinFitParticles(void) const
double Get_ConfidenceLevel(void) const
jerror_t evnt(JEventLoop *locEventLoop, uint64_t eventnumber)
Called every event.
bool Get_KinFitUpdateCovarianceMatricesFlag(void) const
void Set_ControlParams(size_t locGetBatchSize, size_t locNumToAllocateAtOnce, size_t locMaxLocalPoolSize)
shared_ptr< DKinFitParticle > Get_InputKinFitParticle(const shared_ptr< DKinFitParticle > &locKinFitParticle) const
unsigned int Get_NumUnknowns(void) const
unordered_map< const DReaction *, TH1 * > dHistMap_NumParticleCombos
DSourceComboer * dSourceComboer
pair< shared_ptr< const DKinFitChain >, const DKinFitResults * > Fit_Kinematics(const DReactionVertexInfo *locReactionVertexInfo, const DReaction *locReaction, const DParticleCombo *locParticleCombo, DKinFitType locKinFitType, bool locUpdateCovMatricesFlag)
DKinFitUtils_GlueX * dKinFitUtils
virtual void Initialize(JEventLoop *locEventLoop)=0
shared_ptr< const DKinFitChain > Make_KinFitChain(const DReactionVertexInfo *locReactionVertexInfo, const DReaction *locReaction, const DParticleCombo *locParticleCombo, DKinFitType locKinFitType)
vector< DAnalysisAction * > Get_AnalysisActions(void) const
void Recycle(const DType *locResource)
virtual string Get_ActionName(void) const