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