MpdFemtoModelQinvCorrFctn.cxx 7.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274
  1. //
  2. // One-dimensional correlation function for the model estimations
  3. //
  4. // MpdFemtoMaker headers
  5. #include "MpdFemtoModelGausLCMSFreezeOutGenerator.h"
  6. #include "MpdFemtoModelHiddenInfo.h"
  7. #include "MpdFemtoModelManager.h"
  8. #include "MpdFemtoPair.h"
  9. // MpdFemtoMakerUser headers
  10. #include "MpdFemtoModelQinvCorrFctn.h"
  11. // ROOT headers
  12. #include "TH1D.h"
  13. #include "TString.h"
  14. ClassImp(MpdFemtoModelQinvCorrFctn)
  15. //_________________
  16. MpdFemtoModelQinvCorrFctn::MpdFemtoModelQinvCorrFctn() :
  17. MpdFemtoBaseCorrFctn(),
  18. mManager(nullptr),
  19. mNumeratorWeighted(nullptr),
  20. mNumerator(nullptr),
  21. mDenominatorWeighted(nullptr),
  22. mDenominator(nullptr),
  23. mIsUseDenominator(true) {
  24. // Default constructor
  25. mNumeratorWeighted = new TH1D("ModelNumWeighted", "Weighted numerator;q_{inv} (GeV/c);dN/dq_{inv}", 80, 0., 0.8);
  26. mNumerator = new TH1D("ModelNum", "Unweighted numerator;q_{inv} (GeV/c);dN/dq_{inv}", 80, 0., 0.8);
  27. mDenominatorWeighted = new TH1D("ModelDenWeighted", "ModelNumTrueIdeal", 80, 0., 0.8);
  28. mDenominator = new TH1D("ModelDen", "ModelNumFakeIdeal", 80, 0., 0.8);
  29. mNumeratorWeighted->Sumw2();
  30. mNumerator->Sumw2();
  31. mDenominatorWeighted->Sumw2();
  32. mDenominator->Sumw2();
  33. }
  34. //_________________
  35. MpdFemtoModelQinvCorrFctn::MpdFemtoModelQinvCorrFctn(const char *title, const int& aNbins,
  36. const double& aQinvLo, const double& aQinvHi,
  37. const bool useDenominator) :
  38. MpdFemtoBaseCorrFctn(),
  39. mManager(nullptr),
  40. mNumeratorWeighted(nullptr),
  41. mNumerator(nullptr),
  42. mDenominatorWeighted(nullptr),
  43. mDenominator(nullptr),
  44. mIsUseDenominator(useDenominator) {
  45. // Parametrized constructor
  46. TString baseName;
  47. TString baseTitle;
  48. TString appendix = ";q_{inv} (GeV/c);Entries";
  49. baseName = title;
  50. baseName += "_num_wei";
  51. baseTitle = baseName;
  52. baseTitle += appendix;
  53. mNumeratorWeighted = new TH1D(baseName.Data(), baseTitle.Data(),
  54. aNbins, aQinvLo, aQinvHi);
  55. mNumeratorWeighted->Sumw2();
  56. baseName = title;
  57. baseName += "_num";
  58. baseTitle = baseName;
  59. baseTitle += appendix;
  60. mNumerator = new TH1D(baseName.Data(), baseTitle.Data(),
  61. aNbins, aQinvLo, aQinvHi);
  62. mNumerator->Sumw2();
  63. if ( mIsUseDenominator ) {
  64. baseName = title;
  65. baseName += "_den_wei";
  66. baseTitle = baseName;
  67. baseTitle += appendix;
  68. mDenominatorWeighted = new TH1D(baseName.Data(), baseTitle.Data(),
  69. aNbins, aQinvLo, aQinvHi);
  70. mDenominatorWeighted->Sumw2();
  71. baseName = title;
  72. baseName += "_den";
  73. baseTitle = baseName;
  74. baseTitle += appendix;
  75. mDenominator = new TH1D(baseName.Data(), baseTitle.Data(),
  76. aNbins, aQinvLo, aQinvHi);
  77. mDenominator->Sumw2();
  78. }
  79. }
  80. //_________________
  81. MpdFemtoModelQinvCorrFctn::MpdFemtoModelQinvCorrFctn(const MpdFemtoModelQinvCorrFctn& corrFctn) :
  82. MpdFemtoBaseCorrFctn(),
  83. mManager(corrFctn.mManager),
  84. mNumeratorWeighted(nullptr),
  85. mNumerator(nullptr),
  86. mDenominatorWeighted(nullptr),
  87. mDenominator(nullptr),
  88. mIsUseDenominator(corrFctn.mIsUseDenominator) {
  89. // Copy constructor
  90. if (corrFctn.mNumeratorWeighted) {
  91. mNumeratorWeighted = new TH1D(*(corrFctn.mNumeratorWeighted));
  92. }
  93. if (corrFctn.mNumerator) {
  94. mNumerator = new TH1D(*(corrFctn.mNumerator));
  95. }
  96. if ( mIsUseDenominator ) {
  97. if (corrFctn.mDenominatorWeighted) {
  98. mDenominatorWeighted = new TH1D(*(corrFctn.mDenominatorWeighted));
  99. }
  100. if (corrFctn.mDenominator) {
  101. mDenominator = new TH1D(*(corrFctn.mDenominator));
  102. }
  103. }
  104. }
  105. //_________________
  106. MpdFemtoModelQinvCorrFctn::~MpdFemtoModelQinvCorrFctn() {
  107. // Destructor
  108. if (mNumeratorWeighted) {
  109. delete mNumeratorWeighted;
  110. mNumeratorWeighted = nullptr;
  111. }
  112. if (mNumerator) {
  113. delete mNumerator;
  114. mNumerator = nullptr;
  115. }
  116. if ( mIsUseDenominator ) {
  117. if (mDenominatorWeighted) {
  118. delete mDenominatorWeighted;
  119. mDenominatorWeighted = nullptr;
  120. }
  121. if (mDenominator) {
  122. delete mDenominator;
  123. mDenominator = nullptr;
  124. }
  125. }
  126. delete mManager;
  127. }
  128. //_________________
  129. MpdFemtoModelQinvCorrFctn& MpdFemtoModelQinvCorrFctn::operator=(const MpdFemtoModelQinvCorrFctn& corrFctn) {
  130. // Assignment operator
  131. if (this == &corrFctn) {
  132. return *this;
  133. }
  134. if (mNumeratorWeighted) delete mNumeratorWeighted;
  135. mNumeratorWeighted = ((corrFctn.mNumeratorWeighted) ? new TH1D(*corrFctn.mNumeratorWeighted) : nullptr);
  136. if (mNumerator) delete mNumerator;
  137. mNumerator = ((corrFctn.mNumerator) ? new TH1D(*corrFctn.mNumerator) : nullptr);
  138. if ( mIsUseDenominator ) {
  139. if (mDenominatorWeighted) delete mDenominatorWeighted;
  140. mDenominatorWeighted = ((corrFctn.mDenominatorWeighted) ? new TH1D(*corrFctn.mDenominatorWeighted) : nullptr);
  141. if (mDenominator) delete mDenominator;
  142. mDenominator = ((corrFctn.mDenominator) ? new TH1D(*corrFctn.mDenominator) : nullptr);
  143. }
  144. mManager = corrFctn.mManager;
  145. return *this;
  146. }
  147. //_________________
  148. void MpdFemtoModelQinvCorrFctn::connectToManager(MpdFemtoModelManager *manager) {
  149. mManager = manager;
  150. }
  151. //_________________
  152. MpdFemtoString MpdFemtoModelQinvCorrFctn::report() {
  153. // Make a report
  154. MpdFemtoString tStr = "MpdFemtoModelQinvCorrFctn report";
  155. tStr += TString::Format("Number of entries in numerator (w/ weight) :\t%E\n",
  156. mNumeratorWeighted->GetEntries());
  157. tStr += TString::Format("Number of entries in numerator (w/o weight) :\t%E\n",
  158. mNumerator->GetEntries());
  159. if ( mIsUseDenominator ) {
  160. tStr += TString::Format("Number of entries in denominator (w/ weight) :\t%E\n",
  161. mDenominatorWeighted->GetEntries());
  162. tStr += TString::Format("Number of entries in denominator (w/o weight):\t%E\n",
  163. mDenominator->GetEntries());
  164. }
  165. if (mPairCut) {
  166. tStr += "Here is the PairCut specific to this CorrFctn\n";
  167. tStr += mPairCut->report();
  168. } else {
  169. tStr += "No PairCut specific to this CorrFctn\n";
  170. }
  171. return tStr;
  172. }
  173. //_________________
  174. void MpdFemtoModelQinvCorrFctn::addRealPair(MpdFemtoPair* pair) {
  175. // Check if pair passes front-loaded cut if exists
  176. if (mPairCut && !mPairCut->pass(pair)) {
  177. return;
  178. }
  179. // Retrieve femtoscopic weight
  180. double weight = mManager->weight(pair);
  181. mNumeratorWeighted->Fill( fabs(pair->qInv()), weight);
  182. mNumerator->Fill( fabs(pair->qInv()), 1.);
  183. }
  184. //_________________
  185. void MpdFemtoModelQinvCorrFctn::addMixedPair(MpdFemtoPair* pair) {
  186. if ( !mIsUseDenominator ) return;
  187. // Check if pair passes front-loaded cut if exists
  188. if (mPairCut && !mPairCut->pass(pair)) {
  189. return;
  190. }
  191. double weight = mManager->weight(pair);
  192. mDenominatorWeighted->Fill( fabs(pair->qInv()), weight);
  193. mDenominator->Fill( fabs(pair->qInv()), 1.);
  194. }
  195. //_________________
  196. void MpdFemtoModelQinvCorrFctn::eventBegin(const MpdFemtoEvent* /* event */) {
  197. /* empty */
  198. }
  199. //_________________
  200. void MpdFemtoModelQinvCorrFctn::eventEnd(const MpdFemtoEvent* /* event */) {
  201. /* empty */
  202. }
  203. //_________________
  204. void MpdFemtoModelQinvCorrFctn::finish() {
  205. /* empty */
  206. }
  207. //_________________
  208. void MpdFemtoModelQinvCorrFctn::writeOutHistos() {
  209. // Write out data histos
  210. mNumeratorWeighted->Write();
  211. mNumerator->Write();
  212. if ( mIsUseDenominator ) {
  213. mDenominatorWeighted->Write();
  214. mDenominator->Write();
  215. }
  216. }
  217. //_________________
  218. TList* MpdFemtoModelQinvCorrFctn::getOutputList() {
  219. // Prepare the list of objects to be written to the output
  220. TList *tOutputList = new TList();
  221. tOutputList->Add(mNumeratorWeighted);
  222. tOutputList->Add(mNumerator);
  223. if ( mIsUseDenominator ) {
  224. tOutputList->Add(mDenominatorWeighted);
  225. tOutputList->Add(mDenominator);
  226. }
  227. return tOutputList;
  228. }