PWG2/FEMTOSCOPY/AliFemto/AliFemtoQinvCorrFctn.cxx

   1 ///////////////////////////////////////////////////////////////////////////
   2 //                                                                       //
   3 // AliFemtoQinvCorrFctn:                                                 //
   4 // a simple Q-invariant correlation function                             //
   5 //                                                                       //
   6 ///////////////////////////////////////////////////////////////////////////
   7
   8 #include "AliFemtoQinvCorrFctn.h"
   9 //#include "AliFemtoHisto.h"
  10 #include <cstdio>
  11
  12 #ifdef __ROOT__
  13 ClassImp(AliFemtoQinvCorrFctn)
  14 #endif
  15
  16 //____________________________
  17 AliFemtoQinvCorrFctn::AliFemtoQinvCorrFctn(char* title, const int& nbins, const float& QinvLo, const float& QinvHi):
  18   fNumerator(0),
  19   fDenominator(0),
  20   fRatio(0)
  21 {
  22   // set up numerator
  23   //  title = "Num Qinv (MeV/c)";
  24   char tTitNum[100] = "Num";
  25   strcat(tTitNum,title);
  26   fNumerator = new TH1D(tTitNum,title,nbins,QinvLo,QinvHi);
  27   // set up denominator
  28   //title = "Den Qinv (MeV/c)";
  29   char tTitDen[100] = "Den";
  30   strcat(tTitDen,title);
  31   fDenominator = new TH1D(tTitDen,title,nbins,QinvLo,QinvHi);
  32   // set up ratio
  33   //title = "Ratio Qinv (MeV/c)";
  34   char tTitRat[100] = "Rat";
  35   strcat(tTitRat,title);
  36   fRatio = new TH1D(tTitRat,title,nbins,QinvLo,QinvHi);
  37   // this next bit is unfortunately needed so that we can have many histos of same "title"
  38   // it is neccessary if we typedef TH1D to TH1d (which we do)
  39   //fNumerator->SetDirectory(0);
  40   //fDenominator->SetDirectory(0);
  41   //fRatio->SetDirectory(0);
  42
  43   // to enable error bar calculation...
  44   fNumerator->Sumw2();
  45   fDenominator->Sumw2();
  46   fRatio->Sumw2();
  47
  48 }
  49
  50 //____________________________
  51 AliFemtoQinvCorrFctn::AliFemtoQinvCorrFctn(const AliFemtoQinvCorrFctn& aCorrFctn) :
  52   fNumerator(0),
  53   fDenominator(0),
  54   fRatio(0)
  55 {
  56   // copy constructor
  57   fNumerator = new TH1D(*aCorrFctn.fNumerator);
  58   fDenominator = new TH1D(*aCorrFctn.fDenominator);
  59   fRatio = new TH1D(*aCorrFctn.fRatio);
  60 }
  61 //____________________________
  62 AliFemtoQinvCorrFctn::~AliFemtoQinvCorrFctn(){
  63   // destructor
  64   delete fNumerator;
  65   delete fDenominator;
  66   delete fRatio;
  67 }
  68 //_________________________
  69 AliFemtoQinvCorrFctn& AliFemtoQinvCorrFctn::operator=(const AliFemtoQinvCorrFctn& aCorrFctn)
  70 {
  71   // assignment operator
  72   if (this == &aCorrFctn)
  73     return *this;
  74
  75   if (fNumerator) delete fNumerator;
  76   fNumerator = new TH1D(*aCorrFctn.fNumerator);
  77   if (fDenominator) delete fDenominator;
  78   fDenominator = new TH1D(*aCorrFctn.fDenominator);
  79   if (fRatio) delete fRatio;
  80   fRatio = new TH1D(*aCorrFctn.fRatio);
  81
  82   return *this;
  83 }
  84
  85 //_________________________
  86 void AliFemtoQinvCorrFctn::Finish(){
  87   // here is where we should normalize, fit, etc...
  88   // we should NOT Draw() the histos (as I had done it below),
  89   // since we want to insulate ourselves from root at this level
  90   // of the code.  Do it instead at root command line with browser.
  91   //  fNumerator->Draw();
  92   //fDenominator->Draw();
  93   //fRatio->Draw();
  94   fRatio->Divide(fNumerator,fDenominator,1.0,1.0);
  95
  96 }
  97
  98 //____________________________
  99 AliFemtoString AliFemtoQinvCorrFctn::Report(){
 100   // construct report
 101   string stemp = "Qinv Correlation Function Report:\n";
 102   char ctemp[100];
 103   sprintf(ctemp,"Number of entries in numerator:\t%E\n",fNumerator->GetEntries());
 104   stemp += ctemp;
 105   sprintf(ctemp,"Number of entries in denominator:\t%E\n",fDenominator->GetEntries());
 106   stemp += ctemp;
 107   sprintf(ctemp,"Number of entries in ratio:\t%E\n",fRatio->GetEntries());
 108   stemp += ctemp;
 109   //  stemp += mCoulombWeight->Report();
 110   AliFemtoString returnThis = stemp;
 111   return returnThis;
 112 }
 113 //____________________________
 114 void AliFemtoQinvCorrFctn::AddRealPair(AliFemtoPair* pair){
 115   // add true pair
 116   double tQinv = fabs(pair->QInv());   // note - qInv() will be negative for identical pairs...
 117   fNumerator->Fill(tQinv);
 118   //  cout << "AliFemtoQinvCorrFctn::AddRealPair : " << pair->qInv() << " " << tQinv <<
 119   //" " << pair->track1().FourMomentum() << " " << pair->track2().FourMomentum() << endl;
 120 }
 121 //____________________________
 122 void AliFemtoQinvCorrFctn::AddMixedPair(AliFemtoPair* pair){
 123   // add mixed (background) pair
 124   double weight = 1.0;
 125   double tQinv = fabs(pair->QInv());   // note - qInv() will be negative for identical pairs...
 126   fDenominator->Fill(tQinv,weight);
 127 }
 128
 129