added stuff
[u/mrichter/AliRoot.git] / PWG0 / AliCorrectionMatrix2D.cxx
CommitLineData
bf21645b 1/* $Id$ */
2
3// ------------------------------------------------------
4//
5// Class to handle 2d-corrections.
6//
7// ------------------------------------------------------
8//
9
10#include <TH2F.h>
11
12#include <AliLog.h>
13
14#include "AliCorrectionMatrix2D.h"
15
16//____________________________________________________________________
17ClassImp(AliCorrectionMatrix2D)
18
19//____________________________________________________________________
20AliCorrectionMatrix2D::AliCorrectionMatrix2D() :
21 AliCorrectionMatrix()
22{
23 // default constructor
24}
25
26//____________________________________________________________________
27AliCorrectionMatrix2D::AliCorrectionMatrix2D(const AliCorrectionMatrix2D& c)
28 : AliCorrectionMatrix(c)
29{
30 // copy constructor
31 ((AliCorrectionMatrix2D &)c).Copy(*this);
32}
33
34//____________________________________________________________________
61385583 35AliCorrectionMatrix2D &AliCorrectionMatrix2D::operator=(const AliCorrectionMatrix2D &c)
36{
37 // assigment operator
38
39 if (this != &c)
40 ((AliCorrectionMatrix2D &) c).Copy(*this);
41
42 return *this;
43}
44
45//____________________________________________________________________
bf21645b 46AliCorrectionMatrix2D::AliCorrectionMatrix2D(const Char_t* name, const Char_t* title,
47 Int_t nBinX, Float_t Xmin, Float_t Xmax,
48 Int_t nBinY, Float_t Ymin, Float_t Ymax)
49 : AliCorrectionMatrix(name, title)
50{
51 //
52 // constructor
53 //
54
55 fhMeas = new TH2F(Form("meas_%s",name), Form("meas_%s",title), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax);
56 fhGene = new TH2F(Form("gene_%s",name), Form("gene_%s",title), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax);
57 fhCorr = new TH2F(Form("corr_%s",name), Form("corr_%s",title), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax);
58
59 fhMeas->Sumw2();
60 fhGene->Sumw2();
61 fhCorr->Sumw2();
62}
63
64//____________________________________________________________________
65AliCorrectionMatrix2D::AliCorrectionMatrix2D(const Char_t* name, const Char_t* title,
66 Int_t nBinX, Float_t *X, Int_t nBinY, Float_t *Y)
67 : AliCorrectionMatrix(name, title)
68{
69 //
70 // constructor
71 //
72
73 fhMeas = new TH2F(Form("meas_%s",name), Form("meas_%s",title), nBinX, X, nBinY, Y);
74 fhGene = new TH2F(Form("gene_%s",name), Form("gene_%s",title), nBinX, X, nBinY, Y);
75 fhCorr = new TH2F(Form("corr_%s",name), Form("corr_%s",title), nBinX, X, nBinY, Y);
76
77 fhMeas->Sumw2();
78 fhGene->Sumw2();
79 fhCorr->Sumw2();
80}
81
82//____________________________________________________________________
83AliCorrectionMatrix2D::~AliCorrectionMatrix2D()
84{
85 //
86 // destructor
87 //
88
89 // histograms already deleted in base class
90}
91
92TH2F* AliCorrectionMatrix2D::GetGeneratedHistogram() const
93{
94 // return generated histogram casted to correct type
95 return dynamic_cast<TH2F*> (fhGene);
96}
97
98TH2F* AliCorrectionMatrix2D::GetMeasuredHistogram() const
99{
100 // return measured histogram casted to correct type
101 return dynamic_cast<TH2F*> (fhMeas);
102}
103
104//____________________________________________________________________
7584d357 105TH1F* AliCorrectionMatrix2D::Get1DCorrection(Char_t* opt, Float_t min, Float_t max)
bf21645b 106{
107 //
108 // integrate the correction over one variable
109 //
110
111 TH1D* meas1D = 0;
112 TH1D* gene1D = 0;
113
114 if (strcmp(opt,"x")==0) {
7584d357 115 Int_t binMin = GetMeasuredHistogram()->GetYaxis()->FindBin(min);
116 Int_t binMax = GetMeasuredHistogram()->GetYaxis()->FindBin(max);
117
118 if (min==0 && max==0) {
119 meas1D = GetMeasuredHistogram()->ProjectionX();
120 gene1D = GetGeneratedHistogram()->ProjectionX();
121 }
122 else {
123 AliDebug(AliLog::kDebug+1, Form("Getting 1D map. Including y-bins %d to %d \n", binMin, binMax));
124
125 meas1D = GetMeasuredHistogram()->ProjectionX("pm",binMin,binMax);
126 gene1D = GetGeneratedHistogram()->ProjectionX("pg",binMin,binMax);
127 }
bf21645b 128 }
129 if (strcmp(opt,"y")==0) {
7584d357 130 Int_t binMin = GetMeasuredHistogram()->GetXaxis()->FindBin(min);
131 Int_t binMax = GetMeasuredHistogram()->GetXaxis()->FindBin(max);
132
133 if (min==0 && max==0) {
134 meas1D = GetMeasuredHistogram()->ProjectionY();
135 gene1D = GetGeneratedHistogram()->ProjectionY();
136 }
137 else {
138 AliDebug(AliLog::kDebug+1, Form("Getting 1D map. Including x-bins %d to %d \n", binMin, binMax));
139
140 meas1D = GetMeasuredHistogram()->ProjectionY("pm", binMin, binMax);
141 gene1D = GetGeneratedHistogram()->ProjectionY("pg", binMin, binMax);
142 }
bf21645b 143 }
144 gene1D->Sumw2();
145
146 gene1D->SetName(Form("corr_1D_%s",fName.Data()));
147 gene1D->SetTitle(Form("corr_1D_%s",fName.Data()));
148
149 gene1D->Divide(gene1D, meas1D, 1, 1, "B");
150
151 return (TH1F*)gene1D;
152}
153
154//____________________________________________________________________
155void AliCorrectionMatrix2D::FillMeas(Float_t ax, Float_t ay)
156{
157 // add value to measured histogram
158 GetMeasuredHistogram()->Fill(ax, ay);
159}
160
161//____________________________________________________________________
162void AliCorrectionMatrix2D::FillGene(Float_t ax, Float_t ay)
163{
164 // add value to generated histogram
165 GetGeneratedHistogram()->Fill(ax, ay);
166}
167
168//____________________________________________________________________
169Float_t AliCorrectionMatrix2D::GetCorrection(Float_t ax, Float_t ay) const
170{
171 // returns a value of the correction map
172 return fhCorr->GetBinContent(fhCorr->FindBin(ax,ay));
173}
174
175//____________________________________________________________________
176void AliCorrectionMatrix2D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge)
177{
178 // remove edges of correction histogram by removing
179 // - bins with content less than cut
180 // - bins next to bins with zero bin content
181
182 Int_t nBinsX = fhCorr->GetNbinsX();
183 Int_t nBinsY = fhCorr->GetNbinsY();
184
185 // set bin content to zero for bins with content smaller cut
186 for (Int_t bx=0; bx<=nBinsX; bx++) {
187 for (Int_t by=0; by<=nBinsY; by++) {
188 if (fhCorr->GetBinContent(bx,by)>cut) {
189 fhCorr->SetBinContent(bx,by,0);
190 fhCorr->SetBinError(bx,by,0);
191 }
192 }
193 }
194
195 // set bin content to zero for bins next to bins with zero
196 TH2F* tmp = (TH2F*)fhCorr->Clone("tmp");
197 tmp->Reset();
198
199 Bool_t done = kFALSE;
200 Int_t nBinsXCount = 0;
201 Int_t nBinsYCount = 0;
202 while (!done) {
203 if (nBinsXCount<nBinsXedge)
204 for (Int_t bx=0; bx<=nBinsX; bx++) {
205 for (Int_t by=0; by<=nBinsY; by++) {
206 if ((fhCorr->GetBinContent(bx+1,by)==0)||
207 (fhCorr->GetBinContent(bx-1,by)==0))
208 tmp->SetBinContent(bx,by,1);
209
210 }
211 }
212 if (nBinsYCount<nBinsYedge)
213 for (Int_t bx=0; bx<=nBinsX; bx++) {
214 for (Int_t by=0; by<=nBinsY; by++) {
215 if ((fhCorr->GetBinContent(bx,by+1)==0)||
216 (fhCorr->GetBinContent(bx,by-1)==0))
217 tmp->SetBinContent(bx,by,1);
218 }
219 }
220 for (Int_t bx=0; bx<=nBinsX; bx++) {
221 for (Int_t by=0; by<=nBinsY; by++) {
222 if (tmp->GetBinContent(bx,by)==1) {
223 fhCorr->SetBinContent(bx,by,0);
224 fhCorr->SetBinError(bx,by,0);
225 }
226 }
227 }
228 nBinsXCount++;
229 nBinsYCount++;
230 if ((nBinsXCount>=nBinsXedge)&&(nBinsYCount>=nBinsYedge)) done=kTRUE;
231 }
232 tmp->Delete();
233
234}
235