PWG2/FEMTOSCOPY/AliFemtoUser/AliFemtoChi2CorrFctn.cxx

   1 ////////////////////////////////////////////////////////////////////////////////
   2 ///                                                                          ///
   3 /// AliFemtoChi2CorrFctn - A correlation function that saves the correlation ///
   4 /// function as a function of single track quality (chi2/ndof) for its and   ///
   5 /// tpc                                                                      ///
   6 /// Authors: Adam Kisiel kisiel@mps.ohio-state.edu                           ///
   7 ///                                                                          ///
   8 ////////////////////////////////////////////////////////////////////////////////
   9
  10 #include "AliFemtoChi2CorrFctn.h"
  11 //#include "AliFemtoHisto.hh"
  12 #include <cstdio>
  13
  14 #ifdef __ROOT__
  15 ClassImp(AliFemtoChi2CorrFctn)
  16 #endif
  17
  18 //____________________________
  19 AliFemtoChi2CorrFctn::AliFemtoChi2CorrFctn(char* title, const int& nbins, const float& QinvLo, const float& QinvHi):
  20   AliFemtoCorrFctn(),
  21   fChi2ITSSUMNumerator(0),
  22   fChi2ITSSUMDenominator(0),
  23   fChi2TPCSUMNumerator(0),
  24   fChi2TPCSUMDenominator(0),
  25   fChi2ITSONENumerator(0),
  26   fChi2ITSONEDenominator(0),
  27   fChi2TPCONENumerator(0),
  28   fChi2TPCONEDenominator(0),
  29   fSigmaToVertexNumerator(0),
  30   fSigmaToVertexDenominator(0)
  31 {
  32   // set up numerator
  33   char tTitNum[100] = "NumChi2ITSSUM";
  34   strcat(tTitNum,title);
  35   fChi2ITSSUMNumerator = new TH2D(tTitNum,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  36   // set up denominator
  37   char tTitDen[100] = "DenChi2ITSSUM";
  38   strcat(tTitDen,title);
  39   fChi2ITSSUMDenominator = new TH2D(tTitDen,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  40
  41   // set up numerator
  42   char tTit2Num[100] = "NumChi2TPCSUM";
  43   strcat(tTit2Num,title);
  44   fChi2TPCSUMNumerator = new TH2D(tTit2Num,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  45   // set up denominator
  46   char tTit2Den[100] = "DenChi2TPCSUM";
  47   strcat(tTit2Den,title);
  48   fChi2TPCSUMDenominator = new TH2D(tTit2Den,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  49
  50   // to enable error bar calculation...
  51   fChi2ITSSUMNumerator->Sumw2();
  52   fChi2ITSSUMDenominator->Sumw2();
  53
  54   fChi2TPCSUMNumerator->Sumw2();
  55   fChi2TPCSUMDenominator->Sumw2();
  56   // set up numerator
  57   sprintf(tTitNum,"%s%s","NumChi2ITSONE",title);
  58   fChi2ITSONENumerator = new TH2D(tTitNum,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  59   // set up denominator
  60   sprintf(tTitDen,"%s%s", "DenChi2ITSONE", title);
  61   fChi2ITSONEDenominator = new TH2D(tTitDen,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  62
  63   // set up numerator
  64   sprintf(tTit2Num,"%s%s","NumChi2TPCONE",title);
  65   fChi2TPCONENumerator = new TH2D(tTit2Num,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  66   // set up denominator
  67   sprintf(tTit2Den,"%s%s", "DenChi2TPCONE", title);
  68   fChi2TPCONEDenominator = new TH2D(tTit2Den,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  69
  70   // set up numerator
  71   sprintf(tTit2Num,"%s%s","NumSigmaToVertex",title);
  72   fSigmaToVertexNumerator = new TH2D(tTit2Num,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  73   // set up denominator
  74   sprintf(tTit2Den,"%s%s", "DenSigmaToVertex", title);
  75   fSigmaToVertexDenominator = new TH2D(tTit2Den,title,nbins,QinvLo,QinvHi,200,0.0,6.0);
  76
  77   // to enable error bar calculation...
  78   fChi2ITSONENumerator->Sumw2();
  79   fChi2ITSONEDenominator->Sumw2();
  80
  81   fChi2TPCONENumerator->Sumw2();
  82   fChi2TPCONEDenominator->Sumw2();
  83 }
  84
  85 //____________________________
  86 AliFemtoChi2CorrFctn::AliFemtoChi2CorrFctn(const AliFemtoChi2CorrFctn& aCorrFctn) :
  87   AliFemtoCorrFctn(),
  88   fChi2ITSSUMNumerator(0),
  89   fChi2ITSSUMDenominator(0),
  90   fChi2TPCSUMNumerator(0),
  91   fChi2TPCSUMDenominator(0),
  92   fChi2ITSONENumerator(0),
  93   fChi2ITSONEDenominator(0),
  94   fChi2TPCONENumerator(0),
  95   fChi2TPCONEDenominator(0),
  96   fSigmaToVertexNumerator(0),
  97   fSigmaToVertexDenominator(0)
  98 {
  99   // copy constructor
 100   if (aCorrFctn.fChi2ITSSUMNumerator)
 101     fChi2ITSSUMNumerator = new TH2D(*aCorrFctn.fChi2ITSSUMNumerator);
 102   if (aCorrFctn.fChi2ITSSUMDenominator)
 103     fChi2ITSSUMDenominator = new TH2D(*aCorrFctn.fChi2ITSSUMDenominator);
 104   if (aCorrFctn.fChi2TPCSUMNumerator)
 105     fChi2TPCSUMNumerator = new TH2D(*aCorrFctn.fChi2TPCSUMNumerator);
 106   if (aCorrFctn.fChi2TPCSUMDenominator)
 107     fChi2TPCSUMDenominator = new TH2D(*aCorrFctn.fChi2TPCSUMDenominator);
 108   if (aCorrFctn.fChi2ITSONENumerator)
 109     fChi2ITSONENumerator = new TH2D(*aCorrFctn.fChi2ITSONENumerator);
 110   if (aCorrFctn.fChi2ITSONEDenominator)
 111     fChi2ITSONEDenominator = new TH2D(*aCorrFctn.fChi2ITSONEDenominator);
 112   if (aCorrFctn.fChi2TPCONENumerator)
 113     fChi2TPCONENumerator = new TH2D(*aCorrFctn.fChi2TPCONENumerator);
 114   if (aCorrFctn.fChi2TPCONEDenominator)
 115     fChi2TPCONEDenominator = new TH2D(*aCorrFctn.fChi2TPCONEDenominator);
 116   if (aCorrFctn.fSigmaToVertexNumerator)
 117     fSigmaToVertexNumerator = new TH2D(*aCorrFctn.fSigmaToVertexNumerator);
 118   if (aCorrFctn.fSigmaToVertexDenominator)
 119     fSigmaToVertexDenominator = new TH2D(*aCorrFctn.fSigmaToVertexDenominator);
 120 }
 121 //____________________________
 122 AliFemtoChi2CorrFctn::~AliFemtoChi2CorrFctn(){
 123   // destructor
 124   delete fChi2ITSSUMNumerator;
 125   delete fChi2ITSSUMDenominator;
 126   delete fChi2TPCSUMNumerator;
 127   delete fChi2TPCSUMDenominator;
 128   delete fChi2ITSONENumerator;
 129   delete fChi2ITSONEDenominator;
 130   delete fChi2TPCONENumerator;
 131   delete fChi2TPCONEDenominator;
 132   delete fSigmaToVertexNumerator;
 133   delete fSigmaToVertexDenominator;
 134 }
 135 //_________________________
 136 AliFemtoChi2CorrFctn& AliFemtoChi2CorrFctn::operator=(const AliFemtoChi2CorrFctn& aCorrFctn)
 137 {
 138   // assignment operator
 139   if (this == &aCorrFctn)
 140     return *this;
 141
 142   if (aCorrFctn.fChi2ITSSUMNumerator)
 143     fChi2ITSSUMNumerator = new TH2D(*aCorrFctn.fChi2ITSSUMNumerator);
 144   else
 145     fChi2ITSSUMNumerator = 0;
 146   if (aCorrFctn.fChi2ITSSUMDenominator)
 147     fChi2ITSSUMDenominator = new TH2D(*aCorrFctn.fChi2ITSSUMDenominator);
 148   else
 149     fChi2ITSSUMDenominator = 0;
 150   if (aCorrFctn.fChi2TPCSUMNumerator)
 151     fChi2TPCSUMNumerator = new TH2D(*aCorrFctn.fChi2TPCSUMNumerator);
 152   else
 153     fChi2TPCSUMNumerator = 0;
 154   if (aCorrFctn.fChi2TPCSUMDenominator)
 155     fChi2TPCSUMDenominator = new TH2D(*aCorrFctn.fChi2TPCSUMDenominator);
 156   else
 157     fChi2TPCSUMDenominator = 0;
 158   if (aCorrFctn.fChi2ITSONENumerator)
 159     fChi2ITSONENumerator = new TH2D(*aCorrFctn.fChi2ITSONENumerator);
 160   else
 161     fChi2ITSONENumerator = 0;
 162   if (aCorrFctn.fChi2ITSONEDenominator)
 163     fChi2ITSONEDenominator = new TH2D(*aCorrFctn.fChi2ITSONEDenominator);
 164   else
 165     fChi2ITSONEDenominator = 0;
 166   if (aCorrFctn.fChi2TPCONENumerator)
 167     fChi2TPCONENumerator = new TH2D(*aCorrFctn.fChi2TPCONENumerator);
 168   else
 169     fChi2TPCONENumerator = 0;
 170   if (aCorrFctn.fChi2TPCONEDenominator)
 171     fChi2TPCONEDenominator = new TH2D(*aCorrFctn.fChi2TPCONEDenominator);
 172   else
 173     fChi2TPCONEDenominator = 0;
 174   if (aCorrFctn.fSigmaToVertexNumerator)
 175     fSigmaToVertexNumerator = new TH2D(*aCorrFctn.fSigmaToVertexNumerator);
 176   else
 177     fSigmaToVertexNumerator = 0;
 178   if (aCorrFctn.fSigmaToVertexDenominator)
 179     fSigmaToVertexDenominator = new TH2D(*aCorrFctn.fSigmaToVertexDenominator);
 180   else
 181     fSigmaToVertexDenominator = 0;
 182
 183   return *this;
 184 }
 185 //_________________________
 186 void AliFemtoChi2CorrFctn::Finish(){
 187   // here is where we should normalize, fit, etc...
 188   // we should NOT Draw() the histos (as I had done it below),
 189   // since we want to insulate ourselves from root at this level
 190   // of the code.  Do it instead at root command line with browser.
 191   //  mShareNumerator->Draw();
 192   //mShareDenominator->Draw();
 193   //mRatio->Draw();
 194
 195 }
 196
 197 //____________________________
 198 AliFemtoString AliFemtoChi2CorrFctn::Report(){
 199   // create report
 200   string stemp = "ITS and TPC quality Correlation Function Report:\n";
 201   char ctemp[100];
 202   sprintf(ctemp,"Number of entries in numerator:\t%E\n",fChi2ITSSUMNumerator->GetEntries());
 203   stemp += ctemp;
 204   sprintf(ctemp,"Number of entries in denominator:\t%E\n",fChi2ITSSUMDenominator->GetEntries());
 205   stemp += ctemp;
 206   //  stemp += mCoulombWeight->Report();
 207   AliFemtoString returnThis = stemp;
 208   return returnThis;
 209 }
 210 //____________________________
 211 void AliFemtoChi2CorrFctn::AddRealPair( AliFemtoPair* pair){
 212   // add real (effect) pair
 213   double tQinv = fabs(pair->QInv());   // note - qInv() will be negative for identical pairs...
 214
 215   fChi2ITSSUMNumerator->Fill(tQinv,
 216                           (pair->Track1()->Track()->ITSchi2() +
 217                            pair->Track2()->Track()->ITSchi2())/
 218                           (pair->Track1()->Track()->ITSncls() +
 219                            pair->Track2()->Track()->ITSncls()));
 220   fChi2TPCSUMNumerator->Fill(tQinv,
 221                           (pair->Track1()->Track()->TPCchi2() +
 222                            pair->Track2()->Track()->TPCchi2())/
 223                           (pair->Track1()->Track()->TPCncls() +
 224                            pair->Track2()->Track()->TPCncls()));
 225   if ((pair->Track1()->Track()->ITSchi2()/pair->Track1()->Track()->ITSncls()) > (pair->Track2()->Track()->ITSchi2()/pair->Track2()->Track()->ITSncls())) {
 226     fChi2ITSONENumerator->Fill(tQinv,
 227                             (pair->Track1()->Track()->ITSchi2()/
 228                              pair->Track1()->Track()->ITSncls()));
 229   }
 230   else {
 231     fChi2ITSONENumerator->Fill(tQinv,
 232                             (pair->Track2()->Track()->ITSchi2()/
 233                              pair->Track2()->Track()->ITSncls()));
 234   }
 235   if ((pair->Track1()->Track()->TPCchi2()/pair->Track1()->Track()->TPCncls()) > (pair->Track2()->Track()->TPCchi2()/pair->Track2()->Track()->TPCncls())) {
 236     fChi2TPCONENumerator->Fill(tQinv,
 237                             (pair->Track1()->Track()->TPCchi2()/
 238                              pair->Track1()->Track()->TPCncls()));
 239   }
 240   else {
 241     fChi2TPCONENumerator->Fill(tQinv,
 242                             (pair->Track2()->Track()->TPCchi2()/
 243                              pair->Track2()->Track()->TPCncls()));
 244   }
 245   if (pair->Track1()->Track()->SigmaToVertex() > pair->Track2()->Track()->SigmaToVertex()) {
 246     fSigmaToVertexNumerator->Fill(tQinv,
 247                                   pair->Track1()->Track()->SigmaToVertex());
 248   }
 249   else {
 250     fSigmaToVertexNumerator->Fill(tQinv,
 251                                   pair->Track2()->Track()->SigmaToVertex());
 252   }
 253 }
 254 //____________________________
 255 void AliFemtoChi2CorrFctn::AddMixedPair( AliFemtoPair* pair){
 256   // add mixed (background) pair
 257   double tQinv = fabs(pair->QInv());   // note - qInv() will be negative for identical pairs...
 258
 259   fChi2ITSSUMDenominator->Fill(tQinv,
 260                             (pair->Track1()->Track()->ITSchi2() +
 261                              pair->Track2()->Track()->ITSchi2())/
 262                             (pair->Track1()->Track()->ITSncls() +
 263                              pair->Track2()->Track()->ITSncls()));
 264   fChi2TPCSUMDenominator->Fill(tQinv,
 265                             (pair->Track1()->Track()->TPCchi2() +
 266                              pair->Track2()->Track()->TPCchi2())/
 267                             (pair->Track1()->Track()->TPCncls() +
 268                              pair->Track2()->Track()->TPCncls()));
 269   if ((pair->Track1()->Track()->ITSchi2()/pair->Track1()->Track()->ITSncls()) > (pair->Track2()->Track()->ITSchi2()/pair->Track2()->Track()->ITSncls())) {
 270     fChi2ITSONEDenominator->Fill(tQinv,
 271                             (pair->Track1()->Track()->ITSchi2()/
 272                              pair->Track1()->Track()->ITSncls()));
 273   }
 274   else {
 275     fChi2ITSONEDenominator->Fill(tQinv,
 276                             (pair->Track2()->Track()->ITSchi2()/
 277                              pair->Track2()->Track()->ITSncls()));
 278   }
 279   if ((pair->Track1()->Track()->TPCchi2()/pair->Track1()->Track()->TPCncls()) > (pair->Track2()->Track()->TPCchi2()/pair->Track2()->Track()->TPCncls())) {
 280     fChi2TPCONEDenominator->Fill(tQinv,
 281                             (pair->Track1()->Track()->TPCchi2()/
 282                              pair->Track1()->Track()->TPCncls()));
 283   }
 284   else {
 285     fChi2TPCONEDenominator->Fill(tQinv,
 286                             (pair->Track2()->Track()->TPCchi2()/
 287                              pair->Track2()->Track()->TPCncls()));
 288   }
 289   if (pair->Track1()->Track()->SigmaToVertex() > pair->Track2()->Track()->SigmaToVertex()) {
 290     fSigmaToVertexDenominator->Fill(tQinv,
 291                                   pair->Track1()->Track()->SigmaToVertex());
 292   }
 293   else {
 294     fSigmaToVertexDenominator->Fill(tQinv,
 295                                   pair->Track2()->Track()->SigmaToVertex());
 296   }
 297 }
 298
 299
 300 void AliFemtoChi2CorrFctn::WriteHistos()
 301 {
 302   // Write out result histograms
 303   fChi2ITSSUMNumerator->Write();
 304   fChi2ITSSUMDenominator->Write();
 305   fChi2TPCSUMNumerator->Write();
 306   fChi2TPCSUMDenominator->Write();
 307   fChi2ITSONENumerator->Write();
 308   fChi2ITSONEDenominator->Write();
 309   fChi2TPCONENumerator->Write();
 310   fChi2TPCONEDenominator->Write();
 311   fSigmaToVertexNumerator->Write();
 312   fSigmaToVertexDenominator->Write();
 313
 314 }
 315
 316 TList* AliFemtoChi2CorrFctn::GetOutputList()
 317 {
 318   // Prepare the list of objects to be written to the output
 319   TList *tOutputList = new TList();
 320
 321   tOutputList->Add(fChi2ITSSUMNumerator);
 322   tOutputList->Add(fChi2ITSSUMDenominator);
 323   tOutputList->Add(fChi2TPCSUMNumerator);
 324   tOutputList->Add(fChi2TPCSUMDenominator);
 325   tOutputList->Add(fChi2ITSONENumerator);
 326   tOutputList->Add(fChi2ITSONEDenominator);
 327   tOutputList->Add(fChi2TPCONENumerator);
 328   tOutputList->Add(fChi2TPCONEDenominator);
 329   tOutputList->Add(fSigmaToVertexNumerator);
 330   tOutputList->Add(fSigmaToVertexDenominator);
 331
 332   return tOutputList;
 333
 334 }