Interface class for quarkonia efficiencies (Zaida)
[u/mrichter/AliRoot.git] / PWG3 / AliQuarkoniaEfficiency.cxx
1
2 // Root 
3 #include "TAxis.h"
4 #include "TFile.h"
5 #include "TH2.h"
6 #include "TString.h"
7 #include "TDirectory.h"
8
9 // AliRoot includes
10 #include "AliQuarkoniaEfficiency.h"
11 #include "AliLog.h"
12
13
14 AliQuarkoniaEfficiency::AliQuarkoniaEfficiency(Int_t quarkoniaResonance, Int_t decayChannel,
15                                                Int_t simParameterization):
16   fEfficiencyFileName("QuarkoniaEfficiency.root"),
17   fQuarkoniaResonance(quarkoniaResonance),     
18   fDecayChannel(decayChannel),           
19   fParameterization(simParameterization),  
20   fTriggerType(kPairUnlikeApt),
21   fTrigger(kFALSE),
22   fEfficiency(0x0)
23 {
24
25 }
26
27 AliQuarkoniaEfficiency::~AliQuarkoniaEfficiency()
28 {
29   delete fEfficiency;
30 }
31
32 void AliQuarkoniaEfficiency::Init()
33 {
34   switch (fQuarkoniaResonance)  {  
35   case kJpsi:
36     SetTitle("Jpsi");
37     break;
38   case kPsiP:
39     SetTitle("PsiP");
40     break;
41   case kUpsilon:
42     SetTitle("Upsilon");
43     break;
44   case kUpsilonP:
45     SetTitle("UpsilonP");
46     break;
47   case kUpsilonPP:
48     SetTitle("UpsilonPP");
49     break;
50   case kOmega:
51     SetTitle("Omega");
52     break;
53   case kPhi:
54     SetTitle("Phi");
55     break;
56   }
57
58   switch (fDecayChannel) {
59   case kDimuon:
60     SetName("Dimuon");
61     break;
62   case kDielectron:
63     SetName("Dielectron");
64     break;
65   }
66
67   char *param=0;
68   switch (fParameterization){
69   case kFlat:
70     param = "Flat";
71     break;
72   case kCDFscaled:
73     param = "CDFscaled";
74     break;
75   case kCDFscaledPP:
76     param = "CDFscaledPP";
77     break;
78   }
79   
80   char *trig=0;
81   switch (fTriggerType){
82   case kSinglePlusLpt:
83     trig = "SinglePlusLpt";
84     break;
85   case kSinglePlusHpt: 
86     trig = "SinglePlusHpt";
87     break;
88   case kSinglePlusApt:
89     trig = "SinglePlusApt";
90     break;    
91   case kSingleMinusLpt:
92     trig = "SingleMinusLpt";
93     break;    
94   case kSingleMinusHpt:
95     trig = "SingleMinusHpt";
96     break;    
97   case kSingleMinusApt:
98     trig = "SingleMinusApt";
99     break;    
100   case kSingleUndefLpt:
101     trig = "SingleUndefLpt";
102     break;    
103   case kSingleUndefHpt:
104     trig = "SingleUndefHpt";
105     break;    
106   case kSingleUndefApt:
107     trig = "SingleUndefApt";
108     break;    
109   case kPairUnlikeLpt:
110     trig = "PairUnlikeLpt";
111     break;    
112   case kPairUnlikeHpt:
113     trig = "PairUnlikeHpt";
114     break;    
115   case kPairUnlikeApt:
116     trig = "PairUnlikeApt";
117     break;    
118   case kPairLikeLpt:
119     trig = "PairLikeLpt";
120     break;    
121   case kPairLikeHpt:
122     trig = "PairLikeHpt";
123     break;    
124   case kPairLikeApt:
125     trig = "PairLikeApt";
126     break;    
127   }
128
129
130   if(!fEfficiency) delete fEfficiency; 
131
132   TFile efficiencyFile(fEfficiencyFileName);
133   if ( efficiencyFile.IsOpen() ) {
134
135     char quarkoniaDir[15];
136     sprintf(quarkoniaDir,"%s",GetTitle());
137     if (! efficiencyFile.cd(quarkoniaDir) ){
138       AliError(Form("Directory %s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found ",
139                     quarkoniaDir,fEfficiencyFileName.Data(),GetTitle(),GetName() ));
140       return;
141     }
142     
143     char histosDir[30];
144     sprintf(histosDir,"%s/%s_%s_%s",quarkoniaDir,GetTitle(),GetName(),param);
145     if(! efficiencyFile.cd(histosDir) ){
146       AliError(Form("Subdirectory %s/%s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found ",
147                     quarkoniaDir,histosDir,fEfficiencyFileName.Data(),GetTitle(),GetName() ));
148       return;
149     }
150
151     char histoname[50];
152     if(fTrigger) sprintf(histoname,"h%sEfficiencyPtRap_%s",GetTitle(),trig);
153     else sprintf(histoname,"h%sEfficiencyPtRap",GetTitle());
154     char histonameposition[99];
155     sprintf(histonameposition,"%s/%s",histosDir,histoname);
156     fEfficiency = (TH2F*)efficiencyFile.Get(histonameposition);
157
158     if ( !fEfficiency ) {
159       AliError(Form("Histo %s not found in file %s \n Efficiency data for quarkonia %s and channel %s not found",
160                     histoname,GetTitle(), GetName() ));
161     }
162     else {
163       fEfficiency->SetDirectory(0);
164     }
165     efficiencyFile.Close();
166
167   }
168   else {
169     AliError(Form("File %s not found",fEfficiencyFileName.Data()));
170   }
171
172 }
173
174 TH2F*  AliQuarkoniaEfficiency::GetEfficiencyHisto() const
175 {
176   if (fEfficiency) return fEfficiency;
177   else {
178     AliError(Form("Efficiency data for quarkonia %s and channel %s not found",GetTitle(),GetName()));
179     return 0x0;
180   }
181 }
182
183 void  AliQuarkoniaEfficiency::GetEfficiency(Float_t rap, Float_t pT, Double_t &eff, Double_t &error)
184 {
185   Int_t binx=0;
186   Int_t biny=0;
187   
188   if (!fEfficiency) {
189     AliError(Form("Efficiency data for quarkonia %s and channel %s not found",GetTitle(),GetName()));
190   }
191   else {
192     if ( rap < (fEfficiency->GetXaxis())->GetXmin()  ||  
193          rap > (fEfficiency->GetXaxis())->GetXmax()  ||
194          pT  < (fEfficiency->GetYaxis())->GetXmin()  ||  
195          pT  > (fEfficiency->GetYaxis())->GetXmax()   ) {
196       AliInfo("Values out of range");
197       eff   = 0.;
198       error = 0.;
199     }
200     else  { 
201       binx  = fEfficiency->GetXaxis()->FindBin(rap);  
202       biny  = fEfficiency->GetYaxis()->FindBin(pT);
203       eff   = fEfficiency->GetBinContent(binx,biny);
204       error = fEfficiency->GetBinError(binx,biny);
205     }
206   } 
207 }