UParticle.h 8.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234
  1. #ifndef UPARTICLE_H
  2. #define UPARTICLE_H
  3. // C++ headers
  4. #include <limits>
  5. // ROOT headers
  6. #include "TObject.h"
  7. #include "TLorentzVector.h"
  8. #include "TMath.h"
  9. #include "TDatabasePDG.h"
  10. #include "TParticlePDG.h"
  11. #include "TParticle.h"
  12. //_________________
  13. class UParticle : public TObject {
  14. public:
  15. /// Default constructor
  16. UParticle();
  17. /// Constructor that takes parameters
  18. UParticle(Int_t index, Int_t pdg, Int_t status,
  19. Int_t parent, Int_t parentDecay,
  20. Int_t mate, Int_t decay, Int_t child[2],
  21. Double_t px, Double_t py, Double_t pz, Double_t e,
  22. Double_t x, Double_t y, Double_t z, Double_t t,
  23. Double_t weight);
  24. /// Another constructor with parameters
  25. UParticle(Int_t index, Int_t pdg, Int_t status,
  26. Int_t parent, Int_t parentDecay,
  27. Int_t mate, Int_t decay, Int_t child[2],
  28. TLorentzVector mom, TLorentzVector pos,
  29. Double_t weight);
  30. /// Copy constructor that takes UParticle
  31. UParticle(const UParticle& right);
  32. /// Copy constructor that takes TParticle
  33. UParticle(const TParticle& right);
  34. /// Destructor
  35. virtual ~UParticle();
  36. /// Assignement operator
  37. const UParticle& operator=(const UParticle& right);
  38. /// Assignement operator
  39. const UParticle& operator=(const TParticle& right);
  40. /// Comparison operator
  41. const Bool_t operator==(const UParticle& right) const;
  42. void print(Option_t* option = "");
  43. //
  44. // Getters
  45. //
  46. /// Return particle index
  47. Int_t index() const { return (Int_t)fIndex; }
  48. /// Return PDG code
  49. Int_t pdg() const { return fPdg; }
  50. /// Return particle status
  51. Int_t status() const { return (Int_t)fStatus; }
  52. /// Return parent index
  53. Int_t parent() const { return (Int_t)fParent; }
  54. /// Return parent decay index
  55. Int_t parentDecay() const { return (Int_t)fParentDecay; }
  56. /// Return index of the last collision partner
  57. Int_t mate() const { return (Int_t)fMate; }
  58. /// Return decay index (-1 if not decayed)
  59. Int_t decay() const { return (Int_t)fDecay; }
  60. /// Return index of the first child
  61. Int_t firstChild() const { return (Int_t)fChild[0]; }
  62. /// Return index of the second child
  63. Int_t lastChild() const { return (Int_t)fChild[1]; }
  64. /// Return px (GeV/c)
  65. Double_t px() const { return (Double_t)fPx; }
  66. /// Return py (GeV/c)
  67. Double_t py() const { return (Double_t)fPy; }
  68. /// Return pz (GeV/c)
  69. Double_t pz() const { return (Double_t)fPz; }
  70. /// Return p (GeV/c)
  71. Double_t ptot() const
  72. { return TMath::Sqrt( px()*px() + py()*py() + pz()*pz() ); }
  73. /// Return transverse momentum (pT)
  74. Double_t pt() const
  75. { return TMath::Sqrt( px()*px() + py()*py() ); }
  76. /// Return mass (GeV/c^2)
  77. Double_t mass() const
  78. { return TDatabasePDG::Instance()->GetParticle( fPdg )->Mass(); }
  79. /// Return charge
  80. Double_t charge() const
  81. { return TDatabasePDG::Instance()->GetParticle( fPdg )->Charge(); }
  82. /// Return energy of the particle (GeV)
  83. Double_t energy() const
  84. { return TMath::Sqrt( ptot()*ptot() + mass()*mass() ); }
  85. /// Return energy (GeV)
  86. Double_t e() const { return energy(); }
  87. /// Return pseudorapidity
  88. Double_t eta() const { return momentum().Eta(); }
  89. /// Return pseudorapidity
  90. Double_t pseudoRapidity() const { return momentum().Eta(); }
  91. /// Return four-momentum (px,py,pz,E)
  92. TLorentzVector momentum() const
  93. { return TLorentzVector( fPx, fPy, fPz, e() ); }
  94. /// Set four-momentum to the mom vector
  95. void momentum(TLorentzVector& mom) const
  96. { mom.SetXYZM( fPx, fPy, fPz, mass() ); }
  97. /// Return x position (fm)
  98. Double_t x() const { return (Double_t)fX; }
  99. /// Return y position (fm)
  100. Double_t y() const { return (Double_t)fY; }
  101. /// Return z position (fm)
  102. Double_t z() const { return (Double_t)fZ; }
  103. /// Return t position (fm)
  104. Double_t t() const { return (Double_t)fT; }
  105. /// Return four-coordinate (x,y,z,t)
  106. TLorentzVector position() const
  107. { return TLorentzVector( fX, fY, fZ, fT ); }
  108. /// Set four-coordinate to the pos vector
  109. void position(TLorentzVector& pos) const
  110. { pos.SetXYZT( fX, fY, fZ, fT); }
  111. /// Return weight
  112. Double_t weight() const { return fWeight; }
  113. //
  114. // Setters
  115. //
  116. /// Set particle index
  117. void setIndex(const Int_t& index)
  118. { fIndex = ( (index > std::numeric_limits<unsigned short>::max() ) ?
  119. std::numeric_limits<unsigned short>::max() : (UShort_t)index ); }
  120. /// Set PdgId (pdg code)
  121. void setPdg(const Int_t& pdg) {fPdg = pdg;}
  122. /// Set status
  123. void setStatus(const Int_t& status)
  124. { if ( status <= std::numeric_limits<char>::min() ) { fStatus = std::numeric_limits<char>::min(); }
  125. else if ( status >= std::numeric_limits<char>::max() ) { fStatus = std::numeric_limits<char>::max(); }
  126. else { fStatus = (Char_t)status; } }
  127. /// Set parent index
  128. void setParent(const Int_t& parent)
  129. { fParent = ( ( parent > std::numeric_limits<unsigned short>::max() ) ?
  130. std::numeric_limits<unsigned short>::max() : (UShort_t)parent ); }
  131. /// Set parent decay index
  132. void setParentDecay(const Int_t& parentDecay)
  133. { fParentDecay = ( ( parentDecay > std::numeric_limits<unsigned short>::max() ) ?
  134. std::numeric_limits<unsigned short>::max() : (UShort_t)parentDecay ); }
  135. /// Set index of the last collision partner
  136. void setMate(const Int_t& mate)
  137. { fMate = ( (mate > std::numeric_limits<unsigned short>::max() ) ?
  138. std::numeric_limits<unsigned short>::max() : (UShort_t)mate ); }
  139. /// Set decay index (-1 if not decayed)
  140. void setDecay(const Int_t& decay)
  141. { fDecay = ( ( TMath::Abs(decay) > std::numeric_limits<short>::max() ) ?
  142. std::numeric_limits<short>::max() : (Short_t)decay ); }
  143. /// Set 2 childer indeces
  144. void setChild(Int_t child[2])
  145. { setFirstChild( child[0] ); setLastChild( child[1] ); }
  146. /// Set index of the first child
  147. void setFirstChild(const Int_t& child)
  148. { fChild[0] = ( (child > std::numeric_limits<unsigned short>::max() ) ?
  149. std::numeric_limits<unsigned short>::max() : (UShort_t)child ); }
  150. /// Set index of the second child
  151. void setLastChild(const Int_t& child)
  152. { fChild[1] = ( (child > std::numeric_limits<unsigned short>::max() ) ?
  153. std::numeric_limits<unsigned short>::max() : (UShort_t)child ); }
  154. /// Set px (GeV/c)
  155. void setPx(const Double_t& px) {fPx = (Float_t)px; }
  156. /// Set py (GeV/c)
  157. void setPy(const Double_t& py) {fPy = (Float_t)py; }
  158. /// Set pz (GeV/c)
  159. void setPz(const Double_t& pz) {fPz = (Float_t)pz; }
  160. /// Set energy (GeV). IMPORTANT: This is a dummy method.
  161. void setE(const Double_t& e) { /* fE = (Float_t)e; */}
  162. /// Set four-momentum (px,py,pz,E)
  163. void setMomentum(const Double_t& px, const Double_t& py,
  164. const Double_t& pz, const Double_t& e)
  165. { fPx = (Float_t)px; fPy = (Float_t)py; fPz = (Float_t)pz; /* fE = (Float_t)e; */ }
  166. /// Set four-momentum (TLorentzVector)
  167. void setMomentum(const TLorentzVector& mom)
  168. { fPx=(Float_t)mom.Px(); fPy=(Float_t)mom.Py(); fPz=(Float_t)mom.Pz(); /* fE=(Float_t)mom.E(); */ }
  169. /// Set x coordinate (fm)
  170. void setX(const Double_t& x) { fX = (Float_t)x; }
  171. /// Set y coordinate (fm)
  172. void setY(const Double_t& y) { fY = (Float_t)y; }
  173. /// Set z coordinate (fm)
  174. void setZ(const Double_t& z) { fZ = (Float_t)z; }
  175. /// Set t coordinate (fm/c)
  176. void setT(const Double_t& t) { fT = (Float_t)t; }
  177. /// Set four-coordinate (x,y,z,t)
  178. void setPosition(const Double_t& x, const Double_t& y,
  179. const Double_t& z, const Double_t& t)
  180. { fX = (Float_t)x; fY = (Float_t)y; fZ = (Float_t)z; fT = (Float_t)t; }
  181. /// Set four-coordinate (TLorentzVector)
  182. void setPosition(const TLorentzVector& pos)
  183. { fX=(Float_t)pos.X(); fY=(Float_t)pos.Y(); fZ=(Float_t)pos.Z(); fT=(Float_t)pos.T(); }
  184. /// Set weight
  185. void setWeight(const Double_t& weight) { fWeight = (Float_t)weight; }
  186. private:
  187. /// Index of this particle
  188. UShort_t fIndex;
  189. /// PDG code
  190. Int_t fPdg;
  191. /// Status
  192. Char_t fStatus;
  193. /// Index of the parent
  194. UShort_t fParent;
  195. /// Parent decay index
  196. UShort_t fParentDecay;
  197. /// Index of the last collision partner
  198. UShort_t fMate;
  199. /// Decay index (-1 if not decayed)
  200. Short_t fDecay;
  201. /// Index of the first and the last child
  202. UShort_t fChild[2];
  203. /// px (GeV/c)
  204. Float_t fPx;
  205. /// py (GeV/c)
  206. Float_t fPy;
  207. /// pz (GeV/c)
  208. Float_t fPz;
  209. /// x (fm)
  210. Float_t fX;
  211. /// y (fm)
  212. Float_t fY;
  213. /// z (fm)
  214. Float_t fZ;
  215. /// t (fm/c)
  216. Float_t fT;
  217. /// Weight
  218. Float_t fWeight;
  219. ClassDef(UParticle, 3);
  220. };
  221. #endif