MpdFemtoReactionPlaneAnalysis.cxx 8.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297
  1. //
  2. // Class that allows one to perform the azimuthal analysis
  3. //
  4. // C++ headers
  5. #include <iostream>
  6. // MpdFemtoMaker headers
  7. // Base
  8. #include "MpdFemtoBaseTrackCut.h"
  9. #include "MpdFemtoBaseV0Cut.h"
  10. #include "MpdFemtoBaseKinkCut.h"
  11. // Infrastructure
  12. #include "MpdFemtoReactionPlaneAnalysis.h"
  13. #include "MpdFemtoParticleCollection.h"
  14. #include "MpdFemtoPicoEventCollectionVector.h"
  15. #include "MpdFemtoPicoEventCollectionVectorHideAway.h"
  16. #include "PhysicalConstants.h"
  17. // ROOT headers
  18. #include "TMath.h"
  19. #include "TString.h"
  20. ClassImp(MpdFemtoReactionPlaneAnalysis)
  21. //_________________
  22. MpdFemtoReactionPlaneAnalysis::MpdFemtoReactionPlaneAnalysis(unsigned int binsVertex,
  23. double minVertex,
  24. double maxVertex,
  25. unsigned int binsMult,
  26. double minMult,
  27. double maxMult,
  28. unsigned short binsRP,
  29. double minRP,
  30. double maxRP) :
  31. MpdFemtoAnalysis(),
  32. mVertexZBins(binsVertex),
  33. mOverFlowVertexZ(0),
  34. mUnderFlowVertexZ(0),
  35. mMultBins(binsMult),
  36. mOverFlowMult(0),
  37. mUnderFlowMult(0),
  38. mRPBins(binsRP),
  39. mCurrentRP(0) {
  40. // Constructor
  41. mVertexZ[0] = minVertex;
  42. mVertexZ[1] = maxVertex;
  43. mMult[0] = minMult;
  44. mMult[1] = maxMult;
  45. mRP[0] = minRP;
  46. mRP[1] = maxRP;
  47. // Print out warnings, will help user to fix these bugs
  48. if (minVertex >= maxVertex) {
  49. std::cout << "[WARNING] MpdFemtoVertexMultAnalysis - wrong z-vertex positions ("
  50. << minVertex << " >= " << maxVertex << "). No events are expected to pass."
  51. << std::endl;
  52. }
  53. if (minMult >= maxMult) {
  54. std::cout << "[WARNING] MpdFemtoVertexMultAnalysis - wrong multiplicity intervals ("
  55. << minMult << " >= " << maxMult << "). No events are expected to pass."
  56. << std::endl;
  57. }
  58. if (!mCorrFctnCollection) {
  59. mCorrFctnCollection = new MpdFemtoCorrFctnCollection;
  60. }
  61. // If the event collection was already create (it should NOT be) delete
  62. // before we allocate a new one
  63. if (mMixingBuffer) {
  64. delete mMixingBuffer;
  65. mMixingBuffer = nullptr;
  66. }
  67. mPicoEventCollectionVectorHideAway =
  68. new MpdFemtoPicoEventCollectionVectorHideAway(mVertexZBins, mVertexZ[0], mVertexZ[1],
  69. mMultBins, mMult[0], mMult[1],
  70. mRPBins, mRP[0], mRP[1]);
  71. };
  72. //_________________
  73. MpdFemtoReactionPlaneAnalysis::MpdFemtoReactionPlaneAnalysis(const MpdFemtoReactionPlaneAnalysis& a) :
  74. MpdFemtoAnalysis(),
  75. mVertexZBins(a.mVertexZBins),
  76. mOverFlowVertexZ(0),
  77. mUnderFlowVertexZ(0),
  78. mMultBins(a.mMultBins),
  79. mOverFlowMult(0),
  80. mUnderFlowMult(0),
  81. mRPBins(a.mRPBins),
  82. mCurrentRP(0) {
  83. // Copy constructor
  84. mVertexZ[0] = a.mVertexZ[0];
  85. mVertexZ[1] = a.mVertexZ[1];
  86. mMult[0] = a.mMult[0];
  87. mMult[1] = a.mMult[1];
  88. mRP[0] = a.mRP[0];
  89. mRP[1] = a.mRP[1];
  90. if (mMixingBuffer) {
  91. delete mMixingBuffer;
  92. mMixingBuffer = nullptr;
  93. }
  94. if (mPicoEventCollectionVectorHideAway) {
  95. delete mPicoEventCollectionVectorHideAway;
  96. }
  97. mPicoEventCollectionVectorHideAway =
  98. new MpdFemtoPicoEventCollectionVectorHideAway(mVertexZBins, mVertexZ[0], mVertexZ[1],
  99. mMultBins, mMult[0], mMult[1],
  100. mRPBins, mRP[0], mRP[1]);
  101. if (mVerbose) {
  102. std::cout << "MpdFemtoReactionPlaneAnalysis::MpdFemtoReactionPlaneAnalysis(const MpdFemtoReactionPlaneAnalysis& a) - "
  103. << "analysis copied" << std::endl;
  104. }
  105. }
  106. //_________________
  107. MpdFemtoReactionPlaneAnalysis& MpdFemtoReactionPlaneAnalysis::operator=(const MpdFemtoReactionPlaneAnalysis& a) {
  108. // Assignment operator
  109. if (this != &a) {
  110. // Allow parent class to copy the cuts and correlation functions
  111. MpdFemtoAnalysis::operator=(a);
  112. mVertexZBins = a.mVertexZBins;
  113. mVertexZ[0] = a.mVertexZ[0];
  114. mVertexZ[1] = a.mVertexZ[1];
  115. mUnderFlowVertexZ = 0;
  116. mOverFlowVertexZ = 0;
  117. mMultBins = a.mMultBins;
  118. mMult[0] = a.mMult[0];
  119. mMult[1] = a.mMult[1];
  120. mUnderFlowMult = 0;
  121. mOverFlowMult = 0;
  122. mRPBins = a.mRPBins;
  123. mCurrentRP = 0;
  124. if (mMixingBuffer) {
  125. delete mMixingBuffer;
  126. mMixingBuffer = nullptr;
  127. }
  128. delete mPicoEventCollectionVectorHideAway;
  129. mPicoEventCollectionVectorHideAway =
  130. new MpdFemtoPicoEventCollectionVectorHideAway(mVertexZBins, mVertexZ[0], mVertexZ[1],
  131. mMultBins, mMult[0], mMult[1],
  132. mRPBins, mRP[0], mRP[1]);
  133. } // if ( this != &a )
  134. return *this;
  135. }
  136. //_________________
  137. MpdFemtoReactionPlaneAnalysis::~MpdFemtoReactionPlaneAnalysis() {
  138. // Destructor
  139. // Now delete every PicoEvent in the EventMixingBuffer and
  140. // then the Buffer itself
  141. delete mPicoEventCollectionVectorHideAway;
  142. }
  143. //____________________________
  144. MpdFemtoString MpdFemtoReactionPlaneAnalysis::report() {
  145. if (mVerbose) {
  146. std::cout << "MpdFemtoReactionPlaneAnalysis - constructing Report..." << std::endl;
  147. }
  148. TString report("-----------\nHbt MpdFemtoReactionPlaneAnalysis Report:\n");
  149. report += TString::Format("Events are mixed in %d VertexZ bins in the range %E cm to %E cm.\n",
  150. mVertexZBins, mVertexZ[0], mVertexZ[1])
  151. + TString::Format("Events underflowing: %d\n", mUnderFlowVertexZ)
  152. + TString::Format("Events overflowing: %d\n", mOverFlowVertexZ)
  153. + TString::Format("Events are mixed in %d Mult bins in the range %E cm to %E cm.\n",
  154. mMultBins, mMult[0], mMult[1])
  155. + TString::Format("Events underflowing: %d\n", mUnderFlowMult)
  156. + TString::Format("Events overflowing: %d\n", mOverFlowMult)
  157. + TString::Format("Now adding MpdFemtoAnalysis(base) report\n")
  158. + MpdFemtoAnalysis::report();
  159. return MpdFemtoString((const char *) report);
  160. }
  161. //_________________________
  162. void MpdFemtoReactionPlaneAnalysis::processEvent(const MpdFemtoEvent* hbtEvent) {
  163. // Perform event processing in bins of z vertex, multiplicity and reaction plane
  164. // Get right mixing buffer
  165. double vertexZ = hbtEvent->primaryVertex().Z();
  166. double mult = hbtEvent->refMult();
  167. mCurrentRP = hbtEvent->eventPlaneAngle();
  168. if (mCurrentRP < -990) {
  169. // std::cout << "Event planevalue < 990!" << std::endl;
  170. return;
  171. }
  172. if (mCurrentRP < 0.0) {
  173. mCurrentRP += 2 * TMath::Pi();
  174. }
  175. mMixingBuffer = mPicoEventCollectionVectorHideAway->picoEventCollection(vertexZ, mult, mCurrentRP);
  176. if (!mMixingBuffer) {
  177. if (vertexZ < mVertexZ[0]) {
  178. mUnderFlowVertexZ++;
  179. }
  180. if (vertexZ > mVertexZ[1]) {
  181. mOverFlowVertexZ++;
  182. }
  183. if (mult < mMult[0]) {
  184. mUnderFlowMult++;
  185. }
  186. if (mult > mMult[1]) {
  187. mOverFlowMult++;
  188. }
  189. return;
  190. }
  191. // call ProcessEvent() from MpdFemtoAnalysis-base
  192. MpdFemtoAnalysis::processEvent(hbtEvent);
  193. }
  194. //_________________
  195. TList* MpdFemtoReactionPlaneAnalysis::listSettings() {
  196. TList *settings = MpdFemtoAnalysis::listSettings();
  197. settings->AddVector(new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.vertex_z.bins=%d", mVertexZBins)),
  198. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.vertex_z.min=%f", mVertexZ[0])),
  199. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.vertex_z.max=%f", mVertexZ[1])),
  200. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.multiplicity.bins=%d", mMultBins)),
  201. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.multiplicity.min=%f", mMult[0])),
  202. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.multiplicity.max=%f", mMult[1])),
  203. new TObjString(TString::Format("MpdFemtoReactionPlaneAnalysis.eventPlane.bins=%d", mRPBins)),
  204. NULL);
  205. return settings;
  206. }
  207. //_________________________
  208. TList* MpdFemtoReactionPlaneAnalysis::getOutputList() {
  209. // Collect the list of output objects to be written
  210. TList *tOutputList = new TList();
  211. TList *p1Cut = mFirstParticleCut->getOutputList();
  212. TListIter nextp1(p1Cut);
  213. while (TObject * obj = nextp1.Next()) {
  214. tOutputList->Add(obj);
  215. }
  216. delete p1Cut;
  217. if (mSecondParticleCut != mFirstParticleCut) {
  218. TList *p2Cut = mSecondParticleCut->getOutputList();
  219. TIter nextp2(p2Cut);
  220. while (TObject * obj = nextp2()) {
  221. tOutputList->Add(obj);
  222. }
  223. delete p2Cut;
  224. } //if (fSecondParticleCut != fFirstParticleCut)
  225. TList *pairCut = mPairCut->getOutputList();
  226. TIter nextpair(pairCut);
  227. while (TObject * obj = nextpair()) {
  228. tOutputList->Add(obj);
  229. }
  230. delete pairCut;
  231. TList *eventCut = mEventCut->getOutputList();
  232. TIter nextevent(eventCut);
  233. while (TObject * obj = nextevent()) {
  234. tOutputList->Add(obj);
  235. }
  236. delete eventCut;
  237. for (auto &cf : *mCorrFctnCollection) {
  238. TList *tListCf = cf->getOutputList();
  239. TIter nextListCf(tListCf);
  240. while (TObject * obj = nextListCf()) {
  241. tOutputList->Add(obj);
  242. }
  243. delete tListCf;
  244. }
  245. return tOutputList;
  246. }