PWG3/AliQuarkoniaEfficiency.cxx

   1 //===================================================================
   2 //  Class AliQuarkoniaEfficiency
   3 //
   4 //  This class will provide the quarkonia reconstruction efficiency
   5 //  in ALICE without acceptance consideration.
   6 //
   7 //
   8 //  Reconstruction efficiency has been evaluated by means of a flat
   9 //  y and pt distribution of quarkonia in -4 < y < -2.5,
  10 //  0 < pt < 20 GeV/c. Weights have been used to evaluate the
  11 //  reconstruction efficiency in different parameterizations.
  12 //
  13 //  Beware that efficiency histos are preliminary.
  14 //  Just Jpsi, Dimuon, UnlikePair trigger
  15 //  efficiencies should be considered.
  16 //
  17 //
  18 //  Example:
  19 //   Double_t eff,error;
  20 //   AliQuarkoniaEfficiency * JPsiEff = new AliQuarkoniaEfficiency();
  21 //   JPsiEff->Init();
  22 //   JPsiEff->GetEfficiency(Rapidity,Pt,eff,error);
  23 //   printf(" eff = %2.2e  error %2.2e \n",eff,error);
  24 //
  25 //
  26 //  Subatech 2006
  27 //===================================================================
  28
  29 // Root
  30 #include "TFile.h"
  31 #include "TH2.h"
  32
  33 // AliRoot includes
  34 #include "AliQuarkoniaEfficiency.h"
  35 #include "AliLog.h"
  36
  37
  38 //_______________________________________________________________________
  39 AliQuarkoniaEfficiency::AliQuarkoniaEfficiency(Int_t quarkoniaResonance, Int_t decayChannel,
  40                                                Int_t simParameterization):
  41   fEfficiencyFileName("QuarkoniaEfficiency.root"),
  42   fQuarkoniaResonance(quarkoniaResonance),
  43   fDecayChannel(decayChannel),
  44   fParameterization(simParameterization),
  45   fTriggerType(kPairUnlikeApt),
  46   fTrigger(kFALSE),
  47   fEfficiency(0x0)
  48 {
  49   // Constructor
  50 }
  51
  52 //_______________________________________________________________________
  53 AliQuarkoniaEfficiency::~AliQuarkoniaEfficiency()
  54 {
  55   // Destructor
  56   delete fEfficiency;
  57 }
  58
  59 //_______________________________________________________________________
  60 void AliQuarkoniaEfficiency::Init()
  61 {
  62   // Initialize method
  63   switch (fQuarkoniaResonance)  {
  64   case kJpsi:
  65     SetTitle("Jpsi");
  66     break;
  67   case kPsiP:
  68     SetTitle("PsiP");
  69     break;
  70   case kUpsilon:
  71     SetTitle("Upsilon");
  72     break;
  73   case kUpsilonP:
  74     SetTitle("UpsilonP");
  75     break;
  76   case kUpsilonPP:
  77     SetTitle("UpsilonPP");
  78     break;
  79   case kOmega:
  80     SetTitle("Omega");
  81     break;
  82   case kPhi:
  83     SetTitle("Phi");
  84     break;
  85   }
  86
  87   switch (fDecayChannel) {
  88   case kDimuon:
  89     SetName("Dimuon");
  90     break;
  91   case kDielectron:
  92     SetName("Dielectron");
  93     break;
  94   }
  95
  96   char *param=0;
  97   switch (fParameterization){
  98   case kFlat:
  99     param = "Flat";
 100     break;
 101   case kCDFscaled:
 102     param = "CDFscaled";
 103     break;
 104   case kCDFscaledPP:
 105     param = "CDFscaledPP";
 106     break;
 107   }
 108
 109   char *trig=0;
 110   switch (fTriggerType){
 111   case kSinglePlusLpt:
 112     trig = "SinglePlusLpt";
 113     break;
 114   case kSinglePlusHpt:
 115     trig = "SinglePlusHpt";
 116     break;
 117   case kSinglePlusApt:
 118     trig = "SinglePlusApt";
 119     break;
 120   case kSingleMinusLpt:
 121     trig = "SingleMinusLpt";
 122     break;
 123   case kSingleMinusHpt:
 124     trig = "SingleMinusHpt";
 125     break;
 126   case kSingleMinusApt:
 127     trig = "SingleMinusApt";
 128     break;
 129   case kSingleUndefLpt:
 130     trig = "SingleUndefLpt";
 131     break;
 132   case kSingleUndefHpt:
 133     trig = "SingleUndefHpt";
 134     break;
 135   case kSingleUndefApt:
 136     trig = "SingleUndefApt";
 137     break;
 138   case kPairUnlikeLpt:
 139     trig = "PairUnlikeLpt";
 140     break;
 141   case kPairUnlikeHpt:
 142     trig = "PairUnlikeHpt";
 143     break;
 144   case kPairUnlikeApt:
 145     trig = "PairUnlikeApt";
 146     break;
 147   case kPairLikeLpt:
 148     trig = "PairLikeLpt";
 149     break;
 150   case kPairLikeHpt:
 151     trig = "PairLikeHpt";
 152     break;
 153   case kPairLikeApt:
 154     trig = "PairLikeApt";
 155     break;
 156   }
 157
 158
 159   if(!fEfficiency) delete fEfficiency;
 160
 161   TFile efficiencyFile(fEfficiencyFileName);
 162   if ( efficiencyFile.IsOpen() ) {
 163
 164     char quarkoniaDir[15];
 165     sprintf(quarkoniaDir,"%s",GetTitle());
 166     if (! efficiencyFile.cd(quarkoniaDir) ){
 167       AliError(Form("Directory %s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found ",
 168                     quarkoniaDir,fEfficiencyFileName.Data(),GetTitle(),GetName() ));
 169       return;
 170     }
 171
 172     char histosDir[30];
 173     sprintf(histosDir,"%s/%s_%s_%s",quarkoniaDir,GetTitle(),GetName(),param);
 174     if(! efficiencyFile.cd(histosDir) ){
 175       AliError(Form("Subdirectory %s/%s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found ",
 176                     quarkoniaDir,histosDir,fEfficiencyFileName.Data(),GetTitle(),GetName() ));
 177       return;
 178     }
 179
 180     char histoname[50];
 181     if(fTrigger) sprintf(histoname,"h%sEfficiencyPtRap_%s",GetTitle(),trig);
 182     else sprintf(histoname,"h%sEfficiencyPtRap",GetTitle());
 183     char histonameposition[99];
 184     sprintf(histonameposition,"%s/%s",histosDir,histoname);
 185     fEfficiency = (TH2F*)efficiencyFile.Get(histonameposition);
 186
 187     if ( !fEfficiency ) {
 188       AliError(Form("Histo %s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found",
 189                     histoname,GetTitle(), GetName() ));
 190     }
 191     else {
 192       fEfficiency->SetDirectory(0);
 193     }
 194     efficiencyFile.Close();
 195
 196   }
 197   else {
 198     AliError(Form("File %s not found",fEfficiencyFileName.Data()));
 199   }
 200
 201 }
 202
 203 //_______________________________________________________________________
 204 TH2F*  AliQuarkoniaEfficiency::GetEfficiencyHisto() const
 205 {
 206   // Returns the efficiency histogram
 207   if (fEfficiency) return fEfficiency;
 208   else {
 209     AliError(Form("Efficiency data for quarkonia %s and channel %s not found",GetTitle(),GetName()));
 210     return 0x0;
 211   }
 212 }
 213
 214 //_______________________________________________________________________
 215 void  AliQuarkoniaEfficiency::GetEfficiency(Float_t rap, Float_t pT, Double_t &eff, Double_t &error)
 216 {
 217   // Evaluates the efficiency for a given (y,Pt) of the quarkonia
 218   Int_t binx=0;
 219   Int_t biny=0;
 220
 221   if (!fEfficiency) {
 222     AliError(Form("Efficiency data for quarkonia %s and channel %s not found",GetTitle(),GetName()));
 223   }
 224   else {
 225     if ( rap < (fEfficiency->GetXaxis())->GetXmin()  ||
 226          rap > (fEfficiency->GetXaxis())->GetXmax()  ||
 227          pT  < (fEfficiency->GetYaxis())->GetXmin()  ||
 228          pT  > (fEfficiency->GetYaxis())->GetXmax()   ) {
 229       AliInfo("Values out of range");
 230       eff   = 0.;
 231       error = 0.;
 232     }
 233     else  {
 234       binx  = fEfficiency->GetXaxis()->FindBin(rap);
 235       biny  = fEfficiency->GetYaxis()->FindBin(pT);
 236       eff   = fEfficiency->GetBinContent(binx,biny);
 237       error = fEfficiency->GetBinError(binx,biny);
 238     }
 239   }
 240 }