Major dielectron framework update; includes "alignment" to updates in
[u/mrichter/AliRoot.git] / PWG3 / dielectron / AliDielectronCFdraw.cxx
1 /*************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3 *                                                                        *
4 * Author: The ALICE Off-line Project.                                    *
5 * Contributors are mentioned in the code where appropriate.              *
6 *                                                                        *
7 * Permission to use, copy, modify and distribute this software and its   *
8 * documentation strictly for non-commercial purposes is hereby granted   *
9 * without fee, provided that the above copyright notice appears in all   *
10 * copies and that both the copyright notice and this permission notice   *
11 * appear in the supporting documentation. The authors make no claims     *
12 * about the suitability of this software for any purpose. It is          *
13 * provided "as is" without express or implied warranty.                  *
14 **************************************************************************/
15
16 ///////////////////////////////////////////////////////////////////////////
17 //       Dielectron Correction framework draw helper                     //
18 //                                                                       //
19 /*
20
21
22
23
24
25
26
27
28
29 */
30 //                                                                       //
31 ///////////////////////////////////////////////////////////////////////////
32
33 #include <TSeqCollection.h>
34 #include <TObjArray.h>
35 #include <TKey.h>
36 #include <TList.h>
37 #include <TClass.h>
38 #include <TObject.h>
39 #include <TVirtualPS.h>
40 #include <TFile.h>
41 #include <TString.h>
42 #include <TObjString.h>
43 #include <TVectorD.h>
44 #include <TMath.h>
45 #include <TH1.h>
46 #include <TH2.h>
47 #include <TH3.h>
48 #include <TPad.h>
49 #include <TCanvas.h>
50 #include <TLegend.h>
51 #include <AliCFEffGrid.h>
52
53 #include <AliLog.h>
54
55 #include "AliDielectronCFdraw.h"
56
57 ClassImp(AliDielectronCFdraw)
58
59 AliDielectronCFdraw::AliDielectronCFdraw() :
60   TNamed(),
61   fCfContainer(0x0),
62   fEffGrid(0x0),
63   fVdata(0)
64 {
65   //
66   // Ctor
67   //
68 }
69
70 //________________________________________________________________
71 AliDielectronCFdraw::AliDielectronCFdraw(const char*name, const char* title) :
72   TNamed(name,title),
73   fCfContainer(0x0),
74   fEffGrid(0x0),
75   fVdata(0)
76 {
77   //
78   // Named Ctor
79   //
80   
81 }
82
83 //________________________________________________________________
84 AliDielectronCFdraw::AliDielectronCFdraw(AliCFContainer *cont) :
85   TNamed(cont->GetName(), cont->GetTitle()),
86   fCfContainer(cont),
87   fEffGrid(new AliCFEffGrid("eff","eff",*cont)),
88   fVdata(0)
89 {
90   //
91   // directly set the CF container
92   //
93
94 }
95
96 //________________________________________________________________
97 AliDielectronCFdraw::AliDielectronCFdraw(const char* filename) :
98   TNamed(),
99   fCfContainer(0x0),
100   fEffGrid(0x0),
101   fVdata(0)
102 {
103   //
104   // get CF container(s) from file 'filename'
105   //
106   SetCFContainers(filename);
107 }
108
109 //________________________________________________________________
110 void AliDielectronCFdraw::SetCFContainers(const TSeqCollection *arr)
111 {
112   //
113   // Merge CF Container out of several containers
114   //
115
116   TIter next(arr);
117   TObject *o=0x0;
118
119   Int_t nstep=0;
120   while ( (o=next()) ){
121     AliCFContainer *cf=dynamic_cast<AliCFContainer*>(o);
122     if (!cf) continue;
123     nstep+=cf->GetNStep();
124   }
125   if (nstep==0) return;
126   Int_t nbins[1]={1};
127   fCfContainer=new AliCFContainer(GetName(), GetTitle(), nstep, 1, nbins);
128
129   //delete unneeded steps
130 //   for (Int_t istep=0; istep<nstep; ++istep) delete fCfContainer->GetGrid(istep);
131
132   //add step to the new container
133   Int_t istep=0;
134   for (Int_t icf=0; icf<arr->GetEntries(); ++icf){
135     AliCFContainer *cf=dynamic_cast<AliCFContainer*>(arr->At(icf));
136     if (!cf) continue;
137     for (Int_t istepCurr=0; istepCurr<cf->GetNStep(); ++istepCurr){
138       fCfContainer->SetGrid(istep, cf->GetGrid(istepCurr));
139       fCfContainer->SetStepTitle(istep,Form("%s, Pair: %s",cf->GetTitle(),cf->GetStepTitle(istepCurr)));
140       ++istep;
141     }
142   }
143   if (fEffGrid) delete fEffGrid;
144   fEffGrid=new AliCFEffGrid("eff","eff",*fCfContainer);
145 }
146
147 //________________________________________________________________
148 void AliDielectronCFdraw::SetCFContainers(const char* filename)
149 {
150   //
151   // get CF containers from file
152   //
153
154   TFile f(filename);
155   TList *l=f.GetListOfKeys();
156   TIter nextKey(l);
157   TKey *k=0x0;
158   while ( (k=static_cast<TKey*>(nextKey())) ){
159     TObject *o=k->ReadObj();
160     if (o->IsA()->InheritsFrom(TSeqCollection::Class())){
161       TSeqCollection *arr=static_cast<TSeqCollection*>(o);
162       SetCFContainers(arr);
163     } else if (o->IsA()==AliCFContainer::Class()){
164       fCfContainer=static_cast<AliCFContainer*>(o);
165       if (fEffGrid) delete fEffGrid;
166       fEffGrid=new AliCFEffGrid("eff","eff",*fCfContainer);
167     }
168   }
169 }
170
171 //________________________________________________________________
172 void AliDielectronCFdraw::SetRangeUser(Int_t ivar, Double_t min, Double_t max, const char* slices)
173 {
174   //
175   // Set range of cut steps defined in slices
176   // Steps may be separated by one the the characteres ,;:
177   //
178   TObjArray *arr=TString(slices).Tokenize(",:;");
179
180   if (arr->GetEntriesFast()==0){
181     // all slices in case of 0 entries
182     for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
183       fCfContainer->GetGrid(istep)->SetRangeUser(ivar,min,max);
184     }
185   } else {
186     TIter next(arr);
187     TObjString *ostr=0x0;
188     while ( (ostr=static_cast<TObjString*>(next())) ) {
189       Int_t istep=ostr->GetString().Atoi();
190       fCfContainer->GetGrid(istep)->SetRangeUser(ivar,min,max);
191     }
192   }
193   delete arr;
194 }
195
196 //________________________________________________________________
197 void AliDielectronCFdraw::UnsetRangeUser(Int_t ivar, const char* slices)
198 {
199   //
200   // Unset range of cut steps defined in slices
201   // Steps may be separated by one the the characteres ,;:
202   //
203   TObjArray *arr=TString(slices).Tokenize(",:;");
204   
205   if (arr->GetEntriesFast()==0){
206     // all slices in case of 0 entries
207     for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
208       fCfContainer->GetAxis(ivar,istep)->SetRange(0,0);
209     }
210   } else {
211     TIter next(arr);
212     TObjString *ostr=0x0;
213     while ( (ostr=static_cast<TObjString*>(next())) ) {
214       Int_t istep=ostr->GetString().Atoi();
215       fCfContainer->GetAxis(ivar,istep)->SetRange(0,0);
216     }
217   }
218   delete arr;
219 }
220
221 //________________________________________________________________
222 void AliDielectronCFdraw::Draw(const Option_t* varnames, const char* opt, const char* slices)
223 {
224   //
225   // Draw 'variables' of 'slices'
226   // for multidimensional draw variables may be separated by a ':'
227   // slice numbers may be separated by any of ,:;
228   //
229   // variables may be called by either their name or number
230   //
231
232   TObjArray *arrVars=TString(varnames).Tokenize(":");
233   Int_t entries=arrVars->GetEntriesFast();
234   if (entries<1||entries>3){
235     AliError("Wrong number of variables, supported are 1 - 3 dimensions");
236     delete arrVars;
237     return;
238   }
239   
240   TIter next(arrVars);
241   TObjString *ostr=0x0;
242   Int_t ivar[3]={-1,-1,-1};
243   for (Int_t i=entries-1; i>=0; --i){
244     ostr=static_cast<TObjString*>(next());
245     if (ostr->GetString().IsDigit()){
246       ivar[i]=ostr->GetString().Atoi();
247     } else {
248       ivar[i]=fCfContainer->GetVar(ostr->GetName());
249     }
250   }
251
252   Draw(ivar[0],ivar[1],ivar[2],opt,slices);
253   delete arrVars;
254 }
255
256 //________________________________________________________________
257 void AliDielectronCFdraw::Draw(Int_t var, const char* opt, const char* slices)
258 {
259   //
260   // Draw variable var for all slices
261   // slices may be divided by and of ,;:
262   //
263   // if opt contains 'same' all histograms are drawn in the same pad
264   // otherwise the pad will be divided in sub pads and the histograms
265   // are drawn in each sub pad
266   //
267
268   Int_t vars[3]={var,-1,-1};
269   TObjArray *arr=CollectHistosProj(vars,slices);
270   Draw(arr,opt);
271   delete arr; 
272 }
273
274 //________________________________________________________________
275 void AliDielectronCFdraw::Draw(Int_t var0, Int_t var1, const char* opt, const char* slices)
276 {
277   //
278   // Draw 2D case
279   //
280   Int_t vars[3]={var0,var1,-1};
281   TObjArray *arr=CollectHistosProj(vars,slices);
282   Draw(arr,opt);
283   delete arr;
284 }
285
286 //________________________________________________________________
287 void AliDielectronCFdraw::Draw(Int_t var0, Int_t var1, Int_t var2, const char* opt, const char* slices)
288 {
289   //
290   // Draw 3D case
291   //
292   Int_t vars[3]={var0,var1,var2};
293   TObjArray *arr=CollectHistosProj(vars,slices);
294   Draw(arr,opt);
295   delete arr;
296 }
297
298 //________________________________________________________________
299 TObjArray* AliDielectronCFdraw::CollectHistosProj(const Int_t vars[3], const char* slices)
300 {
301   //
302   // Collect histos with up to 3 dimension of the 'slices' separated by one of "':;'"
303   // in a TObjArray and return it
304   // if a dimension is not used it must be set to -1
305   //
306   TObjArray *arr=TString(slices).Tokenize(",:;");
307   TObjArray *arrHists=0x0;
308   if (arr->GetEntriesFast()==0){
309     // all slices in case of 0 entries
310     arrHists=new TObjArray(fCfContainer->GetNStep());
311     for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
312       TH1 *hproj=Project(vars,istep);
313       if (!hproj) continue;
314       hproj->SetName(Form("proj_%02d",istep));
315       hproj->SetTitle(fCfContainer->GetStepTitle(istep));
316       arrHists->Add(hproj);
317     }
318   } else {
319     arrHists=new TObjArray(arr->GetEntriesFast());
320     TIter next(arr);
321     TObjString *ostr=0x0;
322     while ( (ostr=static_cast<TObjString*>(next())) ) {
323       Int_t istep=ostr->GetString().Atoi();
324       TH1 *hproj=Project(vars,istep);
325       if (!hproj) continue;
326       hproj->SetName(Form("proj_%02d",istep));
327       hproj->SetTitle(fCfContainer->GetStepTitle(istep));
328       arrHists->Add(hproj);
329     }
330   }
331   delete arr;
332
333   return arrHists;
334 }
335
336 //________________________________________________________________
337 TH1* AliDielectronCFdraw::Project(const Int_t *vars, Int_t slice)
338 {
339   //
340   // Do an ndim projection
341   //
342   return fCfContainer->Project(slice,vars[0],vars[1],vars[2]);
343 }
344
345 //________________________________________________________________
346 TH1* AliDielectronCFdraw::Project(const Option_t* var, Int_t slice)
347 {
348   //
349   // translate variable names and do projection
350   //
351   TObjArray *arrVars=TString(var).Tokenize(":");
352   Int_t entries=arrVars->GetEntriesFast();
353   if (entries<1||entries>3){
354     AliError("Wrong number of variables, supported are 1 - 3 dimensions");
355     delete arrVars;
356     return 0x0;
357   }
358   
359   TObjString *ostr=0x0;
360   Int_t ivar[3]={-1,-1,-1};
361   for (Int_t i=entries-1; i>=0; --i){
362     ostr=static_cast<TObjString*>(arrVars->At(i));
363     if (ostr->GetString().IsDigit()){
364       ivar[i]=ostr->GetString().Atoi();
365     } else {
366       ivar[i]=fCfContainer->GetVar(ostr->GetName());
367     }
368   }
369   if (ivar[0]==-1) return 0x0;
370   delete arrVars;
371   return fCfContainer->Project(slice,ivar[0],ivar[1],ivar[2]);
372 }
373
374 //________________________________________________________________
375 void AliDielectronCFdraw::DrawEfficiency(const char* varnames, const char* numerators, Int_t denominator, const char* opt)
376 {
377   //
378   // plot efficiencies for variables. Variable may be up to 3 dim, separated by a ':'
379   // you may have several numerators, sparated by one of ',:;'
380   //
381   
382   TObjArray *arrVars=TString(varnames).Tokenize(":");
383   Int_t entries=arrVars->GetEntriesFast();
384   if (entries<1||entries>3){
385     AliError("Wrong number of variables, supported are 1 - 3 dimensions");
386     delete arrVars;
387     return;
388   }
389   
390   TIter next(arrVars);
391   TObjString *ostr=0x0;
392   Int_t ivar[3]={-1,-1,-1};
393   for (Int_t i=0; i<entries; ++i){
394     ostr=static_cast<TObjString*>(next());
395     if (ostr->GetString().IsDigit()){
396       ivar[i]=ostr->GetString().Atoi();
397     } else {
398       ivar[i]=fCfContainer->GetVar(ostr->GetName());
399     }
400   }
401
402   Int_t type=0;
403   TString optStr(opt);
404   if (optStr.Contains("2")) type=1;
405   
406   DrawEfficiency(ivar[2],ivar[1],ivar[0],numerators, denominator,opt,type);
407   delete arrVars;
408 }
409
410 //________________________________________________________________
411 void AliDielectronCFdraw::DrawEfficiency(Int_t var, const char* numerators, Int_t denominator, const char* opt, Int_t type)
412 {
413   //
414   // Draw Efficiencies for all numerators
415   // numerators may be divided by and of ,;:
416   //
417   // if opt contains 'same' all histograms are drawn in the same pad
418   // otherwise the pad will be divided in sub pads and the histograms
419   // are drawn in each sub pad
420   //
421   
422   Int_t vars[3]={var,-1,-1};
423   TObjArray *arr=CollectHistosEff(vars,numerators,denominator,type);
424   TString drawOpt=opt;
425   drawOpt+="eff";
426   Draw(arr,drawOpt);
427   delete arr;
428 }
429
430 //________________________________________________________________
431 void AliDielectronCFdraw::DrawEfficiency(Int_t var0, Int_t var1, const char* numerators, Int_t denominator, const char* opt, Int_t type)
432 {
433   //
434   // Draw 2D case
435   //
436   Int_t vars[3]={var0,var1,-1};
437   TObjArray *arr=CollectHistosEff(vars,numerators,denominator,type);
438   TString drawOpt=opt;
439   drawOpt+="eff";
440   Draw(arr,drawOpt);
441   delete arr;
442 }
443
444 //________________________________________________________________
445 void AliDielectronCFdraw::DrawEfficiency(Int_t var0, Int_t var1, Int_t var2, const char* numerators, Int_t denominator, const char* opt, Int_t type)
446 {
447   //
448   // Draw 3D case
449   //
450   Int_t vars[3]={var0,var1,var2};
451   TObjArray *arr=CollectHistosEff(vars,numerators,denominator,type);
452   TString drawOpt=opt;
453   drawOpt+="eff";
454   Draw(arr,drawOpt);
455   delete arr;
456 }
457
458 //________________________________________________________________
459 TObjArray* AliDielectronCFdraw::CollectHistosEff(const  Int_t vars[3], const char* numerators, Int_t denominator, Int_t type)
460 {
461   //
462   // Collect histos with 'dim'ension of the 'slices' separated by one of "':;'"
463   // in a TObjArray and return it
464   //
465   TObjArray *arr=TString(numerators).Tokenize(",:;");
466   TObjArray *arrHists=0x0;
467
468   if (type==0){
469     if (arr->GetEntriesFast()==0){
470     // all slices in case of 0 entries
471       arrHists=new TObjArray(fCfContainer->GetNStep());
472       fVdata.ResizeTo(arrHists->GetSize());
473       for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
474         fEffGrid->CalculateEfficiency(istep,denominator);
475         TH1 *hproj=ProjectEff(vars);
476         if (!hproj) continue;
477         Float_t eff=fEffGrid->GetAverage();
478         fVdata(istep)=eff;
479         hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
480         hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
481         arrHists->Add(hproj);
482       }
483     } else {
484       arrHists=new TObjArray(arr->GetEntriesFast());
485       TIter next(arr);
486       TObjString *ostr=0x0;
487       fVdata.ResizeTo(arrHists->GetSize());
488       Int_t count=0;
489       while ( (ostr=static_cast<TObjString*>(next())) ) {
490         Int_t istep=ostr->GetString().Atoi();
491         fEffGrid->CalculateEfficiency(istep,denominator);
492         TH1 *hproj=ProjectEff(vars);
493         if (!hproj) continue;
494         Float_t eff=fEffGrid->GetAverage();
495         fVdata(count++)=eff;
496         hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
497         hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
498         arrHists->Add(hproj);
499       }
500     }
501   }
502
503   //second approach
504   if (type==1){
505     TH1 *hDen=Project(vars,denominator);
506     Double_t entriesDen=hDen->GetEffectiveEntries();
507     if (arr->GetEntriesFast()==0){
508     // all slices in case of 0 entries
509       arrHists=new TObjArray(fCfContainer->GetNStep());
510       fVdata.ResizeTo(arrHists->GetSize());
511       for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
512         TH1 *hproj=Project(vars,istep);
513         if (!hproj) continue;
514         Float_t eff=0;
515         if (entriesDen>0) eff=hproj->GetEffectiveEntries()/entriesDen;
516         fVdata(istep)=eff;
517         hproj->Divide(hproj,hDen,1,1,"B");
518         hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
519         hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
520         arrHists->Add(hproj);
521       }
522     } else {
523       arrHists=new TObjArray(arr->GetEntriesFast());
524       fVdata.ResizeTo(arrHists->GetSize());
525       TIter next(arr);
526       TObjString *ostr=0x0;
527       Int_t count=0;
528       while ( (ostr=static_cast<TObjString*>(next())) ) {
529         Int_t istep=ostr->GetString().Atoi();
530         TH1 *hproj=Project(vars,istep);
531         if (!hproj) continue;
532         Float_t eff=0;
533         if (entriesDen>0) eff=hproj->GetEffectiveEntries()/entriesDen;
534         fVdata(count++)=eff;
535         hproj->Divide(hproj,hDen,1,1,"B");
536         hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
537         hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
538         arrHists->Add(hproj);
539       }
540     }
541     delete hDen;
542   }
543   
544
545   delete arr;
546   return arrHists;
547 }
548
549 //________________________________________________________________
550 TH1* AliDielectronCFdraw::ProjectEff(const Int_t vars[3])
551 {
552   //
553   // Do an nim projection
554   //
555   return fEffGrid->Project(vars[0],vars[1],vars[2]);
556 }
557
558 //________________________________________________________________
559 void AliDielectronCFdraw::Draw(const TObjArray *arr, const char* opt)
560 {
561   //
562   // draw all objects in arr
563   //
564   TString optStr(opt);
565   optStr.ToLower();
566   Bool_t drawSame     = optStr.Contains("same");
567   Bool_t drawSamePlus = optStr.Contains("same+");
568   Bool_t drawEff      = optStr.Contains("eff");
569   Bool_t optLeg       = optStr.Contains("leg");
570   Bool_t optScaleMax  = optStr.Contains("max");
571   
572   if (!drawSamePlus) optStr.ReplaceAll("same","");
573   
574   optStr.ReplaceAll("+","");
575   optStr.ReplaceAll("eff","");
576   optStr.ReplaceAll("leg","");
577   optStr.ReplaceAll("max","");
578   
579   if (!gPad) new TCanvas;
580   
581   Int_t nPads = arr->GetEntriesFast();
582   if (nPads==0) return;
583   
584 //   if (nPads==1){
585 //     arr->UncheckedAt(0)->Draw(optStr.Data());
586 //     return;
587 //   }
588   
589   TCanvas *c=gPad->GetCanvas();
590   if (!gVirtualPS&&!drawSamePlus&&!drawSame&&nPads>1) c->Clear();
591   
592   if (!drawSame&&nPads>1){
593     //optimised division
594     Int_t nCols = (Int_t)TMath::Ceil( TMath::Sqrt(nPads) );
595     Int_t nRows = (Int_t)TMath::Ceil( (Double_t)nPads/(Double_t)nCols );
596     c->Divide(nCols,nRows);
597     for (Int_t i=0; i<nPads; ++i){
598       c->cd(i+1);
599       arr->UncheckedAt(i)->Draw(optStr.Data());
600     }
601   } else {
602     TLegend *leg=0;
603     if (drawSamePlus){
604       //find already existing legend to attach entries to it
605       TIter nextPrimitive(gPad->GetListOfPrimitives());
606       TObject *o=0x0;
607       while ((o=nextPrimitive())) if (o->IsA()==TLegend::Class()) leg=(TLegend*)o;
608     }
609     
610     if (optLeg&&!leg) leg=new TLegend(.2,.3,.99,.9);
611     Int_t addColor=0;
612     if (leg) addColor=leg->GetNRows();
613     Int_t offset=20;
614     if (nPads<7) offset=24;
615     for (Int_t i=0; i<nPads; ++i){
616       if (i==1&&!drawSamePlus) optStr+="same";
617       TH1 *hist=static_cast<TH1*>(arr->UncheckedAt(i));
618       hist->SetLineColor(i+1+addColor);
619       hist->SetLineWidth(2);
620       hist->SetMarkerColor(i+1+addColor);
621       hist->SetMarkerStyle(offset+i+addColor);
622       hist->Draw(optStr.Data());
623       
624       if (drawEff&&i==0&&!drawSamePlus) {
625         hist->GetYaxis()->SetRangeUser(0,2);
626         hist->SetYTitle("Rec. Signal / Gen. Signal");
627       }
628       
629       if (leg) leg->AddEntry(hist,hist->GetTitle(),"lp");
630     }
631     if (leg){
632       leg->SetFillColor(10);
633       leg->SetY1(.9-leg->GetNRows()*.05);
634       leg->SetY1NDC(.9-leg->GetNRows()*.05);
635       leg->SetMargin(.1);
636       if (!drawSamePlus) leg->Draw();
637     }
638     if (optScaleMax){
639       TIter nextPrimitive(gPad->GetListOfPrimitives());
640       TObject *o=0x0;
641       TH1 *hfirst=0x0;
642       Float_t max=0;
643       while ((o=nextPrimitive())) if (o->InheritsFrom(TH1::Class())){
644         TH1 *h=(TH1*)o;
645         if (!hfirst) hfirst=h;
646         if (h->GetMaximum()>max) max=h->GetMaximum();
647       }
648       max*=1.1;
649       hfirst->SetMaximum(max);
650     }
651   }
652   
653 }
654
655 //________________________________________________________________
656 Double_t AliDielectronCFdraw::GetAverageEfficiency(Int_t numerator, Int_t denominator, Double_t &effErr)
657 {
658   //
659   // Extract the mean efficiency of the numerator and denominator
660   //
661
662   //use variable 0 as default, since for the average it doesn't matter
663   TH1 *hDen=fCfContainer->Project(denominator,0);
664   Double_t entriesDen=hDen->GetEffectiveEntries();
665   TH1 *hproj=fCfContainer->Project(numerator,0);
666   if (!hproj) return -1.;
667   Double_t entriesNum=hproj->GetEffectiveEntries();
668   if (entriesDen<1||entriesNum<1) return -1;
669   
670   Double_t eff=-1.;
671   if (entriesDen>0) eff=entriesNum/entriesDen;
672   effErr=TMath::Sqrt(1./entriesNum+1./entriesDen)*eff;
673   delete hDen;
674   delete hproj;
675   return eff;
676 }