Logo ROOT   6.10/09
Reference Guide
TProfile2Poly.cxx
Go to the documentation of this file.
1 #include "TProfile2Poly.h"
2 #include "TProfileHelper.h"
3 
4 #include "TMultiGraph.h"
5 #include "TGraph.h"
6 #include "TClass.h"
7 #include "TList.h"
8 #include "TMath.h"
9 
10 #include <cassert>
11 #include <cmath>
12 
14 
15 // -------------- TProfile2PolyBin --------------
16 
18 {
19  fSumw = 0;
20  fSumvw = 0;
21  fSumw2 = 0;
22  fSumwv2 = 0;
23  fError = 0;
24  fAverage = 0;
25  fErrorMode = kERRORMEAN;
26 }
27 
28 TProfile2PolyBin::TProfile2PolyBin(TObject *poly, Int_t bin_number) : TH2PolyBin(poly, bin_number)
29 {
30  fSumw = 0;
31  fSumvw = 0;
32  fSumw2 = 0;
33  fSumwv2 = 0;
34  fError = 0;
35  fAverage = 0;
37 }
38 
40 {
41  this->fSumw += toMerge->fSumw;
42  this->fSumvw += toMerge->fSumvw;
43  this->fSumw2 += toMerge->fSumw2;
44  this->fSumwv2 += toMerge->fSumwv2;
45 }
46 
48 {
49  UpdateAverage();
50  UpdateError();
51  SetChanged(true);
52 }
53 
55 {
56  if (fSumw != 0) fAverage = fSumvw / fSumw;
57 }
58 
60 {
61  Double_t tmp = 0;
62  if (fSumw != 0) tmp = std::sqrt((fSumwv2 / fSumw) - (fAverage * fAverage));
63 
64  fError = tmp;
65 
66  return;
67 
68 }
69 
71 {
72  fSumw = 0;
73  fSumvw = 0;
74  fSumw2 = 0;
75  fSumwv2 = 0;
76  fError = 0;
77  fAverage = 0;
78 }
79 
81 {
82  fSumw += weight;
83  fSumvw += value * weight;
84  fSumw2 += weight * weight;
85  fSumwv2 += weight * value * value;
86  this->Update();
87 }
88 
89 // -------------- TProfile2Poly --------------
90 
91 TProfile2Poly::TProfile2Poly(const char *name, const char *title, Double_t xlow, Double_t xup, Double_t ylow,
92  Double_t yup)
93  : TH2Poly(name, title, xlow, xup, ylow, yup)
94 {
95 }
96 
97 TProfile2Poly::TProfile2Poly(const char *name, const char *title, Int_t nX, Double_t xlow, Double_t xup, Int_t nY,
98  Double_t ylow, Double_t yup)
99  : TH2Poly(name, title, nX, xlow, xup, nY, ylow, yup)
100 {
101 }
102 
104 {
105  if (!poly) return 0;
106 
107  if (fBins == 0) {
108  fBins = new TList();
109  fBins->SetOwner();
110  }
111 
112  fNcells++;
113  Int_t ibin = fNcells - kNOverflow;
114  return new TProfile2PolyBin(poly, ibin);
115 }
116 
118 {
119  return Fill(xcoord, ycoord, value, 1);
120 }
121 
123 {
124  // Find region in which the hit occured
125  Int_t tmp = GetOverflowRegionFromCoordinates(xcoord, ycoord);
126  Int_t overflow_idx = OverflowIdxToArrayIdx(tmp);
127  fOverflowBins[overflow_idx].Fill(value, weight);
128  fOverflowBins[overflow_idx].SetContent(fOverflowBins[overflow_idx].fAverage );
129 
130  // Find the cell to which (x,y) coordinates belong to
131  Int_t n = (Int_t)(floor((xcoord - fXaxis.GetXmin()) / fStepX));
132  Int_t m = (Int_t)(floor((ycoord - fYaxis.GetXmin()) / fStepY));
133 
134  // Make sure the array indices are correct.
135  if (n >= fCellX) n = fCellX - 1;
136  if (m >= fCellY) m = fCellY - 1;
137  if (n < 0) n = 0;
138  if (m < 0) m = 0;
139 
140  // ------------ Update global (per histo) statistics
141  fTsumw += weight;
142  fTsumw2 += weight * weight;
143  fTsumwx += weight * xcoord;
144  fTsumwx2 += weight * xcoord * xcoord;
145  fTsumwy += weight * ycoord;
146  fTsumwy2 += weight * ycoord * ycoord;
147  fTsumwxy += weight * xcoord * ycoord;
148  fTsumwz += weight * value;
149  fTsumwz2 += weight * value * value;
150 
151  // ------------ Update local (per bin) statistics
152  TProfile2PolyBin *bin;
153  TIter next(&fCells[n + fCellX * m]);
154  TObject *obj;
155  while ((obj = next())) {
156  bin = (TProfile2PolyBin *)obj;
157  if (bin->IsInside(xcoord, ycoord)) {
158  fEntries++;
159  bin->Fill(value, weight);
160  bin->Update();
161  bin->SetContent(bin->fAverage);
162  }
163  }
164 
165  return tmp;
166 }
167 
169 {
170  Int_t size = in->GetSize();
171 
172  std::vector<TProfile2Poly *> list;
173  list.reserve(size);
174 
175  for (int i = 0; i < size; i++) {
176  list.push_back((TProfile2Poly *)((TList *)in)->At(i));
177  }
178  return this->Merge(list);
179 }
180 
181 Long64_t TProfile2Poly::Merge(const std::vector<TProfile2Poly *> &list)
182 {
183  if (list.size() == 0) {
184  std::cout << "[FAIL] TProfile2Poly::Merge: No objects to be merged " << std::endl;
185  return -1;
186  }
187 
188  // ------------ Check that bin numbers of TP2P's to be merged are equal
189  std::set<Int_t> numBinUnique;
190  for (const auto &histo : list) {
191  numBinUnique.insert(histo->fBins->GetSize());
192  }
193  if (numBinUnique.size() != 1) {
194  std::cout << "[FAIL] TProfile2Poly::Merge: Bin numbers of TProfile2Polys to be merged differ!" << std::endl;
195  return -1;
196  }
197  Int_t nbins = *numBinUnique.begin();
198 
199  // ------------ Update global (per histo) statistics
200  for (const auto &histo : list) {
201  this->fEntries += histo->fEntries;
202  this->fTsumw += histo->fTsumw;
203  this->fTsumw2 += histo->fTsumw2;
204  this->fTsumwx += histo->fTsumwx;
205  this->fTsumwx2 += histo->fTsumwx2;
206  this->fTsumwy += histo->fTsumwy;
207  this->fTsumwy2 += histo->fTsumwy2;
208  this->fTsumwxy += histo->fTsumwxy;
209  this->fTsumwz += histo->fTsumwz;
210  this->fTsumwz2 += histo->fTsumwz2;
211 
212  // Merge overflow bins
213  for (Int_t i = 0; i < kNOverflow; ++i) {
214  this->fOverflowBins[i].Merge(&histo->fOverflowBins[i]);
215  }
216  }
217 
218  // ------------ Update local (per bin) statistics
219  TProfile2PolyBin *dst = nullptr;
220  TProfile2PolyBin *src = nullptr;
221  for (Int_t i = 0; i < nbins; i++) {
222  dst = (TProfile2PolyBin *)fBins->At(i);
223 
224  for (const auto &e : list) {
225  src = (TProfile2PolyBin *)e->fBins->At(i);
226  dst->Merge(src);
227  }
228 
229  dst->Update();
230  }
231 
232  this->SetContentToAverage();
233  return 1;
234 }
235 
237 {
238  Int_t nbins = fBins->GetSize();
239  for (Int_t i = 0; i < nbins; i++) {
241  bin->Update();
242  bin->SetContent(bin->fAverage);
243  }
244  for (Int_t i = 0; i < kNOverflow; ++i) {
245  TProfile2PolyBin & bin = fOverflowBins[i];
246  bin.Update();
247  bin.SetContent(bin.fAverage);
248  }
249 }
250 
252 {
253  Int_t nbins = fBins->GetSize();
254  for (Int_t i = 0; i < nbins; i++) {
256  bin->Update();
257  bin->SetContent(bin->fError);
258  }
259  for (Int_t i = 0; i < kNOverflow; ++i) {
260  TProfile2PolyBin & bin = fOverflowBins[i];
261  bin.Update();
262  bin.SetContent(bin.fError);
263  }
264 }
265 
267 {
268  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
269  if (bin<0) return fOverflowBins[-bin - 1].GetContent();
270  return ((TProfile2PolyBin*) fBins->At(bin-1))->GetContent();
271 }
272 
273 
275 {
276  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
277  if (bin < 0) return fOverflowBins[-bin - 1].GetEffectiveEntries();
278  return ((TProfile2PolyBin *)fBins->At(bin - 1))->GetEffectiveEntries();
279 }
280 
282 {
283  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
284  if (bin < 0) return fOverflowBins[-bin - 1].GetEntries();
285  return ((TProfile2PolyBin *)fBins->At(bin - 1))->GetEntries();
286 }
287 
289 {
290  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
291  if (bin < 0) return fOverflowBins[-bin - 1].GetEntriesW2();
292  return ((TProfile2PolyBin *)fBins->At(bin - 1))->GetEntriesW2();
293 }
294 
296 {
297  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
298  if (bin < 0) return fOverflowBins[-bin - 1].GetEntriesVW();
299  return ((TProfile2PolyBin *)fBins->At(bin - 1))->GetEntriesVW();
300 }
301 
303 {
304  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
305  if (bin < 0) return fOverflowBins[-bin - 1].GetEntriesWV2();
306  return ((TProfile2PolyBin *)fBins->At(bin - 1))->GetEntriesWV2();
307 }
308 
310 {
311  Double_t tmp = 0;
312  if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
313  if (bin < 0)
314  tmp = fOverflowBins[-bin - 1].GetError();
315  else
316  tmp = ((TProfile2PolyBin *)fBins->At(bin - 1))->GetError();
317 
318  return (fErrorMode == kERRORSPREAD) ? tmp : tmp / std::sqrt(GetBinEffectiveEntries(bin));
319 
320 }
321 
322 ////////////////////////////////////////////////////////////////////////////////
323 /// Fill the array stats from the contents of this profile.
324 /// The array stats must be correctly dimensioned in the calling program.
325 ///
326 /// - stats[0] = sumw
327 /// - stats[1] = sumw2
328 /// - stats[2] = sumwx
329 /// - stats[3] = sumwx2
330 /// - stats[4] = sumwy
331 /// - stats[5] = sumwy2
332 /// - stats[6] = sumwxy
333 /// - stats[7] = sumwz
334 /// - stats[8] = sumwz2
335 ///
336 /// If no axis-subrange is specified (via TAxis::SetRange), the array stats
337 /// is simply a copy of the statistics quantities computed at filling time.
338 /// If a sub-range is specified, the function recomputes these quantities
339 /// from the bin contents in the current axis range.
340 
342 {
343  stats[0] = fTsumw;
344  stats[1] = fTsumw2;
345  stats[2] = fTsumwx;
346  stats[3] = fTsumwx2;
347  stats[4] = fTsumwy;
348  stats[5] = fTsumwy2;
349  stats[6] = fTsumwxy;
350  stats[7] = fTsumwz;
351  stats[8] = fTsumwz2;
352 }
353 
355 {
356  Double_t total = 0;
357  Double_t cont = 0;
358  for (Int_t i = 0; i < kNOverflow; ++i) {
359  cont = GetOverflowContent(i);
360  total += cont;
361  std::cout << "\t" << cont << "\t";
362  if ((i + 1) % 3 == 0) std::cout << std::endl;
363  }
364 
365  std::cout << "Total: " << total << std::endl;
366 }
367 
369 {
370  TIter next(fBins);
371  TObject *obj;
372  TProfile2PolyBin *bin;
373 
374  // Clears bin contents
375  while ((obj = next())) {
376  bin = (TProfile2PolyBin *)obj;
377  bin->ClearContent();
378  bin->ClearStats();
379  }
380  TH2::Reset(opt);
381 }
382 
384 {
385  // The overflow regions are calculated by considering x, y coordinates.
386  // The Middle bin at -5 contains all the TProfile2Poly bins.
387  //
388  // -0 -1 -2
389  // ________
390  // -1: |__|__|__|
391  // -4: |__|__|__|
392  // -7: |__|__|__|
393 
394  Int_t region = 0;
395 
396  if (fNcells <= kNOverflow) return 0;
397 
398  // --- y offset
399  if (y > fYaxis.GetXmax())
400  region += -1;
401  else if (y > fYaxis.GetXmin())
402  region += -4;
403  else
404  region += -7;
405 
406  // --- x offset
407  if (x > fXaxis.GetXmax())
408  region += -2;
409  else if (x > fXaxis.GetXmin())
410  region += -1;
411  else
412  region += 0;
413 
414  return region;
415 }
416 
418 {
419  fErrorMode = type;
420 }
void Fill(Double_t value, Double_t weight)
virtual TProfile2PolyBin * CreateBin(TObject *poly) override
Create appropriate histogram bin.
void SetContentToAverage()
void printOverflowRegions()
friend class TProfile2PolyBin
Definition: TProfile2Poly.h:50
EErrorType fErrorMode
Definition: TProfile2Poly.h:95
void SetContentToError()
virtual Double_t GetEffectiveEntries() const
Number of effective entries of the histogram.
Definition: TH1.cxx:4079
long long Long64_t
Definition: RtypesCore.h:69
Double_t GetEntriesWV2() const
Definition: TProfile2Poly.h:27
EErrorType fErrorMode
Definition: TProfile2Poly.h:38
Double_t fTsumwz2
Definition: TProfile2Poly.h:97
const char Option_t
Definition: RtypesCore.h:62
Double_t fStepX
Definition: TH2Poly.h:140
virtual Double_t GetBinError(Int_t bin) const override
Returns the value of error associated to bin number bin.
TAxis fYaxis
Y axis descriptor.
Definition: TH1.h:81
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
void SetChanged(Bool_t flag)
Definition: TH2Poly.h:44
Double_t GetBinEntriesVW(Int_t bin) const
void Merge(const TProfile2PolyBin *toMerge)
int Int_t
Definition: RtypesCore.h:41
Double_t GetEffectiveEntries() const
Definition: TProfile2Poly.h:23
Bool_t IsInside(Double_t x, Double_t y) const
Return "true" if the point (x,y) is inside the bin.
Definition: TH2Poly.cxx:1509
int nbins[3]
Double_t GetError() const
Definition: TProfile2Poly.h:28
Helper class to represent a bin in the TH2Poly histogram.
Definition: TH2Poly.h:25
EErrorType
Definition: TProfile.h:28
Double_t fTsumwx2
Total Sum of weight*X*X.
Definition: TH1.h:89
Int_t fCellX
Definition: TH2Poly.h:136
double sqrt(double)
Double_t GetXmin() const
Definition: TAxis.h:133
TList * fCells
Definition: TH2Poly.h:139
Double_t x[n]
Definition: legend1.C:17
Long64_t Merge(const std::vector< TProfile2Poly *> &list)
TProfile2PolyBin fOverflowBins[kNOverflow]
Definition: TProfile2Poly.h:94
Double_t fTsumwy2
Definition: TH2.h:35
Double_t GetOverflowContent(Int_t idx)
Definition: TProfile2Poly.h:90
Double_t GetEntries() const
Definition: TProfile2Poly.h:24
void SetContent(Double_t content)
Definition: TH2Poly.h:45
Double_t fTsumwx
Total Sum of weight*X.
Definition: TH1.h:88
A doubly linked list.
Definition: TList.h:43
Double_t GetEntriesW2() const
Definition: TProfile2Poly.h:25
virtual Int_t Fill(Double_t xcoord, Double_t ycoord, Double_t value) override
Increment the bin containing (x,y) by w.
Int_t fCellY
Definition: TH2Poly.h:137
Int_t GetOverflowRegionFromCoordinates(Double_t x, Double_t y)
Double_t GetContent() const
Definition: TH2Poly.h:35
Double_t fTsumwxy
Definition: TH2.h:36
Double_t fTsumwy
Definition: TH2.h:34
Collection abstract base class.
Definition: TCollection.h:42
Double_t fEntries
Number of entries.
Definition: TH1.h:85
TMarker * m
Definition: textangle.C:8
Int_t GetNumberOfBins() const
Definition: TH2Poly.h:110
virtual void GetStats(Double_t *stats) const override
Fill the array stats from the contents of this profile.
double floor(double)
Double_t GetBinEntriesW2(Int_t bin) const
Double_t GetEntriesVW() const
Definition: TProfile2Poly.h:26
virtual TObject * At(Int_t idx) const
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:315
void SetErrorOption(EErrorType type)
TList * fBins
Definition: TH2Poly.h:134
static unsigned int total
Double_t fTsumw2
Total Sum of squares of weights.
Definition: TH1.h:87
virtual Double_t GetBinContent(Int_t bin) const override
Returns the content of the input bin For the overflow/underflow/sea bins: -1 | -2 | -3 ---+----+---- ...
Int_t OverflowIdxToArrayIdx(Int_t val)
#define ClassImp(name)
Definition: Rtypes.h:336
double Double_t
Definition: RtypesCore.h:55
Double_t fTsumw
Total Sum of weights.
Definition: TH1.h:86
int type
Definition: TGX11.cxx:120
Double_t y[n]
Definition: legend1.C:17
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
virtual Double_t GetEntries() const
Return the current number of entries.
Definition: TH1.cxx:4054
void ClearContent()
Definition: TH2Poly.h:32
Mother of all ROOT objects.
Definition: TObject.h:37
virtual void Reset(Option_t *option="") override
Reset this histogram: contents, errors, etc.
Double_t fTsumwz
Definition: TProfile2Poly.h:96
Double_t fStepY
Definition: TH2Poly.h:140
Double_t GetBinEntries(Int_t bin) const
TAxis fXaxis
X axis descriptor.
Definition: TH1.h:80
Double_t GetBinEntriesWV2(Int_t bin) const
Double_t GetBinEffectiveEntries(Int_t bin) const
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
Definition: TH2.cxx:2455
virtual Int_t GetSize() const
Definition: TCollection.h:89
Double_t GetXmax() const
Definition: TAxis.h:134
const Int_t n
Definition: legend1.C:16
2D Histogram with Polygonal Bins
Definition: TH2Poly.h:66
Int_t fNcells
number of bins(1D), cells (2D) +U/Overflows
Definition: TH1.h:79