MpdFemtoPair.h 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701
  1. /**
  2. * \class MpdFemtoPair
  3. * \brief Holds information about pair of particles
  4. *
  5. * The class describes a pair of MpdFemtoParticle instances and allows various
  6. * relative 4-momentum quantities, such as, qInv, qOut, qSide and qLong in different
  7. * frames (CMS - center of mass of the experiment, LCMS - longitudinally co-moving
  8. * system, PF - pair frame), as well as, sum of 4-momenta, such as, invariant mass,
  9. * kT, mT and so on.
  10. * It also provides esitmations for two-particle effects, for instance, track-merging
  11. * and track-splitiing effects.
  12. *
  13. * \author Grigory Nigmatkulov (NRNU MEPhI)
  14. * \date May 18, 2019
  15. * \email nigmatkulov@gmail.com
  16. */
  17. #ifndef MpdFemtoPair_h
  18. #define MpdFemtoPair_h
  19. // C++ headers
  20. #include <utility>
  21. #include <vector>
  22. // MpdFemtoMaker headers
  23. #include "MpdFemtoParticle.h"
  24. #include "MpdFemtoTypes.h"
  25. // ROOT headers
  26. #include "TVector3.h"
  27. #include "TLorentzVector.h"
  28. #include "TMath.h"
  29. //_________________
  30. class MpdFemtoPair {
  31. public:
  32. /// Default constructor
  33. MpdFemtoPair();
  34. /// Constructor with particles
  35. MpdFemtoPair(MpdFemtoParticle*, MpdFemtoParticle*);
  36. /// Copy constructor
  37. MpdFemtoPair(const MpdFemtoPair& copy);
  38. /// Copy constructor
  39. MpdFemtoPair& operator=(const MpdFemtoPair& copy);
  40. /// Default destructor
  41. virtual ~MpdFemtoPair();
  42. //
  43. // Getters
  44. //
  45. /// First particle from a pair
  46. MpdFemtoParticle* track1() const {
  47. return mTrack1;
  48. }
  49. /// Second particle from a pair
  50. MpdFemtoParticle* track2() const {
  51. return mTrack2;
  52. }
  53. /// Relative four-momentum
  54. TLorentzVector fourMomentumDiff() const {
  55. return ( mTrack1->fourMomentum() - mTrack2->fourMomentum());
  56. }
  57. /// Sum of four-momenta
  58. TLorentzVector fourMomentumSum() const {
  59. return ( mTrack1->fourMomentum() + mTrack2->fourMomentum());
  60. }
  61. /// Invariant relative momentum of the pair
  62. double qInv() const {
  63. return (-1.)*fourMomentumDiff().M();
  64. }
  65. /// Pair transverse momentum of the pair
  66. double pT() const {
  67. return p().Perp();
  68. }
  69. /// Half of pair transverse momentum of the pair
  70. double kT() const {
  71. return 0.5 * pT();
  72. }
  73. /// Invariant mass of the pair
  74. double mInv() const {
  75. return fourMomentumSum().M();
  76. }
  77. /// Three-momentum of the pair
  78. TVector3 momentum() const {
  79. return fourMomentumSum().Vect();
  80. }
  81. /// Three-momentum of the pair
  82. TVector3 p() const {
  83. return momentum();
  84. }
  85. /// Px of the pair
  86. double px() const {
  87. return p().X();
  88. }
  89. /// Py of the pair
  90. double py() const {
  91. return p().Y();
  92. }
  93. /// Pz of the pair
  94. double pz() const {
  95. return p().Z();
  96. }
  97. /// Momentum of the pair
  98. double ptot() const {
  99. return p().Mag();
  100. }
  101. /// Squared momentum of the pair
  102. double ptot2() const {
  103. return p().Mag2();
  104. }
  105. /// Energy of the pair
  106. double energy() const {
  107. return fourMomentumSum().Energy();
  108. }
  109. /// Rapidity of the pair
  110. double rap() const {
  111. return fourMomentumSum().Rapidity();
  112. }
  113. /// Rapidity of the pair
  114. double rapidity() const {
  115. return rap();
  116. }
  117. /// Pair emission angle
  118. double emissionAngle() const;
  119. /// Pseudorapidity of the pair
  120. double eta() const {
  121. return fourMomentumSum().Eta();
  122. }
  123. /// Pseudorapidity of the pair
  124. double pseudoRapidity() const {
  125. return eta();
  126. }
  127. /// Azimuthal angle of the pair
  128. double phi() const {
  129. return fourMomentumSum().Phi();
  130. }
  131. /// Delta phi between two particles
  132. double deltaPhi() const {
  133. return ( mTrack1->phi() - mTrack2->phi());
  134. }
  135. /// Delta eta between two particles
  136. double deltaEta() const {
  137. return ( mTrack1->eta() - mTrack2->eta());
  138. }
  139. /// Sqrt(dEta^2+dPhi^2)
  140. double rValue() const {
  141. return TMath::Sqrt(deltaEta() * deltaEta() + deltaPhi() * deltaPhi());
  142. }
  143. /// Bertsch-Pratt side momentum component in Pair Frame
  144. double qSidePf() const;
  145. /// Bertsch-Pratt out momentum component in Pair Frame
  146. double qOutPf() const;
  147. /// Bertsch-Pratt long momentum component in Pair Frame
  148. double qLongPf() const;
  149. /// Bertsch-Pratt side momentum component in Local CMS (longitudinally comoving) frame
  150. double qSideCMS() const;
  151. /// Bertsch-Pratt out momentum component in Local CMS (longitudinally comoving) frame
  152. double qOutCMS() const;
  153. /// Bertsch-Pratt long momentum component in Local CMS (longitudinally comoving) frame
  154. double qLongCMS() const;
  155. /// Bertsch-Pratt side momentum component in Local CMS (longitudinally comoving) frame
  156. /// for non-identical pairs
  157. double dKSide() const {
  158. if (mNonIdParNotCalculated) {
  159. calcNonIdPar();
  160. }
  161. return mDKSide;
  162. }
  163. /// Bertsch-Pratt out momentum component in Local CMS (longitudinally comoving) frame
  164. /// for non-identical pairs
  165. double dKOut() const {
  166. if (mNonIdParNotCalculated) {
  167. calcNonIdPar();
  168. }
  169. return mDKOut;
  170. }
  171. /// Bertsch-Pratt long momentum component in Local CMS (longitudinally comoving) frame
  172. /// for non-identical pairs
  173. double dKLong() const {
  174. if (mNonIdParNotCalculated) {
  175. calcNonIdPar();
  176. }
  177. return mDKLong;
  178. }
  179. /// Bertsch-Pratt side momentum component in a longitudinally boosted frame
  180. /// the argument is the beta of the longitudinal boost (default is 0.0, meaning lab frame)
  181. double qSideBf(double beta = 0.0) const;
  182. /// Bertsch-Pratt side momentum component in a longitudinally boosted frame
  183. /// the argument is the beta of the longitudinal boost (default is 0.0, meaning lab frame)
  184. double qOutBf(double beta = 0.0) const;
  185. /// Bertsch-Pratt side momentum component in a longitudinally boosted frame
  186. /// the argument is the beta of the longitudinal boost (default is 0.0, meaning lab frame)
  187. double qLongBf(double beta = 0.0) const;
  188. /// Yano-Koonin-Podgoretskii Parametrisation in
  189. /// source rest frame (usually lab frame)
  190. void qYKPCMS(double& qP, double& qT, double& q0) const;
  191. /// Yano-Koonin-Podgoretskii Parametrisation in
  192. /// longitudinal co-moving frame
  193. void qYKPLCMS(double& qP, double& qT, double& q0) const;
  194. /// Yano-Koonin-Podgoretskii Parametrisation in
  195. /// pair rest frame
  196. void qYKPPF(double& qP, double& qT, double& q0) const;
  197. /// Track-splitting quantity
  198. double quality() const;
  199. /// Track-splitting quantity
  200. double splittingLevel() const {
  201. return quality();
  202. }
  203. // The following two methods calculate the "nominal" separation of the tracks
  204. // at the inner field cage (EntranceSeparation) and when they exit the TPC,
  205. // which may be at the outer field cage, or at the endcaps.
  206. // "nominal" means that the tracks are assumed to start at (0,0,0). Making this
  207. // assumption is important for the Event Mixing -- it is not a mistake.
  208. /// Distance between particles at TPC exit points
  209. double nominalTpcExitSeparation() const;
  210. /// Distance between particles at TPC entrance points
  211. double nominalTpcEntranceSeparation() const;
  212. /// Average distance between particles within TPC
  213. double nominalTpcAverageSeparation() const;
  214. // Adopted calculation of Entrance/Exit/Average TPC separation to V0 daughters
  215. double tpcExitSeparationTrackV0Pos() const {
  216. return ( mTrack1->nominalTpcExitPoint() - mTrack2->tpcV0PosExitPoint()).Mag();
  217. }
  218. double tpcEntranceSeparationTrackV0Pos() const {
  219. return ( mTrack1->nominalTpcEntrancePoint() - mTrack2->tpcV0PosEntrancePoint()).Mag();
  220. }
  221. double tpcAverageSeparationTrackV0Pos() const;
  222. double tpcExitSeparationTrackV0Neg() const {
  223. return ( mTrack1->nominalTpcExitPoint() - mTrack2->tpcV0NegExitPoint()).Mag();
  224. }
  225. double tpcEntranceSeparationTrackV0Neg() const {
  226. return ( mTrack1->nominalTpcEntrancePoint() - mTrack2->tpcV0NegEntrancePoint()).Mag();
  227. }
  228. double tpcAverageSeparationTrackV0Neg() const;
  229. double tpcExitSeparationV0PosV0Pos() const {
  230. return ( mTrack1->tpcV0PosExitPoint() - mTrack2->tpcV0PosExitPoint()).Mag();
  231. }
  232. double tpcEntranceSeparationV0PosV0Pos() const {
  233. return ( mTrack1->tpcV0PosEntrancePoint() - mTrack2->tpcV0PosEntrancePoint()).Mag();
  234. }
  235. double tpcAverageSeparationV0PosV0Pos() const;
  236. double tpcExitSeparationV0PosV0Neg() const {
  237. return ( mTrack1->tpcV0PosExitPoint() - mTrack2->tpcV0NegExitPoint()).Mag();
  238. }
  239. double tpcEntranceSeparationV0PosV0Neg() const {
  240. return ( mTrack1->tpcV0PosEntrancePoint() - mTrack2->tpcV0NegEntrancePoint()).Mag();
  241. }
  242. double tpcAverageSeparationV0PosV0Neg() const;
  243. double tpcExitSeparationV0NegV0Pos() const {
  244. return ( mTrack1->tpcV0NegExitPoint() - mTrack2->tpcV0PosExitPoint()).Mag();
  245. }
  246. double tpcEntranceSeparationV0NegV0Pos() const {
  247. return ( mTrack1->tpcV0NegEntrancePoint() - mTrack2->tpcV0PosEntrancePoint()).Mag();
  248. }
  249. double tpcAverageSeparationV0NegV0Pos() const;
  250. double tpcExitSeparationV0NegV0Neg() const {
  251. return ( mTrack1->tpcV0NegExitPoint() - mTrack2->tpcV0NegExitPoint()).Mag();
  252. }
  253. double tpcEntranceSeparationV0NegV0Neg() const {
  254. return ( mTrack1->tpcV0NegEntrancePoint() - mTrack2->tpcV0NegEntrancePoint()).Mag();
  255. }
  256. double tpcAverageSeparationV0NegV0Neg() const;
  257. /// Invariant total momentum
  258. double pInv() const;
  259. /// k^* of non-identical particles
  260. double kStar() const {
  261. if (mNonIdParNotCalculated) {
  262. calcNonIdPar();
  263. }
  264. return mKStarCalc;
  265. }
  266. /// k^* of non-identical particles when one of particles was rotated
  267. double kStarFlipped() const;
  268. double cvk() const {
  269. if (mNonIdParNotCalculated) {
  270. calcNonIdPar();
  271. }
  272. return mCVK;
  273. }
  274. double cvkFlipped() const;
  275. /// Invariant momentum of a pair with the first particle
  276. /// flipped in xy-plane (px,py,pz) -> (-px,-py,pz)
  277. double qInvFlippedXY() const;
  278. /// Invariant momentum of a pair with one of the particles
  279. /// flipped in xy-plane (px,py,pz) -> (-px,-py,pz)
  280. double qInvRandomFlippedXY() const;
  281. /// Invariant momentum of a pair with the first particle
  282. /// rotated to opposite-hemisphere (px,py,pz) -> (-px,-py,-pz)
  283. double qInvFlippedXYZ() const;
  284. /// Invariant momentum of a pair with one of the particles
  285. /// rotated to opposite-hemisphere (px,py,pz) -> (-px,-py,-pz)
  286. double qInvRandomFlippedXYZ() const;
  287. /// Opening angle between two particles
  288. double openingAngle() const {
  289. return TMath::RadToDeg() * (mTrack1->p().Angle(mTrack2->p()));
  290. }
  291. /// Side component of k^*
  292. double kStarSide() const {
  293. if (mNonIdParNotCalculated) {
  294. calcNonIdPar();
  295. }
  296. return mDKSide;
  297. }
  298. /// Out component of k^*
  299. double kStarOut() const {
  300. if (mNonIdParNotCalculated) {
  301. calcNonIdPar();
  302. }
  303. return mDKOut;
  304. }
  305. /// Long component of k^*
  306. double kStarLong() const {
  307. if (mNonIdParNotCalculated) {
  308. calcNonIdPar();
  309. }
  310. return mDKLong;
  311. }
  312. double electronPairProbability() const {
  313. return (mTrack1->track()->pidProbElectron()) * (mTrack2->track()->pidProbElectron());
  314. }
  315. double pionPairProbability() const {
  316. return (mTrack1->track()->pidProbPion()) * (mTrack2->track()->pidProbPion());
  317. }
  318. double kaonPairProbability() const {
  319. return (mTrack1->track()->pidProbKaon()) * (mTrack2->track()->pidProbKaon());
  320. }
  321. double protonPairProbability() const {
  322. return (mTrack1->track()->pidProbProton()) * (mTrack2->track()->pidProbProton());
  323. }
  324. double kaonPionPairProbability() const {
  325. return (mTrack1->track()->pidProbKaon()) * (mTrack2->track()->pidProbPion());
  326. }
  327. double dcaInsideTpc() const;
  328. double quality2() const;
  329. double kStarGlobal() const {
  330. if (mNonIdParNotCalculatedGlobal) {
  331. calcNonIdParGlobal();
  332. }
  333. return mKStarCalcGlobal;
  334. }
  335. double cvkGlobal() const {
  336. if (mNonIdParNotCalculatedGlobal) {
  337. calcNonIdParGlobal();
  338. }
  339. return mCVKGlobal;
  340. }
  341. double kStarSideGlobal() const {
  342. if (mNonIdParNotCalculatedGlobal) {
  343. calcNonIdParGlobal();
  344. }
  345. return mDKSideGlobal;
  346. }
  347. double kStarOutGlobal() const {
  348. if (mNonIdParNotCalculatedGlobal) {
  349. calcNonIdParGlobal();
  350. }
  351. return mDKOutGlobal;
  352. }
  353. double kStarLongGlobal() const {
  354. if (mNonIdParNotCalculatedGlobal) {
  355. calcNonIdParGlobal();
  356. }
  357. return mDKLongGlobal;
  358. }
  359. double fractionOfMergedRow() const {
  360. if (mMergingParNotCalculated) {
  361. calcMergingPar();
  362. }
  363. return mFracOfMergedRow;
  364. }
  365. double closestRowAtDCA() const {
  366. if (mMergingParNotCalculated) {
  367. calcMergingPar();
  368. }
  369. return mClosestRowAtDCA;
  370. }
  371. double weightedAvSep() const {
  372. if (mMergingParNotCalculated) {
  373. calcMergingPar();
  374. }
  375. return mWeightedAvSep;
  376. }
  377. double fractionOfMergedRowTrkV0Pos() const {
  378. if (mMergingParNotCalculatedTrkV0Pos) {
  379. calcMergingParFctn(&mMergingParNotCalculatedTrkV0Pos,
  380. mTrack1->z(), mTrack1->u(),
  381. mTrack2->z(), mTrack2->u(),
  382. mTrack1->sect(), mTrack2->sect(),
  383. &mFracOfMergedRowTrkV0Pos,
  384. &mClosestRowAtDCATrkV0Pos);
  385. }
  386. return mFracOfMergedRowTrkV0Pos;
  387. }
  388. double closestRowAtDCATrkV0Pos() const {
  389. if (mMergingParNotCalculatedTrkV0Pos) {
  390. calcMergingParFctn(&mMergingParNotCalculatedTrkV0Pos,
  391. mTrack1->z(), mTrack1->u(),
  392. mTrack2->z(), mTrack2->u(),
  393. mTrack1->sect(), mTrack2->sect(),
  394. &mFracOfMergedRowTrkV0Pos,
  395. &mClosestRowAtDCATrkV0Pos);
  396. }
  397. return mClosestRowAtDCATrkV0Pos;
  398. }
  399. double fractionOfMergedRowTrkV0Neg() const {
  400. if (mMergingParNotCalculatedTrkV0Neg) {
  401. calcMergingParFctn(&mMergingParNotCalculatedTrkV0Neg,
  402. mTrack1->z(), mTrack1->u(),
  403. mTrack2->v0NegZ(), mTrack2->v0NegU(),
  404. mTrack1->sect(), mTrack2->v0NegSect(),
  405. &mFracOfMergedRowTrkV0Neg,
  406. &mClosestRowAtDCATrkV0Neg);
  407. }
  408. return mFracOfMergedRowTrkV0Neg;
  409. }
  410. double closestRowAtDCATrkV0Neg() const {
  411. if (mMergingParNotCalculatedTrkV0Neg) {
  412. calcMergingParFctn(&mMergingParNotCalculatedTrkV0Neg,
  413. mTrack1->z(), mTrack1->u(),
  414. mTrack2->v0NegZ(), mTrack2->v0NegU(),
  415. mTrack1->sect(), mTrack2->v0NegSect(),
  416. &mFracOfMergedRowTrkV0Neg,
  417. &mClosestRowAtDCATrkV0Neg);
  418. }
  419. return mClosestRowAtDCATrkV0Neg;
  420. }
  421. double fractionOfMergedRowV0PosV0Neg() const {
  422. if (mMergingParNotCalculatedV0PosV0Neg) {
  423. calcMergingParFctn(&mMergingParNotCalculatedV0PosV0Neg,
  424. mTrack1->z(), mTrack1->u(),
  425. mTrack2->v0NegZ(), mTrack2->v0NegU(),
  426. mTrack1->sect(), mTrack2->v0NegSect(),
  427. &mFracOfMergedRowV0PosV0Neg,
  428. &mClosestRowAtDCAV0PosV0Neg);
  429. }
  430. return mFracOfMergedRowV0PosV0Neg;
  431. }
  432. double fractionOfMergedRowV0NegV0Pos() const {
  433. if (mMergingParNotCalculatedV0NegV0Pos) {
  434. calcMergingParFctn(&mMergingParNotCalculatedV0NegV0Pos,
  435. mTrack1->v0NegZ(), mTrack1->v0NegU(),
  436. mTrack2->z(), mTrack2->u(),
  437. mTrack1->v0NegSect(), mTrack2->sect(),
  438. &mFracOfMergedRowV0NegV0Pos,
  439. &mClosestRowAtDCAV0NegV0Pos);
  440. }
  441. return mFracOfMergedRowV0NegV0Pos;
  442. }
  443. double fractionOfMergedRowV0PosV0Pos() const {
  444. if (mMergingParNotCalculatedV0PosV0Pos) {
  445. calcMergingParFctn(&mMergingParNotCalculatedV0PosV0Pos,
  446. mTrack1->z(), mTrack1->u(),
  447. mTrack2->z(), mTrack2->u(),
  448. mTrack1->sect(), mTrack2->sect(),
  449. &mFracOfMergedRowV0PosV0Pos,
  450. &mClosestRowAtDCAV0PosV0Pos);
  451. }
  452. return mFracOfMergedRowV0PosV0Pos;
  453. }
  454. double fractionOfMergedRowV0NegV0Neg() const {
  455. if (mMergingParNotCalculatedV0NegV0Neg) {
  456. calcMergingParFctn(&mMergingParNotCalculatedV0NegV0Neg,
  457. mTrack1->v0NegZ(), mTrack1->v0NegU(),
  458. mTrack2->v0NegZ(), mTrack2->v0NegU(),
  459. mTrack1->v0NegSect(), mTrack2->v0NegSect(),
  460. &mFracOfMergedRowV0NegV0Neg,
  461. &mClosestRowAtDCAV0NegV0Neg);
  462. }
  463. return mFracOfMergedRowV0NegV0Neg;
  464. }
  465. //
  466. // Setters
  467. //
  468. /// Set the first particle from a pair
  469. void setTrack1(const MpdFemtoParticle* trkPtr) {
  470. mTrack1 = (MpdFemtoParticle*) trkPtr;
  471. resetParCalculated();
  472. }
  473. /// Set the second particle from a pair
  474. void setTrack2(const MpdFemtoParticle* trkPtr) {
  475. mTrack2 = (MpdFemtoParticle*) trkPtr;
  476. resetParCalculated();
  477. }
  478. /// Set track-merging parameters
  479. void setMergingPar(float aMaxDuInner, float aMaxDzInner,
  480. float aMaxDuOuter, float aMaxDzOuter);
  481. /// Set parameters used for track-merging estimation assuming half field
  482. void setDefaultHalfFieldMergingPar();
  483. /// Set parameters used for track-merging estimation assuming full field
  484. void setDefaultFullFieldMergingPar();
  485. /// Calculate the \Delta\phi^{*} between two particles, which is minal across several points
  486. /// \param p_a momentum of first particle
  487. /// \param charge_a charge of the first particle
  488. /// \param p_b momentum of second particle
  489. /// \param charge_b charge of the second particle
  490. /// \param radius_in_meters Radial distance at which the angle should be taken [Meters]
  491. /// \param magnetic_field_in_tesla Strength and direction of the magnetic field in the detector [Tesla]
  492. ///
  493. static double calculateDPhiStarMin(const TVector3& p_a,
  494. const short& charge_a,
  495. const TVector3& p_b,
  496. const short& charge_b,
  497. const double& rad_step_in_meters,
  498. const double& rad_min_in_meters,
  499. const double& rad_max_in_meters,
  500. const double& magnetic_field);
  501. /// Calculate the \Delta\phi^{*} between two particles at each radial step between
  502. ///
  503. /// \param p_a momentum of first particle
  504. /// \param charge_a charge of the first particle
  505. /// \param p_b momentum of second particle
  506. /// \param charge_b charge of the second particle
  507. /// \param radius_in_meters Radial distance at which the angle should be taken [Meters]
  508. /// \param magnetic_field_in_tesla Strength and direction of the magnetic field in the detector [Tesla]
  509. ///
  510. static std::vector<double> calculateDPhiStarValues(const TVector3& p_a,
  511. const short& charge_a,
  512. const TVector3& p_b,
  513. const short& charge_b,
  514. const double& rad_step_in_meters,
  515. const double& rad_min_in_meters,
  516. const double& rad_max_in_meters,
  517. const double& magnetic_field);
  518. /// Calculate the \Delta\phi^{*} between two particles.
  519. /// \param p_a momentum of first particle
  520. /// \param charge_a charge of the first particle
  521. /// \param p_b momentum of second particle
  522. /// \param charge_b charge of the second particle
  523. /// \param radius_in_meters Radial distance at which the angle should be taken [Meters]
  524. /// \param magnetic_field_in_tesla Strength and direction of the magnetic field in the detector [Tesla]
  525. ///
  526. static double calculateDPhiStar(const TVector3& p_a,
  527. const short& charge_a,
  528. const TVector3& p_b,
  529. const short& charge_b,
  530. const double& radius_in_meters,
  531. const double& magnetic_field);
  532. /// Calculate \Delta\phi between two particles.
  533. /// \param a Momentum of first particle
  534. /// \param b Momentum of second particle
  535. ///
  536. /// The calculation returns $\Delta\phi = \phi2 - \phi1$
  537. ///
  538. static double calculateDPhi(const TVector3& a, const TVector3& b);
  539. /// Calculate \Delta\eta between two particles.
  540. /// \param a Momentum of first particle
  541. /// \param b Momentum of second particle
  542. ///
  543. /// The calculation returns $\Delta\eta = \eta_2 - \eta_1$
  544. ///
  545. static double calculateDEta(const TVector3& a, const TVector3& b);
  546. /// Calculate \Delta\eta between two particles.
  547. /// \param a Momentum of first particle
  548. /// \param b Momentum of second particle
  549. /// \param minRad Radial distance at which the eta value should be taken?
  550. static double calculateDEtaStar(const TVector3& a, const TVector3& b,
  551. const double& radius_in_meters);
  552. private:
  553. /// The first particle from a pair
  554. MpdFemtoParticle* mTrack1;
  555. /// The second particle from a pair
  556. MpdFemtoParticle* mTrack2;
  557. mutable short mNonIdParNotCalculated;
  558. mutable float mDKSide;
  559. mutable float mDKOut;
  560. mutable float mDKLong;
  561. mutable float mCVK;
  562. mutable float mKStarCalc;
  563. void calcNonIdPar() const;
  564. mutable short mNonIdParNotCalculatedGlobal;
  565. mutable float mDKSideGlobal;
  566. mutable float mDKOutGlobal;
  567. mutable float mDKLongGlobal;
  568. mutable float mKStarCalcGlobal;
  569. mutable float mCVKGlobal;
  570. void calcNonIdParGlobal() const;
  571. mutable short mMergingParNotCalculated;
  572. mutable float mWeightedAvSep;
  573. mutable float mFracOfMergedRow;
  574. mutable float mClosestRowAtDCA;
  575. mutable short mMergingParNotCalculatedTrkV0Pos;
  576. mutable float mFracOfMergedRowTrkV0Pos;
  577. mutable float mClosestRowAtDCATrkV0Pos;
  578. mutable short mMergingParNotCalculatedTrkV0Neg;
  579. mutable float mFracOfMergedRowTrkV0Neg;
  580. mutable float mClosestRowAtDCATrkV0Neg;
  581. mutable short mMergingParNotCalculatedV0PosV0Neg;
  582. mutable float mFracOfMergedRowV0PosV0Neg;
  583. mutable float mClosestRowAtDCAV0PosV0Neg;
  584. mutable short mMergingParNotCalculatedV0NegV0Pos;
  585. mutable float mFracOfMergedRowV0NegV0Pos;
  586. mutable float mClosestRowAtDCAV0NegV0Pos;
  587. mutable short mMergingParNotCalculatedV0PosV0Pos;
  588. mutable float mFracOfMergedRowV0PosV0Pos;
  589. mutable float mClosestRowAtDCAV0PosV0Pos;
  590. mutable short mMergingParNotCalculatedV0NegV0Neg;
  591. mutable float mFracOfMergedRowV0NegV0Neg;
  592. mutable float mClosestRowAtDCAV0NegV0Neg;
  593. static float mMaxDuInner;
  594. static float mMaxDzInner;
  595. static float mMaxDuOuter;
  596. static float mMaxDzOuter;
  597. static float mTpcRadiusMin; //[cm]
  598. static float mTpcRadiusMax; //[cm]
  599. void calcMergingPar() const;
  600. /// Calculate track-merging parameters
  601. void calcMergingParFctn(short* tmpMergingParNotCalculatedFctn,
  602. float* tmpZ1, float* tmpU1,
  603. float* tmpZ2, float* tmpU2,
  604. int *tmpSect1, int *tmpSect2,
  605. float* tmpFracOfMergedRow,
  606. float* tmpClosestRowAtDCA) const;
  607. /// Reset calculated parameters
  608. void resetParCalculated() {
  609. mNonIdParNotCalculated = 1;
  610. mNonIdParNotCalculatedGlobal = 1;
  611. mMergingParNotCalculated = 1;
  612. mMergingParNotCalculatedTrkV0Pos = 1;
  613. mMergingParNotCalculatedTrkV0Neg = 1;
  614. mMergingParNotCalculatedV0PosV0Pos = 1;
  615. mMergingParNotCalculatedV0NegV0Pos = 1;
  616. mMergingParNotCalculatedV0PosV0Neg = 1;
  617. mMergingParNotCalculatedV0NegV0Neg = 1;
  618. }
  619. ClassDef(MpdFemtoPair, 0)
  620. };
  621. #endif // #define MpdFemtoPair_h