Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DReactionStep.h
Go to the documentation of this file.
1 #ifndef _DReactionStep_
2 #define _DReactionStep_
3 
4 #include <vector>
5 #include <memory>
6 #include <string>
7 #include <iostream>
8 #include <algorithm>
9 #include <stdlib.h>
10 #include <functional>
11 #include <numeric>
12 
13 #include "particleType.h"
14 
15 using namespace std;
16 
17 //--------------------------------------------------------------------------------------
18 // At one point, the Apple compiler did not support std::accumulate so the block of code
19 // below was added. Recent versions of Xcode include a compiler that does support this
20 // so the code below is now disabled. If you are compiling on Mac OS X with an older
21 // compiler that complains that accumulate is not defined, then you'll need to comment
22 // out the next line and uncomment the one following it. 5/16/2018 DL
23 #if 0
24 //#ifdef __APPLE__
25 namespace std{
26 template <class InputIterator, class T>
27  T accumulate (InputIterator first, InputIterator last, T init)
28 {
29  while (first!=last) {
30  init = init + *first;
31  ++first;
32  }
33  return init;
34 }
35 template <class InputIterator, class T, class BinaryOperation>
36  T accumulate (InputIterator first, InputIterator last, T init, BinaryOperation binary_op)
37 {
38  while (first!=last) {
39  binary_op(init,*first);
40  ++first;
41  }
42  return init;
43 }
44 };
45 #endif
46 //--------------------------------------------------------------------------------------
47 
48 namespace DAnalysis
49 {
50 
52 {
53  public:
54 
55  //CONSTRUCTORS
56 
57  //if fixed-target production or rescattering:
58  DReactionStep(Particle_t locScatteringPID, Particle_t locTargetPID, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID = Unknown,
59  bool locInclusiveFlag = false, bool locBeamMissingFlag = false);
60  //if 2 beams: collider experiment
61  DReactionStep(pair<Particle_t, Particle_t> locBeamPIDs, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID = Unknown,
62  bool locInclusiveFlag = false, bool locFirstBeamMissingFlag = false, bool locSecondBeamMissingFlag = false);
63  //decaying particle
64  DReactionStep(Particle_t locDecayingPID, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID = Unknown, bool locInclusiveFlag = false);
65  // default
66  DReactionStep(void); //DEPRECATED
67 
68  void Reset(void);
69 
70  // OPERATORS
71  DReactionStep& operator=(const DReactionStep& locSourceData);
72  bool operator<(const DReactionStep& locStep) const{return *dReactionStepInfo < *(locStep.dReactionStepInfo);} //ignores dKinFitConstrainInitMassFlag!!!!
73 
74  // MANUALLY SET PIDs: //DEPRECATED
75  void Set_InitialParticleID(Particle_t locPID, bool locIsMissingFlag = false);
76  void Set_TargetParticleID(Particle_t locPID){dReactionStepInfo->dTargetPID = locPID;}
77  void Add_FinalParticleID(Particle_t locPID, bool locIsMissingFlag = false);
78  void Set_KinFitConstrainInitMassFlag(bool locFlag){dReactionStepInfo->dKinFitConstrainInitMassFlag = locFlag;}
79 
80  // GET INITIAL, TARGET, AND MISSING FINAL PIDs:
81  Particle_t Get_InitialPID(void) const{return dReactionStepInfo->dInitialPID;}
82  Particle_t Get_SecondBeamPID(void) const{return dReactionStepInfo->dSecondBeamPID;}
83  Particle_t Get_TargetPID(void) const{return dReactionStepInfo->dTargetPID;}
84  Particle_t Get_MissingPID(void) const;
85 
86  // GET FINAL PARTICLE PIDs:
87  Particle_t Get_FinalPID(size_t locIndex) const{return dReactionStepInfo->dFinalPIDs[locIndex];}
88  Particle_t Get_PID(int locParticlceIndex) const;
89  size_t Get_NumFinalPIDs(void) const{return dReactionStepInfo->dFinalPIDs.size();}
90  vector<Particle_t> Get_FinalPIDs(bool locIncludeMissingFlag = true, Charge_t locCharge = d_AllCharges, bool locIncludeDuplicatesFlag = true) const;
91 
92  //GET CONTROL FLAGS
93  int Get_MissingParticleIndex(void) const{return dReactionStepInfo->dMissingParticleIndex;}
94  bool Get_IsInclusiveFlag(void) const{return (dReactionStepInfo->dMissingParticleIndex == DReactionStep::Get_ParticleIndex_Inclusive());}
95  bool Get_IsBeamMissingFlag(void) const{return (dReactionStepInfo->dMissingParticleIndex == DReactionStep::Get_ParticleIndex_Initial());}
96  bool Get_IsSecondBeamMissingFlag(void) const{return (dReactionStepInfo->dMissingParticleIndex == DReactionStep::Get_ParticleIndex_SecondBeam());}
97  bool Get_KinFitConstrainInitMassFlag(void) const{return dReactionStepInfo->dKinFitConstrainInitMassFlag;}
98 
99  //definitions of negative values for any particle index //in order such that operator< returns order expected for string (e.g. gp->...)
100  constexpr static int Get_ParticleIndex_None(void){return -5;}
101  constexpr static int Get_ParticleIndex_Inclusive(void){return -4;}
102  constexpr static int Get_ParticleIndex_Initial(void){return -3;}
103  constexpr static int Get_ParticleIndex_SecondBeam(void){return -2;}
104  constexpr static int Get_ParticleIndex_Target(void){return -1;}
105 
106  private:
107 
108  int Prepare_InfoArguments(vector<Particle_t>& locFinalPIDs, Particle_t locMissingFinalPID, bool locInclusiveFlag, bool locBeamMissingFlag, bool locSecondBeamMissingFlag) const;
109  static void Check_IsResonance(Particle_t locPID);
110 
112  {
113  //CONSTRUCTORS
114  DReactionStepInfo(void) = default;
115  DReactionStepInfo(Particle_t locInitialPID, Particle_t locSecondBeamPID, Particle_t locTargetPID, vector<Particle_t> locFinalPIDs,
116  int locMissingParticleIndex = -1);
117 
118  bool operator<(const DReactionStepInfo& locInfo) const; //ignores dKinFitConstrainInitMassFlag!!! //not ideal, but in a hurry: for DParticleComboCreator::dComboStepMap
119 
120  // PID MEMBERS:
121  Particle_t dInitialPID = Unknown; //e.g. lambda, gamma
122  Particle_t dSecondBeamPID = Unknown; //second beam, Unknown if not present
123  Particle_t dTargetPID = Unknown; //unknown if none
124  vector<Particle_t> dFinalPIDs; //if inclusive, there is no indication of it here!
125 
126  // CONTROL MEMBERS:
127  int dMissingParticleIndex = DReactionStep::Get_ParticleIndex_None(); //final state particle at this index is missing
128  bool dKinFitConstrainInitMassFlag = true; //default true, is ignored when not applicable (e.g. init is non-decaying (beam) particle)
129  };
130 
131  //memory of object in shared_ptr is managed automatically: deleted automatically when no references are left
132  //This is done because sometimes a new object is needed (e.g. DParticleComboStep) for which this info hasn't changed (from DReactionStep)
133  //Thus, just share this between the two objects, instead of doubling the memory usage
134  //By inheriting this class, you also get to share the same interface
135  shared_ptr<DReactionStepInfo> dReactionStepInfo;
136 };
137 
138 /********************************************************* NAMESPACE-SCOPE NON-INLINE FUNCTION DECLARATIONS **********************************************************/
139 
140 string Get_InitialParticlesName(const DReactionStep* locStep, bool locTLatexFlag);
141 bool Are_ParticlesIdentical(const DReactionStep* locStep1, const DReactionStep* locStep2, bool locExceptMissingUnknownInInputFlag);
142 vector<string> Get_FinalParticleNames(const DReactionStep* locStep, bool locIncludeMissingFlag, bool locTLatexFlag);
143 
144 /************************************************************** DREACTIONSTEP CONSTRUCTORS & OPERATORS ***************************************************************/
145 
146 //if fixed-target production or rescattering:
147 inline DReactionStep::DReactionStep(Particle_t locScatteringPID, Particle_t locTargetPID, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID,
148  bool locInclusiveFlag, bool locBeamMissingFlag)
149 {
150  //Prepare arguments
151  int locMissingParticleIndex = Prepare_InfoArguments(locNonMissingFinalPIDs, locMissingFinalPID, locInclusiveFlag, locBeamMissingFlag, false);
152  dReactionStepInfo = std::make_shared<DReactionStepInfo>(locScatteringPID, Unknown, locTargetPID, locNonMissingFinalPIDs, locMissingParticleIndex);
153 }
154 
155 //if 2 beams: collider experiment
156 inline DReactionStep::DReactionStep(pair<Particle_t, Particle_t> locBeamPIDs, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID,
157  bool locInclusiveFlag, bool locFirstBeamMissingFlag, bool locSecondBeamMissingFlag)
158 {
159  //Prepare arguments
160  int locMissingParticleIndex = Prepare_InfoArguments(locNonMissingFinalPIDs, locMissingFinalPID, locInclusiveFlag, locFirstBeamMissingFlag, locSecondBeamMissingFlag);
161  dReactionStepInfo = std::make_shared<DReactionStepInfo>(locBeamPIDs.first, locBeamPIDs.second, Unknown, locNonMissingFinalPIDs, locMissingParticleIndex);
162 }
163 
164 //decaying particle
165 inline DReactionStep::DReactionStep(Particle_t locDecayingPID, vector<Particle_t> locNonMissingFinalPIDs, Particle_t locMissingFinalPID, bool locInclusiveFlag)
166 {
167  //Prepare arguments
168  int locMissingParticleIndex = Prepare_InfoArguments(locNonMissingFinalPIDs, locMissingFinalPID, locInclusiveFlag, false, false);
169  dReactionStepInfo = std::make_shared<DReactionStepInfo>(locDecayingPID, Unknown, Unknown, locNonMissingFinalPIDs, locMissingParticleIndex);
170 }
171 
172 // default
173 inline DReactionStep::DReactionStep(void) : dReactionStepInfo(std::make_shared<DReactionStepInfo>())
174 {
175 }
176 
177 //operator=
179 {
180  //Replace current data with a new, independent copy of the input data: tracked separately from input so it can be modified
181  dReactionStepInfo = std::make_shared<DReactionStepInfo>(*(locSourceData.dReactionStepInfo));
182  return *this;
183 }
184 
185 inline void DReactionStep::Reset(void)
186 {
187  dReactionStepInfo = std::make_shared<DReactionStepInfo>();
188 }
189 
190 /*************************************************************** DREACTIONSTEPINFO FUNCTIONS *****************************************************************/
191 
192 inline DReactionStep::DReactionStepInfo::DReactionStepInfo(Particle_t locInitialPID, Particle_t locSecondBeamPID, Particle_t locTargetPID,
193  vector<Particle_t> locFinalPIDs, int locMissingParticleIndex) :
194  dInitialPID(locInitialPID), dSecondBeamPID(locSecondBeamPID), dTargetPID(locTargetPID), dFinalPIDs(locFinalPIDs),
195  dMissingParticleIndex(locMissingParticleIndex)
196 {
200  for(const auto& locPID : dFinalPIDs)
201  Check_IsResonance(locPID);
202 }
203 
205 {
206  if(dInitialPID != locInfo.dInitialPID)
207  return (dInitialPID < locInfo.dInitialPID);
208 
209  if(dSecondBeamPID != locInfo.dSecondBeamPID)
210  return (dSecondBeamPID < locInfo.dSecondBeamPID);
211 
212  if(dTargetPID != locInfo.dTargetPID)
213  return (dTargetPID < locInfo.dTargetPID);
214 
215  if(dMissingParticleIndex != locInfo.dMissingParticleIndex)
216  return (dMissingParticleIndex < locInfo.dMissingParticleIndex);
217 
218  return (dFinalPIDs < locInfo.dFinalPIDs);
219 }
220 
221 /************************************************************** DREACTIONSTEP INLINE FUNCTIONS ***************************************************************/
222 
224 {
225  if(IsResonance(locPID))
226  {
227  cout << "ERROR: CANNOT SET RESONANCE PID. ABORTING." << endl;
228  abort();
229  }
230 }
231 
232 inline void DReactionStep::Set_InitialParticleID(Particle_t locPID, bool locIsMissingFlag)
233 {
234  Check_IsResonance(locPID);
235  dReactionStepInfo->dInitialPID = locPID;
236  if(locIsMissingFlag)
238 }
239 
241 {
242  return (dReactionStepInfo->dMissingParticleIndex < 0) ? Unknown : dReactionStepInfo->dFinalPIDs[dReactionStepInfo->dMissingParticleIndex];
243 }
244 
245 inline Particle_t DReactionStep::Get_PID(int locParticlceIndex) const
246 {
247  if(locParticlceIndex >= 0)
248  return Get_FinalPID(locParticlceIndex);
249 
250  switch (locParticlceIndex)
251  {
252  case Get_ParticleIndex_Initial(): return dReactionStepInfo->dInitialPID;
253  case Get_ParticleIndex_SecondBeam(): return dReactionStepInfo->dSecondBeamPID;
254  case Get_ParticleIndex_Target(): return dReactionStepInfo->dTargetPID;
255  default: return Unknown;
256  }
257 }
258 
259 /****************************************************** NAMESPACE-SCOPE INLINE FUNCTIONS *******************************************************/
260 
261 inline size_t Get_NumFinalPIDs(const DReactionStep* locStep, Particle_t locInputPID, bool locIncludeMissingFlag)
262 {
263  auto locFinalPIDs = locStep->Get_FinalPIDs(locIncludeMissingFlag);
264  auto locComparator = [&locInputPID](size_t locTotal, Particle_t locPID) -> size_t {return ((locPID == locInputPID) ? locTotal + 1 : locTotal);};
265  return std::accumulate(locFinalPIDs.begin(), locFinalPIDs.end(), size_t(0), locComparator);
266 }
267 
268 inline string Get_FinalParticlesName(const DReactionStep* locStep, bool locIncludeMissingFlag, bool locTLatexFlag)
269 {
270  auto locNames = Get_FinalParticleNames(locStep, locIncludeMissingFlag, locTLatexFlag);
271  if(locTLatexFlag)
272  return std::accumulate(locNames.begin(), locNames.end(), string(""));
273  else
274  {
275  auto locRetriever = [](const string& locTotal, const string& locString) -> string {return locTotal + string("_") + locString;};
276  return std::accumulate(std::next(locNames.begin()), locNames.end(), locNames.front(), locRetriever);
277  }
278 }
279 
280 inline string Get_StepName(const DReactionStep* locStep, bool locIncludeMissingFlag, bool locTLatexFlag)
281 {
282  auto locStepName = Get_InitialParticlesName(locStep, locTLatexFlag);
283  locStepName += locTLatexFlag ? "#rightarrow" : "__";
284  locStepName += Get_FinalParticlesName(locStep, locIncludeMissingFlag, locTLatexFlag);
285  return locStepName;
286 }
287 
288 inline void Print_ReactionStep(const DReactionStep* locReactionStep)
289 {
290  cout << DAnalysis::Get_StepName(locReactionStep, true, false) << endl;
291 }
292 
293 inline int Get_ParticleIndex(const DReactionStep* locStep, Particle_t locInputPID, size_t locInstance)
294 {
295  //locInstance starts from 1!!!
296  auto locFinalPIDs = locStep->Get_FinalPIDs(false); //exclude missing
297  size_t locCount = 0;
298  auto Instance_Finder = [locInputPID, locInstance, &locCount](Particle_t locPID) -> bool
299  {return (locPID != locInputPID) ? false : (++locCount == locInstance) ? true : false;};
300  auto locIterator = std::find_if(locFinalPIDs.begin(), locFinalPIDs.end(), Instance_Finder);
301  return (locIterator == locFinalPIDs.end()) ? DReactionStep::Get_ParticleIndex_None() : std::distance(locFinalPIDs.begin(), locIterator);
302 }
303 
305 {
306  auto locMissingParticleIndex = locStep->Get_MissingParticleIndex();
307  return (locStep->Get_IsInclusiveFlag() || (locMissingParticleIndex >= 0));
308 }
309 
310 inline bool Check_ChannelEquality(const DReactionStep* lhs, const DReactionStep* rhs, bool locSameOrderFlag = true, bool locRightSubsetOfLeftFlag = false)
311 {
312  if((lhs->Get_InitialPID() != rhs->Get_InitialPID()) || (lhs->Get_SecondBeamPID() != rhs->Get_SecondBeamPID()) || (lhs->Get_TargetPID() != rhs->Get_TargetPID()))
313  return false;
314 
315  if(locSameOrderFlag)
316  return ((lhs->Get_MissingParticleIndex() == rhs->Get_MissingParticleIndex()) && (lhs->Get_FinalPIDs() == rhs->Get_FinalPIDs()));
317 
318  //MISSING PARTICLE INDEX
319  if(lhs->Get_MissingPID() != rhs->Get_MissingPID())
320  return false;
321  auto locMissingIndex_lhs = lhs->Get_MissingParticleIndex();
322  auto locMissingIndex_rhs = rhs->Get_MissingParticleIndex();
323 
324  //if it's not, then ... one must be inclusive, the other none
325  if((lhs->Get_MissingPID() == Unknown) && (locMissingIndex_lhs != locMissingIndex_rhs))
326  {
327  if(!locRightSubsetOfLeftFlag)
328  return false;
329  if((locMissingIndex_lhs != DReactionStep::Get_ParticleIndex_Inclusive()) || (locMissingIndex_rhs != DReactionStep::Get_ParticleIndex_None()))
330  return false; //e.g. missing beam vs. inclusive
331  }
332 
333  // GET FINAL PARTICLE PIDs:
334  auto locPIDS_lhs = lhs->Get_FinalPIDs();
335  auto locPIDS_rhs = rhs->Get_FinalPIDs();
336  std::sort(locPIDS_lhs.begin(), locPIDS_lhs.end());
337  std::sort(locPIDS_rhs.begin(), locPIDS_rhs.end());
338  if(locRightSubsetOfLeftFlag)
339  {
340  decltype(locPIDS_lhs) locPIDDifference{};
341  locPIDDifference.reserve(locPIDS_rhs.size());
342  std::set_difference(locPIDS_rhs.begin(), locPIDS_rhs.end(), locPIDS_lhs.begin(), locPIDS_lhs.end(), std::back_inserter(locPIDDifference)); //right - left better be empty
343  return locPIDDifference.empty();
344  }
345  if(locPIDS_lhs.size() != locPIDS_rhs.size())
346  return false;
347 
348  return std::is_permutation(locPIDS_lhs.begin(), locPIDS_lhs.end(), locPIDS_rhs.begin());
349 }
350 
351 } //end DAnalysis namespace
352 
353 #endif // _DReactionStep_
354 
int Get_ParticleIndex(const DReactionStep *locStep, Particle_t locInputPID, size_t locInstance)
static constexpr int Get_ParticleIndex_Initial(void)
bool operator<(const DSourceComboUse &lhs, const DSourceComboUse &rhs)
Definition: DSourceCombo.h:95
bool Check_ChannelEquality(const DReaction *lhs, const DReaction *rhs, bool locSameOrderFlag=true, bool locRightSubsetOfLeftFlag=false)
Definition: DReaction.h:185
bool Get_KinFitConstrainInitMassFlag(void) const
Definition: DReactionStep.h:97
bool operator<(const DReactionStepInfo &locInfo) const
bool Are_ParticlesIdentical(const DReactionStep *locStep1, const DReactionStep *locStep2, bool locExceptMissingUnknownInInputFlag)
char string[256]
static unsigned short int IsResonance(Particle_t p)
Definition: particleType.h:802
void Set_KinFitConstrainInitMassFlag(bool locFlag)
Definition: DReactionStep.h:78
Particle_t Get_TargetPID(void) const
Definition: DReactionStep.h:83
string Get_StepName(const DReactionStep *locStep, bool locIncludeMissingFlag, bool locTLatexFlag)
Charge_t
bool Get_HasMissingParticle_FinalState(const DReactionStep *locStep)
vector< string > Get_FinalParticleNames(const DReactionStep *locStep, bool locIncludeMissingFlag, bool locTLatexFlag)
DReactionStep & operator=(const DReactionStep &locSourceData)
Particle_t Get_SecondBeamPID(void) const
Definition: DReactionStep.h:82
void Print_ReactionStep(const DReactionStep *locReactionStep)
Particle_t Get_InitialPID(void) const
Definition: DReactionStep.h:81
Particle_t Get_PID(int locParticlceIndex) const
shared_ptr< DReactionStepInfo > dReactionStepInfo
Particle_t Get_FinalPID(size_t locIndex) const
Definition: DReactionStep.h:87
static constexpr int Get_ParticleIndex_Inclusive(void)
static constexpr int Get_ParticleIndex_None(void)
string Get_FinalParticlesName(const DReactionStep *locStep, bool locIncludeMissingFlag, bool locTLatexFlag)
int Get_MissingParticleIndex(void) const
Definition: DReactionStep.h:93
bool operator<(const DReactionStep &locStep) const
Definition: DReactionStep.h:72
static constexpr int Get_ParticleIndex_Target(void)
static void Check_IsResonance(Particle_t locPID)
size_t Get_NumFinalPIDs(void) const
Definition: DReactionStep.h:89
void Set_InitialParticleID(Particle_t locPID, bool locIsMissingFlag=false)
vector< Particle_t > Get_FinalPIDs(bool locIncludeMissingFlag=true, Charge_t locCharge=d_AllCharges, bool locIncludeDuplicatesFlag=true) const
static constexpr int Get_ParticleIndex_SecondBeam(void)
bool Get_IsInclusiveFlag(void) const
Definition: DReactionStep.h:94
bool Get_IsSecondBeamMissingFlag(void) const
Definition: DReactionStep.h:96
string Get_InitialParticlesName(const DReactionStep *locStep, bool locTLatexFlag)
bool Get_IsBeamMissingFlag(void) const
Definition: DReactionStep.h:95
void Set_TargetParticleID(Particle_t locPID)
Definition: DReactionStep.h:76
Particle_t Get_MissingPID(void) const
Particle_t
Definition: particleType.h:12