MpdFemtoQinvCorrFctnKt.cxx 6.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232
  1. //
  2. // One-dimensional correlation function with kT-binning
  3. //
  4. // MpdFemtoMakerUser headers
  5. #include "MpdFemtoQinvCorrFctnKt.h"
  6. // ROOT headers
  7. #include "TString.h"
  8. #include "TH1D.h"
  9. // C++ headers
  10. #include <iostream>
  11. ClassImp(MpdFemtoQinvCorrFctnKt)
  12. //_________________
  13. MpdFemtoQinvCorrFctnKt::MpdFemtoQinvCorrFctnKt(const char* title,
  14. const int &nQbins, const double &qLo, const double &qHi,
  15. const int& nKtBins, const double& ktLo, const double& ktHi) :
  16. MpdFemtoBaseCorrFctn() {
  17. // Constructor
  18. // Define number of kT bins and kT step
  19. setKtRange(nKtBins, ktLo, ktHi);
  20. TString baseName(title);
  21. TString baseTitle(title);
  22. TString appendix = ";q_{inv} (GeV/c);Entries";
  23. // Histogram loop
  24. for (int iKtBin = 0; iKtBin < mNKtBins; iKtBin++) {
  25. baseName = title;
  26. baseName += "_num_";
  27. baseName += iKtBin;
  28. baseTitle = baseName;
  29. baseTitle += " ";
  30. baseTitle += Form("%3.2f", (mKtRange[0] + mKtStep * iKtBin));
  31. baseTitle += "#leq k_{T} (GeV/c) #leq ";
  32. baseTitle += Form("%3.2f", (mKtRange[0] + mKtStep * (iKtBin + 1)));
  33. baseTitle += appendix;
  34. TH1D *histoNum = new TH1D(baseName.Data(), baseTitle.Data(),
  35. nQbins, qLo, qHi);
  36. histoNum->Sumw2();
  37. mNumerator.push_back(histoNum);
  38. baseName = title;
  39. baseName += "_den_";
  40. baseName += iKtBin;
  41. baseTitle = baseName;
  42. baseTitle += " ";
  43. baseTitle += Form("%3.2f", (mKtRange[0] + mKtStep * iKtBin));
  44. baseTitle += "#leq k_{T} (GeV/c) #leq ";
  45. baseTitle += Form("%3.2f", (mKtRange[0] + mKtStep * (iKtBin + 1)));
  46. baseTitle += appendix;
  47. TH1D *histoDen = new TH1D(baseName.Data(), baseTitle.Data(),
  48. nQbins, qLo, qHi);
  49. histoDen->Sumw2();
  50. mDenominator.push_back(histoDen);
  51. } // for(int i=0; i<mNumberKt; i++)
  52. }
  53. //_________________
  54. MpdFemtoQinvCorrFctnKt::MpdFemtoQinvCorrFctnKt(const MpdFemtoQinvCorrFctnKt& c) :
  55. MpdFemtoBaseCorrFctn(c),
  56. mNKtBins(c.mNKtBins),
  57. mKtStep(c.mKtStep) {
  58. // Copy constructor
  59. mKtRange[0] = c.mKtRange[0];
  60. mKtRange[1] = c.mKtRange[1];
  61. for (int iCF = 0; iCF < mNKtBins; iCF++) {
  62. mNumerator.push_back((c.mNumerator.at(iCF)) ? new TH1D(*c.mNumerator.at(iCF)) : nullptr);
  63. mDenominator.push_back((c.mDenominator.at(iCF)) ? new TH1D(*c.mDenominator.at(iCF)) : nullptr);
  64. } // for ( int iBin=0; iBin<mNumberKt; iBin++ )
  65. }
  66. //_________________
  67. MpdFemtoQinvCorrFctnKt& MpdFemtoQinvCorrFctnKt::operator=(const MpdFemtoQinvCorrFctnKt& c) {
  68. // Assignment operator
  69. if (this != &c) {
  70. mNKtBins = c.mNKtBins;
  71. mKtRange[0] = c.mKtRange[0];
  72. mKtRange[1] = c.mKtRange[1];
  73. mKtStep = c.mKtStep;
  74. for (int iCF = 0; iCF < mNKtBins; iCF++) {
  75. mNumerator.push_back((c.mNumerator.at(iCF)) ? new TH1D(*c.mNumerator.at(iCF)) : nullptr);
  76. mDenominator.push_back((c.mDenominator.at(iCF)) ? new TH1D(*c.mDenominator.at(iCF)) : nullptr);
  77. } // for ( int iBin=0; iBin<mNumberKt; iBin++ )
  78. } // if (this != &c)
  79. return *this;
  80. }
  81. //_________________
  82. MpdFemtoQinvCorrFctnKt::~MpdFemtoQinvCorrFctnKt() {
  83. // Destructor
  84. if (!mNumerator.empty()) {
  85. mNumerator.clear();
  86. }
  87. if (!mDenominator.empty()) {
  88. mDenominator.clear();
  89. }
  90. }
  91. //_________________
  92. void MpdFemtoQinvCorrFctnKt::writeOutHistos() {
  93. // Write histograms to the file
  94. if (!mNumerator.empty()) {
  95. for (unsigned int i = 0; i < mNumerator.size(); i++) {
  96. mNumerator.at(i)->Write();
  97. }
  98. }
  99. if (!mDenominator.empty()) {
  100. for (unsigned int i = 0; i < mDenominator.size(); i++) {
  101. mDenominator.at(i)->Write();
  102. }
  103. }
  104. }
  105. //_________________
  106. TList* MpdFemtoQinvCorrFctnKt::getOutputList() {
  107. // Prepare the list of objects to be written to the output
  108. TList *outputList = new TList();
  109. for (int i = 0; i < mNKtBins; i++) {
  110. outputList->Add(mNumerator.at(i));
  111. outputList->Add(mDenominator.at(i));
  112. }
  113. return outputList;
  114. }
  115. //_________________
  116. void MpdFemtoQinvCorrFctnKt::eventBegin(const MpdFemtoEvent* /* event */) {
  117. /* empty */
  118. }
  119. //_________________
  120. void MpdFemtoQinvCorrFctnKt::eventEnd(const MpdFemtoEvent* /* event */) {
  121. /* empty */
  122. }
  123. //_________________
  124. void MpdFemtoQinvCorrFctnKt::finish() {
  125. /* empty */
  126. }
  127. //_________________
  128. void MpdFemtoQinvCorrFctnKt::setKtRange(const int& nbins, const double& ktLo, const double& ktHi) {
  129. mNKtBins = nbins;
  130. mKtRange[0] = ktLo;
  131. mKtRange[1] = ktHi;
  132. mKtStep = (mKtRange[1] - mKtRange[0]) / mNKtBins;
  133. }
  134. //_________________
  135. void MpdFemtoQinvCorrFctnKt::addRealPair(MpdFemtoPair* pair) {
  136. // Check if pair passes front-loaded cut if exists
  137. if (mPairCut && !mPairCut->pass(pair)) {
  138. return;
  139. }
  140. int mIndexKt = (int) (((pair->kT() - mKtRange[0]) / mKtStep));
  141. if ((mIndexKt >= 0) && (mIndexKt < mNKtBins)) {
  142. if (mNumerator.at(mIndexKt)) {
  143. mNumerator.at(mIndexKt)->Fill(fabs(pair->qInv()));
  144. } else {
  145. std::cout << "[ERROR] void MpdFemtoQinvCorrFctnKt::addRealPair(MpdFemtoPair* pair) - "
  146. << "No correlation function with index: " << mIndexKt << " was found"
  147. << std::endl;
  148. }
  149. } // if ( ( mIndexKt>=0 ) && ( mIndexKt<mNumberKt ) )
  150. }
  151. //_________________
  152. void MpdFemtoQinvCorrFctnKt::addMixedPair(MpdFemtoPair* pair) {
  153. // Check if pair passes front-loaded cut if exists
  154. if (mPairCut && !mPairCut->pass(pair)) {
  155. return;
  156. }
  157. int mIndexKt = (int) (((fabs(pair->kT()) - mKtRange[0]) / mKtStep));
  158. if ((mIndexKt >= 0) && (mIndexKt < mNKtBins)) {
  159. if (mDenominator.at(mIndexKt)) {
  160. mDenominator.at(mIndexKt)->Fill(fabs(pair->qInv()));
  161. } else {
  162. std::cout << "[ERROR] void MpdFemtoQinvCorrFctnKt::addMixedPair(MpdFemtoPair* pair) - "
  163. << "No correlation function with index: " << mIndexKt << " was found"
  164. << std::endl;
  165. }
  166. } // if ( ( mIndexKt>=0 ) && ( mIndexKt<mNumberKt ) )
  167. }
  168. //_________________
  169. MpdFemtoString MpdFemtoQinvCorrFctnKt::report() {
  170. // Make a report
  171. MpdFemtoString tStr = "MpdFemtoQinvCorrFctnKt report";
  172. int mNumeratorEntries = 0;
  173. int mDenominatorEntries = 0;
  174. for (int i = 0; i < mNKtBins; i++) {
  175. mNumeratorEntries += (int) mNumerator.at(i)->GetEntries();
  176. mDenominatorEntries += (int) mDenominator.at(i)->GetEntries();
  177. } // for (int i=0; i<mNumberKt; i++)
  178. tStr += TString::Format("\nTotal number of pairs in the numerator :\t%d\n", mNumeratorEntries);
  179. tStr += TString::Format("Total number of pairs in the denominator:\t%d\n", mDenominatorEntries);
  180. for (int i = 0; i < mNKtBins; i++) {
  181. tStr += TString::Format("Total number of pairs in %d-th numerator :\t%E\n",
  182. i, mNumerator.at(i)->GetEntries());
  183. tStr += TString::Format("Total number of pairs in %d-th denominator :\t%E\n",
  184. i, mDenominator.at(i)->GetEntries());
  185. }
  186. if (mPairCut) {
  187. tStr += "Here is the PairCut specific to this CorrFctn\n";
  188. tStr += mPairCut->report();
  189. } else {
  190. tStr += "No PairCut specific to this CorrFctn\n";
  191. }
  192. return tStr;
  193. }