Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TH2Poly.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id$
2// TH2Poly v2.1
3// Author: Olivier Couet, Deniz Gunceler, Danilo Piparo
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13#include "TH2Poly.h"
14#include "TMultiGraph.h"
15#include "TGraph.h"
16#include "Riostream.h"
17#include "TList.h"
18#include "TMath.h"
19#include <cassert>
20
22
23/** \class TH2Poly
24 \ingroup Hist
252D Histogram with Polygonal Bins
26
27## Overview
28`TH2Poly` is a 2D Histogram class (TH2) allowing to define polygonal
29bins of arbitrary shape.
30
31Each bin in the `TH2Poly` histogram is a `TH2PolyBin` object.
32`TH2PolyBin` is a very simple class containing the vertices (stored
33as `TGraph`s or `TMultiGraph`s ) and contents of the polygonal
34bin as well as several related functions.
35
36Essentially, a `TH2Poly` is a TList of `TH2PolyBin` objects
37with methods to manipulate them.
38
39Bins are defined using one of the `AddBin()` methods. The bin definition
40should be done before filling.
41
42The histogram can be filled with `Fill(Double_t x, Double_t y, Double_t w)
43`. `w` is the weight.
44If no weight is specified, it is assumed to be 1.
45
46Not all histogram's area need to be binned. Filling an area without bins,
47will falls into the overflows. Adding a bin is not retroactive; it doesn't
48affect previous fillings. A `Fill()` call, that
49was previously ignored due to the lack of a bin at the specified location, is
50not reconsidered when that location is binned later.
51
52If there are two overlapping bins, the first one in the list will be incremented
53by `Fill()`.
54
55The histogram may automatically extends its limits if a bin outside the
56histogram limits is added. This is done when the default constructor (with no
57arguments) is used. It generates a histogram with no limits along the X and Y
58axis. Adding bins to it will extend it up to a proper size.
59
60`TH2Poly` implements a partitioning algorithm to speed up bins' filling
61(see the "Partitioning Algorithm" section for details).
62The partitioning algorithm divides the histogram into regions called cells.
63The bins that each cell intersects are recorded in an array of `TList`s.
64When a coordinate in the histogram is to be filled; the method (quickly) finds
65which cell the coordinate belongs. It then only loops over the bins
66intersecting that cell to find the bin the input coordinate corresponds to.
67The partitioning of the histogram is updated continuously as each bin is added.
68The default number of cells on each axis is 25. This number could be set to
69another value in the constructor or adjusted later by calling the
70`ChangePartition(Int_t, Int_t)` method. The partitioning algorithm is
71considerably faster than the brute force algorithm (i.e. checking if each bin
72contains the input coordinates), especially if the histogram is to be filled
73many times.
74
75The following very simple macro shows how to build and fill a `TH2Poly`:
76~~~ {.cpp}
77{
78 TH2Poly *h2p = new TH2Poly();
79
80 Double_t x1[] = {0, 5, 6};
81 Double_t y1[] = {0, 0, 5};
82 Double_t x2[] = {0, -1, -1, 0};
83 Double_t y2[] = {0, 0, -1, 3};
84 Double_t x3[] = {4, 3, 0, 1, 2.4};
85 Double_t y3[] = {4, 3.7, 1, 3.7, 2.5};
86
87 h2p->AddBin(3, x1, y1);
88 h2p->AddBin(4, x2, y2);
89 h2p->AddBin(5, x3, y3);
90
91 h2p->Fill(0.1, 0.01, 3);
92 h2p->Fill(-0.5, -0.5, 7);
93 h2p->Fill(-0.7, -0.5, 1);
94 h2p->Fill(1, 3, 1.5);
96~~~
97
98More examples can be found in th2polyBoxes.C, th2polyEurope.C, th2polyHoneycomb.C
99and th2polyUSA.C.
100
101## Partitioning Algorithm
102The partitioning algorithm forms an essential part of the `TH2Poly`
103class. It is implemented to speed up the filling of bins.
104
105With the brute force approach, the filling is done in the following way: An
106iterator loops over all bins in the `TH2Poly` and invokes the
107method `IsInside()` for each of them.
108This method checks if the input location is in that bin. If the filling
109coordinate is inside, the bin is filled. Looping over all the bin is
110very slow.
111
112The alternative is to divide the histogram into virtual rectangular regions
113called "cells". Each cell stores the pointers of the bins intersecting it.
114When a coordinate is to be filled, the method finds which cell the coordinate
115falls into. Since the cells are rectangular, this can be done very quickly.
116It then only loops over the bins associated with that cell and calls `IsInside()`
117only on that bins. This reduces considerably the number of bins on which `IsInside()`
118is called and therefore speed up by a huge factor the filling compare to the brute force
119approach where `IsInside()` is called for all bins.
120
121The addition of bins to the appropriate cells is done when the bin is added
122to the histogram. To do this, `AddBin()` calls the
123`AddBinToPartition()` method.
124This method adds the input bin to the partitioning matrix.
125
126The number of partition cells per axis can be specified in the constructor.
127If it is not specified, the default value of 25 along each axis will be
128assigned. This value was chosen because it is small enough to avoid slowing
129down AddBin(), while being large enough to enhance Fill() by a considerable
130amount. Regardless of how it is initialized at construction time, it can be
131changed later with the `ChangePartition()` method.
132`ChangePartition()` deletes the
133old partition matrix and generates a new one with the specified number of cells
134on each axis.
135
136The optimum number of partition cells per axis changes with the number of
137times `Fill()` will be called. Although partitioning greatly speeds up
138filling, it also adds a constant time delay into the code. When `Fill()`
139is to be called many times, it is more efficient to divide the histogram into
140a large number cells. However, if the histogram is to be filled only a few
141times, it is better to divide into a small number of cells.
142*/
143
144////////////////////////////////////////////////////////////////////////////////
145/// Default Constructor. No boundaries specified.
146
148{
149 Initialize(0., 0., 0., 0., 25, 25);
150 SetName("NoName");
151 SetTitle("NoTitle");
152 SetFloat();
153}
154
155////////////////////////////////////////////////////////////////////////////////
156/// Constructor with specified name and boundaries,
157/// but no partition cell number.
158
159TH2Poly::TH2Poly(const char *name,const char *title, Double_t xlow,Double_t xup
160 , Double_t ylow,Double_t yup)
161{
162 Initialize(xlow, xup, ylow, yup, 25, 25);
163 SetName(name);
164 SetTitle(title);
166}
167
168////////////////////////////////////////////////////////////////////////////////
169/// Constructor with specified name and boundaries and partition cell number.
170
171TH2Poly::TH2Poly(const char *name,const char *title,
172 Int_t nX, Double_t xlow, Double_t xup,
173 Int_t nY, Double_t ylow, Double_t yup)
174{
175 Initialize(xlow, xup, ylow, yup, nX, nY);
176 SetName(name);
177 SetTitle(title);
179}
180
181////////////////////////////////////////////////////////////////////////////////
182/// Destructor.
183
185{
186 delete[] fCells;
187 delete[] fIsEmpty;
188 delete[] fCompletelyInside;
189 // delete at the end the bin List since it owns the objects
190 delete fBins;
191}
192
193////////////////////////////////////////////////////////////////////////////////
194/// Create appropriate histogram bin.
195/// e.g. TH2Poly creates TH2PolyBin,
196/// TProfile2Poly creates TProfile2PolyBin
197/// This is done so that TH2Poly::AddBin does not have to be duplicated,
198/// but only create needs to be reimplemented for additional histogram types
199
201{
202 if (!poly) return 0;
203
204 if (fBins == 0) {
205 fBins = new TList();
206 fBins->SetOwner();
207 }
208
209 fNcells++;
210 Int_t ibin = fNcells - kNOverflow;
211 // if structure fsumw2 is created extend it
212 if (fSumw2.fN) fSumw2.Set(fNcells);
213 return new TH2PolyBin(poly, ibin);
214}
215
216////////////////////////////////////////////////////////////////////////////////
217/// Adds a new bin to the histogram. It can be any object having the method
218/// IsInside(). It returns the bin number in the histogram. It returns 0 if
219/// it failed to add. To allow the histogram limits to expand when a bin
220/// outside the limits is added, call SetFloat() before adding the bin.
221
223{
224 auto *bin = CreateBin(poly);
225 Int_t ibin = fNcells-kNOverflow;
226 if(!bin) return 0;
227
228 // If the bin lies outside histogram boundaries, then extends the boundaries.
229 // Also changes the partition information accordingly
230 Bool_t flag = kFALSE;
231 if (fFloat) {
232 if (fXaxis.GetXmin() > bin->GetXMin()) {
233 fXaxis.Set(100, bin->GetXMin(), fXaxis.GetXmax());
234 flag = kTRUE;
235 }
236 if (fXaxis.GetXmax() < bin->GetXMax()) {
237 fXaxis.Set(100, fXaxis.GetXmin(), bin->GetXMax());
238 flag = kTRUE;
239 }
240 if (fYaxis.GetXmin() > bin->GetYMin()) {
241 fYaxis.Set(100, bin->GetYMin(), fYaxis.GetXmax());
242 flag = kTRUE;
243 }
244 if (fYaxis.GetXmax() < bin->GetYMax()) {
245 fYaxis.Set(100, fYaxis.GetXmin(), bin->GetYMax());
246 flag = kTRUE;
247 }
248 if (flag) ChangePartition(fCellX, fCellY);
249 } else {
250 /*Implement polygon clipping code here*/
251 }
252
253 fBins->Add((TObject*) bin);
255
256 // Adds the bin to the partition matrix
258
259 return ibin;
260}
261
262////////////////////////////////////////////////////////////////////////////////
263/// Adds a new bin to the histogram. The number of vertices and their (x,y)
264/// coordinates are required as input. It returns the bin number in the
265/// histogram.
266
268{
269 TGraph *g = new TGraph(n, x, y);
270 Int_t bin = AddBin(g);
271 return bin;
272}
273
274////////////////////////////////////////////////////////////////////////////////
275/// Add a new bin to the histogram. The bin shape is a rectangle.
276/// It returns the bin number of the bin in the histogram.
277
279{
280 Double_t x[] = {x1, x1, x2, x2, x1};
281 Double_t y[] = {y1, y2, y2, y1, y1};
282 TGraph *g = new TGraph(5, x, y);
283 Int_t bin = AddBin(g);
284 return bin;
285}
286
287////////////////////////////////////////////////////////////////////////////////
288/// Performs the operation: this = this + c1*h1.
289
291{
292 Int_t bin;
293
294 TH2Poly *h1p = (TH2Poly *)h1;
295
296 // Check if number of bins is the same.
297 if (h1p->GetNumberOfBins() != GetNumberOfBins()) {
298 Error("Add", "Attempt to add histograms with different number of bins");
299 return kFALSE;
300 }
301
302 // Check if the bins are the same.
303 TList *h1pBins = h1p->GetBins();
304 TH2PolyBin *thisBin, *h1pBin;
305 for (bin = 1; bin <= GetNumberOfBins(); bin++) {
306 thisBin = (TH2PolyBin *)fBins->At(bin - 1);
307 h1pBin = (TH2PolyBin *)h1pBins->At(bin - 1);
308 if (thisBin->GetXMin() != h1pBin->GetXMin() ||
309 thisBin->GetXMax() != h1pBin->GetXMax() ||
310 thisBin->GetYMin() != h1pBin->GetYMin() ||
311 thisBin->GetYMax() != h1pBin->GetYMax()) {
312 Error("Add", "Attempt to add histograms with different bin limits");
313 return kFALSE;
314 }
315 }
316
317
318 // Create Sumw2 if h1p has Sumw2 set
319 if (fSumw2.fN == 0 && h1p->GetSumw2N() != 0) Sumw2();
320
321 // statistics can be preserved only in case of positive coefficients
322 // otherwise with negative c1 (histogram subtraction) one risks to get negative variances
323 Bool_t resetStats = (c1 < 0);
324 Double_t s1[kNstat] = {0};
325 Double_t s2[kNstat] = {0};
326 if (!resetStats) {
327 // need to initialize to zero s1 and s2 since
328 // GetStats fills only used elements depending on dimension and type
329 GetStats(s1);
330 h1->GetStats(s2);
331 }
332 // get number of entries now because afterwards UpdateBinContent will change it
333 Double_t entries = TMath::Abs( GetEntries() + c1 * h1->GetEntries() );
334
335
336 // Perform the Add.
337 Double_t factor = 1;
338 if (h1p->GetNormFactor() != 0)
339 factor = h1p->GetNormFactor() / h1p->GetSumOfWeights();
340 for (bin = 0; bin < fNcells; bin++) {
341 Double_t y = this->RetrieveBinContent(bin) + c1 * h1p->RetrieveBinContent(bin);
342 UpdateBinContent(bin, y);
343 if (fSumw2.fN) {
344 Double_t esq = factor * factor * h1p->GetBinErrorSqUnchecked(bin);
345 fSumw2.fArray[bin] += c1 * c1 * factor * factor * esq;
346 }
347 }
348
349 // update statistics (do here to avoid changes by SetBinContent)
350 if (resetStats) {
351 // statistics need to be reset in case coefficient are negative
352 ResetStats();
353 } else {
354 for (Int_t i = 0; i < kNstat; i++) {
355 if (i == 1) s1[i] += c1 * c1 * s2[i];
356 else s1[i] += c1 * s2[i];
357 }
358 PutStats(s1);
359 SetEntries(entries);
360 }
361 return kTRUE;
362}
363
364////////////////////////////////////////////////////////////////////////////////
365/// Adds the input bin into the partition cell matrix. This method is called
366/// in AddBin() and ChangePartition().
367
369{
370 // Cell Info
371 Int_t nl, nr, mb, mt; // Max/min indices of the cells that contain the bin
372 Double_t xclipl, xclipr, yclipb, yclipt; // x and y coordinates of a cell
373 Double_t binXmax, binXmin, binYmax, binYmin; // The max/min bin coordinates
374
375 binXmax = bin->GetXMax();
376 binXmin = bin->GetXMin();
377 binYmax = bin->GetYMax();
378 binYmin = bin->GetYMin();
379 nl = (Int_t)(floor((binXmin - fXaxis.GetXmin())/fStepX));
380 nr = (Int_t)(floor((binXmax - fXaxis.GetXmin())/fStepX));
381 mb = (Int_t)(floor((binYmin - fYaxis.GetXmin())/fStepY));
382 mt = (Int_t)(floor((binYmax - fYaxis.GetXmin())/fStepY));
383
384 // Make sure the array indices are correct.
385 if (nr>=fCellX) nr = fCellX-1;
386 if (mt>=fCellY) mt = fCellY-1;
387 if (nl<0) nl = 0;
388 if (mb<0) mb = 0;
389
390 // number of cells in the grid
391 //N.B. not to be confused with fNcells (the number of bins) !
393
394 // Loop over all cells
395 for (int i = nl; i <= nr; i++) {
396 xclipl = fXaxis.GetXmin() + i*fStepX;
397 xclipr = xclipl + fStepX;
398 for (int j = mb; j <= mt; j++) {
399 yclipb = fYaxis.GetXmin() + j*fStepY;
400 yclipt = yclipb + fStepY;
401
402 // If the bin is completely inside the cell,
403 // add that bin to the cell then return
404 if ((binXmin >= xclipl) && (binXmax <= xclipr) &&
405 (binYmax <= yclipt) && (binYmin >= yclipb)){
406 fCells[i + j*fCellX].Add((TObject*) bin);
407 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
408 return;
409 }
410
411 // If any of the sides of the cell intersect with any side of the bin,
412 // add that bin then continue
413 if (IsIntersecting(bin, xclipl, xclipr, yclipb, yclipt)) {
414 fCells[i + j*fCellX].Add((TObject*) bin);
415 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
416 continue;
417 }
418 // If a corner of the cell is inside the bin and since there is no
419 // intersection, then that cell completely inside the bin.
420 if((bin->IsInside(xclipl,yclipb)) || (bin->IsInside(xclipl,yclipt))){
421 fCells[i + j*fCellX].Add((TObject*) bin);
422 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
424 continue;
425 }
426 if((bin->IsInside(xclipr,yclipb)) || (bin->IsInside(xclipr,yclipt))){
427 fCells[i + j*fCellX].Add((TObject*) bin);
428 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
430 continue;
431 }
432 }
433 }
434}
435
436////////////////////////////////////////////////////////////////////////////////
437/// Changes the number of partition cells in the histogram.
438/// Deletes the old partition and constructs a new one.
439
441{
442 fCellX = n; // Set the number of cells
443 fCellY = m; // Set the number of cells
444
445 delete [] fCells; // Deletes the old partition
446
447 // number of cells in the grid
448 //N.B. not to be confused with fNcells (the number of bins) !
450 fCells = new TList [fNCells]; // Sets an empty partition
451
454
455 delete [] fIsEmpty;
456 delete [] fCompletelyInside;
457 fIsEmpty = new Bool_t [fNCells];
459
460 // Initializes the flags
461 for (int i = 0; i<fNCells; i++) {
462 fIsEmpty[i] = kTRUE;
464 }
465
466 // TList iterator
467 TIter next(fBins);
468 TObject *obj;
469
470 while((obj = next())){ // Loop over bins and add them to the partition
472 }
473}
474
475////////////////////////////////////////////////////////////////////////////////
476/// Make a complete copy of the underlying object. If 'newname' is set,
477/// the copy's name will be set to that name.
478
479TObject* TH2Poly::Clone(const char* newname) const
480{
481 // TH1::Clone relies on ::Copy to implemented by the derived class.
482 // Until this is implemented, revert to the much slower default version
483 // (and possibly non-thread safe).
484
485 return TNamed::Clone(newname);
486}
487
488////////////////////////////////////////////////////////////////////////////////
489/// Clears the contents of all bins in the histogram.
490
492{
493 TIter next(fBins);
494 TObject *obj;
495 TH2PolyBin *bin;
496
497 // Clears the bin contents
498 while ((obj = next())) {
499 bin = (TH2PolyBin*) obj;
500 bin->ClearContent();
501 }
502
503 // Clears the statistics
504 fTsumw = 0;
505 fTsumw2 = 0;
506 fTsumwx = 0;
507 fTsumwx2 = 0;
508 fTsumwy = 0;
509 fTsumwy2 = 0;
510 fEntries = 0;
511}
512
513////////////////////////////////////////////////////////////////////////////////
514/// Reset this histogram: contents, errors, etc.
515
517{
518 TIter next(fBins);
519 TObject *obj;
520 TH2PolyBin *bin;
521
522 // Clears the bin contents
523 while ((obj = next())) {
524 bin = (TH2PolyBin*) obj;
525 bin->ClearContent();
526 }
527
528 TH2::Reset(opt);
529}
530
531////////////////////////////////////////////////////////////////////////////////
532/// Returns the bin number of the bin at the given coordinate. -1 to -9 are
533/// the overflow and underflow bins. overflow bin -5 is the unbinned areas in
534/// the histogram (also called "the sea"). The third parameter can be left
535/// blank.
536/// The overflow/underflow bins are:
537///~~~ {.cpp}
538/// -1 | -2 | -3
539/// -------------
540/// -4 | -5 | -6
541/// -------------
542/// -7 | -8 | -9
543///~~~
544/// where -5 means is the "sea" bin (i.e. unbinned areas)
545
547{
548
549 // Checks for overflow/underflow
550 Int_t overflow = 0;
551 if (y > fYaxis.GetXmax()) overflow += -1;
552 else if (y > fYaxis.GetXmin()) overflow += -4;
553 else overflow += -7;
554 if (x > fXaxis.GetXmax()) overflow += -2;
555 else if (x > fXaxis.GetXmin()) overflow += -1;
556 if (overflow != -5) return overflow;
557
558 // Finds the cell (x,y) coordinates belong to
561
562 // Make sure the array indices are correct.
563 if (n>=fCellX) n = fCellX-1;
564 if (m>=fCellY) m = fCellY-1;
565 if (n<0) n = 0;
566 if (m<0) m = 0;
567
568 if (fIsEmpty[n+fCellX*m]) return -5;
569
570 TH2PolyBin *bin;
571
572 TIter next(&fCells[n+fCellX*m]);
573 TObject *obj;
574
575 // Search for the bin in the cell
576 while ((obj=next())) {
577 bin = (TH2PolyBin*)obj;
578 if (bin->IsInside(x,y)) return bin->GetBinNumber();
579 }
580
581 // If the search has not returned a bin, the point must be on "the sea"
582 return -5;
583}
584
585////////////////////////////////////////////////////////////////////////////////
586/// Increment the bin containing (x,y) by 1.
587/// Uses the partitioning algorithm.
588
590{
591 return Fill(x, y, 1.0);
592}
593
594////////////////////////////////////////////////////////////////////////////////
595/// Increment the bin containing (x,y) by w.
596/// Uses the partitioning algorithm.
597
599{
600 // see GetBinCOntent for definition of overflow bins
601 // in case of weighted events store weight square in fSumw2.fArray
602 // but with this indexing:
603 // fSumw2.fArray[0:kNOverflow-1] : sum of weight squares for the overflow bins
604 // fSumw2.fArray[kNOverflow:fNcells] : sum of weight squares for the standard bins
605 // where fNcells = kNOverflow + Number of bins. kNOverflow=9
606
607 if (fNcells <= kNOverflow) return 0;
608
609 // create sum of weight square array if weights are different than 1
610 if (!fSumw2.fN && w != 1.0 && !TestBit(TH1::kIsNotW) ) Sumw2();
611
612 Int_t overflow = 0;
613 if (y > fYaxis.GetXmax()) overflow += -1;
614 else if (y > fYaxis.GetXmin()) overflow += -4;
615 else overflow += -7;
616 if (x > fXaxis.GetXmax()) overflow += -2;
617 else if(x > fXaxis.GetXmin()) overflow += -1;
618 if (overflow != -5) {
619 fOverflow[-overflow - 1]+= w;
620 if (fSumw2.fN) fSumw2.fArray[-overflow - 1] += w*w;
621 return overflow;
622 }
623
624 // Finds the cell (x,y) coordinates belong to
627
628 // Make sure the array indices are correct.
629 if (n>=fCellX) n = fCellX-1;
630 if (m>=fCellY) m = fCellY-1;
631 if (n<0) n = 0;
632 if (m<0) m = 0;
633
634 if (fIsEmpty[n+fCellX*m]) {
635 fOverflow[4]+= w;
636 if (fSumw2.fN) fSumw2.fArray[4] += w*w;
637 return -5;
638 }
639
640 TH2PolyBin *bin;
641 Int_t bi;
642
643 TIter next(&fCells[n+fCellX*m]);
644 TObject *obj;
645
646 while ((obj=next())) {
647 bin = (TH2PolyBin*)obj;
648 // needs to account offset in array for overflow bins
649 bi = bin->GetBinNumber()-1+kNOverflow;
650 if (bin->IsInside(x,y)) {
651 bin->Fill(w);
652
653 // Statistics
654 fTsumw = fTsumw + w;
655 fTsumw2 = fTsumw2 + w*w;
656 fTsumwx = fTsumwx + w*x;
657 fTsumwx2 = fTsumwx2 + w*x*x;
658 fTsumwy = fTsumwy + w*y;
659 fTsumwy2 = fTsumwy2 + w*y*y;
660 if (fSumw2.fN) {
661 assert(bi < fSumw2.fN);
662 fSumw2.fArray[bi] += w*w;
663 }
664 fEntries++;
665
667
668 return bin->GetBinNumber();
669 }
670 }
671
672 fOverflow[4]+= w;
673 if (fSumw2.fN) fSumw2.fArray[4] += w*w;
674 return -5;
675}
676
677////////////////////////////////////////////////////////////////////////////////
678/// Increment the bin named "name" by w.
679
681{
682 TString sname(name);
683
684 TIter next(fBins);
685 TObject *obj;
686 TH2PolyBin *bin;
687
688 while ((obj = next())) {
689 bin = (TH2PolyBin*) obj;
690 if (!sname.CompareTo(bin->GetPolygon()->GetName())) {
691 bin->Fill(w);
692 fEntries++;
694 return bin->GetBinNumber();
695 }
696 }
697
698 return 0;
699}
700
701////////////////////////////////////////////////////////////////////////////////
702/// Fills a 2-D histogram with an array of values and weights.
703///
704/// \param [in] ntimes: number of entries in arrays x and w
705/// (array size must be ntimes*stride)
706/// \param [in] x: array of x values to be histogrammed
707/// \param [in] y: array of y values to be histogrammed
708/// \param [in] w: array of weights
709/// \param [in] stride: step size through arrays x, y and w
710
711void TH2Poly::FillN(Int_t ntimes, const Double_t* x, const Double_t* y,
712 const Double_t* w, Int_t stride)
713{
714 for (int i = 0; i < ntimes; i += stride) {
715 Fill(x[i], y[i], w[i]);
716 }
717}
718
719////////////////////////////////////////////////////////////////////////////////
720/// Returns the integral of bin contents.
721/// By default the integral is computed as the sum of bin contents.
722/// If option "width" or "area" is specified, the integral is the sum of
723/// the bin contents multiplied by the area of the bin.
724
726{
727 TString opt = option;
728 opt.ToLower();
729
730 Double_t w;
731 Double_t integral = 0.;
732
733 TIter next(fBins);
734 TObject *obj;
735 TH2PolyBin *bin;
736 if ((opt.Contains("width")) || (opt.Contains("area"))) {
737 while ((obj = next())) {
738 bin = (TH2PolyBin *)obj;
739 w = bin->GetArea();
740 integral += w * (bin->GetContent());
741 }
742 } else {
743 // need to recompute integral in case somebith called SetBinContent and
744 // cannot use fTsumw since it is not updated in that case
745 while ((obj = next())) {
746 bin = (TH2PolyBin *)obj;
747 integral += (bin->GetContent());
748 }
749 }
750 return integral;
751}
752
753////////////////////////////////////////////////////////////////////////////////
754/// Returns the content of the input bin
755/// For the overflow/underflow/sea bins:
756///~~~ {.cpp}
757/// -1 | -2 | -3
758/// ---+----+----
759/// -4 | -5 | -6
760/// ---+----+----
761/// -7 | -8 | -9
762///~~~
763/// where -5 is the "sea" bin (i.e. unbinned areas)
764
766{
767 if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
768 if (bin<0) return fOverflow[-bin - 1];
769 return ((TH2PolyBin*) fBins->At(bin-1))->GetContent();
770}
771
772////////////////////////////////////////////////////////////////////////////////
773/// Returns the value of error associated to bin number bin.
774/// If the sum of squares of weights has been defined (via Sumw2),
775/// this function returns the sqrt(sum of w2).
776/// otherwise it returns the sqrt(contents) for this bin.
777/// Bins are in range [1:nbins] and for bin < 0 in range [-9:-1] it returns errors for overflow bins.
778/// See also TH2Poly::GetBinContent
779
781{
782 if (bin == 0 || bin > GetNumberOfBins() || bin < - kNOverflow) return 0;
783 if (fBuffer) ((TH1*)this)->BufferEmpty();
784 // in case of weighted events the sum of the weights are stored in a different way than
785 // a normal histogram
786 // fSumw2.fArray[0:kNOverflow-1] : sum of weight squares for the overflow bins (
787 // fSumw2.fArray[kNOverflow:fNcells] : sum of weight squares for the standard bins
788 // fNcells = kNOverflow (9) + Number of bins
789 if (fSumw2.fN) {
790 Int_t binIndex = (bin > 0) ? bin+kNOverflow-1 : -(bin+1);
791 Double_t err2 = fSumw2.fArray[binIndex];
792 return TMath::Sqrt(err2);
793 }
794 Double_t error2 = TMath::Abs(GetBinContent(bin));
795 return TMath::Sqrt(error2);
796}
797
798////////////////////////////////////////////////////////////////////////////////
799/// Set the bin Error.
800/// Re-implementation for TH2Poly given the different bin indexing in the
801/// stored squared error array.
802/// See also notes in TH1::SetBinError
803///
804/// Bins are in range [1:nbins] and for bin < 0 in the range [-9:-1] the errors is set for the overflow bins
805
806
808{
809 if (bin == 0 || bin > GetNumberOfBins() || bin < - kNOverflow) return;
810 if (!fSumw2.fN) Sumw2();
812 // see comment in GetBinError for special convention of bin index in fSumw2 array
813 Int_t binIndex = (bin > 0) ? bin+kNOverflow-1 : -(bin+1);
814 fSumw2.fArray[binIndex] = error * error;
815}
816
817
818
819////////////////////////////////////////////////////////////////////////////////
820/// Returns the bin name.
821
822const char *TH2Poly::GetBinName(Int_t bin) const
823{
824 if (bin > GetNumberOfBins()) return "";
825 if (bin < 0) return "";
826 return ((TH2PolyBin*) fBins->At(bin-1))->GetPolygon()->GetName();
827}
828
829////////////////////////////////////////////////////////////////////////////////
830/// Returns the bin title.
831
832const char *TH2Poly::GetBinTitle(Int_t bin) const
833{
834 if (bin > GetNumberOfBins()) return "";
835 if (bin < 0) return "";
836 return ((TH2PolyBin*) fBins->At(bin-1))->GetPolygon()->GetTitle();
837}
838
839////////////////////////////////////////////////////////////////////////////////
840/// Returns the maximum value of the histogram.
841
843{
844 if (fNcells <= kNOverflow) return 0;
845 if (fMaximum != -1111) return fMaximum;
846
847 TH2PolyBin *b;
848
849 TIter next(fBins);
850 TObject *obj;
851 Double_t max,c;
852
853 max = ((TH2PolyBin*) next())->GetContent();
854
855 while ((obj=next())) {
856 b = (TH2PolyBin*)obj;
857 c = b->GetContent();
858 if (c>max) max = c;
859 }
860 return max;
861}
862
863////////////////////////////////////////////////////////////////////////////////
864/// Returns the maximum value of the histogram that is less than maxval.
865
867{
868 if (fNcells <= kNOverflow) return 0;
869 if (fMaximum != -1111) return fMaximum;
870
871 TH2PolyBin *b;
872
873 TIter next(fBins);
874 TObject *obj;
875 Double_t max,c;
876
877 max = ((TH2PolyBin*) next())->GetContent();
878
879 while ((obj=next())) {
880 b = (TH2PolyBin*)obj;
881 c = b->GetContent();
882 if (c>max && c<maxval) max=c;
883 }
884 return max;
885}
886
887////////////////////////////////////////////////////////////////////////////////
888/// Returns the minimum value of the histogram.
889
891{
892 if (fNcells <= kNOverflow) return 0;
893 if (fMinimum != -1111) return fMinimum;
894
895 TH2PolyBin *b;
896
897 TIter next(fBins);
898 TObject *obj;
899 Double_t min,c;
900
901 min = ((TH2PolyBin*) next())->GetContent();
902
903 while ((obj=next())) {
904 b = (TH2PolyBin*)obj;
905 c = b->GetContent();
906 if (c<min) min=c;
907 }
908 return min;
909}
910
911////////////////////////////////////////////////////////////////////////////////
912/// Returns the minimum value of the histogram that is greater than minval.
913
915{
916 if (fNcells <= kNOverflow) return 0;
917 if (fMinimum != -1111) return fMinimum;
918
919 TH2PolyBin *b;
920
921 TIter next(fBins);
922 TObject *obj;
923 Double_t min,c;
924
925 min = ((TH2PolyBin*) next())->GetContent();
926
927 while ((obj=next())) {
928 b = (TH2PolyBin*)obj;
929 c = b->GetContent();
930 if (c<min && c>minval) min=c;
931 }
932 return min;
933}
934
935////////////////////////////////////////////////////////////////////////////////
936/// Bins the histogram using a honeycomb structure
937
939 Int_t k, Int_t s)
940{
941 // Add the bins
942 Double_t numberOfHexagonsInTheRow;
943 Double_t x[6], y[6];
944 Double_t xloop, yloop, xtemp;
945 xloop = xstart; yloop = ystart + a/2.0;
946 for (int sCounter = 0; sCounter < s; sCounter++) {
947
948 xtemp = xloop; // Resets the temp variable
949
950 // Determine the number of hexagons in that row
951 if(sCounter%2 == 0){numberOfHexagonsInTheRow = k;}
952 else{numberOfHexagonsInTheRow = k - 1;}
953
954 for (int kCounter = 0; kCounter < numberOfHexagonsInTheRow; kCounter++) {
955
956 // Go around the hexagon
957 x[0] = xtemp;
958 y[0] = yloop;
959 x[1] = x[0];
960 y[1] = y[0] + a;
961 x[2] = x[1] + a*TMath::Sqrt(3)/2.0;
962 y[2] = y[1] + a/2.0;
963 x[3] = x[2] + a*TMath::Sqrt(3)/2.0;
964 y[3] = y[1];
965 x[4] = x[3];
966 y[4] = y[0];
967 x[5] = x[2];
968 y[5] = y[4] - a/2.0;
969
970 this->AddBin(6, x, y);
971
972 // Go right
973 xtemp += a*TMath::Sqrt(3);
974 }
975
976 // Increment the starting position
977 if (sCounter%2 == 0) xloop += a*TMath::Sqrt(3)/2.0;
978 else xloop -= a*TMath::Sqrt(3)/2.0;
979 yloop += 1.5*a;
980 }
981}
982
983////////////////////////////////////////////////////////////////////////////////
984/// Initializes the TH2Poly object. This method is called by the constructor.
985
987 Double_t ylow, Double_t yup, Int_t n, Int_t m)
988{
989 Int_t i;
990 fDimension = 2; //The dimension of the histogram
991
992 fBins = 0;
994
995 // Sets the boundaries of the histogram
996 fXaxis.Set(100, xlow, xup);
997 fYaxis.Set(100, ylow, yup);
998
999 for (i=0; i<9; i++) fOverflow[i] = 0.;
1000
1001 // Statistics
1002 fEntries = 0; // The total number of entries
1003 fTsumw = 0.; // Total amount of content in the histogram
1004 fTsumw2 = 0.; // Sum square of the weights
1005 fTsumwx = 0.; // Weighted sum of x coordinates
1006 fTsumwx2 = 0.; // Weighted sum of the squares of x coordinates
1007 fTsumwy2 = 0.; // Weighted sum of the squares of y coordinates
1008 fTsumwy = 0.; // Weighted sum of y coordinates
1009
1010 fCellX = n; // Set the number of cells to default
1011 fCellY = m; // Set the number of cells to default
1012
1013 // number of cells in the grid
1014 //N.B. not to be confused with fNcells (the number of bins) !
1016 fCells = new TList [fNCells]; // Sets an empty partition
1017 fStepX = (fXaxis.GetXmax() - fXaxis.GetXmin())/fCellX; // Cell width
1018 fStepY = (fYaxis.GetXmax() - fYaxis.GetXmin())/fCellY; // Cell height
1019
1020 fIsEmpty = new Bool_t [fNCells]; // Empty partition
1021 fCompletelyInside = new Bool_t [fNCells]; // Cell is completely inside bin
1022
1023 for (i = 0; i<fNCells; i++) { // Initializes the flags
1024 fIsEmpty[i] = kTRUE;
1026 }
1027
1028 // 3D Painter flags
1031}
1032
1033////////////////////////////////////////////////////////////////////////////////
1034/// Returns kTRUE if the input bin is intersecting with the
1035/// input rectangle (xclipl, xclipr, yclipb, yclipt)
1036
1038 Double_t xclipl, Double_t xclipr,
1039 Double_t yclipb, Double_t yclipt)
1040{
1041 Int_t gn;
1042 Double_t *gx;
1043 Double_t *gy;
1044 Bool_t inter = kFALSE;
1045 TObject *poly = bin->GetPolygon();
1046
1047 if (poly->IsA() == TGraph::Class()) {
1048 TGraph *g = (TGraph*)poly;
1049 gx = g->GetX();
1050 gy = g->GetY();
1051 gn = g->GetN();
1052 inter = IsIntersectingPolygon(gn, gx, gy, xclipl, xclipr, yclipb, yclipt);
1053 }
1054
1055 if (poly->IsA() == TMultiGraph::Class()) {
1056 TMultiGraph *mg = (TMultiGraph*)poly;
1057 TList *gl = mg->GetListOfGraphs();
1058 if (!gl) return inter;
1059 TGraph *g;
1060 TIter next(gl);
1061 while ((g = (TGraph*) next())) {
1062 gx = g->GetX();
1063 gy = g->GetY();
1064 gn = g->GetN();
1065 inter = IsIntersectingPolygon(gn, gx, gy, xclipl, xclipr,
1066 yclipb, yclipt);
1067 if (inter) return inter;
1068 }
1069 }
1070
1071 return inter;
1072}
1073
1074////////////////////////////////////////////////////////////////////////////////
1075/// Returns kTRUE if the input polygon (bn, x, y) is intersecting with the
1076/// input rectangle (xclipl, xclipr, yclipb, yclipt)
1077
1079 Double_t xclipl, Double_t xclipr,
1080 Double_t yclipb, Double_t yclipt)
1081{
1082 Bool_t p0R, p0L, p0T, p0B, p0xM, p0yM, p1R, p1L, p1T;
1083 Bool_t p1B, p1xM, p1yM, p0In, p1In;
1084
1085 for (int counter = 0; counter < (bn-1); counter++) {
1086 // If both are on the same side, return kFALSE
1087 p0L = x[counter] <= xclipl; // Point 0 is on the left
1088 p1L = x[counter + 1] <= xclipl; // Point 1 is on the left
1089 if (p0L && p1L) continue;
1090 p0R = x[counter] >= xclipr; // Point 0 is on the right
1091 p1R = x[counter + 1] >= xclipr; // Point 1 is on the right
1092 if (p0R && p1R) continue;
1093 p0T = y[counter] >= yclipt; // Point 0 is at the top
1094 p1T = y[counter + 1] >= yclipt; // Point 1 is at the top
1095 if (p0T && p1T) continue;
1096 p0B = y[counter] <= yclipb; // Point 0 is at the bottom
1097 p1B = y[counter + 1] <= yclipb; // Point 1 is at the bottom
1098 if (p0B && p1B) continue;
1099
1100 // Checks to see if any are inside
1101 p0xM = !p0R && !p0L; // Point 0 is inside along x
1102 p0yM = !p0T && !p0B; // Point 1 is inside along x
1103 p1xM = !p1R && !p1L; // Point 0 is inside along y
1104 p1yM = !p1T && !p1B; // Point 1 is inside along y
1105 p0In = p0xM && p0yM; // Point 0 is inside
1106 p1In = p1xM && p1yM; // Point 1 is inside
1107 if (p0In) {
1108 if (p1In) continue;
1109 return kTRUE;
1110 } else {
1111 if (p1In) return kTRUE;
1112 }
1113
1114 // We know by now that the points are not in the same side and not inside.
1115
1116 // Checks to see if they are opposite
1117
1118 if (p0xM && p1xM) return kTRUE;
1119 if (p0yM && p1yM) return kTRUE;
1120
1121 // We now know that the points are in different x and y indices
1122
1123 Double_t xcoord[3], ycoord[3];
1124 xcoord[0] = x[counter];
1125 xcoord[1] = x[counter + 1];
1126 ycoord[0] = y[counter];
1127 ycoord[1] = y[counter + 1];
1128
1129 if (p0L) {
1130 if(p1T){
1131 xcoord[2] = xclipl;
1132 ycoord[2] = yclipb;
1133 if((TMath::IsInside(xclipl, yclipt, 3, xcoord, ycoord)) ||
1134 (TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord))) continue;
1135 else return kTRUE;
1136 } else if (p1B) {
1137 xcoord[2] = xclipl;
1138 ycoord[2] = yclipt;
1139 if((TMath::IsInside(xclipl, yclipb, 3, xcoord, ycoord)) ||
1140 (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord))) continue;
1141 else return kTRUE;
1142 } else { // p1yM
1143 if (p0T) {
1144 xcoord[2] = xclipl;
1145 ycoord[2] = yclipb;
1146 if (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord)) continue;
1147 else return kTRUE;
1148 }
1149 if (p0B) {
1150 xcoord[2] = xclipl;
1151 ycoord[2] = yclipt;
1152 if (TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord)) continue;
1153 else return kTRUE;
1154 }
1155 }
1156 } else if (p0R) {
1157 if (p1T) {
1158 xcoord[2] = xclipl;
1159 ycoord[2] = yclipb;
1160 if ((TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord)) ||
1161 (TMath::IsInside(xclipl, yclipt, 3, xcoord, ycoord))) continue;
1162 else return kTRUE;
1163 } else if (p1B) {
1164 xcoord[2] = xclipl;
1165 ycoord[2] = yclipt;
1166 if ((TMath::IsInside(xclipl, yclipb, 3, xcoord, ycoord)) ||
1167 (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord))) continue;
1168 else return kTRUE;
1169 } else{ // p1yM
1170 if (p0T) {
1171 xcoord[2] = xclipr;
1172 ycoord[2] = yclipb;
1173 if (TMath::IsInside(xclipl, yclipt, 3, xcoord, ycoord)) continue;
1174 else return kTRUE;
1175 }
1176 if (p0B) {
1177 xcoord[2] = xclipr;
1178 ycoord[2] = yclipt;
1179 if (TMath::IsInside(xclipl, yclipb, 3, xcoord, ycoord)) continue;
1180 else return kTRUE;
1181 }
1182 }
1183 }
1184 }
1185 return kFALSE;
1186}
1187
1188////////////////////////////////////////////////////////////////////////////////
1189/// Merge TH2Polys
1190/// Given the special nature of the TH2Poly, the merge is implemented in
1191/// terms of subsequent TH2Poly::Add calls.
1193{
1194 for (auto h2pAsObj : *coll) {
1195 if (!Add((TH1*)h2pAsObj, 1.)) {
1196 Warning("Merge", "An issue was encountered during the merge operation.");
1197 return 0L;
1198 }
1199 }
1200 return GetEntries();
1201}
1202
1203////////////////////////////////////////////////////////////////////////////////
1204/// Save primitive as a C++ statement(s) on output stream out
1205
1206void TH2Poly::SavePrimitive(std::ostream &out, Option_t *option)
1207{
1208 out <<" "<<std::endl;
1209 out <<" "<< ClassName() <<" *";
1210
1211 //histogram pointer has by default the histogram name.
1212 //however, in case histogram has no directory, it is safer to add a
1213 //incremental suffix
1214 static Int_t hcounter = 0;
1215 TString histName = GetName();
1216 if (!fDirectory && !histName.Contains("Graph")) {
1217 hcounter++;
1218 histName += "__";
1219 histName += hcounter;
1220 }
1221 const char *hname = histName.Data();
1222
1223 //Construct the class initialization
1224 out << hname << " = new " << ClassName() << "(\"" << hname << "\", \""
1225 << GetTitle() << "\", " << fCellX << ", " << fXaxis.GetXmin()
1226 << ", " << fXaxis.GetXmax()
1227 << ", " << fCellY << ", " << fYaxis.GetXmin() << ", "
1228 << fYaxis.GetXmax() << ");" << std::endl;
1229
1230 // Save Bins
1231 TIter next(fBins);
1232 TObject *obj;
1233 TH2PolyBin *th2pBin;
1234
1235 while((obj = next())){
1236 th2pBin = (TH2PolyBin*) obj;
1237 th2pBin->GetPolygon()->SavePrimitive(out,
1238 TString::Format("th2poly%s",histName.Data()));
1239 }
1240
1241 // save bin contents
1242 out<<" "<<std::endl;
1243 Int_t bin;
1244 for (bin=1;bin<=GetNumberOfBins();bin++) {
1245 Double_t bc = GetBinContent(bin);
1246 if (bc) {
1247 out<<" "<<hname<<"->SetBinContent("<<bin<<","<<bc<<");"<<std::endl;
1248 }
1249 }
1250
1251 // save bin errors
1252 if (fSumw2.fN) {
1253 for (bin=1;bin<=GetNumberOfBins();bin++) {
1254 Double_t be = GetBinError(bin);
1255 if (be) {
1256 out<<" "<<hname<<"->SetBinError("<<bin<<","<<be<<");"<<std::endl;
1257 }
1258 }
1259 }
1260 TH1::SavePrimitiveHelp(out, hname, option);
1261}
1262
1263////////////////////////////////////////////////////////////////////////////////
1264/// Multiply this histogram by a constant c1.
1265
1267{
1268 for( int i = 0; i < this->GetNumberOfBins(); i++ ) {
1269 this->SetBinContent(i+1, c1*this->GetBinContent(i+1));
1270 }
1271 for( int i = 0; i < kNOverflow; i++ ) {
1272 this->SetBinContent(-i-1, c1*this->GetBinContent(-i-1) );
1273 }
1274}
1275
1276////////////////////////////////////////////////////////////////////////////////
1277/// Sets the contents of the input bin to the input content
1278/// Negative values between -1 and -9 are for the overflows and the sea
1279
1281{
1282 if (bin > GetNumberOfBins() || bin == 0 || bin < -9 ) return;
1283 if (bin > 0) {
1284 ((TH2PolyBin*) fBins->At(bin-1))->SetContent(content);
1285 }
1286 else
1287 fOverflow[-bin - 1] = content;
1288
1290 fEntries++;
1291}
1292
1293////////////////////////////////////////////////////////////////////////////////
1294/// When set to kTRUE, allows the histogram to expand if a bin outside the
1295/// limits is added.
1296
1298{
1299 fFloat = flag;
1300}
1301
1302////////////////////////////////////////////////////////////////////////////////
1303/// Return "true" if the point (x,y) is inside the bin of binnr.
1304
1306{
1307 if (!fBins) return false;
1308 TH2PolyBin* bin = (TH2PolyBin*)fBins->At(binnr);
1309 if (!bin) return false;
1310 return bin->IsInside(x,y);
1311}
1312
1313void TH2Poly::GetStats(Double_t *stats) const
1314{
1315 stats[0] = fTsumw;
1316 stats[1] = fTsumw2;
1317 stats[2] = fTsumwx;
1318 stats[3] = fTsumwx2;
1319 stats[4] = fTsumwy;
1320 stats[5] = fTsumwy2;
1321 stats[6] = fTsumwxy;
1322}
1323
1324/** \class TH2PolyBin
1325 \ingroup Hist
1326Helper class to represent a bin in the TH2Poly histogram
1327*/
1328
1329////////////////////////////////////////////////////////////////////////////////
1330/// Default constructor.
1331
1333{
1334 fPoly = 0;
1335 fContent = 0.;
1336 fNumber = 0;
1337 fXmax = -1111;
1338 fXmin = -1111;
1339 fYmax = -1111;
1340 fYmin = -1111;
1341 fArea = 0;
1343}
1344
1345////////////////////////////////////////////////////////////////////////////////
1346/// Normal constructor.
1347
1349{
1350 fContent = 0.;
1351 fNumber = bin_number;
1352 fArea = 0.;
1353 fPoly = poly;
1354 fXmax = -1111;
1355 fXmin = -1111;
1356 fYmax = -1111;
1357 fYmin = -1111;
1359}
1360
1361////////////////////////////////////////////////////////////////////////////////
1362/// Destructor.
1363
1365{
1366 if (fPoly) delete fPoly;
1367}
1368
1369////////////////////////////////////////////////////////////////////////////////
1370/// Returns the area of the bin.
1371
1373{
1374 Int_t bn;
1375
1376 if (fArea == 0) {
1377 if (fPoly->IsA() == TGraph::Class()) {
1378 TGraph *g = (TGraph*)fPoly;
1379 bn = g->GetN();
1380 fArea = g->Integral(0,bn-1);
1381 }
1382
1383 if (fPoly->IsA() == TMultiGraph::Class()) {
1385 TList *gl = mg->GetListOfGraphs();
1386 if (!gl) return fArea;
1387 TGraph *g;
1388 TIter next(gl);
1389 while ((g = (TGraph*) next())) {
1390 bn = g->GetN();
1391 fArea = fArea + g->Integral(0,bn-1);
1392 }
1393 }
1394 }
1395
1396 return fArea;
1397}
1398
1399////////////////////////////////////////////////////////////////////////////////
1400/// Returns the maximum value for the x coordinates of the bin.
1401
1403{
1404 if (fXmax != -1111) return fXmax;
1405
1406 Int_t bn,i;
1407 Double_t *bx;
1408
1409 if (fPoly->IsA() == TGraph::Class()) {
1410 TGraph *g = (TGraph*)fPoly;
1411 bx = g->GetX();
1412 bn = g->GetN();
1413 fXmax = bx[0];
1414 for (i=1; i<bn; i++) {if (fXmax < bx[i]) fXmax = bx[i];}
1415 }
1416
1417 if (fPoly->IsA() == TMultiGraph::Class()) {
1419 TList *gl = mg->GetListOfGraphs();
1420 if (!gl) return fXmax;
1421 TGraph *g;
1422 TIter next(gl);
1423 Bool_t first = kTRUE;
1424 while ((g = (TGraph*) next())) {
1425 bx = g->GetX();
1426 bn = g->GetN();
1427 if (first) {fXmax = bx[0]; first = kFALSE;}
1428 for (i=0; i<bn; i++) {if (fXmax < bx[i]) fXmax = bx[i];}
1429 }
1430 }
1431
1432 return fXmax;
1433}
1434
1435////////////////////////////////////////////////////////////////////////////////
1436/// Returns the minimum value for the x coordinates of the bin.
1437
1439{
1440 if (fXmin != -1111) return fXmin;
1441
1442 Int_t bn,i;
1443 Double_t *bx;
1444
1445 if (fPoly->IsA() == TGraph::Class()) {
1446 TGraph *g = (TGraph*)fPoly;
1447 bx = g->GetX();
1448 bn = g->GetN();
1449 fXmin = bx[0];
1450 for (i=1; i<bn; i++) {if (fXmin > bx[i]) fXmin = bx[i];}
1451 }
1452
1453 if (fPoly->IsA() == TMultiGraph::Class()) {
1455 TList *gl = mg->GetListOfGraphs();
1456 if (!gl) return fXmin;
1457 TGraph *g;
1458 TIter next(gl);
1459 Bool_t first = kTRUE;
1460 while ((g = (TGraph*) next())) {
1461 bx = g->GetX();
1462 bn = g->GetN();
1463 if (first) {fXmin = bx[0]; first = kFALSE;}
1464 for (i=0; i<bn; i++) {if (fXmin > bx[i]) fXmin = bx[i];}
1465 }
1466 }
1467
1468 return fXmin;
1469}
1470
1471////////////////////////////////////////////////////////////////////////////////
1472/// Returns the maximum value for the y coordinates of the bin.
1473
1475{
1476 if (fYmax != -1111) return fYmax;
1477
1478 Int_t bn,i;
1479 Double_t *by;
1480
1481 if (fPoly->IsA() == TGraph::Class()) {
1482 TGraph *g = (TGraph*)fPoly;
1483 by = g->GetY();
1484 bn = g->GetN();
1485 fYmax = by[0];
1486 for (i=1; i<bn; i++) {if (fYmax < by[i]) fYmax = by[i];}
1487 }
1488
1489 if (fPoly->IsA() == TMultiGraph::Class()) {
1491 TList *gl = mg->GetListOfGraphs();
1492 if (!gl) return fYmax;
1493 TGraph *g;
1494 TIter next(gl);
1495 Bool_t first = kTRUE;
1496 while ((g = (TGraph*) next())) {
1497 by = g->GetY();
1498 bn = g->GetN();
1499 if (first) {fYmax = by[0]; first = kFALSE;}
1500 for (i=0; i<bn; i++) {if (fYmax < by[i]) fYmax = by[i];}
1501 }
1502 }
1503
1504 return fYmax;
1505}
1506
1507////////////////////////////////////////////////////////////////////////////////
1508/// Returns the minimum value for the y coordinates of the bin.
1509
1511{
1512 if (fYmin != -1111) return fYmin;
1513
1514 Int_t bn,i;
1515 Double_t *by;
1516
1517 if (fPoly->IsA() == TGraph::Class()) {
1518 TGraph *g = (TGraph*)fPoly;
1519 by = g->GetY();
1520 bn = g->GetN();
1521 fYmin = by[0];
1522 for (i=1; i<bn; i++) {if (fYmin > by[i]) fYmin = by[i];}
1523 }
1524
1525 if (fPoly->IsA() == TMultiGraph::Class()) {
1527 TList *gl = mg->GetListOfGraphs();
1528 if (!gl) return fYmin;
1529 TGraph *g;
1530 TIter next(gl);
1531 Bool_t first = kTRUE;
1532 while ((g = (TGraph*) next())) {
1533 by = g->GetY();
1534 bn = g->GetN();
1535 if (first) {fYmin = by[0]; first = kFALSE;}
1536 for (i=0; i<bn; i++) {if (fYmin > by[i]) fYmin = by[i];}
1537 }
1538 }
1539
1540 return fYmin;
1541}
1542
1543////////////////////////////////////////////////////////////////////////////////
1544/// Return "true" if the point (x,y) is inside the bin.
1545
1547{
1548 Int_t in=0;
1549
1550 if (fPoly->IsA() == TGraph::Class()) {
1551 TGraph *g = (TGraph*)fPoly;
1552 in = g->IsInside(x, y);
1553 }
1554
1555 if (fPoly->IsA() == TMultiGraph::Class()) {
1557 in = mg->IsInside(x, y);
1558 }
1559
1560 return in;
1561}
1562
1563////////////////////////////////////////////////////////////////////////
1564/// RE-implement dummy functions to avoid users calling the
1565/// corresponding implemntations in TH1 or TH2
1566//////////////////////////////////////////////////////////////////////////
1567
1568////////////////////////////////////////////////////////////////////////////////
1569/// Performs the operation: this = this + c1*f1. NOT IMPLEMENTED for TH2Poly
1571{
1572 Error("Add","Not implement for TH2Poly");
1573 return kFALSE;
1574}
1575
1576////////////////////////////////////////////////////////////////////////////////
1577/// Replace contents of this histogram by the addition of h1 and h2. NOT IMPLEMENTED for TH2Poly
1579{
1580 Error("Add","Not implement for TH2Poly");
1581 return kFALSE;
1582}
1583
1584////////////////////////////////////////////////////////////////////////////////
1585/// Performs the operation: this = this / c1*f1. NOT IMPLEMENTED for TH2Poly
1587{
1588 Error("Divide","Not implement for TH2Poly");
1589 return kFALSE;
1590}
1591
1592////////////////////////////////////////////////////////////////////////////////
1593/// NOT IMPLEMENTED for TH2Poly
1595{
1596 Error("Multiply","Not implement for TH2Poly");
1597 return kFALSE;
1598}
1599////////////////////////////////////////////////////////////////////////////////
1600/// NOT IMPLEMENTED for TH2Poly
1602{
1603 Error("ComputeIntegral","Not implement for TH2Poly");
1604 return TMath::QuietNaN();
1605}
1606////////////////////////////////////////////////////////////////////////////////
1607/// NOT IMPLEMENTED for TH2Poly
1609{
1610 Error("FFT","Not implement for TH2Poly");
1611 return nullptr;
1612}
1613////////////////////////////////////////////////////////////////////////////////
1614/// NOT IMPLEMENTED for TH2Poly
1616{
1617 Error("GetAsymmetry","Not implement for TH2Poly");
1618 return nullptr;
1619}
1620////////////////////////////////////////////////////////////////////////////////
1621/// NOT IMPLEMENTED for TH2Poly
1623{
1624 Error("Interpolate","Not implement for TH2Poly");
1625 return TMath::QuietNaN();
1626}
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define g(i)
Definition RSha256.hxx:105
#define a(i)
Definition RSha256.hxx:99
#define s1(x)
Definition RSha256.hxx:91
static const double x2[5]
static const double x1[5]
int Int_t
Definition RtypesCore.h:45
const Bool_t kFALSE
Definition RtypesCore.h:92
double Double_t
Definition RtypesCore.h:59
long long Long64_t
Definition RtypesCore.h:73
const Bool_t kTRUE
Definition RtypesCore.h:91
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:364
@ kCounter
Definition TDataType.h:34
char name[80]
Definition TGX11.cxx:110
double floor(double)
Double_t * fArray
Definition TArrayD.h:30
void Set(Int_t n)
Set size of this array to n doubles.
Definition TArrayD.cxx:106
Int_t fN
Definition TArray.h:38
Double_t GetXmax() const
Definition TAxis.h:134
virtual void Set(Int_t nbins, Double_t xmin, Double_t xmax)
Initialize axis with fix bins.
Definition TAxis.cxx:731
Double_t GetXmin() const
Definition TAxis.h:133
Collection abstract base class.
Definition TCollection.h:63
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
1-Dim function class
Definition TF1.h:213
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
Double_t * fBuffer
[fBufferSize] entry buffer
Definition TH1.h:107
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition TH1.cxx:6678
Int_t fNcells
number of bins(1D), cells (2D) +U/Overflows
Definition TH1.h:88
virtual void GetStats(Double_t *stats) const
fill the array stats from the contents of this histogram The array stats must be correctly dimensione...
Definition TH1.cxx:7726
Double_t fTsumw
Total Sum of weights.
Definition TH1.h:95
Double_t fTsumw2
Total Sum of squares of weights.
Definition TH1.h:96
Double_t fTsumwx2
Total Sum of weight*X*X.
Definition TH1.h:98
virtual Double_t GetNormFactor() const
Definition TH1.h:300
@ kIsNotW
Histogram is forced to be not weighted even when the histogram is filled with weighted different than...
Definition TH1.h:171
TDirectory * fDirectory
!Pointer to directory holding this histogram
Definition TH1.h:108
Double_t fMaximum
Maximum value for plotting.
Definition TH1.h:99
Int_t fDimension
!Histogram dimension (1, 2 or 3 dim)
Definition TH1.h:109
virtual void SetContent(const Double_t *content)
Replace bin contents by the contents of array content.
Definition TH1.cxx:8246
virtual void SavePrimitiveHelp(std::ostream &out, const char *hname, Option_t *option="")
Helper function for the SavePrimitive functions from TH1 or classes derived from TH1,...
Definition TH1.cxx:7264
virtual Double_t GetBinErrorSqUnchecked(Int_t bin) const
Definition TH1.h:443
Double_t fMinimum
Minimum value for plotting.
Definition TH1.h:100
@ kNstat
Definition TH1.h:183
virtual Double_t GetEntries() const
Return the current number of entries.
Definition TH1.cxx:4386
virtual void ResetStats()
Reset the statistics including the number of entries and replace with values calculated from bin cont...
Definition TH1.cxx:7795
virtual void SetBinErrorOption(EBinErrorOpt type)
Definition TH1.h:376
Double_t fEntries
Number of entries.
Definition TH1.h:94
virtual void SetName(const char *name)
Change the name of this histogram.
Definition TH1.cxx:8800
@ kNormal
errors with Normal (Wald) approximation: errorUp=errorLow= sqrt(N)
Definition TH1.h:64
TAxis fXaxis
X axis descriptor.
Definition TH1.h:89
TArrayD fSumw2
Array of sum of squares of weights.
Definition TH1.h:103
virtual Int_t GetSumw2N() const
Definition TH1.h:314
TAxis fYaxis
Y axis descriptor.
Definition TH1.h:90
virtual Double_t GetSumOfWeights() const
Return the sum of weights excluding under/overflows.
Definition TH1.cxx:7810
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition TH1.cxx:8860
virtual void SetEntries(Double_t n)
Definition TH1.h:385
Double_t fTsumwx
Total Sum of weight*X.
Definition TH1.h:97
Helper class to represent a bin in the TH2Poly histogram.
Definition TH2Poly.h:25
Double_t GetXMin()
Returns the minimum value for the x coordinates of the bin.
Definition TH2Poly.cxx:1438
Double_t GetYMax()
Returns the maximum value for the y coordinates of the bin.
Definition TH2Poly.cxx:1474
Double_t GetArea()
Returns the area of the bin.
Definition TH2Poly.cxx:1372
void ClearContent()
Definition TH2Poly.h:32
void Fill(Double_t w)
Definition TH2Poly.h:33
Double_t GetYMin()
Returns the minimum value for the y coordinates of the bin.
Definition TH2Poly.cxx:1510
Double_t fArea
Bin area.
Definition TH2Poly.h:51
Double_t fContent
Bin content.
Definition TH2Poly.h:52
Bool_t IsInside(Double_t x, Double_t y) const
Return "true" if the point (x,y) is inside the bin.
Definition TH2Poly.cxx:1546
Double_t fXmax
X maximum value.
Definition TH2Poly.h:55
Int_t fNumber
Bin number of the bin in TH2Poly.
Definition TH2Poly.h:49
TH2PolyBin()
Default constructor.
Definition TH2Poly.cxx:1332
Double_t fYmax
Y maximum value.
Definition TH2Poly.h:56
Double_t GetXMax()
Returns the maximum value for the x coordinates of the bin.
Definition TH2Poly.cxx:1402
Double_t GetContent() const
Definition TH2Poly.h:35
Double_t fYmin
Y minimum value.
Definition TH2Poly.h:54
void SetChanged(Bool_t flag)
Definition TH2Poly.h:44
Int_t GetBinNumber() const
Definition TH2Poly.h:37
Double_t fXmin
X minimum value.
Definition TH2Poly.h:53
TObject * GetPolygon() const
Definition TH2Poly.h:38
virtual ~TH2PolyBin()
Destructor.
Definition TH2Poly.cxx:1364
TObject * fPoly
Object holding the polygon definition.
Definition TH2Poly.h:50
2D Histogram with Polygonal Bins
Definition TH2Poly.h:66
Double_t Integral(Option_t *option="") const
Returns the integral of bin contents.
Definition TH2Poly.cxx:725
virtual void UpdateBinContent(Int_t bin, Double_t content)
Raw update of bin content on internal data structure see convention for numbering bins in TH1::GetBin...
Definition TH2Poly.h:177
TList * GetBins()
Returns the TList of all bins in the histogram.
Definition TH2Poly.h:98
void ClearBinContents()
Clears the contents of all bins in the histogram.
Definition TH2Poly.cxx:491
Double_t fOverflow[kNOverflow]
Overflow bins.
Definition TH2Poly.h:156
Bool_t fFloat
When set to kTRUE, allows the histogram to expand if a bin outside the limits is added.
Definition TH2Poly.h:164
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error.
Definition TH2Poly.cxx:807
Bool_t IsIntersecting(TH2PolyBin *bin, Double_t xclipl, Double_t xclipr, Double_t yclipb, Double_t yclipt)
Returns kTRUE if the input bin is intersecting with the input rectangle (xclipl, xclipr,...
Definition TH2Poly.cxx:1037
void Honeycomb(Double_t xstart, Double_t ystart, Double_t a, Int_t k, Int_t s)
Bins the histogram using a honeycomb structure.
Definition TH2Poly.cxx:938
void SetFloat(Bool_t flag=true)
When set to kTRUE, allows the histogram to expand if a bin outside the limits is added.
Definition TH2Poly.cxx:1297
virtual void GetStats(Double_t *stats) const
Fill the array stats from the contents of this histogram The array stats must be correctly dimensione...
Definition TH2Poly.cxx:1313
virtual TH2PolyBin * CreateBin(TObject *poly)
Create appropriate histogram bin.
Definition TH2Poly.cxx:200
virtual Double_t RetrieveBinContent(Int_t bin) const
Raw retrieval of bin content on internal data structure see convention for numbering bins in TH1::Get...
Definition TH2Poly.h:174
Bool_t * fIsEmpty
[fNCells] The array that returns true if the cell at the given coordinate is empty
Definition TH2Poly.h:162
TList * fBins
List of bins. The list owns the contained objects.
Definition TH2Poly.h:167
void AddBinToPartition(TH2PolyBin *bin)
Adds the input bin into the partition cell matrix.
Definition TH2Poly.cxx:368
@ kNOverflow
Definition TH2Poly.h:154
Int_t GetNumberOfBins() const
Definition TH2Poly.h:110
void SetBinContentChanged(Bool_t flag)
Definition TH2Poly.h:119
virtual ~TH2Poly()
Destructor.
Definition TH2Poly.cxx:184
virtual Int_t Fill(Double_t x, Double_t y)
Increment the bin containing (x,y) by 1.
Definition TH2Poly.cxx:589
Bool_t * fCompletelyInside
[fNCells] The array that returns true if the cell at the given coordinate is completely inside a bin
Definition TH2Poly.h:163
TH2Poly()
Default Constructor. No boundaries specified.
Definition TH2Poly.cxx:147
Bool_t IsInsideBin(Int_t binnr, Double_t x, Double_t y)
Return "true" if the point (x,y) is inside the bin of binnr.
Definition TH2Poly.cxx:1305
void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
Definition TH2Poly.cxx:1206
TObject * Clone(const char *newname="") const
Make a complete copy of the underlying object.
Definition TH2Poly.cxx:479
virtual void Reset(Option_t *option)
Reset this histogram: contents, errors, etc.
Definition TH2Poly.cxx:516
const char * GetBinTitle(Int_t bin) const
Returns the bin title.
Definition TH2Poly.cxx:832
void SetNewBinAdded(Bool_t flag)
Definition TH2Poly.h:121
void ChangePartition(Int_t n, Int_t m)
Changes the number of partition cells in the histogram.
Definition TH2Poly.cxx:440
virtual Double_t Interpolate(Double_t, Double_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1622
Double_t GetMinimum() const
Returns the minimum value of the histogram.
Definition TH2Poly.cxx:890
const char * GetBinName(Int_t bin) const
Returns the bin name.
Definition TH2Poly.cxx:822
Double_t fStepX
Definition TH2Poly.h:161
virtual Bool_t Multiply(TF1 *, Double_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1594
Int_t fNCells
Number of partition cells: fCellX*fCellY.
Definition TH2Poly.h:159
virtual Double_t GetBinContent(Int_t bin) const
Returns the content of the input bin For the overflow/underflow/sea bins:
Definition TH2Poly.cxx:765
void Initialize(Double_t xlow, Double_t xup, Double_t ylow, Double_t yup, Int_t n, Int_t m)
Initializes the TH2Poly object. This method is called by the constructor.
Definition TH2Poly.cxx:986
Int_t fCellX
Number of partition cells in the x-direction of the histogram.
Definition TH2Poly.h:157
virtual Bool_t Add(const TH1 *h1, Double_t c1)
Performs the operation: this = this + c1*h1.
Definition TH2Poly.cxx:290
virtual Int_t AddBin(TObject *poly)
Adds a new bin to the histogram.
Definition TH2Poly.cxx:222
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
Definition TH2Poly.cxx:1266
Long64_t Merge(TCollection *)
Merge TH2Polys Given the special nature of the TH2Poly, the merge is implemented in terms of subseque...
Definition TH2Poly.cxx:1192
Bool_t IsIntersectingPolygon(Int_t bn, Double_t *x, Double_t *y, Double_t xclipl, Double_t xclipr, Double_t yclipb, Double_t yclipt)
Returns kTRUE if the input polygon (bn, x, y) is intersecting with the input rectangle (xclipl,...
Definition TH2Poly.cxx:1078
virtual TH1 * FFT(TH1 *, Option_t *)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1608
virtual Double_t GetBinError(Int_t bin) const
Returns the value of error associated to bin number bin.
Definition TH2Poly.cxx:780
virtual TH1 * GetAsymmetry(TH1 *, Double_t, Double_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1615
virtual void SetBinContent(Int_t bin, Double_t content)
Sets the contents of the input bin to the input content Negative values between -1 and -9 are for the...
Definition TH2Poly.cxx:1280
Double_t GetMaximum() const
Returns the maximum value of the histogram.
Definition TH2Poly.cxx:842
Int_t FindBin(Double_t x, Double_t y, Double_t z=0)
Returns the bin number of the bin at the given coordinate.
Definition TH2Poly.cxx:546
virtual Double_t ComputeIntegral(Bool_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1601
Double_t fStepY
Dimensions of a partition cell.
Definition TH2Poly.h:161
void FillN(Int_t ntimes, const Double_t *x, const Double_t *y, const Double_t *w, Int_t stride=1)
Fills a 2-D histogram with an array of values and weights.
Definition TH2Poly.cxx:711
TList * fCells
[fNCells] The array of TLists that store the bins that intersect with each cell. List do not own the ...
Definition TH2Poly.h:160
virtual Bool_t Divide(TF1 *, Double_t)
Performs the operation: this = this / c1*f1. NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1586
Int_t fCellY
Number of partition cells in the y-direction of the histogram.
Definition TH2Poly.h:158
virtual void PutStats(Double_t *stats)
Replace current statistics with the values in array stats.
Definition TH2.cxx:2376
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
Definition TH2.cxx:2491
Double_t fTsumwxy
Definition TH2.h:36
Double_t fTsumwy2
Definition TH2.h:35
Double_t fTsumwy
Definition TH2.h:34
A doubly linked list.
Definition TList.h:44
virtual void Add(TObject *obj)
Definition TList.h:87
virtual TObject * At(Int_t idx) const
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:357
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition TMultiGraph.h:36
TList * GetListOfGraphs() const
Definition TMultiGraph.h:70
virtual Int_t IsInside(Double_t x, Double_t y) const
Return 1 if the point (x,y) is inside one of the graphs 0 otherwise.
virtual const char * GetTitle() const
Returns title of object.
Definition TNamed.h:48
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition TNamed.cxx:74
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
Mother of all ROOT objects.
Definition TObject.h:37
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:359
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:187
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:130
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:879
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition TObject.cxx:666
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:893
Basic string class.
Definition TString.h:136
void ToLower()
Change string to lower-case.
Definition TString.cxx:1145
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition TString.cxx:438
const char * Data() const
Definition TString.h:369
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2331
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:624
return c1
Definition legend1.C:41
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
TH1F * h1
Definition legend1.C:5
Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y)
Function which returns kTRUE if point xp,yp lies inside the polygon defined by the np points in array...
Definition TMath.h:1209
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754
Definition TMath.h:901
Double_t Sqrt(Double_t x)
Definition TMath.h:691
Short_t Abs(Short_t d)
Definition TMathBase.h:120
Definition first.py:1
auto * m
Definition textangle.C:8