56ee055ad1d736aed231286a691532be0c95c728
[u/mrichter/AliRoot.git] / PWG0 / AliCorrectionMatrix3D.cxx
1 /* $Id$ */
2
3 // ------------------------------------------------------
4 //
5 // Class to handle 3d-corrections.
6 //
7 // ------------------------------------------------------
8 //
9
10 #include <TH3F.h>
11 #include <TH1F.h>
12
13 #include <AliLog.h>
14
15 #include "AliCorrectionMatrix3D.h"
16 #include "AliPWG0Helper.h"
17
18 //____________________________________________________________________
19 ClassImp(AliCorrectionMatrix3D)
20
21 //____________________________________________________________________
22 AliCorrectionMatrix3D::AliCorrectionMatrix3D() :
23   AliCorrectionMatrix()
24 {
25   // default constructor
26 }
27
28 //____________________________________________________________________
29 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const AliCorrectionMatrix3D& c)
30   : AliCorrectionMatrix(c)
31 {
32   // copy constructor
33   ((AliCorrectionMatrix3D &)c).Copy(*this);
34 }
35
36 //____________________________________________________________________
37 AliCorrectionMatrix3D &AliCorrectionMatrix3D::operator=(const AliCorrectionMatrix3D &c)
38 {
39   // assigment operator
40
41   if (this != &c)
42     ((AliCorrectionMatrix3D &) c).Copy(*this);
43
44   return *this;
45 }
46
47 //____________________________________________________________________
48 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
49               Int_t nBinX, Float_t Xmin, Float_t Xmax,
50               Int_t nBinY, Float_t Ymin, Float_t Ymax,
51               Int_t nBinZ, Float_t Zmin, Float_t Zmax)
52   : AliCorrectionMatrix(name, title)
53 {
54   //
55   // constructor
56   //
57
58   Float_t* binLimitsX = new Float_t[nBinX+1];
59   for (Int_t i=0; i<=nBinX; ++i)
60     binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;
61
62   Float_t* binLimitsY = new Float_t[nBinY+1];
63   for (Int_t i=0; i<=nBinY; ++i)
64     binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
65
66   Float_t* binLimitsZ = new Float_t[nBinZ+1];
67   for (Int_t i=0; i<=nBinZ; ++i)
68     binLimitsZ[i] = Zmin + (Zmax - Zmin) / nBinZ * i;
69
70   CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
71
72   delete[] binLimitsX;
73   delete[] binLimitsY;
74   delete[] binLimitsZ;
75 }
76
77 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
78       Int_t nBinX, Float_t Xmin, Float_t Xmax,
79       Int_t nBinY, Float_t Ymin, Float_t Ymax,
80       Int_t nBinZ, const Float_t* zbins)
81   : AliCorrectionMatrix(name, title)
82 {
83   // constructor with variable bin sizes
84
85   Float_t* binLimitsX = new Float_t[nBinX+1];
86   for (Int_t i=0; i<=nBinX; ++i)
87     binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;
88
89   Float_t* binLimitsY = new Float_t[nBinY+1];
90   for (Int_t i=0; i<=nBinY; ++i)
91     binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
92
93   CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);
94
95   delete[] binLimitsX;
96   delete[] binLimitsY;
97 }
98
99 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
100       Int_t nBinX, const Float_t* xbins,
101       Int_t nBinY, Float_t Ymin, Float_t Ymax,
102       Int_t nBinZ, const Float_t* zbins)
103   : AliCorrectionMatrix(name, title)
104 {
105   // constructor with variable bin sizes
106
107   Float_t* binLimitsY = new Float_t[nBinY+1];
108   for (Int_t i=0; i<=nBinY; ++i)
109     binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
110
111   CreateHists(nBinX, xbins, nBinY, binLimitsY, nBinZ, zbins);
112
113   delete[] binLimitsY;
114 }
115
116 //____________________________________________________________________
117 void AliCorrectionMatrix3D::CreateHists(Int_t nBinX, const Float_t* binLimitsX,
118       Int_t nBinY, const Float_t* binLimitsY,
119       Int_t nBinZ, const Float_t* binLimitsZ)
120 {
121   // create the histograms
122
123   // do not add this hists to the directory
124   Bool_t oldStatus = TH1::AddDirectoryStatus();
125   TH1::AddDirectory(kFALSE);
126
127   fhMeas  = new TH3F("measured",   Form("%s measured",GetTitle()),   nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
128   fhGene  = new TH3F("generated",  Form("%s generated",GetTitle()),  nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
129   fhCorr  = new TH3F("correction", Form("%s correction",GetTitle()), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
130
131   fhMeas->Sumw2();
132   fhGene->Sumw2();
133   fhCorr->Sumw2();
134
135   TH1::AddDirectory(oldStatus);
136 }
137
138
139 //____________________________________________________________________
140 AliCorrectionMatrix3D::~AliCorrectionMatrix3D()
141 {
142   //
143   // destructor
144   //
145
146   // histograms already deleted in base class
147 }
148
149 //____________________________________________________________________
150 TH3F* AliCorrectionMatrix3D::GetGeneratedHistogram()
151 {
152   // return generated histogram casted to correct type
153   return dynamic_cast<TH3F*> (fhGene);
154 }
155
156 //____________________________________________________________________
157 TH3F* AliCorrectionMatrix3D::GetMeasuredHistogram()
158 {
159   // return measured histogram casted to correct type
160   return dynamic_cast<TH3F*> (fhMeas);
161 }
162
163 //____________________________________________________________________
164 TH3F* AliCorrectionMatrix3D::GetCorrectionHistogram()
165 {
166   // return correction histogram casted to correct type
167   return dynamic_cast<TH3F*> (fhCorr);
168 }
169
170 //____________________________________________________________________
171 void AliCorrectionMatrix3D::FillMeas(Float_t ax, Float_t ay, Float_t az)
172 {
173   // add value to measured histogram
174   GetMeasuredHistogram()->Fill(ax, ay, az);
175 }
176
177 //____________________________________________________________________
178 void AliCorrectionMatrix3D::FillGene(Float_t ax, Float_t ay, Float_t az)
179 {
180   // add value to generated histogram
181   GetGeneratedHistogram()->Fill(ax, ay, az);
182 }
183
184 //____________________________________________________________________
185 Float_t AliCorrectionMatrix3D::GetCorrection(Float_t ax, Float_t ay, Float_t az) const
186 {
187   // returns a value of the correction map
188   return fhCorr->GetBinContent(fhCorr->FindBin(ax, ay, az));
189 }
190
191 //____________________________________________________________________
192 //void AliCorrectionMatrix3D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge, Int_t nBinsZedge)
193 void AliCorrectionMatrix3D::RemoveEdges(Float_t, Int_t, Int_t, Int_t)
194 {
195   // so what do we do here...
196 }
197
198 //____________________________________________________________________
199 void AliCorrectionMatrix3D::SaveHistograms()
200 {
201   //
202   // saves the histograms
203   //
204
205   AliCorrectionMatrix::SaveHistograms();
206
207   if (GetGeneratedHistogram() && GetMeasuredHistogram())
208   {
209     gDirectory->cd(GetName());
210     
211     AliPWG0Helper::CreateDividedProjections(GetGeneratedHistogram(), GetMeasuredHistogram(), 0, kFALSE, kTRUE);
212     
213     gDirectory->cd("..");
214   }
215 }
216
217 //____________________________________________________________________
218 Int_t AliCorrectionMatrix3D::CheckEmptyBins(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax, Bool_t quiet)
219 {
220   //
221   // counts the number of empty Bins in a given region
222   //
223
224   TH3F* hist = GetGeneratedHistogram();
225   if (!hist)
226     return -1;
227
228   Int_t emptyBins = 0;
229   for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
230     for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
231       for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
232         if (hist->GetBinContent(x, y, z) == 0)
233         {
234           if (!quiet)
235             printf("Empty bin in %s at vtx = %f, eta = %f, pt = %f\n", GetName(), hist->GetXaxis()->GetBinCenter(x), hist->GetYaxis()->GetBinCenter(y), hist->GetZaxis()->GetBinCenter(z));
236           ++emptyBins;
237         }
238
239   return emptyBins;
240 }
241
242 //____________________________________________________________________
243 TH1F* AliCorrectionMatrix3D::PlotBinErrors(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax)
244 {
245   //
246   // makes a 1d plots of the relative bin errors
247   //
248
249   TH3F* hist = GetCorrectionHistogram();
250   if (!hist)
251     return 0;
252     
253   TH1F* target = new TH1F("relerrors", "relerrors", 100, 0, 10);
254
255   for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
256     for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
257       for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
258         if (hist->GetBinContent(x, y, z) != 0)
259           target->Fill(100 * hist->GetBinError(x, y, z) / hist->GetBinContent(x, y, z));
260
261   return target;
262 }
263