qaReader_phqmd.cxx 2.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137
  1. #include <qaReader_phqmd.h>
  2. ClassImp(qaReader_phqmd);
  3. qaReader_phqmd::qaReader_phqmd(/* args */) : is_init(false), fCurrentEvent(-1), fEvent(nullptr), fParticle(nullptr)
  4. {
  5. }
  6. qaReader_phqmd::~qaReader_phqmd()
  7. {
  8. }
  9. Bool_t qaReader_phqmd::ChainCheck()
  10. {
  11. if (!is_init)
  12. {
  13. return false;
  14. }
  15. if (fCurrentEvent == -1)
  16. {
  17. return false;
  18. }
  19. if (!fChainPHSD->GetEntry(fCurrentEvent))
  20. {
  21. return false;
  22. }
  23. if (!fChainFrag->GetEntry(fCurrentEvent))
  24. {
  25. return false;
  26. }
  27. return true;
  28. }
  29. void qaReader_phqmd::SetChain(const TString &inputFileName)
  30. {
  31. fChainPHSD = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainPHSDName);
  32. fChainPHSD->SetBranchAddress("event", &fEvent);
  33. fChainFrag = (TChain *)qaUtility::GetInstance()->initChain(inputFileName, fChainFragName);
  34. fChainFrag->SetBranchAddress("eventFrag", &fEventFrag);
  35. is_init = kTRUE;
  36. }
  37. qaEvent *qaReader_phqmd::ReadEvent(Long64_t iev)
  38. {
  39. fCurrentEvent = iev;
  40. if (!ChainCheck())
  41. {
  42. return nullptr;
  43. }
  44. qaEvent *event = new qaEvent();
  45. event->SetB(fEvent->GetB());
  46. // event->SetPhiRP(0.); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
  47. event->SetPhiRP(TMath::Pi()); // in PHQMD axis orientation is rotated by 180 deg. In HSD/PHSD it's ok.
  48. event->SetNparticles(fEvent->GetNpart() + fEventFrag->GetNfragMST()); // number of particles from PHSD and number of MST fragments
  49. return event;
  50. }
  51. qaParticle *qaReader_phqmd::ReadParticle(Int_t ipart)
  52. {
  53. if (!ChainCheck())
  54. {
  55. return nullptr;
  56. }
  57. if (ipart >= fEvent->GetNpart() + fEventFrag->GetNfragMST())
  58. {
  59. return nullptr;
  60. }
  61. if (ipart < 0)
  62. {
  63. return nullptr;
  64. }
  65. bool is_particle = false;
  66. bool is_fragment = false;
  67. if (ipart < fEvent->GetNpart())
  68. {
  69. is_particle = true;
  70. }
  71. if (ipart >= fEvent->GetNpart())
  72. {
  73. is_fragment = true;
  74. }
  75. if (!is_particle && !is_fragment)
  76. {
  77. return nullptr;
  78. }
  79. if (is_particle && is_fragment)
  80. {
  81. return nullptr;
  82. }
  83. qaParticle *particle = new qaParticle();
  84. if (is_particle)
  85. {
  86. fParticle = fEvent->GetParticle(ipart);
  87. if (fParticle->IsInMST()) // Check whether the particle from PHSD was used to create MST fragment
  88. {
  89. return nullptr;
  90. }
  91. particle->SetEnergy(fParticle->E());
  92. particle->SetPdg(fParticle->GetPdg());
  93. particle->SetPxPyPz(fParticle->Px(), fParticle->Py(), fParticle->Pz());
  94. particle->SetTime(0.);
  95. particle->SetXYZ(0., 0., 0.);
  96. particle->SetCharge(fParticle->GetCharge());
  97. }
  98. if (is_fragment)
  99. {
  100. fFragmentMST = fEventFrag->GetFragmentMST(ipart - fEvent->GetNpart());
  101. particle->SetEnergy(fFragmentMST->E());
  102. particle->SetPdg(fFragmentMST->GetPdg());
  103. particle->SetPxPyPz(fFragmentMST->Px(), fFragmentMST->Py(), fFragmentMST->Pz());
  104. particle->SetTime(0.);
  105. particle->SetXYZ(0., 0., 0.);
  106. particle->SetCharge(fFragmentMST->GetZ());
  107. }
  108. return particle;
  109. }