]>
Commit | Line | Data |
---|---|---|
bf21645b | 1 | /* $Id$ */ |
2 | ||
3 | // ------------------------------------------------------ | |
4 | // | |
5 | // Class to handle 2d-corrections. | |
6 | // | |
7 | // ------------------------------------------------------ | |
8 | // | |
9 | ||
10 | #include <TH2F.h> | |
69b09e3b | 11 | #include <TMath.h> |
bf21645b | 12 | |
13 | #include <AliLog.h> | |
14 | ||
15 | #include "AliCorrectionMatrix2D.h" | |
16 | ||
17 | //____________________________________________________________________ | |
18 | ClassImp(AliCorrectionMatrix2D) | |
19 | ||
20 | //____________________________________________________________________ | |
21 | AliCorrectionMatrix2D::AliCorrectionMatrix2D() : | |
22 | AliCorrectionMatrix() | |
23 | { | |
24 | // default constructor | |
25 | } | |
26 | ||
27 | //____________________________________________________________________ | |
28 | AliCorrectionMatrix2D::AliCorrectionMatrix2D(const AliCorrectionMatrix2D& c) | |
29 | : AliCorrectionMatrix(c) | |
30 | { | |
31 | // copy constructor | |
32 | ((AliCorrectionMatrix2D &)c).Copy(*this); | |
33 | } | |
34 | ||
61385583 | 35 | //____________________________________________________________________ |
36 | AliCorrectionMatrix2D &AliCorrectionMatrix2D::operator=(const AliCorrectionMatrix2D &c) | |
37 | { | |
38 | // assigment operator | |
39 | ||
40 | if (this != &c) | |
41 | ((AliCorrectionMatrix2D &) c).Copy(*this); | |
42 | ||
43 | return *this; | |
44 | } | |
45 | ||
bf21645b | 46 | //____________________________________________________________________ |
47 | AliCorrectionMatrix2D::AliCorrectionMatrix2D(const Char_t* name, const Char_t* title, | |
48 | Int_t nBinX, Float_t Xmin, Float_t Xmax, | |
49 | Int_t nBinY, Float_t Ymin, Float_t Ymax) | |
50 | : AliCorrectionMatrix(name, title) | |
51 | { | |
52 | // | |
53 | // constructor | |
54 | // | |
55 | ||
29771dc8 | 56 | // do not add this hists to the directory |
57 | Bool_t oldStatus = TH1::AddDirectoryStatus(); | |
58 | TH1::AddDirectory(kFALSE); | |
59 | ||
82ac9fc4 | 60 | fhMeas = new TH2F("measured", Form("%s measured", GetTitle()), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax); |
61 | fhGene = new TH2F("generated", Form("%s generated", GetTitle()), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax); | |
62 | fhCorr = new TH2F("correction", Form("%s correction", GetTitle()), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax); | |
29771dc8 | 63 | |
64 | TH1::AddDirectory(oldStatus); | |
bf21645b | 65 | |
66 | fhMeas->Sumw2(); | |
67 | fhGene->Sumw2(); | |
68 | fhCorr->Sumw2(); | |
69 | } | |
70 | ||
71 | //____________________________________________________________________ | |
72 | AliCorrectionMatrix2D::AliCorrectionMatrix2D(const Char_t* name, const Char_t* title, | |
73 | Int_t nBinX, Float_t *X, Int_t nBinY, Float_t *Y) | |
74 | : AliCorrectionMatrix(name, title) | |
75 | { | |
76 | // | |
77 | // constructor | |
78 | // | |
79 | ||
29771dc8 | 80 | // do not add this hists to the directory |
81 | Bool_t oldStatus = TH1::AddDirectoryStatus(); | |
82 | TH1::AddDirectory(kFALSE); | |
bf21645b | 83 | |
29771dc8 | 84 | fhMeas = new TH2F("measured", Form("%s measured",title), nBinX, X, nBinY, Y); |
85 | fhGene = new TH2F("generated", Form("%s generated",title), nBinX, X, nBinY, Y); | |
86 | fhCorr = new TH2F("correction", Form("%s correction",title), nBinX, X, nBinY, Y); | |
87 | ||
88 | TH1::AddDirectory(oldStatus); | |
89 | ||
90 | fhMeas->Sumw2(); | |
bf21645b | 91 | fhGene->Sumw2(); |
92 | fhCorr->Sumw2(); | |
93 | } | |
94 | ||
95 | //____________________________________________________________________ | |
96 | AliCorrectionMatrix2D::~AliCorrectionMatrix2D() | |
97 | { | |
98 | // | |
99 | // destructor | |
100 | // | |
101 | ||
102 | // histograms already deleted in base class | |
103 | } | |
104 | ||
5a6310fe | 105 | TH2* AliCorrectionMatrix2D::GetGeneratedHistogram() const |
bf21645b | 106 | { |
107 | // return generated histogram casted to correct type | |
5a6310fe | 108 | return dynamic_cast<TH2*> (fhGene); |
bf21645b | 109 | } |
110 | ||
5a6310fe | 111 | TH2* AliCorrectionMatrix2D::GetMeasuredHistogram() const |
bf21645b | 112 | { |
113 | // return measured histogram casted to correct type | |
5a6310fe | 114 | return dynamic_cast<TH2*> (fhMeas); |
bf21645b | 115 | } |
116 | ||
117 | //____________________________________________________________________ | |
69b09e3b | 118 | TH1* AliCorrectionMatrix2D::Get1DCorrectionHistogram(const Char_t* opt, Float_t min, Float_t max, Bool_t binomialErrors) |
bf21645b | 119 | { |
120 | // | |
121 | // integrate the correction over one variable | |
122 | // | |
123 | ||
124 | TH1D* meas1D = 0; | |
125 | TH1D* gene1D = 0; | |
126 | ||
127 | if (strcmp(opt,"x")==0) { | |
06e4b91b | 128 | Int_t binMin = fhMeas->GetYaxis()->FindBin(min); |
129 | Int_t binMax = fhGene->GetYaxis()->FindBin(max); | |
7584d357 | 130 | |
06e4b91b | 131 | if (min>=max) { |
132 | meas1D = ((TH2F*)fhMeas)->ProjectionX(); | |
133 | gene1D = ((TH2F*)fhGene)->ProjectionX(); | |
7584d357 | 134 | } |
135 | else { | |
dd367a14 | 136 | Printf("Getting 1D map. Including y-bins %d to %d", binMin, binMax); |
7584d357 | 137 | |
5a6310fe | 138 | meas1D = ((TH2F*)fhMeas)->ProjectionX(Form("%s_x_pm", GetName()),binMin,binMax); |
139 | gene1D = ((TH2F*)fhGene)->ProjectionX(Form("%s_x_pg", GetName()),binMin,binMax); | |
7584d357 | 140 | } |
bf21645b | 141 | } |
5a6310fe | 142 | else if (strcmp(opt,"y")==0) { |
06e4b91b | 143 | Int_t binMin = fhMeas->GetXaxis()->FindBin(min); |
144 | Int_t binMax = fhMeas->GetXaxis()->FindBin(max); | |
dd367a14 | 145 | |
06e4b91b | 146 | if (min>=max) { |
147 | meas1D = ((TH2F*)fhMeas)->ProjectionY(); | |
148 | gene1D = ((TH2F*)fhGene)->ProjectionY(); | |
7584d357 | 149 | } |
150 | else { | |
5a6310fe | 151 | Printf("Getting 1D map. Including x-bins %d to %d \n", binMin, binMax); |
dd367a14 | 152 | |
5a6310fe | 153 | meas1D = ((TH2F*)fhMeas)->ProjectionY(Form("%s_y_pm", GetName()), binMin, binMax); |
154 | gene1D = ((TH2F*)fhGene)->ProjectionY(Form("%s_y_pg", GetName()), binMin, binMax); | |
7584d357 | 155 | } |
bf21645b | 156 | } |
5a6310fe | 157 | else { |
158 | Printf("ERROR: Invalid option"); | |
159 | return 0; | |
160 | } | |
161 | ||
69b09e3b | 162 | if (!binomialErrors) |
163 | { | |
164 | // set the errors on gene manually, and clear the ones on meas. | |
165 | gene1D->Sumw2(); | |
166 | for (Int_t bin=0; bin <= gene1D->GetNbinsX()+1; bin++) | |
167 | { | |
168 | gene1D->SetBinError(bin, TMath::Sqrt(gene1D->GetBinContent(bin))); | |
169 | meas1D->SetBinError(bin, 0); | |
170 | } | |
171 | } | |
172 | ||
bf21645b | 173 | gene1D->SetName(Form("corr_1D_%s",fName.Data())); |
174 | gene1D->SetTitle(Form("corr_1D_%s",fName.Data())); | |
51f6de65 | 175 | |
176 | TH1* divided = (TH1*) gene1D->Clone(Form("corr_1D_%s",fName.Data())); | |
177 | divided->Reset(); | |
69b09e3b | 178 | |
179 | divided->Divide(gene1D, meas1D, 1, 1, (binomialErrors) ? "B" : ""); | |
bf21645b | 180 | |
51f6de65 | 181 | Printf("%p %p", gene1D, meas1D); |
bf21645b | 182 | |
51f6de65 | 183 | return (TH1F*)divided; |
bf21645b | 184 | } |
185 | ||
186 | //____________________________________________________________________ | |
187 | void AliCorrectionMatrix2D::FillMeas(Float_t ax, Float_t ay) | |
188 | { | |
189 | // add value to measured histogram | |
06e4b91b | 190 | ((TH2F*)fhMeas)->Fill(ax, ay); |
bf21645b | 191 | } |
192 | ||
193 | //____________________________________________________________________ | |
194 | void AliCorrectionMatrix2D::FillGene(Float_t ax, Float_t ay) | |
195 | { | |
196 | // add value to generated histogram | |
06e4b91b | 197 | ((TH2F*)fhGene)->Fill(ax, ay); |
bf21645b | 198 | } |
199 | ||
200 | //____________________________________________________________________ | |
201 | Float_t AliCorrectionMatrix2D::GetCorrection(Float_t ax, Float_t ay) const | |
202 | { | |
203 | // returns a value of the correction map | |
204 | return fhCorr->GetBinContent(fhCorr->FindBin(ax,ay)); | |
205 | } | |
206 | ||
207 | //____________________________________________________________________ | |
208 | void AliCorrectionMatrix2D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge) | |
209 | { | |
210 | // remove edges of correction histogram by removing | |
211 | // - bins with content less than cut | |
212 | // - bins next to bins with zero bin content | |
213 | ||
214 | Int_t nBinsX = fhCorr->GetNbinsX(); | |
215 | Int_t nBinsY = fhCorr->GetNbinsY(); | |
216 | ||
217 | // set bin content to zero for bins with content smaller cut | |
218 | for (Int_t bx=0; bx<=nBinsX; bx++) { | |
219 | for (Int_t by=0; by<=nBinsY; by++) { | |
220 | if (fhCorr->GetBinContent(bx,by)>cut) { | |
221 | fhCorr->SetBinContent(bx,by,0); | |
222 | fhCorr->SetBinError(bx,by,0); | |
223 | } | |
224 | } | |
225 | } | |
226 | ||
227 | // set bin content to zero for bins next to bins with zero | |
228 | TH2F* tmp = (TH2F*)fhCorr->Clone("tmp"); | |
229 | tmp->Reset(); | |
230 | ||
231 | Bool_t done = kFALSE; | |
232 | Int_t nBinsXCount = 0; | |
233 | Int_t nBinsYCount = 0; | |
234 | while (!done) { | |
235 | if (nBinsXCount<nBinsXedge) | |
236 | for (Int_t bx=0; bx<=nBinsX; bx++) { | |
237 | for (Int_t by=0; by<=nBinsY; by++) { | |
238 | if ((fhCorr->GetBinContent(bx+1,by)==0)|| | |
239 | (fhCorr->GetBinContent(bx-1,by)==0)) | |
240 | tmp->SetBinContent(bx,by,1); | |
241 | ||
242 | } | |
243 | } | |
244 | if (nBinsYCount<nBinsYedge) | |
245 | for (Int_t bx=0; bx<=nBinsX; bx++) { | |
246 | for (Int_t by=0; by<=nBinsY; by++) { | |
247 | if ((fhCorr->GetBinContent(bx,by+1)==0)|| | |
248 | (fhCorr->GetBinContent(bx,by-1)==0)) | |
249 | tmp->SetBinContent(bx,by,1); | |
250 | } | |
251 | } | |
252 | for (Int_t bx=0; bx<=nBinsX; bx++) { | |
253 | for (Int_t by=0; by<=nBinsY; by++) { | |
254 | if (tmp->GetBinContent(bx,by)==1) { | |
255 | fhCorr->SetBinContent(bx,by,0); | |
256 | fhCorr->SetBinError(bx,by,0); | |
257 | } | |
258 | } | |
259 | } | |
260 | nBinsXCount++; | |
261 | nBinsYCount++; | |
262 | if ((nBinsXCount>=nBinsXedge)&&(nBinsYCount>=nBinsYedge)) done=kTRUE; | |
263 | } | |
264 | tmp->Delete(); | |
265 | ||
266 | } | |
267 | ||
745d6088 | 268 | //____________________________________________________________________ |
269 | void AliCorrectionMatrix2D::Rebin(Int_t x, Int_t y) | |
270 | { | |
271 | // rebins the histograms, recalculates the correction | |
272 | ||
273 | GetGeneratedHistogram()->Rebin2D(x, y); | |
274 | GetMeasuredHistogram()->Rebin2D(x, y); | |
275 | GetCorrectionHistogram()->Rebin2D(x, y); | |
276 | Divide(); | |
277 | } |