Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DSourceComboVertexer.h
Go to the documentation of this file.
1 #ifndef DSourceComboVertexer_h
2 #define DSourceComboVertexer_h
3 
4 #include <set>
5 #include <unordered_map>
6 #include <utility>
7 #include <memory>
8 #include <vector>
9 #include <algorithm>
10 #include <iostream>
11 #include <sstream>
12 #include <map>
13 
14 #include "JANA/JEventLoop.h"
15 
16 #include "particleType.h"
18 #include "PID/DNeutralShower.h"
19 #include "ANALYSIS/DSourceCombo.h"
20 #include "ANALYSIS/DReaction.h"
24 
25 #include "PID/DVertex.h"
26 using namespace std;
27 
28 namespace DAnalysis
29 {
30 
31 class DSourceComboer;
34 
36 {
37  public:
38 
39  //CONSTRUCTORS
40  DSourceComboVertexer(void) = delete;
41  DSourceComboVertexer(JEventLoop* locEventLoop, DSourceComboer* locSourceComboer, DSourceComboP4Handler* locSourceComboP4Handler);
42  void Reset(void);
43 void Set_Vertex(const DVertex* locVertex){dVertex = locVertex;} //COMPARE
44  //SETUP
45  void Set_SourceComboTimeHandler(const DSourceComboTimeHandler* locSourceComboTimeHandler){dSourceComboTimeHandler = locSourceComboTimeHandler;}
46  void Set_DebugLevel(int locDebugLevel){dDebugLevel = locDebugLevel;}
47 
48  //COMPUTE
49  void Calc_VertexTimeOffsets_WithCharged(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locReactionChargedCombo);
50  void Calc_VertexTimeOffsets_WithPhotons(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locReactionChargedCombo, const DSourceCombo* locReactionFullCombo);
51  void Calc_VertexTimeOffsets_WithBeam(const DReactionVertexInfo* locReactionVertexInfo, const DSourceComboUse& locReactionFullComboUse, const DSourceCombo* locReactionFullCombo, const DKinematicData* locBeamParticle);
52 
53  bool Get_VertexDeterminableWithCharged(const DReactionStepVertexInfo* locStepVertexInfo) const;
54  bool Get_VertexDeterminableWithPhotons(const DReactionStepVertexInfo* locStepVertexInfo) const;
55 
56  //GET IS-KNOWN
57  bool Get_IsVertexKnown(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const;
58  bool Get_IsVertexKnown_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const;
59  bool Get_IsTimeOffsetKnown(bool locIsPrimaryProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle) const;
60 
61  //GET VERTEX
62  DVector3 Get_Vertex_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const;
63  DVector3 Get_Vertex(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const;
64  DVector3 Get_Vertex(bool locIsProductionVertex, const vector<const DKinematicData*>& locVertexParticles) const;
65  DVector3 Get_Vertex(const DReactionStepVertexInfo* locStepVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle, bool locComboIsFullyCharged) const; //Only one that will handle dangling correctly!
66  DVector3 Get_PrimaryVertex(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle) const;
67 
68  //GET TIME OFFSET
69  double Get_TimeOffset(bool locIsPrimaryProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle) const;
70  double Get_TimeOffset(const DReactionVertexInfo* locReactionVertexInfo, const DReactionStepVertexInfo* locStepVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle) const; //Only one that will handle dangling correctly!
71 
72  //GET CONSTRAINING PARTICLES (for vertex)
73  vector<const DKinematicData*> Get_ConstrainingParticles(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const;
74  vector<const DKinematicData*> Get_ConstrainingParticles_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const;
75 
76  //GET VERTEX-Z BINS
77  signed char Get_VertexZBin(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locPrimaryVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const;
78  signed char Get_VertexZBin_NoBeam(bool locIsProductionVertex, const DSourceCombo* locPrimaryVertexCombo, bool locIsCombo2ndVertex) const;
79  signed char Get_VertexZBin(const DReactionStepVertexInfo* locStepVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle, bool locComboIsFullyCharged) const; //Only one that will handle dangling correctly!
80  vector<signed char> Get_VertexZBins(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle, bool locComboIsFullyCharged) const; //This will call the above
81 
82  private:
83  vector<const DKinematicData*>::const_iterator Get_ThetaNearest90Iterator(const vector<const DKinematicData*>& locParticles);
84  vector<const DKinematicData*> Get_FullConstrainDecayingParticles(const DReactionStepVertexInfo* locStepVertexInfo, const map<pair<int, int>, const DKinematicData*>& locReconDecayParticleMap);
85 
86  DVector3 Calc_Vertex(bool locIsProductionVertexFlag, const vector<pair<Particle_t, const JObject*>>& locChargedSourceParticles, const vector<const DKinematicData*>& locDecayingParticles, vector<const DKinematicData*>& locVertexParticles);
87  void Calc_TimeOffsets(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locChargedReactionCombo, const DSourceCombo* locFullReactionCombo = nullptr);
88 
89  void Construct_DecayingParticle_InvariantMass(const DReactionStepVertexInfo* locReactionStepVertexInfo, const DSourceCombo* locVertexCombo, DVector3 locVertex, map<pair<int, int>, const DKinematicData*>& locReconDecayParticleMap);
90  void Construct_DecayingParticle_MissingMass(const DReactionStepVertexInfo* locReactionStepVertexInfo, const DSourceComboUse& locReactionFullComboUse, const DSourceCombo* locReactionFullCombo, const DSourceCombo* locFullVertexCombo, const DKinematicData* locBeamParticle, DVector3 locVertex, int locRFBunch, double locRFVertexTime, map<pair<int, int>, const DKinematicData*>& locReconDecayParticleMap);
91 
92  //HANDLERS/ETC.
95  const DSourceComboTimeHandler* dSourceComboTimeHandler = nullptr;
97  int dDebugLevel = 0;
98 
99  //EXPERIMENT INFORMATION
101 const DVertex* dVertex; //COMPARE
102  double dMinThetaForVertex = 30.0;
103 
104  //DETERMINABILITY
105  unordered_map<const DReactionStepVertexInfo*, bool> dVertexDeterminableWithChargedMap; //excludes dangling vertex infos!! //only includes primary combos at each vertex
106  unordered_map<const DReactionStepVertexInfo*, bool> dVertexDeterminableWithPhotonsMap; //excludes determinable-by-charged & dangling vertex infos!! //only includes primary combos at each vertex
107 
108  //TIME OFFSETS
109  //time offsets & (sometimes) vertices depend on the ENTIRE reaction combo, not just the downstream ones! //time offset is from the RF time
110  //bool: is the PRIMARY vertex a production vertex //why is bool used throughout here?
111  //because the vertexing is different whether it's a production vertex or not, and a given combo of particles can be used either way
112  map<tuple<bool, const DSourceCombo*, const DKinematicData*>, unordered_map<const DSourceCombo*, double>> dTimeOffsets; //first combo: primary reaction combo //kinematics: beam particle
113 
114  //VERTEX-CONSTRAINING PARTICLES
115  //kinematic data: beam particle
116  //first bool: is vertex combo production-vertex flag. However, this flag is true for EVERY vertex that the beam is needed to find
117  //second bool: is vertex combo a charged combo re-used in place of a has-neutrals combo (e.g. Xi0 -> pi0, Lambda: The lambda combo is used for Xi0 during charged-only stage)
118  //false if actual slot (e.g. Lambda decay), true if re-used (e.g. Xi0 decay) //only true in charged-only stage!
119  //when looping in dependency-order over vertex-infos, it will always be false first, then true //only true in charged-only stage!
120  //worst-case?: g, p -> K0, Sigma+, Sigma+ -> pi0, (p)
121  //the pi+, pi- combo is used for the K0 decay, the K0 is used to find the production vertex
122  //AND the production vertex is used for the Sigma+ decay vertex (is dangling, indeterminate): 3 vertices!
123  //however, the 3rd vertex result is the same as the 2nd since copied: call with bool = true
124  //If beam is NOT needed to find the vertex then the primary reaction combo (first combo) is nullptr!!! (also not needed)
125  map<tuple<bool, const DSourceCombo*, const DSourceCombo*, const DKinematicData*, bool>, vector<const DKinematicData*>> dConstrainingParticlesByCombo; //first combo: primary reaction combo
126 
127  //VERTICES
128  //Note that this only includes the particles used to find the vertex (+ rarely an extra or 2), not necessarily ALL of those at the vertex
129  //bool: is production vertex
130  map<pair<bool, vector<const DKinematicData*>>, DVector3> dVertexMap; //vector: from dConstrainingParticlesByCombo
131 
132  //Reconstructed Decaying Particles
133  map<tuple<Particle_t, bool, const DSourceCombo*, const DKinematicData*>, const DKinematicData*> dReconDecayParticles_FromProducts; //pid, is prod vertex, full vertex combo, beam particle
134  map<tuple<Particle_t, const DSourceCombo*, bool, const DSourceCombo*, const DKinematicData*>, const DKinematicData*> dReconDecayParticles_FromMissing; //decay pid, full reaction combo, is prod vertex, full vertex combo, beam particle (this order)
135 
136  //RESOURCE POOL
138 };
139 
140 inline void DSourceComboVertexer::Reset(void)
141 {
142  dConstrainingParticlesByCombo.clear();
143  dVertexMap.clear();
144  dTimeOffsets.clear();
145 
146  for(const auto& locParticlePair : dReconDecayParticles_FromProducts)
147  dResourcePool_KinematicData.Recycle(locParticlePair.second);
148  for(const auto& locParticlePair : dReconDecayParticles_FromMissing)
149  dResourcePool_KinematicData.Recycle(locParticlePair.second);
150 
151  dReconDecayParticles_FromProducts.clear();
152  dReconDecayParticles_FromMissing.clear();
153 
154  //undeterminable vertices
155  dVertexMap.emplace(std::make_pair(false, vector<const DKinematicData*>()), dTargetCenter);
156  dVertexMap.emplace(std::make_pair(true, vector<const DKinematicData*>()), dTargetCenter);
157 }
158 
159 inline bool DSourceComboVertexer::Get_IsTimeOffsetKnown(bool locIsPrimaryProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle) const
160 {
161  //the data member MAY be dependent on the beam particle, but it may not
162  //so, first search with the beam particle; if not found then search without it
163  if(locBeamParticle == nullptr)
164  {
165  auto locIterator = dTimeOffsets.find(std::make_tuple(locIsPrimaryProductionVertex, locReactionCombo, nullptr));
166  if(locIterator == dTimeOffsets.end())
167  return false;
168  auto& locComboMap = locIterator->second;
169  auto locComboIterator = locComboMap.find(locVertexCombo);
170  return (locComboIterator != locComboMap.end());
171  }
172 
173  auto locIterator = dTimeOffsets.find(std::make_tuple(locIsPrimaryProductionVertex, locReactionCombo, locBeamParticle));
174  if(locIterator == dTimeOffsets.end())
175  return Get_TimeOffset(locIsPrimaryProductionVertex, locReactionCombo, locVertexCombo, nullptr); //try without beam
176 
177  auto& locComboMap = locIterator->second;
178  auto locComboIterator = locComboMap.find(locVertexCombo);
179  if(locComboIterator == locComboMap.end())
180  return Get_TimeOffset(locIsPrimaryProductionVertex, locReactionCombo, locVertexCombo, nullptr); //try without beam
181 
182  return true;
183 }
184 
185 inline double DSourceComboVertexer::Get_TimeOffset(bool locIsPrimaryProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle) const
186 {
187  //the data member MAY be dependent on the beam particle, but it may not
188  //so, first search with the beam particle; if not found then search without it
189  if(locBeamParticle == nullptr)
190  {
191  auto locIterator = dTimeOffsets.find(std::make_tuple(locIsPrimaryProductionVertex, locReactionCombo, nullptr));
192  if(locIterator == dTimeOffsets.end())
193  return 0.0;
194  auto& locComboMap = locIterator->second;
195  auto locComboIterator = locComboMap.find(locVertexCombo);
196  return ((locComboIterator != locComboMap.end()) ? locComboIterator->second : 0.0);
197  }
198 
199  auto locIterator = dTimeOffsets.find(std::make_tuple(locIsPrimaryProductionVertex, locReactionCombo, locBeamParticle));
200  if(locIterator == dTimeOffsets.end())
201  return Get_TimeOffset(locIsPrimaryProductionVertex, locReactionCombo, locVertexCombo, nullptr); //try without beam
202 
203  auto& locComboMap = locIterator->second;
204  auto locComboIterator = locComboMap.find(locVertexCombo);
205  if(locComboIterator == locComboMap.end())
206  return Get_TimeOffset(locIsPrimaryProductionVertex, locReactionCombo, locVertexCombo, nullptr); //try without beam
207 
208  return locComboIterator->second;
209 }
210 
211 inline vector<const DKinematicData*> DSourceComboVertexer::Get_ConstrainingParticles(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const
212 {
213  //the data member MAY be dependent on the beam particle, but it may not
214  //so, first search with the beam particle; if not found then search without it
215 
216 //cout << "Query tuple: " << locIsProductionVertex << ", " << locReactionCombo << ", " << locVertexCombo << ", " << locBeamParticle << ", " << locIsCombo2ndVertex << endl;
217 
218  if(locBeamParticle == nullptr)
219  {
220  auto locIterator = dConstrainingParticlesByCombo.find(std::make_tuple(locIsProductionVertex, (const DSourceCombo*)nullptr, locVertexCombo, (const DKinematicData*)nullptr, locIsCombo2ndVertex));
221  if(locIterator != dConstrainingParticlesByCombo.end())
222  return locIterator->second;
223  if(!locIsCombo2ndVertex)
224  return {};
225 
226  //We MAY have the wrong flag for is-prod-vertex, due to dangling vertex on 2nd combo (e.g. worst-case scenario detailed above). try the other one
227  locIterator = dConstrainingParticlesByCombo.find(std::make_tuple(!locIsProductionVertex, (const DSourceCombo*)nullptr, locVertexCombo, (const DKinematicData*)nullptr, locIsCombo2ndVertex));
228  if(locIterator != dConstrainingParticlesByCombo.end())
229  return locIterator->second;
230  return {};
231  }
232 
233  auto locIterator = dConstrainingParticlesByCombo.find(std::make_tuple(true, locReactionCombo, locVertexCombo, locBeamParticle, locIsCombo2ndVertex));
234  if(locIterator == dConstrainingParticlesByCombo.end()) //see if beam is not needed
235  locIterator = dConstrainingParticlesByCombo.find(std::make_tuple(locIsProductionVertex, (const DSourceCombo*)nullptr, locVertexCombo, (const DKinematicData*)nullptr, locIsCombo2ndVertex));
236  if(locIterator != dConstrainingParticlesByCombo.end())
237  return locIterator->second;
238  return {};
239 }
240 
241 inline DVector3 DSourceComboVertexer::Get_Vertex(bool locIsProductionVertex, const vector<const DKinematicData*>& locVertexParticles) const
242 {
243  auto locIterator = dVertexMap.find(std::make_pair(locIsProductionVertex, locVertexParticles));
244  if(locIterator != dVertexMap.end())
245  return locIterator->second;
246  return dTargetCenter;
247 }
248 
249 inline signed char DSourceComboVertexer::Get_VertexZBin_NoBeam(bool locIsProductionVertex, const DSourceCombo* locPrimaryVertexCombo, bool locIsCombo2ndVertex) const
250 {
251  return Get_VertexZBin(locIsProductionVertex, nullptr, locPrimaryVertexCombo, nullptr, locIsCombo2ndVertex);
252 }
253 
254 inline vector<const DKinematicData*> DSourceComboVertexer::Get_ConstrainingParticles_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const
255 {
256  return Get_ConstrainingParticles(locIsProductionVertex, nullptr, locVertexCombo, nullptr, locIsCombo2ndVertex);
257 }
258 
259 inline bool DSourceComboVertexer::Get_IsVertexKnown(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const
260 {
261  return !Get_ConstrainingParticles(locIsProductionVertex, locReactionCombo, locVertexCombo, locBeamParticle, locIsCombo2ndVertex).empty();
262 }
263 
264 inline bool DSourceComboVertexer::Get_IsVertexKnown_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const
265 {
266  return !Get_ConstrainingParticles_NoBeam(locIsProductionVertex, locVertexCombo, locIsCombo2ndVertex).empty();
267 }
268 
269 //2 cases: you KNOW beam photon MUST be nullptr, and one where it MAY be nullptr
270 inline DVector3 DSourceComboVertexer::Get_Vertex_NoBeam(bool locIsProductionVertex, const DSourceCombo* locVertexCombo, bool locIsCombo2ndVertex) const
271 {
272  return Get_Vertex(locIsProductionVertex, Get_ConstrainingParticles(locIsProductionVertex, nullptr, locVertexCombo, nullptr, locIsCombo2ndVertex));
273 }
274 
275 inline DVector3 DSourceComboVertexer::Get_Vertex(bool locIsProductionVertex, const DSourceCombo* locReactionCombo, const DSourceCombo* locVertexCombo, const DKinematicData* locBeamParticle, bool locIsCombo2ndVertex) const
276 {
277  //bool: for the vertex we want, not the primary
278  auto locConstrainingParticles = Get_ConstrainingParticles(true, locReactionCombo, locVertexCombo, locBeamParticle, locIsCombo2ndVertex);
279  if(locConstrainingParticles.empty())
280  locConstrainingParticles = Get_ConstrainingParticles(locIsProductionVertex, nullptr, locVertexCombo, nullptr, locIsCombo2ndVertex);
281  return Get_Vertex(locIsProductionVertex, locConstrainingParticles);
282 }
283 
284 inline DVector3 DSourceComboVertexer::Get_PrimaryVertex(const DReactionVertexInfo* locReactionVertexInfo, const DSourceCombo* locReactionCombo, const DKinematicData* locBeamParticle) const
285 {
286  auto locStepVertexInfo = locReactionVertexInfo->Get_StepVertexInfo(0);
287  auto locIsProductionVertex = locStepVertexInfo->Get_ProductionVertexFlag();
288  auto locComboIsFullyCharged = locReactionCombo->Get_SourceParticles(true, d_Neutral).empty();
289  auto locIsCombo2ndVertex = (locComboIsFullyCharged && locStepVertexInfo->Get_FullConstrainParticles(false, d_FinalState, d_Charged, false).empty());
290  return Get_Vertex(locIsProductionVertex, locReactionCombo, locReactionCombo, locBeamParticle, locIsCombo2ndVertex);
291 }
292 
293 inline vector<const DKinematicData*>::const_iterator DSourceComboVertexer::Get_ThetaNearest90Iterator(const vector<const DKinematicData*>& locParticles)
294 {
295  //true if first less than second
296  auto Get_Nearer90Theta = [](const DKinematicData* lhs, const DKinematicData* rhs) -> bool
297  {return fabs(rhs->momentum().Theta() - 0.5*TMath::Pi()) < fabs(lhs->momentum().Theta() - 0.5*TMath::Pi());};
298  return std::max_element(locParticles.begin(), locParticles.end(), Get_Nearer90Theta);
299 }
300 
302 {
303  auto locIterator = dVertexDeterminableWithChargedMap.find(locStepVertexInfo);
304  if(locIterator == dVertexDeterminableWithChargedMap.end())
305  return false;
306  return locIterator->second;
307 }
308 
310 {
311  auto locIterator = dVertexDeterminableWithPhotonsMap.find(locStepVertexInfo);
312  if(locIterator == dVertexDeterminableWithPhotonsMap.end())
313  return false;
314  return locIterator->second;
315 }
316 
317 } //end DAnalysis namespace
318 
319 #endif // DSourceComboVertexer_h
void Set_SourceComboTimeHandler(const DSourceComboTimeHandler *locSourceComboTimeHandler)
bool Get_VertexDeterminableWithPhotons(const DReactionStepVertexInfo *locStepVertexInfo) const
signed char Get_VertexZBin_NoBeam(bool locIsProductionVertex, const DSourceCombo *locPrimaryVertexCombo, bool locIsCombo2ndVertex) const
unordered_map< const DReactionStepVertexInfo *, bool > dVertexDeterminableWithPhotonsMap
TVector3 DVector3
Definition: DVector3.h:14
void Set_DebugLevel(int locDebugLevel)
const DReactionStepVertexInfo * Get_StepVertexInfo(size_t locStepIndex) const
map< tuple< bool, const DSourceCombo *, const DSourceCombo *, const DKinematicData *, bool >, vector< const DKinematicData * > > dConstrainingParticlesByCombo
vector< const DKinematicData * > Get_ConstrainingParticles_NoBeam(bool locIsProductionVertex, const DSourceCombo *locVertexCombo, bool locIsCombo2ndVertex) const
bool Get_VertexDeterminableWithCharged(const DReactionStepVertexInfo *locStepVertexInfo) const
DVector3 Get_Vertex(bool locIsProductionVertex, const DSourceCombo *locReactionCombo, const DSourceCombo *locVertexCombo, const DKinematicData *locBeamParticle, bool locIsCombo2ndVertex) const
DVector3 Get_PrimaryVertex(const DReactionVertexInfo *locReactionVertexInfo, const DSourceCombo *locReactionCombo, const DKinematicData *locBeamParticle) const
tuple< Particle_t, signed char, const DSourceComboInfo *, bool, Particle_t > DSourceComboUse
Definition: DSourceCombo.h:34
vector< pair< Particle_t, const JObject * > > Get_SourceParticles(bool locEntireChainFlag=false, Charge_t locCharge=d_AllCharges) const
Definition: DSourceCombo.h:290
map< pair< bool, vector< const DKinematicData * > >, DVector3 > dVertexMap
const DAnalysisUtilities * dAnalysisUtilities
map< tuple< Particle_t, const DSourceCombo *, bool, const DSourceCombo *, const DKinematicData * >, const DKinematicData * > dReconDecayParticles_FromMissing
signed char Get_VertexZBin(bool locIsProductionVertex, const DSourceCombo *locReactionCombo, const DSourceCombo *locPrimaryVertexCombo, const DKinematicData *locBeamParticle, bool locIsCombo2ndVertex) const
unordered_map< const DReactionStepVertexInfo *, bool > dVertexDeterminableWithChargedMap
void Set_Vertex(const DVertex *locVertex)
map< tuple< Particle_t, bool, const DSourceCombo *, const DKinematicData * >, const DKinematicData * > dReconDecayParticles_FromProducts
vector< const DKinematicData * > Get_ConstrainingParticles(bool locIsProductionVertex, const DSourceCombo *locReactionCombo, const DSourceCombo *locVertexCombo, const DKinematicData *locBeamParticle, bool locIsCombo2ndVertex) const
DSourceComboP4Handler * dSourceComboP4Handler
const DVector3 & momentum(void) const
DVector3 Get_Vertex_NoBeam(bool locIsProductionVertex, const DSourceCombo *locVertexCombo, bool locIsCombo2ndVertex) const
bool Get_IsVertexKnown(bool locIsProductionVertex, const DSourceCombo *locReactionCombo, const DSourceCombo *locVertexCombo, const DKinematicData *locBeamParticle, bool locIsCombo2ndVertex) const
bool Get_IsVertexKnown_NoBeam(bool locIsProductionVertex, const DSourceCombo *locVertexCombo, bool locIsCombo2ndVertex) const
map< tuple< bool, const DSourceCombo *, const DKinematicData * >, unordered_map< const DSourceCombo *, double > > dTimeOffsets
DResourcePool< DKinematicData > dResourcePool_KinematicData
vector< const DKinematicData * >::const_iterator Get_ThetaNearest90Iterator(const vector< const DKinematicData * > &locParticles)
void Recycle(const DType *locResource)