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 "TInterpreter.h"
17#include "Riostream.h"
18#include "TList.h"
19#include "TMath.h"
20#include <cassert>
21
22
23/** \class TH2Poly
24 \ingroup Histograms
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 auto 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~~~
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/// Copy constructor
183TH2Poly::TH2Poly(const TH2Poly & rhs) : TH2(), fCells(nullptr),
184fIsEmpty(nullptr), fCompletelyInside(nullptr), fBins(nullptr) {
185 rhs.Copy(*this);
186}
187
188////////////////////////////////////////////////////////////////////////////////
189/// Destructor.
190
192{
193 delete[] fCells;
194 delete[] fIsEmpty;
195 delete[] fCompletelyInside;
196 // delete at the end the bin List since it owns the objects
197 delete fBins;
198}
199
200/////////////////////////////////////////////////////////////////////////////////
201/// Assignment operator
203 if (this != &rhs)
204 rhs.Copy(*this);
205 return *this;
206}
207
208////////////////////////////////////////////////////////////////////////////////
209/// Copy function for TH2Poly
210
212{
213 // copy first TH2 information
215 auto & newth2p = dynamic_cast<TH2Poly &>(newobj);
216 newth2p.SetName(GetName());
217 newth2p.SetTitle(GetTitle());
218
219 newth2p.fCellX = fCellX;
220 newth2p.fCellY = fCellY;
221 newth2p.fNCells = fNCells;
222 newth2p.fStepX = fStepX;
223 newth2p.fStepY = fStepY;
224
225 // deallocate previous arrays, if existing
226 delete[] newth2p.fCells;
227 delete[] newth2p.fIsEmpty;
228 delete[] newth2p.fCompletelyInside;
229 // allocate arrays
230 newth2p.fCells = new TList [fNCells];
231 newth2p.fIsEmpty = new Bool_t [fNCells]; // Empty partition
232 newth2p.fCompletelyInside = new Bool_t [fNCells]; // Cell is completely inside bin
233 // Initializes the flags
234 for (int i = 0; i<fNCells; i++) {
235 newth2p.fIsEmpty[i] = fIsEmpty[i];
236 newth2p.fCompletelyInside[i] = fCompletelyInside[i];
237 }
238 // need to use Clone to copy the contained bin list
239 delete newth2p.fBins; // in case there was something before there
240 if (!fBins) {
241 newth2p.fBins = nullptr;
242 }
243 else {
244 newth2p.fBins = dynamic_cast<TList *>(fBins->Clone());
245 if (!newth2p.fBins)
246 Error("Copy","Error cloning the TH2Poly bin list");
247 else {
248 // add bins in the fCells partition. We need to add the TH2PolyBin objects
249 // of the new copied histograms. For this we call AddBinToPartition
250 // we could probably optimize this by implementing a copy of the partition
251 for (auto bin : *(newth2p.fBins)) {
252 newth2p.AddBinToPartition(dynamic_cast<TH2PolyBin*>(bin));
253 }
254 }
255 }
256 // copy overflow contents
257 for(int i = 0; i < kNOverflow; i++ ) {
258 newth2p.fOverflow[i] = fOverflow[i];
259 }
260 // copy other data members
261 newth2p.fFloat = fFloat;
262 newth2p.fNewBinAdded = fNewBinAdded;
263 newth2p.fBinContentChanged = fBinContentChanged;
264}
265
266
267////////////////////////////////////////////////////////////////////////////////
268/// Create appropriate histogram bin.
269/// e.g. TH2Poly creates TH2PolyBin,
270/// TProfile2Poly creates TProfile2PolyBin
271/// This is done so that TH2Poly::AddBin does not have to be duplicated,
272/// but only create needs to be reimplemented for additional histogram types
273
275{
276 if (!poly) return nullptr;
277
278 if (fBins == nullptr) {
279 fBins = new TList();
280 fBins->SetOwner();
281 }
282
283 fNcells++;
285 // if structure fsumw2 is created extend it
286 if (fSumw2.fN) fSumw2.Set(fNcells);
287 return new TH2PolyBin(poly, ibin);
288}
289
290////////////////////////////////////////////////////////////////////////////////
291/// Adds a new bin to the histogram. It can be any object having the method
292/// IsInside(). It returns the bin number in the histogram. It returns 0 if
293/// it failed to add. To allow the histogram limits to expand when a bin
294/// outside the limits is added, call SetFloat() before adding the bin.
295
297{
298 auto *bin = CreateBin(poly);
300 if(!bin) return 0;
301
302 // If the bin lies outside histogram boundaries, then extends the boundaries.
303 // Also changes the partition information accordingly
305 if (fFloat) {
306 if (fXaxis.GetXmin() > bin->GetXMin()) {
307 fXaxis.Set(100, bin->GetXMin(), fXaxis.GetXmax());
308 flag = kTRUE;
309 }
310 if (fXaxis.GetXmax() < bin->GetXMax()) {
311 fXaxis.Set(100, fXaxis.GetXmin(), bin->GetXMax());
312 flag = kTRUE;
313 }
314 if (fYaxis.GetXmin() > bin->GetYMin()) {
315 fYaxis.Set(100, bin->GetYMin(), fYaxis.GetXmax());
316 flag = kTRUE;
317 }
318 if (fYaxis.GetXmax() < bin->GetYMax()) {
319 fYaxis.Set(100, fYaxis.GetXmin(), bin->GetYMax());
320 flag = kTRUE;
321 }
323 } else {
324 /*Implement polygon clipping code here*/
325 }
326
327 fBins->Add((TObject*) bin);
329
330 // Adds the bin to the partition matrix
332
333 return ibin;
334}
335
336////////////////////////////////////////////////////////////////////////////////
337/// Adds a new bin to the histogram. The number of vertices and their (x,y)
338/// coordinates are required as input. It returns the bin number in the
339/// histogram.
340
342{
343 TGraph *g = new TGraph(n, x, y);
344 Int_t bin = AddBin(g);
345 return bin;
346}
347
348////////////////////////////////////////////////////////////////////////////////
349/// Add a new bin to the histogram. The bin shape is a rectangle.
350/// It returns the bin number of the bin in the histogram.
351
353{
354 Double_t x[] = {x1, x1, x2, x2, x1};
355 Double_t y[] = {y1, y2, y2, y1, y1};
356 TGraph *g = new TGraph(5, x, y);
357 Int_t bin = AddBin(g);
358 return bin;
359}
360
361////////////////////////////////////////////////////////////////////////////////
362/// Performs the operation: this = this + c1*h1.
363
365{
366 Int_t bin;
367
368 TH2Poly *h1p = (TH2Poly *)h1;
369
370 // Check if number of bins is the same.
371 if (h1p->GetNumberOfBins() != GetNumberOfBins()) {
372 Error("Add", "Attempt to add histograms with different number of bins");
373 return kFALSE;
374 }
375
376 // Check if the bins are the same.
377 TList *h1pBins = h1p->GetBins();
379 for (bin = 1; bin <= GetNumberOfBins(); bin++) {
380 thisBin = (TH2PolyBin *)fBins->At(bin - 1);
381 h1pBin = (TH2PolyBin *)h1pBins->At(bin - 1);
382 if (thisBin->GetXMin() != h1pBin->GetXMin() ||
383 thisBin->GetXMax() != h1pBin->GetXMax() ||
384 thisBin->GetYMin() != h1pBin->GetYMin() ||
385 thisBin->GetYMax() != h1pBin->GetYMax()) {
386 Error("Add", "Attempt to add histograms with different bin limits");
387 return kFALSE;
388 }
389 }
390
391
392 // Create Sumw2 if h1p has Sumw2 set
393 if (fSumw2.fN == 0 && h1p->GetSumw2N() != 0) Sumw2();
394
395 // statistics can be preserved only in case of positive coefficients
396 // otherwise with negative c1 (histogram subtraction) one risks to get negative variances
397 Bool_t resetStats = (c1 < 0);
398 Double_t s1[kNstat] = {0};
399 Double_t s2[kNstat] = {0};
400 if (!resetStats) {
401 // need to initialize to zero s1 and s2 since
402 // GetStats fills only used elements depending on dimension and type
403 GetStats(s1);
404 h1->GetStats(s2);
405 }
406 // get number of entries now because afterwards UpdateBinContent will change it
407 Double_t entries = TMath::Abs( GetEntries() + c1 * h1->GetEntries() );
408
409
410 // Perform the Add.
411 Double_t factor = 1;
412 if (h1p->GetNormFactor() != 0)
413 factor = h1p->GetNormFactor() / h1p->GetSumOfWeights();
414 for (bin = 0; bin < fNcells; bin++) {
415 Double_t y = this->RetrieveBinContent(bin) + c1 * h1p->RetrieveBinContent(bin);
416 UpdateBinContent(bin, y);
417 if (fSumw2.fN) {
418 Double_t esq = factor * factor * h1p->GetBinErrorSqUnchecked(bin);
419 fSumw2.fArray[bin] += c1 * c1 * factor * factor * esq;
420 }
421 }
422
423 // update statistics (do here to avoid changes by SetBinContent)
424 if (resetStats) {
425 // statistics need to be reset in case coefficient are negative
426 ResetStats();
427 } else {
428 for (Int_t i = 0; i < kNstat; i++) {
429 if (i == 1) s1[i] += c1 * c1 * s2[i];
430 else s1[i] += c1 * s2[i];
431 }
432 PutStats(s1);
433 SetEntries(entries);
434 }
435 return kTRUE;
436}
437
438////////////////////////////////////////////////////////////////////////////////
439/// Adds the input bin into the partition cell matrix. This method is called
440/// in AddBin() and ChangePartition().
441
443{
444 // Cell Info
445 Int_t nl, nr, mb, mt; // Max/min indices of the cells that contain the bin
446 Double_t xclipl, xclipr, yclipb, yclipt; // x and y coordinates of a cell
447 Double_t binXmax, binXmin, binYmax, binYmin; // The max/min bin coordinates
448
449 binXmax = bin->GetXMax();
450 binXmin = bin->GetXMin();
451 binYmax = bin->GetYMax();
452 binYmin = bin->GetYMin();
453 nl = (Int_t)(floor((binXmin - fXaxis.GetXmin())/fStepX));
454 nr = (Int_t)(floor((binXmax - fXaxis.GetXmin())/fStepX));
455 mb = (Int_t)(floor((binYmin - fYaxis.GetXmin())/fStepY));
456 mt = (Int_t)(floor((binYmax - fYaxis.GetXmin())/fStepY));
457
458 // Make sure the array indices are correct.
459 if (nr>=fCellX) nr = fCellX-1;
460 if (mt>=fCellY) mt = fCellY-1;
461 if (nl<0) nl = 0;
462 if (mb<0) mb = 0;
463
464 // number of cells in the grid
465 //N.B. not to be confused with fNcells (the number of bins) !
467
468 // Loop over all cells
469 for (int i = nl; i <= nr; i++) {
470 xclipl = fXaxis.GetXmin() + i*fStepX;
471 xclipr = xclipl + fStepX;
472 for (int j = mb; j <= mt; j++) {
474 yclipt = yclipb + fStepY;
475
476 // If the bin is completely inside the cell,
477 // add that bin to the cell then return
478 if ((binXmin >= xclipl) && (binXmax <= xclipr) &&
479 (binYmax <= yclipt) && (binYmin >= yclipb)){
480 fCells[i + j*fCellX].Add((TObject*) bin);
481 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
482 return;
483 }
484
485 // If any of the sides of the cell intersect with any side of the bin,
486 // add that bin then continue
488 fCells[i + j*fCellX].Add((TObject*) bin);
489 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
490 continue;
491 }
492 // If a corner of the cell is inside the bin and since there is no
493 // intersection, then that cell completely inside the bin.
494 if((bin->IsInside(xclipl,yclipb)) || (bin->IsInside(xclipl,yclipt))){
495 fCells[i + j*fCellX].Add((TObject*) bin);
496 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
498 continue;
499 }
500 if((bin->IsInside(xclipr,yclipb)) || (bin->IsInside(xclipr,yclipt))){
501 fCells[i + j*fCellX].Add((TObject*) bin);
502 fIsEmpty[i + j*fCellX] = kFALSE; // Makes the cell non-empty
504 continue;
505 }
506 }
507 }
508}
509
510////////////////////////////////////////////////////////////////////////////////
511/// Changes the number of partition cells in the histogram.
512/// Deletes the old partition and constructs a new one.
513
515{
516 fCellX = n; // Set the number of cells
517 fCellY = m; // Set the number of cells
518
519 delete [] fCells; // Deletes the old partition
520
521 // number of cells in the grid
522 //N.B. not to be confused with fNcells (the number of bins) !
524 fCells = new TList [fNCells]; // Sets an empty partition
525
528
529 delete [] fIsEmpty;
530 delete [] fCompletelyInside;
531 fIsEmpty = new Bool_t [fNCells];
533
534 // Initializes the flags
535 for (int i = 0; i<fNCells; i++) {
536 fIsEmpty[i] = kTRUE;
538 }
539
540 // TList iterator
541 TIter next(fBins);
542 TObject *obj;
543
544 while((obj = next())){ // Loop over bins and add them to the partition
546 }
547}
548
549////////////////////////////////////////////////////////////////////////////////
550/// Clears the contents of all bins in the histogram.
551
553{
554 TIter next(fBins);
555 TObject *obj;
556 TH2PolyBin *bin;
557
558 // Clears the bin contents
559 while ((obj = next())) {
560 bin = (TH2PolyBin*) obj;
561 bin->ClearContent();
562 }
563
564 // Clears the statistics
565 fTsumw = 0;
566 fTsumw2 = 0;
567 fTsumwx = 0;
568 fTsumwx2 = 0;
569 fTsumwy = 0;
570 fTsumwy2 = 0;
571 fEntries = 0;
572}
573
574////////////////////////////////////////////////////////////////////////////////
575/// Reset this histogram: contents, errors, etc.
576
578{
579 TIter next(fBins);
580 TObject *obj;
581 TH2PolyBin *bin;
582
583 // Clears the bin contents
584 while ((obj = next())) {
585 bin = (TH2PolyBin*) obj;
586 bin->ClearContent();
587 }
588
589 TH2::Reset(opt);
590}
591
592////////////////////////////////////////////////////////////////////////////////
593/// Returns the bin number of the bin at the given coordinate. -1 to -9 are
594/// the overflow and underflow bins. overflow bin -5 is the unbinned areas in
595/// the histogram (also called "the sea"). The third parameter can be left
596/// blank.
597/// The overflow/underflow bins are:
598///~~~ {.cpp}
599/// -1 | -2 | -3
600/// -------------
601/// -4 | -5 | -6
602/// -------------
603/// -7 | -8 | -9
604///~~~
605/// where -5 means is the "sea" bin (i.e. unbinned areas)
606
608{
609
610 // Checks for overflow/underflow
611 Int_t overflow = 0;
612 if (y > fYaxis.GetXmax()) overflow += -1;
613 else if (y > fYaxis.GetXmin()) overflow += -4;
614 else overflow += -7;
615 if (x > fXaxis.GetXmax()) overflow += -2;
616 else if (x > fXaxis.GetXmin()) overflow += -1;
617 if (overflow != -5) return overflow;
618
619 // Finds the cell (x,y) coordinates belong to
620 Int_t n = (Int_t)(floor((x-fXaxis.GetXmin())/fStepX));
621 Int_t m = (Int_t)(floor((y-fYaxis.GetXmin())/fStepY));
622
623 // Make sure the array indices are correct.
624 if (n>=fCellX) n = fCellX-1;
625 if (m>=fCellY) m = fCellY-1;
626 if (n<0) n = 0;
627 if (m<0) m = 0;
628
629 if (fIsEmpty[n+fCellX*m]) return -5;
630
631 TH2PolyBin *bin;
632
633 TIter next(&fCells[n+fCellX*m]);
634 TObject *obj;
635
636 // Search for the bin in the cell
637 while ((obj=next())) {
638 bin = (TH2PolyBin*)obj;
639 if (bin->IsInside(x,y)) return bin->GetBinNumber();
640 }
641
642 // If the search has not returned a bin, the point must be on "the sea"
643 return -5;
644}
645
646////////////////////////////////////////////////////////////////////////////////
647/// Increment the bin containing (x,y) by 1.
648/// Uses the partitioning algorithm.
649
651{
652 return Fill(x, y, 1.0);
653}
654
655////////////////////////////////////////////////////////////////////////////////
656/// Increment the bin containing (x,y) by w.
657/// Uses the partitioning algorithm.
658
660{
661 // see GetBinCOntent for definition of overflow bins
662 // in case of weighted events store weight square in fSumw2.fArray
663 // but with this indexing:
664 // fSumw2.fArray[0:kNOverflow-1] : sum of weight squares for the overflow bins
665 // fSumw2.fArray[kNOverflow:fNcells] : sum of weight squares for the standard bins
666 // where fNcells = kNOverflow + Number of bins. kNOverflow=9
667
668 if (fNcells <= kNOverflow) return 0;
669
670 // create sum of weight square array if weights are different than 1
671 if (!fSumw2.fN && w != 1.0 && !TestBit(TH1::kIsNotW) ) Sumw2();
672
673 Int_t overflow = 0;
674 if (y > fYaxis.GetXmax()) overflow += -1;
675 else if (y > fYaxis.GetXmin()) overflow += -4;
676 else overflow += -7;
677 if (x > fXaxis.GetXmax()) overflow += -2;
678 else if(x > fXaxis.GetXmin()) overflow += -1;
679 if (overflow != -5) {
680 fOverflow[-overflow - 1]+= w;
681 if (fSumw2.fN) fSumw2.fArray[-overflow - 1] += w*w;
682 return overflow;
683 }
684
685 // Finds the cell (x,y) coordinates belong to
686 Int_t n = (Int_t)(floor((x-fXaxis.GetXmin())/fStepX));
687 Int_t m = (Int_t)(floor((y-fYaxis.GetXmin())/fStepY));
688
689 // Make sure the array indices are correct.
690 if (n>=fCellX) n = fCellX-1;
691 if (m>=fCellY) m = fCellY-1;
692 if (n<0) n = 0;
693 if (m<0) m = 0;
694
695 if (fIsEmpty[n+fCellX*m]) {
696 fOverflow[4]+= w;
697 if (fSumw2.fN) fSumw2.fArray[4] += w*w;
698 return -5;
699 }
700
701 TH2PolyBin *bin;
702 Int_t bi;
703
704 TIter next(&fCells[n+fCellX*m]);
705 TObject *obj;
706
707 while ((obj=next())) {
708 bin = (TH2PolyBin*)obj;
709 // needs to account offset in array for overflow bins
710 bi = bin->GetBinNumber()-1+kNOverflow;
711 if (bin->IsInside(x,y)) {
712 bin->Fill(w);
713
714 // Statistics
715 fTsumw = fTsumw + w;
716 fTsumw2 = fTsumw2 + w*w;
717 fTsumwx = fTsumwx + w*x;
718 fTsumwx2 = fTsumwx2 + w*x*x;
719 fTsumwy = fTsumwy + w*y;
720 fTsumwy2 = fTsumwy2 + w*y*y;
721 if (fSumw2.fN) {
722 assert(bi < fSumw2.fN);
723 fSumw2.fArray[bi] += w*w;
724 }
725 fEntries++;
726
728
729 return bin->GetBinNumber();
730 }
731 }
732
733 fOverflow[4]+= w;
734 if (fSumw2.fN) fSumw2.fArray[4] += w*w;
735 return -5;
736}
737
738////////////////////////////////////////////////////////////////////////////////
739/// Increment the bin named "name" by w.
740
742{
744
745 TIter next(fBins);
746 TObject *obj;
747 TH2PolyBin *bin;
748
749 while ((obj = next())) {
750 bin = (TH2PolyBin*) obj;
751 if (!sname.CompareTo(bin->GetPolygon()->GetName())) {
752 bin->Fill(w);
753 fEntries++;
755 return bin->GetBinNumber();
756 }
757 }
758
759 return 0;
760}
761
762////////////////////////////////////////////////////////////////////////////////
763/// Fills a 2-D histogram with an array of values and weights.
764///
765/// \param [in] ntimes: number of entries in arrays x and w
766/// (array size must be ntimes*stride)
767/// \param [in] x: array of x values to be histogrammed
768/// \param [in] y: array of y values to be histogrammed
769/// \param [in] w: array of weights
770/// \param [in] stride: step size through arrays x, y and w
771
773 const Double_t* w, Int_t stride)
774{
775 for (int i = 0; i < ntimes; i += stride) {
776 Fill(x[i], y[i], w[i]);
777 }
778}
779
780////////////////////////////////////////////////////////////////////////////////
781/// Returns the integral of bin contents.
782/// By default the integral is computed as the sum of bin contents.
783/// If option "width" or "area" is specified, the integral is the sum of
784/// the bin contents multiplied by the area of the bin.
785
787{
788 TString opt = option;
789 opt.ToLower();
790
791 Double_t w;
792 Double_t integral = 0.;
793
794 TIter next(fBins);
795 TObject *obj;
796 TH2PolyBin *bin;
797 if ((opt.Contains("width")) || (opt.Contains("area"))) {
798 while ((obj = next())) {
799 bin = (TH2PolyBin *)obj;
800 w = bin->GetArea();
801 integral += w * (bin->GetContent());
802 }
803 } else {
804 // need to recompute integral in case SetBinContent was called.
805 // fTsumw cannot be used since it is not updated in that case
806 while ((obj = next())) {
807 bin = (TH2PolyBin *)obj;
808 integral += (bin->GetContent());
809 }
810 }
811 return integral;
812}
813
814////////////////////////////////////////////////////////////////////////////////
815/// Returns the content of the input bin
816/// Bin numbers are from [1,nbins] and
817/// for the overflow/underflow/sea bins the range is [-9,-1]:
818///~~~ {.cpp}
819/// -1 | -2 | -3
820/// ---+----+----
821/// -4 | -5 | -6
822/// ---+----+----
823/// -7 | -8 | -9
824///~~~
825/// where -5 is the "sea" bin (i.e. unbinned areas)
826
828{
829 if (bin > GetNumberOfBins() || bin == 0 || bin < -kNOverflow) return 0;
830 if (bin<0) return fOverflow[-bin - 1];
831 return ((TH2PolyBin*) fBins->At(bin-1))->GetContent();
832}
833
834////////////////////////////////////////////////////////////////////////////////
835/// Returns the value of error associated to bin number bin.
836/// If the sum of squares of weights has been defined (via Sumw2),
837/// this function returns the sqrt(sum of w2).
838/// otherwise it returns the sqrt(contents) for this bin.
839/// Bins are in range [1:nbins] and for bin < 0 in range [-9:-1] it returns errors for overflow bins.
840/// See also TH2Poly::GetBinContent
841
843{
844 if (bin == 0 || bin > GetNumberOfBins() || bin < - kNOverflow) return 0;
845 if (fBuffer) ((TH1*)this)->BufferEmpty();
846 // in case of weighted events the sum of the weights are stored in a different way than
847 // a normal histogram
848 // fSumw2.fArray[0:kNOverflow-1] : sum of weight squares for the overflow bins (
849 // fSumw2.fArray[kNOverflow:fNcells] : sum of weight squares for the standard bins
850 // fNcells = kNOverflow (9) + Number of bins
851 if (fSumw2.fN) {
852 Int_t binIndex = (bin > 0) ? bin+kNOverflow-1 : -(bin+1);
854 return TMath::Sqrt(err2);
855 }
857 return TMath::Sqrt(error2);
858}
859
860////////////////////////////////////////////////////////////////////////////////
861/// Return the number of bins :
862/// it should be the size of the bin list
864 Int_t nbins = fNcells-kNOverflow;
865 if (nbins != (fBins ? fBins->GetSize() : 0))
866 Fatal("GetNumberOfBins","Object has an invalid number of bins");
867 return nbins;
868}
869
870////////////////////////////////////////////////////////////////////////////////
871/// Set the bin Error.
872/// Re-implementation for TH2Poly given the different bin indexing in the
873/// stored squared error array.
874/// See also notes in TH1::SetBinError
875///
876/// Bins are in range [1:nbins] and for bin < 0 in the range [-9:-1] the errors is set for the overflow bins
877
878
880{
881 if (bin == 0 || bin > GetNumberOfBins() || bin < - kNOverflow) return;
882 if (!fSumw2.fN) Sumw2();
884 // see comment in GetBinError for special convention of bin index in fSumw2 array
885 Int_t binIndex = (bin > 0) ? bin+kNOverflow-1 : -(bin+1);
886 fSumw2.fArray[binIndex] = error * error;
887}
888
889
890
891////////////////////////////////////////////////////////////////////////////////
892/// Returns the bin name.
893
894const char *TH2Poly::GetBinName(Int_t bin) const
895{
896 if (bin > GetNumberOfBins()) return "";
897 if (bin <= 0) return "";
898 return ((TH2PolyBin*) fBins->At(bin-1))->GetPolygon()->GetName();
899}
900
901////////////////////////////////////////////////////////////////////////////////
902/// Returns the bin title.
903
904const char *TH2Poly::GetBinTitle(Int_t bin) const
905{
906 if (bin > GetNumberOfBins()) return "";
907 if (bin <= 0) return "";
908 return ((TH2PolyBin*) fBins->At(bin-1))->GetPolygon()->GetTitle();
909}
910
911////////////////////////////////////////////////////////////////////////////////
912/// Returns the maximum value of the histogram.
913
915{
916 if (fNcells <= kNOverflow) return 0;
917 if (fMaximum != -1111) return fMaximum;
918
919 TH2PolyBin *b;
920
921 TIter next(fBins);
922 TObject *obj;
923 Double_t max,c;
924
925 max = ((TH2PolyBin*) next())->GetContent();
926
927 while ((obj=next())) {
928 b = (TH2PolyBin*)obj;
929 c = b->GetContent();
930 if (c>max) max = c;
931 }
932 return max;
933}
934
935////////////////////////////////////////////////////////////////////////////////
936/// Returns the maximum value of the histogram that is less than maxval.
937
939{
940 if (fNcells <= kNOverflow) return 0;
941 if (fMaximum != -1111) return fMaximum;
942
943 TH2PolyBin *b;
944
945 TIter next(fBins);
946 TObject *obj;
947 Double_t max,c;
948
949 max = ((TH2PolyBin*) next())->GetContent();
950
951 while ((obj=next())) {
952 b = (TH2PolyBin*)obj;
953 c = b->GetContent();
954 if (c>max && c<maxval) max=c;
955 }
956 return max;
957}
958
959////////////////////////////////////////////////////////////////////////////////
960/// Returns the minimum value of the histogram.
961
963{
964 if (fNcells <= kNOverflow) return 0;
965 if (fMinimum != -1111) return fMinimum;
966
967 TH2PolyBin *b;
968
969 TIter next(fBins);
970 TObject *obj;
971 Double_t min,c;
972
973 min = ((TH2PolyBin*) next())->GetContent();
974
975 while ((obj=next())) {
976 b = (TH2PolyBin*)obj;
977 c = b->GetContent();
978 if (c<min) min=c;
979 }
980 return min;
981}
982
983////////////////////////////////////////////////////////////////////////////////
984/// Returns the minimum value of the histogram that is greater than minval.
985
987{
988 if (fNcells <= kNOverflow) return 0;
989 if (fMinimum != -1111) return fMinimum;
990
991 TH2PolyBin *b;
992
993 TIter next(fBins);
994 TObject *obj;
995 Double_t min,c;
996
997 min = ((TH2PolyBin*) next())->GetContent();
998
999 while ((obj=next())) {
1000 b = (TH2PolyBin*)obj;
1001 c = b->GetContent();
1002 if (c<min && c>minval) min=c;
1003 }
1004 return min;
1005}
1006
1007////////////////////////////////////////////////////////////////////////////////
1008/// Bins the histogram using a honeycomb structure
1009/// If the option "v" is specified, the hexagons are drawn "vertically" (default).
1010/// If the option "h" is selected they are drawn "horizontally".
1011
1013{
1014 TString opt = option;
1015 opt.ToLower();
1017 Double_t x[6], y[6];
1019
1020 // Add the bins
1021 if (opt.Contains("v")) {
1022 yloop = ystart + a / 2.0;
1023 for (int sCounter = 0; sCounter < s; sCounter++) {
1024 xtemp = xloop;
1025 // Determine the number of hexagons in that row
1026 if (sCounter%2 == 0)
1028 else
1030
1032 // Go around the hexagon
1033 x[0] = xtemp;
1034 y[0] = yloop;
1035 x[1] = x[0];
1036 y[1] = y[0] + a;
1037 x[2] = x[1] + a * TMath::Sqrt(3) / 2.0;
1038 y[2] = y[1] + a / 2.0;
1039 x[3] = x[2] + a * TMath::Sqrt(3) / 2.0;
1040 y[3] = y[1];
1041 x[4] = x[3];
1042 y[4] = y[0];
1043 x[5] = x[2];
1044 y[5] = y[4] - a/2.0;
1045 this->AddBin(6, x, y);
1046 // Go right
1047 xtemp += a * TMath::Sqrt(3);
1048 }
1049 // Increment the starting position
1050 if (sCounter%2 == 0)
1051 xloop += a * TMath::Sqrt(3) / 2.0;
1052 else
1053 xloop -= a * TMath::Sqrt(3) / 2.0;
1054 yloop += 1.5 * a;
1055 }
1056 } else if (opt.Contains("h")) {
1057 yloop = ystart + a*TMath::Sqrt(3)/2.0;
1058 for (int sCounter = 0; sCounter < s; sCounter++) {
1059 ytemp = yloop;
1060 // Determine the number of hexagons in that row
1061 if (sCounter%2 == 0)
1063 else
1066 // Go around the hexagon
1067 x[0] = xloop;
1068 y[0] = ytemp;
1069 x[1] = x[0] + a/2.0;
1070 y[1] = y[0] + a * TMath::Sqrt(3) / 2.0;
1071 x[2] = x[1] + a;
1072 y[2] = y[1];
1073 x[3] = x[2] + a/2.0;
1074 y[3] = y[0];
1075 x[4] = x[2];
1076 y[4] = y[3] - a * TMath::Sqrt(3) / 2.0;
1077 x[5] = x[1];
1078 y[5] = y[4];
1079 this->AddBin(6, x, y);
1080 // Go up
1081 ytemp += a * TMath::Sqrt(3);
1082 }
1083 // Increment the starting position
1084 if (sCounter%2 == 0)
1085 yloop += a * TMath::Sqrt(3) / 2.0;
1086 else
1087 yloop -= a * TMath::Sqrt(3) / 2.0;
1088 xloop += 1.5 * a;
1089 }
1090 } else {
1091 Error("Honeycomb", "Unknown option");
1092 }
1093
1094}
1095
1096////////////////////////////////////////////////////////////////////////////////
1097/// Initializes the TH2Poly object. This method is called by the constructor.
1098
1100 Double_t ylow, Double_t yup, Int_t n, Int_t m)
1101{
1102 Int_t i;
1103 fDimension = 2; //The dimension of the histogram
1104
1105 fBins = nullptr;
1107
1108 // Sets the boundaries of the histogram
1109 fXaxis.Set(100, xlow, xup);
1110 fYaxis.Set(100, ylow, yup);
1111
1112 for (i=0; i<9; i++) fOverflow[i] = 0.;
1113
1114 // Statistics
1115 fEntries = 0; // The total number of entries
1116 fTsumw = 0.; // Total amount of content in the histogram
1117 fTsumw2 = 0.; // Sum square of the weights
1118 fTsumwx = 0.; // Weighted sum of x coordinates
1119 fTsumwx2 = 0.; // Weighted sum of the squares of x coordinates
1120 fTsumwy2 = 0.; // Weighted sum of the squares of y coordinates
1121 fTsumwy = 0.; // Weighted sum of y coordinates
1122
1123 fCellX = n; // Set the number of cells to default
1124 fCellY = m; // Set the number of cells to default
1125
1126 // number of cells in the grid
1127 //N.B. not to be confused with fNcells (the number of bins) !
1129 fCells = new TList [fNCells]; // Sets an empty partition
1130 fStepX = (fXaxis.GetXmax() - fXaxis.GetXmin())/fCellX; // Cell width
1131 fStepY = (fYaxis.GetXmax() - fYaxis.GetXmin())/fCellY; // Cell height
1132
1133 fIsEmpty = new Bool_t [fNCells]; // Empty partition
1134 fCompletelyInside = new Bool_t [fNCells]; // Cell is completely inside bin
1135
1136 for (i = 0; i<fNCells; i++) { // Initializes the flags
1137 fIsEmpty[i] = kTRUE;
1139 }
1140
1141 // 3D Painter flags
1144}
1145
1146////////////////////////////////////////////////////////////////////////////////
1147/// Returns kTRUE if the input bin is intersecting with the
1148/// input rectangle (xclipl, xclipr, yclipb, yclipt)
1149
1153{
1154 Int_t gn;
1155 Double_t *gx;
1156 Double_t *gy;
1158 TObject *poly = bin->GetPolygon();
1159
1160 if (poly->IsA() == TGraph::Class()) {
1161 TGraph *g = (TGraph*)poly;
1162 gx = g->GetX();
1163 gy = g->GetY();
1164 gn = g->GetN();
1166 }
1167
1168 if (poly->IsA() == TMultiGraph::Class()) {
1170 TList *gl = mg->GetListOfGraphs();
1171 if (!gl) return inter;
1172 TGraph *g;
1173 TIter next(gl);
1174 while ((g = (TGraph*) next())) {
1175 gx = g->GetX();
1176 gy = g->GetY();
1177 gn = g->GetN();
1179 yclipb, yclipt);
1180 if (inter) return inter;
1181 }
1182 }
1183
1184 return inter;
1185}
1186
1187////////////////////////////////////////////////////////////////////////////////
1188/// Returns kTRUE if the input polygon (bn, x, y) is intersecting with the
1189/// input rectangle (xclipl, xclipr, yclipb, yclipt)
1190
1194{
1195 Bool_t p0R, p0L, p0T, p0B, p0xM, p0yM, p1R, p1L, p1T;
1197
1198 for (int counter = 0; counter < (bn-1); counter++) {
1199 // If both are on the same side, return kFALSE
1200 p0L = x[counter] <= xclipl; // Point 0 is on the left
1201 p1L = x[counter + 1] <= xclipl; // Point 1 is on the left
1202 if (p0L && p1L) continue;
1203 p0R = x[counter] >= xclipr; // Point 0 is on the right
1204 p1R = x[counter + 1] >= xclipr; // Point 1 is on the right
1205 if (p0R && p1R) continue;
1206 p0T = y[counter] >= yclipt; // Point 0 is at the top
1207 p1T = y[counter + 1] >= yclipt; // Point 1 is at the top
1208 if (p0T && p1T) continue;
1209 p0B = y[counter] <= yclipb; // Point 0 is at the bottom
1210 p1B = y[counter + 1] <= yclipb; // Point 1 is at the bottom
1211 if (p0B && p1B) continue;
1212
1213 // Checks to see if any are inside
1214 p0xM = !p0R && !p0L; // Point 0 is inside along x
1215 p0yM = !p0T && !p0B; // Point 1 is inside along x
1216 p1xM = !p1R && !p1L; // Point 0 is inside along y
1217 p1yM = !p1T && !p1B; // Point 1 is inside along y
1218 p0In = p0xM && p0yM; // Point 0 is inside
1219 p1In = p1xM && p1yM; // Point 1 is inside
1220 if (p0In) {
1221 if (p1In) continue;
1222 return kTRUE;
1223 } else {
1224 if (p1In) return kTRUE;
1225 }
1226
1227 // We know by now that the points are not in the same side and not inside.
1228
1229 // Checks to see if they are opposite
1230
1231 if (p0xM && p1xM) return kTRUE;
1232 if (p0yM && p1yM) return kTRUE;
1233
1234 // We now know that the points are in different x and y indices
1235
1236 Double_t xcoord[3], ycoord[3];
1237 xcoord[0] = x[counter];
1238 xcoord[1] = x[counter + 1];
1239 ycoord[0] = y[counter];
1240 ycoord[1] = y[counter + 1];
1241
1242 if (p0L) {
1243 if(p1T){
1244 xcoord[2] = xclipl;
1245 ycoord[2] = yclipb;
1247 (TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord))) continue;
1248 else return kTRUE;
1249 } else if (p1B) {
1250 xcoord[2] = xclipl;
1251 ycoord[2] = yclipt;
1253 (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord))) continue;
1254 else return kTRUE;
1255 } else { // p1yM
1256 if (p0T) {
1257 xcoord[2] = xclipl;
1258 ycoord[2] = yclipb;
1259 if (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord)) continue;
1260 else return kTRUE;
1261 }
1262 if (p0B) {
1263 xcoord[2] = xclipl;
1264 ycoord[2] = yclipt;
1265 if (TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord)) continue;
1266 else return kTRUE;
1267 }
1268 }
1269 } else if (p0R) {
1270 if (p1T) {
1271 xcoord[2] = xclipl;
1272 ycoord[2] = yclipb;
1273 if ((TMath::IsInside(xclipr, yclipb, 3, xcoord, ycoord)) ||
1274 (TMath::IsInside(xclipl, yclipt, 3, xcoord, ycoord))) continue;
1275 else return kTRUE;
1276 } else if (p1B) {
1277 xcoord[2] = xclipl;
1278 ycoord[2] = yclipt;
1279 if ((TMath::IsInside(xclipl, yclipb, 3, xcoord, ycoord)) ||
1280 (TMath::IsInside(xclipr, yclipt, 3, xcoord, ycoord))) continue;
1281 else return kTRUE;
1282 } else{ // p1yM
1283 if (p0T) {
1284 xcoord[2] = xclipr;
1285 ycoord[2] = yclipb;
1286 if (TMath::IsInside(xclipl, yclipt, 3, xcoord, ycoord)) continue;
1287 else return kTRUE;
1288 }
1289 if (p0B) {
1290 xcoord[2] = xclipr;
1291 ycoord[2] = yclipt;
1292 if (TMath::IsInside(xclipl, yclipb, 3, xcoord, ycoord)) continue;
1293 else return kTRUE;
1294 }
1295 }
1296 }
1297 }
1298 return kFALSE;
1299}
1300
1301////////////////////////////////////////////////////////////////////////////////
1302/// Merge TH2Polys
1303/// Given the special nature of the TH2Poly, the merge is implemented in
1304/// terms of subsequent TH2Poly::Add calls.
1306{
1307 for (auto h2pAsObj : *coll) {
1308 if (!Add((TH1*)h2pAsObj, 1.)) {
1309 Warning("Merge", "An issue was encountered during the merge operation.");
1310 return 0L;
1311 }
1312 }
1313 return GetEntries();
1314}
1315
1316////////////////////////////////////////////////////////////////////////////////
1317/// Save primitive as a C++ statement(s) on output stream out
1318
1319void TH2Poly::SavePrimitive(std::ostream &out, Option_t *option)
1320{
1321 // histogram pointer has by default the histogram name.
1322 // however, in case histogram has no directory, it is safer to add a
1323 // incremental suffix
1325
1326 out << " \n";
1327
1328 // Construct the class initialization
1329 out << " " << ClassName() << " *" << hname << " = new " << ClassName() << "(\"" << hname << "\", \""
1330 << TString(GetTitle()).ReplaceSpecialCppChars() << "\", " << fCellX << ", " << fXaxis.GetXmin() << ", "
1331 << fXaxis.GetXmax() << ", " << fCellY << ", " << fYaxis.GetXmin() << ", " << fYaxis.GetXmax() << ");\n";
1332
1333 // Save Bins
1334 TIter next(fBins);
1335 while (auto obj = next()) {
1336 auto th2pBin = (TH2PolyBin *)obj;
1337 th2pBin->GetPolygon()->SavePrimitive(out, TString::Format("th2poly%s", hname.Data()));
1338 }
1339
1340 // save bin contents
1341 out << " \n";
1342 for (Int_t bin = 1; bin <= GetNumberOfBins(); bin++) {
1343 Double_t bc = GetBinContent(bin);
1344 if (bc)
1345 out << " " << hname << "->SetBinContent(" << bin << "," << bc << ");\n";
1346 }
1347
1348 // save bin errors
1349 if (fSumw2.fN) {
1350 for (Int_t bin = 1; bin <= GetNumberOfBins(); bin++) {
1351 Double_t be = GetBinError(bin);
1352 if (be)
1353 out << " " << hname << "->SetBinError(" << bin << "," << be << ");\n";
1354 }
1355 }
1356
1358}
1359
1360////////////////////////////////////////////////////////////////////////////////
1361/// Multiply this histogram by a constant c1.
1362
1364{
1365 for( int i = 0; i < this->GetNumberOfBins(); i++ ) {
1366 this->SetBinContent(i+1, c1*this->GetBinContent(i+1));
1367 }
1368 for( int i = 0; i < kNOverflow; i++ ) {
1369 this->SetBinContent(-i-1, c1*this->GetBinContent(-i-1) );
1370 }
1371}
1372
1373////////////////////////////////////////////////////////////////////////////////
1374/// Sets the contents of the input bin to the input content
1375/// Negative values between -1 and -9 are for the overflows and the sea
1376
1378{
1379 if (bin > GetNumberOfBins() || bin == 0 || bin < -9 ) return;
1380 if (bin > 0) {
1381 ((TH2PolyBin*) fBins->At(bin-1))->SetContent(content);
1382 }
1383 else
1384 fOverflow[-bin - 1] = content;
1385
1387 fEntries++;
1388}
1389
1390////////////////////////////////////////////////////////////////////////////////
1391/// When set to kTRUE, allows the histogram to expand if a bin outside the
1392/// limits is added.
1393
1395{
1396 fFloat = flag;
1397}
1398
1399////////////////////////////////////////////////////////////////////////////////
1400/// Return "true" if the point (x,y) is inside the bin of binnr.
1401
1403{
1404 if (!fBins) return false;
1405 TH2PolyBin* bin = (TH2PolyBin*)fBins->At(binnr);
1406 if (!bin) return false;
1407 return bin->IsInside(x,y);
1408}
1409
1411{
1412 stats[0] = fTsumw;
1413 stats[1] = fTsumw2;
1414 stats[2] = fTsumwx;
1415 stats[3] = fTsumwx2;
1416 stats[4] = fTsumwy;
1417 stats[5] = fTsumwy2;
1418 stats[6] = fTsumwxy;
1419}
1420
1421/** \class TH2PolyBin
1422 \ingroup Histograms
1423Helper class to represent a bin in the TH2Poly histogram
1424*/
1425
1426////////////////////////////////////////////////////////////////////////////////
1427/// Default constructor.
1428
1430{
1431 fPoly = nullptr;
1432 fContent = 0.;
1433 fNumber = 0;
1434 fXmax = -1111;
1435 fXmin = -1111;
1436 fYmax = -1111;
1437 fYmin = -1111;
1438 fArea = 0;
1440}
1441
1442////////////////////////////////////////////////////////////////////////////////
1443/// Normal constructor.
1444
1446{
1447 fContent = 0.;
1449 fArea = 0.;
1450 fPoly = poly;
1451 fXmax = -1111;
1452 fXmin = -1111;
1453 fYmax = -1111;
1454 fYmin = -1111;
1456}
1457
1458////////////////////////////////////////////////////////////////////////////////
1459/// Destructor.
1460
1462{
1463 if (fPoly) delete fPoly;
1464}
1465
1466////////////////////////////////////////////////////////////////////////////////
1467/// Returns the area of the bin.
1468
1470{
1471 Int_t bn;
1472
1473 if (fArea == 0) {
1474 if (fPoly->IsA() == TGraph::Class()) {
1475 TGraph *g = (TGraph*)fPoly;
1476 bn = g->GetN();
1477 fArea = g->Integral(0,bn-1);
1478 }
1479
1480 if (fPoly->IsA() == TMultiGraph::Class()) {
1482 TList *gl = mg->GetListOfGraphs();
1483 if (!gl) return fArea;
1484 TGraph *g;
1485 TIter next(gl);
1486 while ((g = (TGraph*) next())) {
1487 bn = g->GetN();
1488 fArea = fArea + g->Integral(0,bn-1);
1489 }
1490 }
1491 }
1492
1493 return fArea;
1494}
1495
1496////////////////////////////////////////////////////////////////////////////////
1497/// Returns the maximum value for the x coordinates of the bin.
1498
1500{
1501 if (fXmax != -1111) return fXmax;
1502
1503 Int_t bn,i;
1504 Double_t *bx;
1505
1506 if (fPoly->IsA() == TGraph::Class()) {
1507 TGraph *g = (TGraph*)fPoly;
1508 bx = g->GetX();
1509 bn = g->GetN();
1510 fXmax = bx[0];
1511 for (i=1; i<bn; i++) {if (fXmax < bx[i]) fXmax = bx[i];}
1512 }
1513
1514 if (fPoly->IsA() == TMultiGraph::Class()) {
1516 TList *gl = mg->GetListOfGraphs();
1517 if (!gl) return fXmax;
1518 TGraph *g;
1519 TIter next(gl);
1520 Bool_t first = kTRUE;
1521 while ((g = (TGraph*) next())) {
1522 bx = g->GetX();
1523 bn = g->GetN();
1524 if (first) {fXmax = bx[0]; first = kFALSE;}
1525 for (i=0; i<bn; i++) {if (fXmax < bx[i]) fXmax = bx[i];}
1526 }
1527 }
1528
1529 return fXmax;
1530}
1531
1532////////////////////////////////////////////////////////////////////////////////
1533/// Returns the minimum value for the x coordinates of the bin.
1534
1536{
1537 if (fXmin != -1111) return fXmin;
1538
1539 Int_t bn,i;
1540 Double_t *bx;
1541
1542 if (fPoly->IsA() == TGraph::Class()) {
1543 TGraph *g = (TGraph*)fPoly;
1544 bx = g->GetX();
1545 bn = g->GetN();
1546 fXmin = bx[0];
1547 for (i=1; i<bn; i++) {if (fXmin > bx[i]) fXmin = bx[i];}
1548 }
1549
1550 if (fPoly->IsA() == TMultiGraph::Class()) {
1552 TList *gl = mg->GetListOfGraphs();
1553 if (!gl) return fXmin;
1554 TGraph *g;
1555 TIter next(gl);
1556 Bool_t first = kTRUE;
1557 while ((g = (TGraph*) next())) {
1558 bx = g->GetX();
1559 bn = g->GetN();
1560 if (first) {fXmin = bx[0]; first = kFALSE;}
1561 for (i=0; i<bn; i++) {if (fXmin > bx[i]) fXmin = bx[i];}
1562 }
1563 }
1564
1565 return fXmin;
1566}
1567
1568////////////////////////////////////////////////////////////////////////////////
1569/// Returns the maximum value for the y coordinates of the bin.
1570
1572{
1573 if (fYmax != -1111) return fYmax;
1574
1575 Int_t bn,i;
1576 Double_t *by;
1577
1578 if (fPoly->IsA() == TGraph::Class()) {
1579 TGraph *g = (TGraph*)fPoly;
1580 by = g->GetY();
1581 bn = g->GetN();
1582 fYmax = by[0];
1583 for (i=1; i<bn; i++) {if (fYmax < by[i]) fYmax = by[i];}
1584 }
1585
1586 if (fPoly->IsA() == TMultiGraph::Class()) {
1588 TList *gl = mg->GetListOfGraphs();
1589 if (!gl) return fYmax;
1590 TGraph *g;
1591 TIter next(gl);
1592 Bool_t first = kTRUE;
1593 while ((g = (TGraph*) next())) {
1594 by = g->GetY();
1595 bn = g->GetN();
1596 if (first) {fYmax = by[0]; first = kFALSE;}
1597 for (i=0; i<bn; i++) {if (fYmax < by[i]) fYmax = by[i];}
1598 }
1599 }
1600
1601 return fYmax;
1602}
1603
1604////////////////////////////////////////////////////////////////////////////////
1605/// Returns the minimum value for the y coordinates of the bin.
1606
1608{
1609 if (fYmin != -1111) return fYmin;
1610
1611 Int_t bn,i;
1612 Double_t *by;
1613
1614 if (fPoly->IsA() == TGraph::Class()) {
1615 TGraph *g = (TGraph*)fPoly;
1616 by = g->GetY();
1617 bn = g->GetN();
1618 fYmin = by[0];
1619 for (i=1; i<bn; i++) {if (fYmin > by[i]) fYmin = by[i];}
1620 }
1621
1622 if (fPoly->IsA() == TMultiGraph::Class()) {
1624 TList *gl = mg->GetListOfGraphs();
1625 if (!gl) return fYmin;
1626 TGraph *g;
1627 TIter next(gl);
1628 Bool_t first = kTRUE;
1629 while ((g = (TGraph*) next())) {
1630 by = g->GetY();
1631 bn = g->GetN();
1632 if (first) {fYmin = by[0]; first = kFALSE;}
1633 for (i=0; i<bn; i++) {if (fYmin > by[i]) fYmin = by[i];}
1634 }
1635 }
1636
1637 return fYmin;
1638}
1639
1640////////////////////////////////////////////////////////////////////////////////
1641/// Return "true" if the point (x,y) is inside the bin.
1642
1644{
1645 Int_t in=0;
1646
1647 if (fPoly->IsA() == TGraph::Class()) {
1648 TGraph *g = (TGraph*)fPoly;
1649 in = g->IsInside(x, y);
1650 }
1651
1652 if (fPoly->IsA() == TMultiGraph::Class()) {
1654 in = mg->IsInside(x, y);
1655 }
1656
1657 return in;
1658}
1659
1660////////////////////////////////////////////////////////////////////////////////
1661// RE-implement dummy functions to avoid users calling the
1662// corresponding implementations in TH1 or TH2
1663////////////////////////////////////////////////////////////////////////////////
1664
1665////////////////////////////////////////////////////////////////////////////////
1666/// NOT IMPLEMENTED for TH2Poly
1668{
1669 Error("Add", "Not implemented for TH2Poly");
1670 return kFALSE;
1671}
1672
1673////////////////////////////////////////////////////////////////////////////////
1674/// NOT IMPLEMENTED for TH2Poly
1676{
1677 Error("Add", "Not implemented for TH2Poly");
1678 return kFALSE;
1679}
1680
1681////////////////////////////////////////////////////////////////////////////////
1682/// NOT IMPLEMENTED for TH2Poly
1684{
1685 Error("Divide", "Not implemented for TH2Poly");
1686 return kFALSE;
1687}
1688
1689////////////////////////////////////////////////////////////////////////////////
1690/// NOT IMPLEMENTED for TH2Poly
1692{
1693 Error("Multiply", "Not implemented for TH2Poly");
1694 return kFALSE;
1695}
1696////////////////////////////////////////////////////////////////////////////////
1697/// NOT IMPLEMENTED for TH2Poly
1699{
1700 Error("ComputeIntegral", "Not implemented for TH2Poly");
1701 return TMath::QuietNaN();
1702}
1703////////////////////////////////////////////////////////////////////////////////
1704/// NOT IMPLEMENTED for TH2Poly
1706{
1707 Error("FFT", "Not implemented for TH2Poly");
1708 return nullptr;
1709}
1710////////////////////////////////////////////////////////////////////////////////
1711/// NOT IMPLEMENTED for TH2Poly
1713{
1714 Error("GetAsymmetry", "Not implemented for TH2Poly");
1715 return nullptr;
1716}
1717////////////////////////////////////////////////////////////////////////////////
1718/// NOT IMPLEMENTED for TH2Poly
1720{
1721 Error("Interpolate", "Not implemented for TH2Poly");
1722 return TMath::QuietNaN();
1723}
1724////////////////////////////////////////////////////////////////////////////////
1725/// NOT IMPLEMENTED for TH2Poly
1727{
1728 Error("AddBinContent", "Not implemented for TH2Poly");
1729}
1730////////////////////////////////////////////////////////////////////////////////
1731/// NOT IMPLEMENTED for TH2Poly
1733{
1734 Error("AddBinContent", "Not implemented for TH2Poly");
1735}
#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
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
@ kCounter
Definition TDataType.h:34
Option_t Option_t option
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char y2
Option_t Option_t TPoint TPoint const char y1
char name[80]
Definition TGX11.cxx:110
Double_t * fArray
Definition TArrayD.h:30
void Set(Int_t n) override
Set size of this array to n doubles.
Definition TArrayD.cxx:105
Int_t fN
Definition TArray.h:38
Double_t GetXmax() const
Definition TAxis.h:142
virtual void Set(Int_t nbins, Double_t xmin, Double_t xmax)
Initialize axis with fix bins.
Definition TAxis.cxx:783
Double_t GetXmin() const
Definition TAxis.h:141
Collection abstract base class.
Definition TCollection.h:65
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
TObject * Clone(const char *newname="") const override
Make a clone of an collection using the Streamer facility.
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
1-Dim function class
Definition TF1.h:182
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
static TClass * Class()
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:109
Double_t * fBuffer
[fBufferSize] entry buffer
Definition TH1.h:169
Int_t fNcells
Number of bins(1D), cells (2D) +U/Overflows.
Definition TH1.h:150
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:7856
void SetTitle(const char *title) override
Change/set the title.
Definition TH1.cxx:6753
Double_t fTsumw
Total Sum of weights.
Definition TH1.h:157
Double_t fTsumw2
Total Sum of squares of weights.
Definition TH1.h:158
Double_t fTsumwx2
Total Sum of weight*X*X.
Definition TH1.h:160
@ kIsNotW
Histogram is forced to be not weighted even when the histogram is filled with weighted.
Definition TH1.h:411
Double_t fMaximum
Maximum value for plotting.
Definition TH1.h:161
TString ProvideSaveName(Option_t *option, Bool_t testfdir=kFALSE)
Provide variable name for histogram for saving as primitive Histogram pointer has by default the hist...
Definition TH1.cxx:7277
Int_t fDimension
! Histogram dimension (1, 2 or 3 dim)
Definition TH1.h:171
virtual void SetContent(const Double_t *content)
Replace bin contents by the contents of array content.
Definition TH1.cxx:8426
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:7404
Double_t fMinimum
Minimum value for plotting.
Definition TH1.h:162
virtual Double_t GetEntries() const
Return the current number of entries.
Definition TH1.cxx:4400
void SetName(const char *name) override
Change the name of this histogram.
Definition TH1.cxx:8988
virtual void ResetStats()
Reset the statistics including the number of entries and replace with values calculated from bin cont...
Definition TH1.cxx:7925
virtual void SetBinErrorOption(EBinErrorOpt type)
Definition TH1.h:630
@ kNstat
Size of statistics data (up to TProfile3D)
Definition TH1.h:423
Double_t fEntries
Number of entries.
Definition TH1.h:156
@ kNormal
Errors with Normal (Wald) approximation: errorUp=errorLow= sqrt(N)
Definition TH1.h:115
TAxis fXaxis
X axis descriptor.
Definition TH1.h:151
TArrayD fSumw2
Array of sum of squares of weights.
Definition TH1.h:165
TAxis fYaxis
Y axis descriptor.
Definition TH1.h:152
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition TH1.cxx:9048
virtual void SetEntries(Double_t n)
Definition TH1.h:640
Double_t fTsumwx
Total Sum of weight*X.
Definition TH1.h:159
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:1535
~TH2PolyBin() override
Destructor.
Definition TH2Poly.cxx:1461
Double_t GetYMax()
Returns the maximum value for the y coordinates of the bin.
Definition TH2Poly.cxx:1571
Double_t GetArea()
Returns the area of the bin.
Definition TH2Poly.cxx:1469
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:1607
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:1643
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:1429
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:1499
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
TObject * fPoly
Object holding the polygon definition.
Definition TH2Poly.h:50
2D Histogram with Polygonal Bins
Definition TH2Poly.h:66
Bool_t Multiply(TF1 *, Double_t) override
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1691
void UpdateBinContent(Int_t bin, Double_t content) override
Raw update of bin content on internal data structure see convention for numbering bins in TH1::GetBin...
Definition TH2Poly.h:182
TH2Poly & operator=(const TH2Poly &rhs)
Assignment operator.
Definition TH2Poly.cxx:202
void ClearBinContents()
Clears the contents of all bins in the histogram.
Definition TH2Poly.cxx:552
Double_t fOverflow[kNOverflow]
Overflow bins.
Definition TH2Poly.h:161
@ kNOverflow
Definition TH2Poly.h:159
Bool_t fFloat
When set to kTRUE, allows the histogram to expand if a bin outside the limits is added.
Definition TH2Poly.h:169
Double_t ComputeIntegral(Bool_t) override
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1698
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:1150
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:1394
Double_t Integral(Option_t *option="") const override
Returns the integral of bin contents.
Definition TH2Poly.cxx:786
virtual TH2PolyBin * CreateBin(TObject *poly)
Create appropriate histogram bin.
Definition TH2Poly.cxx:274
Double_t GetBinContent(Int_t bin) const override
Returns the content of the input bin Bin numbers are from [1,nbins] and for the overflow/underflow/se...
Definition TH2Poly.cxx:827
Bool_t * fIsEmpty
[fNCells] The array that returns true if the cell at the given coordinate is empty
Definition TH2Poly.h:167
Double_t RetrieveBinContent(Int_t bin) const override
Raw retrieval of bin content on internal data structure see convention for numbering bins in TH1::Get...
Definition TH2Poly.h:179
TList * fBins
List of bins. The list owns the contained objects.
Definition TH2Poly.h:172
void AddBinToPartition(TH2PolyBin *bin)
Adds the input bin into the partition cell matrix.
Definition TH2Poly.cxx:442
Bool_t fBinContentChanged
!For the 3D Painter
Definition TH2Poly.h:171
Int_t Fill(Double_t x, Double_t y) override
Increment the bin containing (x,y) by 1.
Definition TH2Poly.cxx:650
Int_t GetNumberOfBins() const
Return the number of bins : it should be the size of the bin list.
Definition TH2Poly.cxx:863
void SetBinContentChanged(Bool_t flag)
Definition TH2Poly.h:121
void GetStats(Double_t *stats) const override
Fill the array stats from the contents of this histogram The array stats must be correctly dimensione...
Definition TH2Poly.cxx:1410
Bool_t * fCompletelyInside
[fNCells] The array that returns true if the cell at the given coordinate is completely inside a bin
Definition TH2Poly.h:168
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:1402
Bool_t fNewBinAdded
!For the 3D Painter
Definition TH2Poly.h:170
void Copy(TObject &newth2p) const override
Copy function for TH2Poly.
Definition TH2Poly.cxx:211
void AddBinContent(Int_t) override
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1726
void SetBinError(Int_t bin, Double_t error) override
Set the bin Error.
Definition TH2Poly.cxx:879
const char * GetBinTitle(Int_t bin) const
Returns the bin title.
Definition TH2Poly.cxx:904
void SetNewBinAdded(Bool_t flag)
Definition TH2Poly.h:123
void Scale(Double_t c1=1, Option_t *option="") override
Multiply this histogram by a constant c1.
Definition TH2Poly.cxx:1363
void ChangePartition(Int_t n, Int_t m)
Changes the number of partition cells in the histogram.
Definition TH2Poly.cxx:514
virtual Double_t Interpolate(Double_t, Double_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1719
Double_t GetMinimum() const
Returns the minimum value of the histogram.
Definition TH2Poly.cxx:962
const char * GetBinName(Int_t bin) const
Returns the bin name.
Definition TH2Poly.cxx:894
Double_t fStepX
Definition TH2Poly.h:166
void FillN(Int_t ntimes, const Double_t *x, const Double_t *y, const Double_t *w, Int_t stride=1) override
Fills a 2-D histogram with an array of values and weights.
Definition TH2Poly.cxx:772
Int_t fNCells
Number of partition cells: fCellX*fCellY.
Definition TH2Poly.h:164
Int_t FindBin(Double_t x, Double_t y, Double_t z=0) override
Returns the bin number of the bin at the given coordinate.
Definition TH2Poly.cxx:607
Bool_t Divide(TF1 *, Double_t) override
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1683
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:1099
Int_t fCellX
Number of partition cells in the x-direction of the histogram.
Definition TH2Poly.h:162
~TH2Poly() override
Destructor.
Definition TH2Poly.cxx:191
virtual Int_t AddBin(TObject *poly)
Adds a new bin to the histogram.
Definition TH2Poly.cxx:296
Double_t GetBinError(Int_t bin) const override
Returns the value of error associated to bin number bin.
Definition TH2Poly.cxx:842
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:1191
Bool_t Add(const TH1 *h1, Double_t c1) override
Performs the operation: this = this + c1*h1.
Definition TH2Poly.cxx:364
virtual TH1 * GetAsymmetry(TH1 *, Double_t, Double_t)
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1712
void Honeycomb(Double_t xstart, Double_t ystart, Double_t a, Int_t k, Int_t s, Option_t *option="v")
Bins the histogram using a honeycomb structure If the option "v" is specified, the hexagons are drawn...
Definition TH2Poly.cxx:1012
Long64_t Merge(TCollection *) override
Merge TH2Polys Given the special nature of the TH2Poly, the merge is implemented in terms of subseque...
Definition TH2Poly.cxx:1305
void SetBinContent(Int_t bin, Double_t content) override
Sets the contents of the input bin to the input content Negative values between -1 and -9 are for the...
Definition TH2Poly.cxx:1377
Double_t GetMaximum() const
Returns the maximum value of the histogram.
Definition TH2Poly.cxx:914
void Reset(Option_t *option) override
Reset this histogram: contents, errors, etc.
Definition TH2Poly.cxx:577
Double_t fStepY
Dimensions of a partition cell.
Definition TH2Poly.h:166
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
Definition TH2Poly.cxx:1319
TH1 * FFT(TH1 *, Option_t *) override
NOT IMPLEMENTED for TH2Poly.
Definition TH2Poly.cxx:1705
TList * fCells
[fNCells] The array of TLists that store the bins that intersect with each cell. List do not own the ...
Definition TH2Poly.h:165
Int_t fCellY
Number of partition cells in the y-direction of the histogram.
Definition TH2Poly.h:163
Service class for 2-D histogram classes.
Definition TH2.h:30
void Copy(TObject &hnew) const override
Copy.
Definition TH2.cxx:343
Double_t fTsumwxy
Total Sum of weight*X*Y.
Definition TH2.h:36
void Reset(Option_t *option="") override
Reset this histogram: contents, errors, etc.
Definition TH2.cxx:2560
Double_t fTsumwy2
Total Sum of weight*Y*Y.
Definition TH2.h:35
Double_t fTsumwy
Total Sum of weight*Y.
Definition TH2.h:34
void PutStats(Double_t *stats) override
Replace current statistics with the values in array stats.
Definition TH2.cxx:2445
A doubly linked list.
Definition TList.h:38
void Add(TObject *obj) override
Definition TList.h:81
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:354
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition TMultiGraph.h:34
TList * GetListOfGraphs() const
Definition TMultiGraph.h:68
static TClass * Class()
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.
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:457
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:226
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition TObject.cxx:1099
virtual TClass * IsA() const
Definition TObject.h:246
Basic string class.
Definition TString.h:138
void ToLower()
Change string to lower-case.
Definition TString.cxx:1189
TString & ReplaceSpecialCppChars()
Find special characters which are typically used in printf() calls and replace them by appropriate es...
Definition TString.cxx:1121
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:2384
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:640
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
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:1320
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754.
Definition TMath.h:913
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:673
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124
TMarker m
Definition textangle.C:8