Added debug flags
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliFMDCorrections.cxx
1 #include "AliFMDCorrections.h"
2 #include <AliESDFMD.h>
3 #include <TAxis.h>
4 #include <TList.h>
5 #include <TMath.h>
6 #include "AliFMDAnaParameters.h"
7 #include "AliLog.h"
8 #include <TH2D.h>
9
10 ClassImp(AliFMDCorrections)
11 #if 0
12 ; // For Emacs
13 #endif 
14
15 //____________________________________________________________________
16 AliFMDCorrections::AliFMDCorrections()
17   : TNamed(), 
18     fRingHistos(),
19     fMultCut(0.3),
20     fDebug(0)
21 {}
22
23 //____________________________________________________________________
24 AliFMDCorrections::AliFMDCorrections(const char* title)
25   : TNamed("fmdCorrections", title), 
26     fRingHistos(), 
27     fMultCut(0.3),
28     fDebug(0)
29 {
30   fRingHistos.SetName(GetName());
31   fRingHistos.Add(new RingHistos(1, 'I'));
32   fRingHistos.Add(new RingHistos(2, 'I'));
33   fRingHistos.Add(new RingHistos(2, 'O'));
34   fRingHistos.Add(new RingHistos(3, 'I'));
35   fRingHistos.Add(new RingHistos(3, 'O'));
36 }
37
38 //____________________________________________________________________
39 AliFMDCorrections::AliFMDCorrections(const AliFMDCorrections& o)
40   : TNamed(o), 
41     fRingHistos(), 
42     fMultCut(o.fMultCut),
43     fDebug(o.fDebug)
44 {
45   TIter    next(&o.fRingHistos);
46   TObject* obj = 0;
47   while ((obj = next())) fRingHistos.Add(obj);
48 }
49
50 //____________________________________________________________________
51 AliFMDCorrections::~AliFMDCorrections()
52 {
53   fRingHistos.Delete();
54 }
55
56 //____________________________________________________________________
57 AliFMDCorrections&
58 AliFMDCorrections::operator=(const AliFMDCorrections& o)
59 {
60   TNamed::operator=(o);
61
62   fMultCut = o.fMultCut;
63   fDebug   = o.fDebug;
64   fRingHistos.Delete();
65   TIter    next(&o.fRingHistos);
66   TObject* obj = 0;
67   while ((obj = next())) fRingHistos.Add(obj);
68   
69   return *this;
70 }
71
72 //____________________________________________________________________
73 AliFMDCorrections::RingHistos*
74 AliFMDCorrections::GetRingHistos(UShort_t d, Char_t r) const
75 {
76   Int_t idx = -1;
77   switch (d) { 
78   case 1: idx = 0; break;
79   case 2: idx = 1 + (r == 'I' || r == 'i' ? 0 : 1); break;
80   case 3: idx = 3 + (r == 'I' || r == 'i' ? 0 : 1); break;
81   }
82   if (idx < 0 || idx >= fRingHistos.GetEntries()) return 0;
83   
84   return static_cast<RingHistos*>(fRingHistos.At(idx));
85 }
86     
87 //____________________________________________________________________
88 Bool_t
89 AliFMDCorrections::Correct(AliForwardUtil::Histos& hists,
90                            Int_t                   vtxbin)
91 {
92   AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
93
94   for (UShort_t d=1; d<=3; d++) { 
95     UShort_t nr = (d == 1 ? 1 : 2);
96     for (UShort_t q=0; q<nr; q++) { 
97       Char_t      r = (q == 0 ? 'I' : 'O');
98       TH2D*       h = hists.Get(d,r);
99       RingHistos* rh= GetRingHistos(d,r);
100       TH2F*       bg= pars->GetBackgroundCorrection(d, r, vtxbin);
101       TH2F*       ef= pars->GetEventSelectionEfficiency("INEL",vtxbin,r);
102       if (!bg) { 
103         AliWarning(Form("No secondary correction for FMDM%d%c in vertex bin %d",
104                         d, r, vtxbin));
105         continue;
106       }
107       if (!ef) { 
108         AliWarning(Form("No event %s selection efficiency in vertex bin %d",
109                         "INEL", vtxbin));
110         continue;
111       }
112
113       // Divide by primary/total ratio
114       h->Divide(bg);
115       
116       // Divide by the event selection efficiency 
117       h->Divide(ef);
118
119       if(!pars->SharingEffPresent()) { 
120         AliWarning("No sharing efficiencies");
121         continue;
122       }
123
124       TH1F* sf = pars->GetSharingEfficiencyTrVtx(d,r,vtxbin); 
125       
126       for (Int_t ieta = 1; ieta <= h->GetNbinsX(); ieta++) {
127         Float_t c  = sf->GetBinContent(ieta);
128         Float_t ec = sf->GetBinError(ieta);
129         
130         if (c == 0) continue;
131
132         for (Int_t iphi = 1; iphi <= h->GetNbinsY(); iphi++) { 
133           Double_t m  = h->GetBinContent(ieta, iphi) / c;
134           Double_t em = h->GetBinError(ieta, iphi);
135           
136           Double_t e  = TMath::Sqrt(em * em + (m * ec) * (m * ec));
137           
138           h->SetBinContent(ieta,iphi,m);
139           h->SetBinError(ieta,iphi,e);
140         }
141       }
142
143       rh->fDensity->Add(h);
144     }
145   }
146   
147   return kTRUE;
148 }
149
150 //____________________________________________________________________
151 void
152 AliFMDCorrections::ScaleHistograms(TList* dir, Int_t nEvents)
153 {
154   if (nEvents <= 0) return;
155   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
156   if (!d) return;
157
158   TIter    next(&fRingHistos);
159   RingHistos* o = 0;
160   while ((o = static_cast<RingHistos*>(next())))
161     o->ScaleHistograms(d, nEvents);
162 }
163 //____________________________________________________________________
164 void
165 AliFMDCorrections::DefineOutput(TList* dir)
166 {
167   TList* d = new TList;
168   d->SetName(GetName());
169   dir->Add(d);
170   TIter    next(&fRingHistos);
171   RingHistos* o = 0;
172   while ((o = static_cast<RingHistos*>(next()))) {
173     o->Output(d);
174   }
175 }
176
177 //====================================================================
178 AliFMDCorrections::RingHistos::RingHistos()
179   : AliForwardUtil::RingHistos(), 
180     fDensity(0)
181 {}
182
183 //____________________________________________________________________
184 AliFMDCorrections::RingHistos::RingHistos(UShort_t d, Char_t r)
185   : AliForwardUtil::RingHistos(d,r), 
186     fDensity(0)
187 {
188   fDensity = new TH2D(Form("FMD%d%c_Primary_Density", d, r), 
189                       Form("in FMD%d%c", d, r), 
190                       200, -4, 6, (r == 'I' || r == 'i' ? 20 : 40), 
191                       0, 2*TMath::Pi());
192   fDensity->SetDirectory(0);
193   fDensity->SetXTitle("#eta");
194   fDensity->SetYTitle("#phi [radians]");
195   fDensity->SetZTitle("Primary N_{ch} density");
196 }
197 //____________________________________________________________________
198 AliFMDCorrections::RingHistos::RingHistos(const RingHistos& o)
199   : AliForwardUtil::RingHistos(o), 
200     fDensity(o.fDensity)
201 {}
202
203 //____________________________________________________________________
204 AliFMDCorrections::RingHistos&
205 AliFMDCorrections::RingHistos::operator=(const RingHistos& o)
206 {
207   AliForwardUtil::RingHistos::operator=(o);
208   
209   if (fDensity) delete fDensity;
210   
211   fDensity = static_cast<TH2D*>(o.fDensity->Clone());
212   
213   return *this;
214 }
215 //____________________________________________________________________
216 AliFMDCorrections::RingHistos::~RingHistos()
217 {
218   if (fDensity) delete fDensity;
219 }
220
221 //____________________________________________________________________
222 void
223 AliFMDCorrections::RingHistos::Output(TList* dir)
224 {
225   TList* d = DefineOutputList(dir);
226   d->Add(fDensity);
227 }
228
229 //____________________________________________________________________
230 void
231 AliFMDCorrections::RingHistos::ScaleHistograms(TList* dir, Int_t nEvents)
232
233   TList* l = GetOutputList(dir);
234   if (!l) return; 
235
236   TH1* density = GetOutputHist(l,Form("%s_Primary_Density", fName.Data()));
237   if (density) density->Scale(1./nEvents);
238 }
239
240 //____________________________________________________________________
241 //
242 // EOF
243 //
244           
245
246