MpdFemtoAverageSeparation.cxx 4.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192
  1. //
  2. // One-dimensional distr. of Average Separation
  3. //
  4. // MpdFemtoMaker headers
  5. #include "MpdFemtoPair.h"
  6. // MpdFemtoMakerUser headers
  7. #include "MpdFemtoAverageSeparation.h"
  8. // ROOT headers
  9. #include "TMath.h"
  10. #include "TString.h"
  11. ClassImp(MpdFemtoAverageSeparation)
  12. //_________________
  13. MpdFemtoAverageSeparation::MpdFemtoAverageSeparation(const char* title,
  14. const int& nBinsAvSep, const double& AvSepLo,
  15. const double& AvSepHi) :
  16. MpdFemtoBaseCorrFctn() {
  17. // Parametrized constructor
  18. // Set initial parameters
  19. setHistoParameters(nBinsAvSep, AvSepLo, AvSepHi);
  20. // Define string parameters
  21. TString baseName(title);
  22. TString baseTitle(title);
  23. TString appendix = "; AvSep#Delta#phi^{*}_{min} (rad);#Delta#eta";
  24. // Set numerator
  25. baseName = title;
  26. baseName += "_num_";
  27. baseTitle += "_num_";
  28. TH1F *hNumer = new TH1F(baseName.Data(), baseTitle.Data(), nBinsAvSep,
  29. AvSepLo, AvSepHi);
  30. hNumer->Sumw2();
  31. mNumerator.push_back(hNumer);
  32. // Denominator
  33. baseName = title;
  34. baseName += "_den_";
  35. baseTitle += "_den_";
  36. TH1F *hDenom = new TH1F(baseName.Data(), baseTitle.Data(),nBinsAvSep,
  37. AvSepLo, AvSepHi);
  38. hDenom->Sumw2();
  39. mDenominator.push_back(hDenom);
  40. }
  41. //____________________________
  42. MpdFemtoAverageSeparation::~MpdFemtoAverageSeparation() {
  43. // Destructor
  44. if (!mNumerator.empty()) {
  45. mNumerator.clear();
  46. }
  47. if (!mDenominator.empty()) {
  48. mDenominator.clear();
  49. }
  50. }
  51. //_________________
  52. void MpdFemtoAverageSeparation::eventBegin(const MpdFemtoEvent* event) {
  53. // setBField(event->bField());
  54. }
  55. //_________________
  56. void MpdFemtoAverageSeparation::eventEnd(const MpdFemtoEvent* /* event */) {
  57. /* empty */
  58. }
  59. //_________________________
  60. void MpdFemtoAverageSeparation::writeOutHistos() {
  61. // Write out all histograms to file
  62. for (unsigned int i = 0; i < mNumerator.size(); i++) {
  63. mNumerator.at(i)->Write();
  64. }
  65. for (unsigned int i = 0; i < mDenominator.size(); i++) {
  66. mDenominator.at(i)->Write();
  67. }
  68. }
  69. //______________________________
  70. TList* MpdFemtoAverageSeparation::getOutputList() {
  71. // Prepare the list of objects to be written to the output
  72. TList *tOutputList = new TList();
  73. // for (int i = 0; i < mNKtBins; i++) {
  74. tOutputList->Add(mNumerator.at(0));
  75. tOutputList->Add(mDenominator.at(0));
  76. // } // for (int i=0; i<mNKtBins; i++)
  77. return tOutputList;
  78. }
  79. //_________________________
  80. void MpdFemtoAverageSeparation::finish() {
  81. /* empty */
  82. }
  83. //____________________________
  84. MpdFemtoString MpdFemtoAverageSeparation::report() {
  85. // Construct the report
  86. TString report("Average Separation report:\n");
  87. // for (int iKt = 0; iKt < mNKtBins; iKt++) {
  88. int iKt =0;
  89. report += TString::Format("Number of entries in %d-th numerator :\t%E\n",
  90. iKt, mNumerator.at(iKt)->GetEntries());
  91. report += TString::Format("Number of entries in %d-th denominator :\t%E\n",
  92. iKt, mDenominator.at(iKt)->GetEntries());
  93. // } // for (int iKt=0; iKt<mNKtBins; iKt++)
  94. if (mPairCut) {
  95. report += "Report from the front-loaded PairCut specific to this CorrFctn\n";
  96. report += mPairCut->report();
  97. } else {
  98. report += "No front-loaded PairCut specific to this CorrFctn\n";
  99. }
  100. return MpdFemtoString((const char *) report);
  101. }
  102. //-------------------------------------------
  103. void MpdFemtoAverageSeparation::addRealPair(MpdFemtoPair* pair) {
  104. // Perform operations on real pairs
  105. // Check front-loaded pair cut
  106. if (mPairCut && !mPairCut->pass(pair)) {
  107. return;
  108. }
  109. if (pair->qInv() < 0.15) {
  110. double avsep = pair->nominalTpcAverageSeparation();
  111. // Check that histrogram exists
  112. if (mNumerator.at(0)) {
  113. mNumerator.at(0)->Fill(avsep, 1.);
  114. } // if ( mNumerator.at(mIndexKt) )
  115. } // if ( (mIndexKt>=0) && (mIndexKt<mNKtBins) )
  116. }
  117. //_________________
  118. void MpdFemtoAverageSeparation::addMixedPair(MpdFemtoPair* pair) {
  119. // Check front-loaded pair cut
  120. if (mPairCut && !mPairCut->pass(pair)) {
  121. return;
  122. }
  123. if (pair->qInv() < 0.15) {
  124. double avsep = pair->nominalTpcAverageSeparation();
  125. // Check that histrogram exists
  126. if (mDenominator.at(0)) {
  127. mDenominator.at(0)->Fill(avsep, 1.);
  128. } // if ( mDenominator.at(mIndexKt) )
  129. } // if ( (mIndexKt>=0) && (mIndexKt<mNKtBins) )
  130. }
  131. void MpdFemtoAverageSeparation::setHistoParameters(const int& nBinsAvSep, const double& AvSepLo, const double& AvSepHi) {
  132. if (AvSepLo > AvSepHi) {
  133. std::cout << "[WARNING] void MpdFemtoAverageSeparation::void setHistoParameters - "
  134. << " AvSepLo > AvSepHi. Resetting to defaults." << std::endl;
  135. mAvSepRange[0] = 0;
  136. mAvSepRange[1] = 100.;
  137. }
  138. else {
  139. mAvSepRange[0] = AvSepLo;
  140. mAvSepRange[1] = AvSepHi;
  141. }
  142. if (nBinsAvSep < 1) {
  143. std::cout << "[WARNING] void MpdFemtoAverageSeparation::void setHistoParameters - "
  144. << " nBinsAvSep<1. Resetting to defaults." << std::endl;
  145. mAvSepBins = 200;
  146. } else {
  147. mAvSepBins = nBinsAvSep;
  148. }
  149. }