Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGraph2D.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id: TGraph2D.cxx,v 1.00
2// Author: Olivier Couet
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "TROOT.h"
13#include "TBuffer.h"
14#include "TMath.h"
15#include "TH2.h"
16#include "TF2.h"
17#include "TList.h"
18#include "TGraph2D.h"
19#include "TGraphDelaunay.h"
20#include "TGraphDelaunay2D.h"
21#include "TVirtualPad.h"
22#include "TVirtualFitter.h"
23#include "TVirtualHistPainter.h"
24#include "TPluginManager.h"
25#include "TSystem.h"
26#include "strtok.h"
27#include "snprintf.h"
28
29#include <cstdlib>
30#include <cassert>
31#include <iostream>
32#include <fstream>
33
34#include "HFitInterface.h"
35#include "Fit/DataRange.h"
37
38
39
40/** \class TGraph2D
41 \ingroup Graphs
42Graphics object made of three arrays X, Y and Z with the same number of points each.
43
44- [Creating a TGraph2D](\ref G2D00)
45- [Drawing options](\ref G2D01)
46- [Examples](\ref G2D02)
47 - [SURF1 Example](\ref G2D021)
48 - [Fitting Example](\ref G2D022)
49 - [PCOL Example](\ref G2D023)
50- [Definition of the Delaunay triangulation (After B. Delaunay)](\ref G2D03)
51
52
53\anchor G2D00
54## Creating a TGraph2D
55
56This class has different constructors:
57- With an array's dimension and three arrays x, y, and z:
58~~~ {.cpp}
59 auto g = new TGraph2D(n, x, y, z);
60~~~
61 x, y, z arrays can be doubles, floats, or ints.
62- With an array's dimension only:
63~~~ {.cpp}
64 auto g = new TGraph2D(n);
65~~~
66 The internal arrays are then filled with `SetPoint()`. The following line
67 fills the internal arrays at the position `i` with the values
68 `x`, `y`, `z`.
69~~~ {.cpp}
70 g->SetPoint(i, x, y, z);
71~~~
72- Without parameters:
73~~~ {.cpp}
74 auto g = new TGraph2D();
75~~~
76 again `SetPoint()` must be used to fill the internal arrays.
77- From a file:
78~~~ {.cpp}
79 auto g = new TGraph2D("graph.dat");
80~~~
81 Arrays are read from the ASCII file "graph.dat" according to a specifies
82 format. The default format is `%%lg %%lg %%lg`
83
84Note that in any of these three cases, `SetPoint()` can be used to change a data
85point or add a new one. If the data point index (`i`) is greater than the
86current size of the internal arrays, they are automatically extended.
87
88Like TGraph some TGraph2D constructors do not have the TGraph2D title and name as parameters.
89For these constructors TGraph2D has the default title and name "Graph2D". To change the
90default title and name `SetTitle` and `SetName` should be called on the TGraph2D after its
91creation.
92
93\anchor G2D01
94## Drawing options
95
96Specific drawing options can be used to paint a TGraph2D:
97
98| Option | Description |
99|----------|-------------------------------------------------------------------|
100| "TRI" | The Delaunay triangles are drawn using filled area. An hidden surface drawing technique is used. The surface is painted with the current fill area color. The edges of each triangles are painted with the current line color. |
101| "TRIW" | The Delaunay triangles are drawn as wire frame. |
102| "TRI1" | The Delaunay triangles are painted with color levels. The edges of each triangles are painted with the current line color. |
103| "TRI2" | The Delaunay triangles are painted with color levels. |
104| "P" | Draw a marker at each vertex. |
105| "P0" | Draw a circle at each vertex. Each circle background is white. |
106| "PCOL" | Draw a marker at each vertex. The color of each marker is defined according to its Z position. |
107| "LINE" | Draw a 3D polyline. |
108| "CONT5" | Draw a contour plot using Delaunay triangles.|
109
110The Delaunay triangulation algorithm assumes that each (x, y) coordinate corresponds to a unique z value,
111meaning duplicate (x, y) points are not allowed. Consequently, when using drawing options that rely on this
112algorithm (e.g., TRI, SURF, etc.), a warning may appear instructing you to remove duplicates
113(see RemoveDuplicates()).
114
115A TGraph2D can be also drawn with any options valid to draw a 2D histogram
116(like `COL`, `SURF`, `LEGO`, `CONT` etc..).
117
118When a TGraph2D is drawn with one of the 2D histogram drawing option,
119an intermediate 2D histogram is filled using the Delaunay triangles
120to interpolate the data set. The 2D histogram has equidistant bins along the X
121and Y directions. The number of bins along each direction can be change using
122`SetNpx()` and `SetNpy()`. Each bin is filled with the Z
123value found via a linear interpolation on the plane defined by the triangle above
124the (X,Y) coordinates of the bin center.
125
126The existing (X,Y,Z) points can be randomly scattered.
127The Delaunay triangles are build in the (X,Y) plane. These 2D triangles are then
128used to define flat planes in (X,Y,Z) over which the interpolation is done to fill
129the 2D histogram. The 3D triangles int takes build a 3D surface in
130the form of tessellating triangles at various angles. The triangles found can be
131drawn in 3D with one of the TGraph2D specific drawing options.
132
133The histogram generated by the Delaunay interpolation can be accessed using the
134`GetHistogram()` method.
135
136The axis settings (title, ranges etc ...) can be changed accessing the axis via
137the GetXaxis GetYaxis and GetZaxis methods. They access the histogram axis created
138at drawing time only. Therefore they should called after the TGraph2D is drawn:
139
140~~~ {.cpp}
141 auto g = new TGraph2D();
142
143 [...]
144
145 g->Draw("tri1");
146 gPad->Update();
147 g->GetXaxis()->SetTitle("X axis title");
148~~~
149
150\anchor G2D02
151## Examples
152
153\anchor G2D021
154### SURF1 Example
155
156Begin_Macro(source)
157{
158 auto c = new TCanvas("c","Graph2D example",0,0,600,400);
159 Double_t x, y, z, P = 6.;
160 Int_t np = 200;
161 auto dt = new TGraph2D();
162 dt->SetTitle("Graph title; X axis title; Y axis title; Z axis title");
163 auto r = new TRandom();
164 for (Int_t N=0; N<np; N++) {
165 x = 2*P*(r->Rndm(N))-P;
166 y = 2*P*(r->Rndm(N))-P;
167 z = (sin(x)/x)*(sin(y)/y)+0.2;
168 dt->SetPoint(N,x,y,z);
169 }
170 gStyle->SetPalette(1);
171 dt->Draw("surf1");
172}
173End_Macro
174
175\anchor G2D022
176### Fitting Example
177
1782D graphs can be fitted as shown by the following example:
179
180Begin_Macro(source)
181../../../tutorials/math/fit/graph2dfit.C
182End_Macro
183
184\anchor G2D023
185### PCOL Example
186
187Example showing the PCOL option.
188
189Begin_Macro(source)
190{
191 auto c = new TCanvas("c","Graph2D example",0,0,600,400);
192 Double_t P = 5.;
193 Int_t npx = 20 ;
194 Int_t npy = 20 ;
195 Double_t x = -P;
196 Double_t y = -P;
197 Double_t z;
198 Int_t k = 0;
199 Double_t dx = (2*P)/npx;
200 Double_t dy = (2*P)/npy;
201 auto dt = new TGraph2D(npx*npy);
202 dt->SetNpy(41);
203 dt->SetNpx(40);
204 for (Int_t i=0; i<npx; i++) {
205 for (Int_t j=0; j<npy; j++) {
206 z = sin(sqrt(x*x+y*y))+1;
207 dt->SetPoint(k,x,y,z);
208 k++;
209 y = y+dy;
210 }
211 x = x+dx;
212 y = -P;
213 }
214 gStyle->SetPalette(1);
215 dt->SetMarkerStyle(20);
216 dt->Draw("pcol");
217}
218End_Macro
219
220\anchor G2D03
221## Definition of the Delaunay triangulation (After B. Delaunay)
222
223For a set S of points in the Euclidean plane, the unique triangulation DT(S)
224of S such that no point in S is inside the circumcircle of any triangle in
225DT(S). DT(S) is the dual of the Voronoi diagram of S.
226If n is the number of points in S, the Voronoi diagram of S is the partitioning
227of the plane containing S points into n convex polygons such that each polygon
228contains exactly one point and every point in a given polygon is closer to its
229central point than to any other. A Voronoi diagram is sometimes also known as
230a Dirichlet tessellation.
231
232\image html tgraph2d_delaunay.png
233
234[This applet](http://www.cs.cornell.edu/Info/People/chew/Delaunay.html)
235gives a nice practical view of Delaunay triangulation and Voronoi diagram.
236*/
237
238
239////////////////////////////////////////////////////////////////////////////////
240/// Graph2D default constructor
241
243 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(0)
244{
245 fSize = 0;
246 fMargin = 0.;
247 fNpx = 40;
248 fNpy = 40;
249 fDirectory = nullptr;
250 fHistogram = nullptr;
251 fDelaunay = nullptr;
252 fMaximum = -1111;
253 fMinimum = -1111;
254 fX = nullptr;
255 fY = nullptr;
256 fZ = nullptr;
257 fZout = 0;
258 fMaxIter = 100000;
259 fPainter = nullptr;
260 fFunctions = new TList;
262}
263
264
265////////////////////////////////////////////////////////////////////////////////
266/// Graph2D constructor with three vectors of ints as input.
267
269 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(n)
270{
271 Build(n);
272
273 // Copy the input vectors into local arrays
274 for (Int_t i = 0; i < fNpoints; ++i) {
275 fX[i] = (Double_t)x[i];
276 fY[i] = (Double_t)y[i];
277 fZ[i] = (Double_t)z[i];
278 }
279}
280
281
282////////////////////////////////////////////////////////////////////////////////
283/// Graph2D constructor with three vectors of floats as input.
284
286 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(n)
287{
288 Build(n);
289
290 // Copy the input vectors into local arrays
291 for (Int_t i = 0; i < fNpoints; ++i) {
292 fX[i] = x[i];
293 fY[i] = y[i];
294 fZ[i] = z[i];
295 }
296}
297
298
299////////////////////////////////////////////////////////////////////////////////
300/// Graph2D constructor with three vectors of doubles as input.
301
303 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(n)
304{
305 Build(n);
306
307 // Copy the input vectors into local arrays
308 for (Int_t i = 0; i < fNpoints; ++i) {
309 fX[i] = x[i];
310 fY[i] = y[i];
311 fZ[i] = z[i];
312 }
313}
314
315
316////////////////////////////////////////////////////////////////////////////////
317/// Graph2D constructor with a TH2 (h2) as input.
318/// Only the h2's bins within the X and Y axis ranges are used.
319/// Empty bins, recognized when both content and errors are zero, are excluded.
320
322 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(0)
323{
324 Build(h2->GetNbinsX()*h2->GetNbinsY());
325
326 TString gname = "Graph2D_from_" + TString(h2->GetName());
327 SetName(gname);
328 // need to call later because sets title in ref histogram
329 SetTitle(h2->GetTitle());
330
331 TAxis *xaxis = h2->GetXaxis();
332 TAxis *yaxis = h2->GetYaxis();
333 Int_t xfirst = xaxis->GetFirst();
334 Int_t xlast = xaxis->GetLast();
335 Int_t yfirst = yaxis->GetFirst();
336 Int_t ylast = yaxis->GetLast();
337
338 Double_t x, y, z;
339 Int_t k = 0;
340
341 for (Int_t i = xfirst; i <= xlast; i++) {
342 for (Int_t j = yfirst; j <= ylast; j++) {
343 x = xaxis->GetBinCenter(i);
344 y = yaxis->GetBinCenter(j);
345 z = h2->GetBinContent(i, j);
346 Double_t ez = h2->GetBinError(i, j);
347 if (z != 0. || ez != 0) {
348 SetPoint(k, x, y, z);
349 k++;
350 }
351 }
352 }
353}
354
355
356////////////////////////////////////////////////////////////////////////////////
357/// Graph2D constructor with name, title and three vectors of doubles as input.
358/// name : name of 2D graph (avoid blanks)
359/// title : 2D graph title
360/// if title is of the form "stringt;stringx;stringy;stringz"
361/// the 2D graph title is set to stringt, the x axis title to stringx,
362/// the y axis title to stringy,etc
363
364TGraph2D::TGraph2D(const char *name, const char *title,
366 : TNamed(name, title), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(n)
367{
368 Build(n);
369
370 // Copy the input vectors into local arrays
371 for (Int_t i = 0; i < fNpoints; ++i) {
372 fX[i] = x[i];
373 fY[i] = y[i];
374 fZ[i] = z[i];
375 }
376}
377
378
379////////////////////////////////////////////////////////////////////////////////
380/// Graph2D constructor. The arrays fX, fY and fZ should be filled via
381/// calls to SetPoint
382
384 : TNamed("Graph2D", "Graph2D"), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(n)
385{
386 Build(n);
387 for (Int_t i = 0; i < fNpoints; i++) {
388 fX[i] = 0.;
389 fY[i] = 0.;
390 fZ[i] = 0.;
391 }
392}
393
394
395////////////////////////////////////////////////////////////////////////////////
396/// Graph2D constructor reading input from filename
397/// filename is assumed to contain at least three columns of numbers.
398/// For files separated by a specific delimiter different from ' ' and '\\t' (e.g. ';' in csv files)
399/// you can avoid using %*s to bypass this delimiter by explicitly specify the "option" argument,
400/// e.g. option=" \\t,;" for columns of figures separated by any of these characters (' ', '\\t', ',', ';')
401/// used once (e.g. "1;1") or in a combined way (" 1;,;; 1").
402/// Note in that case, the instantiation is about 2 times slower.
403
405 : TNamed("Graph2D", filename), TAttLine(1, 1, 1), TAttFill(0, 1001), fNpoints(0)
406{
407 Double_t x, y, z;
410
411 std::ifstream infile(fname.Data());
412 if (!infile.good()) {
413 MakeZombie();
414 Error("TGraph2D", "Cannot open file: %s, TGraph2D is Zombie", filename);
415 return;
416 } else {
417 Build(100);
418 }
419 std::string line;
420 Int_t np = 0;
421
422 if (strcmp(option, "") == 0) { // No delimiters specified (standard constructor).
423
424 while (std::getline(infile, line, '\n')) {
425 if (3 != sscanf(line.c_str(), format, &x, &y, &z)) {
426 continue; // skip empty and ill-formed lines
427 }
428 SetPoint(np, x, y, z);
429 np++;
430 }
431
432 } else { // A delimiter has been specified in "option"
433
434 // Checking format and creating its boolean equivalent
436 format_.ReplaceAll(" ", "") ;
437 format_.ReplaceAll("\t", "") ;
438 format_.ReplaceAll("lg", "") ;
439 format_.ReplaceAll("s", "") ;
440 format_.ReplaceAll("%*", "0") ;
441 format_.ReplaceAll("%", "1") ;
442 if (!format_.IsDigit()) {
443 Error("TGraph2D", "Incorrect input format! Allowed format tags are {\"%%lg\",\"%%*lg\" or \"%%*s\"}");
444 return;
445 }
446 Int_t ntokens = format_.Length() ;
447 if (ntokens < 3) {
448 Error("TGraph2D", "Incorrect input format! Only %d tag(s) in format whereas 3 \"%%lg\" tags are expected!", ntokens);
449 return;
450 }
453 for (Int_t idx = 0; idx < ntokens; idx++) {
454 isTokenToBeSaved[idx] = TString::Format("%c", format_[idx]).Atoi() ; //atoi(&format_[idx]) does not work for some reason...
455 if (isTokenToBeSaved[idx] == 1) {
457 }
458 }
459 if (ntokens >= 3 && ntokensToBeSaved != 3) { //first condition not to repeat the previous error message
460 Error("TGraph2D", "Incorrect input format! There are %d \"%%lg\" tag(s) in format whereas 3 and only 3 are expected!", ntokensToBeSaved);
461 delete [] isTokenToBeSaved ;
462 return;
463 }
464
465 // Initializing loop variables
466 Bool_t isLineToBeSkipped = kFALSE ; //empty and ill-formed lines
467 char * token = nullptr ;
468 TString token_str = "" ;
469 Int_t token_idx = 0 ;
470 Double_t * value = new Double_t [3] ; //x,y,z buffers
471 Int_t value_idx = 0 ;
472
473 // Looping
474 char *rest;
475 while (std::getline(infile, line, '\n')) {
476 if (!line.empty()) {
477 if (line[line.size() - 1] == char(13)) { // removing DOS CR character
478 line.erase(line.end() - 1, line.end()) ;
479 }
480 token = R__STRTOK_R(const_cast<char*>(line.c_str()), option, &rest);
481 while (token != nullptr && value_idx < 3) {
483 token_str = TString(token) ;
484 token_str.ReplaceAll("\t", "") ;
485 if (!token_str.IsFloat()) {
487 break ;
488 } else {
489 value[value_idx] = token_str.Atof() ;
490 value_idx++ ;
491 }
492 }
493 token = R__STRTOK_R(nullptr, option, &rest); // next token
494 token_idx++ ;
495 }
496 if (!isLineToBeSkipped && value_idx == 3) {
497 x = value[0] ;
498 y = value[1] ;
499 z = value[2] ;
500 SetPoint(np, x, y, z) ;
501 np++ ;
502 }
503 }
505 token = nullptr ;
506 token_idx = 0 ;
507 value_idx = 0 ;
508 }
509
510 // Cleaning
511 delete [] isTokenToBeSaved ;
512 delete [] value ;
513 delete token ;
514 }
515 infile.close();
516}
517
518
519////////////////////////////////////////////////////////////////////////////////
520/// Graph2D copy constructor.
521/// copy everything apart from the list of contained functions
522
525 fX(nullptr), fY(nullptr), fZ(nullptr),
526 fHistogram(nullptr), fDirectory(nullptr), fPainter(nullptr)
527{
528 fFunctions = new TList(); // do not copy the functions
529
530 // use operator=
531 (*this) = g;
532
533 // append TGraph2D to gdirectory
536 if (fDirectory) {
537 // append without replacing existing objects
538 fDirectory->Append(this);
539 }
540 }
541}
542
543
544////////////////////////////////////////////////////////////////////////////////
545/// TGraph2D destructor.
546
548{
549 Clear();
550}
551
552
553////////////////////////////////////////////////////////////////////////////////
554/// Graph2D operator "="
555
557{
558 if (this == &g) return *this;
559
560 // delete before existing contained objects
561 if (fX) delete [] fX;
562 if (fY) delete [] fY;
563 if (fZ) delete [] fZ;
564 if (fHistogram && !fUserHisto) {
565 delete fHistogram;
566 fHistogram = nullptr;
567 fDelaunay = nullptr;
568 }
569 // copy everything except the function list
570 fNpoints = g.fNpoints;
571 fNpx = g.fNpx;
572 fNpy = g.fNpy;
573 fMaxIter = g.fMaxIter;
574 fSize = fNpoints; // force size to be the same of npoints
575 fX = (fSize > 0) ? new Double_t[fSize] : nullptr;
576 fY = (fSize > 0) ? new Double_t[fSize] : nullptr;
577 fZ = (fSize > 0) ? new Double_t[fSize] : nullptr;
578 fMinimum = g.fMinimum;
579 fMaximum = g.fMaximum;
580 fMargin = g.fMargin;
581 fZout = g.fZout;
582 fUserHisto = g.fUserHisto;
583 if (g.fHistogram)
584 fHistogram = (fUserHisto ) ? g.fHistogram : new TH2D(*g.fHistogram);
585
586 // copy the points
587 for (Int_t n = 0; n < fSize; n++) {
588 fX[n] = g.fX[n];
589 fY[n] = g.fY[n];
590 fZ[n] = g.fZ[n];
591 }
592
593 return *this;
594}
595
596////////////////////////////////////////////////////////////////////////////////
597/// Creates the 2D graph basic data structure
598
600{
601 if (n <= 0) {
602 Error("TGraph2D", "Invalid number of points (%d)", n);
603 return;
604 }
605
606 fSize = n;
607 fMargin = 0.;
608 fNpx = 40;
609 fNpy = 40;
610 fDirectory = nullptr;
611 fHistogram = nullptr;
612 fDelaunay = nullptr;
613 fMaximum = -1111;
614 fMinimum = -1111;
615 fX = new Double_t[fSize];
616 fY = new Double_t[fSize];
617 fZ = new Double_t[fSize];
618 fZout = 0;
619 fMaxIter = 100000;
620 fFunctions = new TList;
621 fPainter = nullptr;
623
626 if (fDirectory) {
627 fDirectory->Append(this, kTRUE);
628 }
629 }
630}
631
632////////////////////////////////////////////////////////////////////////////////
633/// Performs the operation: `z = z + c1*f(x,y,z)`
634/// Errors are not recalculated.
635///
636/// \param f may be a 2-D function TF2 or 3-d function TF3
637/// \param c1 a scaling factor, 1 by default
638
640{
641 //if (fHistogram) SetBit(kResetHisto);
642
643 for (Int_t i = 0; i < fNpoints; i++) {
644 fZ[i] += c1*f->Eval(fX[i], fY[i], fZ[i]);
645 }
646 if (gPad) gPad->Modified();
647}
648
649////////////////////////////////////////////////////////////////////////////////
650/// Apply function f to all the data points
651/// f may be a 2-D function TF2 or 3-d function TF3
652/// The Z values of the 2D graph are replaced by the new values computed
653/// using the function
654
656{
657 //if (fHistogram) SetBit(kResetHisto);
658
659 for (Int_t i = 0; i < fNpoints; i++) {
660 fZ[i] = f->Eval(fX[i], fY[i], fZ[i]);
661 }
662 if (gPad) gPad->Modified();
663}
664
665////////////////////////////////////////////////////////////////////////////////
666/// Browse
667
669{
670 Draw("p0");
671 gPad->Update();
672}
673
674
675////////////////////////////////////////////////////////////////////////////////
676/// Free all memory allocated by this object.
677
678void TGraph2D::Clear(Option_t * /*option = "" */)
679{
680 if (fX) delete [] fX;
681 fX = nullptr;
682 if (fY) delete [] fY;
683 fY = nullptr;
684 if (fZ) delete [] fZ;
685 fZ = nullptr;
686 fSize = fNpoints = 0;
687 if (fHistogram && !fUserHisto) {
688 delete fHistogram;
689 fHistogram = nullptr;
690 fDelaunay = nullptr;
691 }
692 if (fFunctions) {
695 delete fFunctions;
696 fFunctions = nullptr;
697 }
698 if (fDirectory) {
699 fDirectory->Remove(this);
700 fDirectory = nullptr;
701 }
702}
703
704
705////////////////////////////////////////////////////////////////////////////////
706/// Perform the automatic addition of the graph to the given directory
707///
708/// Note this function is called in place when the semantic requires
709/// this object to be added to a directory (I.e. when being read from
710/// a TKey or being Cloned)
711
713{
715 if (addStatus) {
716 SetDirectory(dir);
717 if (dir) {
719 }
720 }
721}
722
723
724////////////////////////////////////////////////////////////////////////////////
725/// Computes distance from point px,py to a graph
726
733
734
735////////////////////////////////////////////////////////////////////////////////
736/// Specific drawing options can be used to paint a TGraph2D:
737///
738/// - "TRI" : The Delaunay triangles are drawn using filled area.
739/// An hidden surface drawing technique is used. The surface is
740/// painted with the current fill area color. The edges of each
741/// triangles are painted with the current line color.
742/// - "TRIW" : The Delaunay triangles are drawn as wire frame
743/// - "TRI1" : The Delaunay triangles are painted with color levels. The edges
744/// of each triangles are painted with the current line color.
745/// - "TRI2" : the Delaunay triangles are painted with color levels.
746/// - "P" : Draw a marker at each vertex
747/// - "P0" : Draw a circle at each vertex. Each circle background is white.
748/// - "PCOL" : Draw a marker at each vertex. The color of each marker is
749/// defined according to its Z position.
750/// - "CONT" : Draw contours
751/// - "LINE" : Draw a 3D polyline
752///
753/// A TGraph2D can be also drawn with ANY options valid to draw a 2D histogram.
754///
755/// When a TGraph2D is drawn with one of the 2D histogram drawing option,
756/// a intermediate 2D histogram is filled using the Delaunay triangles
757/// technique to interpolate the data set.
758
760{
761 TString opt = option;
762 opt.ToLower();
763 if (gPad) {
764 if (!gPad->IsEditable()) gROOT->MakeDefCanvas();
765 if (!opt.Contains("same")) {
766 //the following statement is necessary in case one attempts to draw
767 //a temporary histogram already in the current pad
768 if (TestBit(kCanDelete)) gPad->GetListOfPrimitives()->Remove(this);
769 gPad->Clear();
770 }
771 }
772 AppendPad(opt.Data());
773}
774
775
776////////////////////////////////////////////////////////////////////////////////
777/// Executes action corresponding to one event
778
780{
781 if (fHistogram) fHistogram->ExecuteEvent(event, px, py);
782}
783
784
785////////////////////////////////////////////////////////////////////////////////
786/// search object named name in the list of functions
787
789{
790 return fFunctions ? fFunctions->FindObject(name) : nullptr;
791}
792
793
794////////////////////////////////////////////////////////////////////////////////
795/// search object obj in the list of functions
796
798{
799 return fFunctions ? fFunctions->FindObject(obj) : nullptr;
800}
801
802
803////////////////////////////////////////////////////////////////////////////////
804/// Fits this graph with function with name fname
805/// Predefined functions such as gaus, expo and poln are automatically
806/// created by ROOT.
807/// fname can also be a formula, accepted by the linear fitter (linear parts divided
808/// by "++" sign), for example "x++sin(y)" for fitting "[0]*x+[1]*sin(y)"
809
811{
812
813 char *linear;
814 linear = (char*)strstr(fname, "++");
815
816 if (linear) {
817 TF2 f2(fname, fname);
818 return Fit(&f2, option, "");
819 }
820 TF2 * f2 = (TF2*)gROOT->GetFunction(fname);
821 if (!f2) {
822 Printf("Unknown function: %s", fname);
823 return -1;
824 }
825 return Fit(f2, option, "");
826
827}
828
829
830////////////////////////////////////////////////////////////////////////////////
831/// Fits this 2D graph with function f2
832///
833/// f2 is an already predefined function created by TF2.
834///
835/// See TGraph::Fit for the available fitting options and fitting notes
836///
838{
839 // internal graph2D fitting methods
841 Option_t *goption = "";
843
844 // create range and minimizer options with default values
848}
849
850
851////////////////////////////////////////////////////////////////////////////////
852/// Display a GUI panel with all graph fit options.
853///
854/// See class TFitEditor for example
855
857{
858 if (!gPad)
859 gROOT->MakeDefCanvas();
860
861 if (!gPad) {
862 Error("FitPanel", "Unable to create a default canvas");
863 return;
864 }
865
866 // use plugin manager to create instance of TFitEditor
867 TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
868 if (handler && handler->LoadPlugin() != -1) {
869 if (handler->ExecPlugin(2, gPad, this) == 0)
870 Error("FitPanel", "Unable to crate the FitPanel");
871 } else
872 Error("FitPanel", "Unable to find the FitPanel plug-in");
873
874}
875
876
877////////////////////////////////////////////////////////////////////////////////
878/// Get x axis of the graph.
879
881{
882 TH1 *h = ((TGraph2D*)this)->GetHistogram("empty");
883 if (!h) return nullptr;
884 return h->GetXaxis();
885}
886
887
888////////////////////////////////////////////////////////////////////////////////
889/// Get y axis of the graph.
890
892{
893 TH1 *h = ((TGraph2D*)this)->GetHistogram("empty");
894 if (!h) return nullptr;
895 return h->GetYaxis();
896}
897
898
899////////////////////////////////////////////////////////////////////////////////
900/// Get z axis of the graph.
901
903{
904 TH1 *h = ((TGraph2D*)this)->GetHistogram("empty");
905 if (!h) return nullptr;
906 return h->GetZaxis();
907}
908
909
910////////////////////////////////////////////////////////////////////////////////
911/// Returns the X and Y graphs building a contour. A contour level may
912/// consist in several parts not connected to each other. This function
913/// returns them in a graphs' list.
914
916{
917 if (fNpoints <= 0) {
918 Error("GetContourList", "Empty TGraph2D");
919 return nullptr;
920 }
921
922 if (!fHistogram) GetHistogram("empty");
923
925
926 return fPainter->GetContourList(contour);
927}
928
929
930////////////////////////////////////////////////////////////////////////////////
931/// This function is called by Graph2DFitChisquare.
932/// It always returns a negative value. Real implementation in TGraph2DErrors
933
935{
936 return -1;
937}
938
939
940////////////////////////////////////////////////////////////////////////////////
941/// This function is called by Graph2DFitChisquare.
942/// It always returns a negative value. Real implementation in TGraph2DErrors
943
945{
946 return -1;
947}
948
949
950////////////////////////////////////////////////////////////////////////////////
951/// This function is called by Graph2DFitChisquare.
952/// It always returns a negative value. Real implementation in TGraph2DErrors
953
955{
956 return -1;
957}
958
959
960////////////////////////////////////////////////////////////////////////////////
961/// Add a TGraphDelaunay in the list of the fHistogram's functions
962
964{
965
967
968 if (oldInterp) {
969 TGraphDelaunay *dt = new TGraphDelaunay(this);
970 dt->SetMaxIter(fMaxIter);
971 dt->SetMarginBinsContent(fZout);
972 fDelaunay = dt;
974 if (!hl->FindObject("TGraphDelaunay")) hl->Add(fDelaunay);
975 } else {
977 dt->SetMarginBinsContent(fZout);
978 fDelaunay = dt;
980 if (!hl->FindObject("TGraphDelaunay2D")) hl->Add(fDelaunay);
981 }
982}
983
984////////////////////////////////////////////////////////////////////////////////
985/// Return pointer to function with name.
986///
987/// Functions such as TGraph2D::Fit store the fitted function in the list of
988/// functions of this graph.
989
990TF2 *TGraph2D::GetFunction(const char *name) const
991{
992 return dynamic_cast<TF2*>(FindObject(name));
993}
994
995////////////////////////////////////////////////////////////////////////////////
996/// By default returns a pointer to the Delaunay histogram. If fHistogram
997/// doesn't exist, books the 2D histogram fHistogram with a margin around
998/// the hull. Calls TGraphDelaunay::Interpolate at each bin centre to build up
999/// an interpolated 2D histogram.
1000///
1001/// If the "empty" option is selected, returns an empty histogram booked with
1002/// the limits of fX, fY and fZ. This option is used when the data set is
1003/// drawn with markers only. In that particular case there is no need to
1004/// find the Delaunay triangles.
1005///
1006/// By default use the new interpolation routine based on Triangles
1007/// If the option "old" the old interpolation is used
1008
1010{
1011 // for an empty graph create histogram in [0,1][0,1]
1012 if (fNpoints <= 0) {
1013 if (!fHistogram) {
1014 // do not add the histogram to gDirectory
1015 TDirectory::TContext ctx(nullptr);
1016 fHistogram = new TH2D(GetName(), GetTitle(), fNpx , 0., 1., fNpy, 0., 1.);
1018 }
1019 return fHistogram;
1020 }
1021
1022 TString opt = option;
1023 opt.ToLower();
1024 Bool_t empty = opt.Contains("empty");
1025 Bool_t oldInterp = opt.Contains("old");
1026
1027 if (fHistogram) {
1028 if (!empty && fHistogram->GetEntries() == 0) {
1029 if (!fUserHisto) {
1030 delete fHistogram;
1031 fHistogram = nullptr;
1032 fDelaunay = nullptr;
1033 }
1034 } else if (fHistogram->GetEntries() == 0)
1035 {; }
1036 // check case if interpolation type has changed
1037 else if ( (TestBit(kOldInterpolation) && !oldInterp) || ( !TestBit(kOldInterpolation) && oldInterp ) ) {
1038 delete fHistogram;
1039 fHistogram = nullptr;
1040 fDelaunay = nullptr;
1041 }
1042 // normal case return existing histogram
1043 else {
1044 return fHistogram;
1045 }
1046 }
1047
1049
1050 // Book fHistogram if needed. It is not added in the current directory
1051 if (!fUserHisto) {
1056 hxmin = xmin - fMargin * (xmax - xmin);
1057 hymin = ymin - fMargin * (ymax - ymin);
1058 hxmax = xmax + fMargin * (xmax - xmin);
1059 hymax = ymax + fMargin * (ymax - ymin);
1060 Double_t epsilon = 1e-9;
1061 if (TMath::AreEqualRel(hxmax,hxmin,epsilon)) {
1062 if (TMath::Abs(hxmin) < epsilon) {
1063 hxmin = -0.001;
1064 hxmax = 0.001;
1065 } else {
1066 hxmin = hxmin-TMath::Abs(hxmin)*(epsilon/2.);
1067 hxmax = hxmax+TMath::Abs(hxmax)*(epsilon/2.);
1068 }
1069 }
1070 if (TMath::AreEqualRel(hymax, hymin, epsilon)) {
1071 if (TMath::Abs(hymin) < epsilon) {
1072 hymin = -0.001;
1073 hymax = 0.001;
1074 } else {
1075 hymin = hymin-TMath::Abs(hymin)*(epsilon/2.);
1076 hymax = hymax+TMath::Abs(hymax)*(epsilon/2.);
1077 }
1078 }
1079 if (fHistogram) {
1082 } else {
1083 TDirectory::TContext ctx(nullptr); // to avoid adding fHistogram to gDirectory
1084 fHistogram = new TH2D(GetName(), GetTitle(),
1085 fNpx , hxmin, hxmax,
1086 fNpy, hymin, hymax);
1088 }
1091 } else {
1096 }
1097
1098 // Option "empty" is selected. An empty histogram is returned.
1100 if (empty) {
1101 if (fMinimum != -1111) {
1102 hzmin = fMinimum;
1103 } else {
1104 hzmin = GetZminE();
1105 }
1106 if (fMaximum != -1111) {
1107 hzmax = fMaximum;
1108 } else {
1109 hzmax = GetZmaxE();
1110 }
1111 if (hzmin == hzmax) {
1112 Double_t hz = hzmin;
1113 if (hz==0) {
1114 hzmin = -0.01;
1115 hzmax = 0.01;
1116 } else {
1117 hzmin = hz - 0.01 * TMath::Abs(hz);
1118 hzmax = hz + 0.01 * TMath::Abs(hz);
1119 }
1120 }
1123 return fHistogram;
1124 }
1125
1126 Double_t dx = (hxmax - hxmin) / fNpx;
1127 Double_t dy = (hymax - hymin) / fNpy;
1128
1129 Double_t x, y, z;
1130
1131 for (Int_t ix = 1; ix <= fNpx; ix++) {
1132 x = hxmin + (ix - 0.5) * dx;
1133 for (Int_t iy = 1; iy <= fNpy; iy++) {
1134 y = hymin + (iy - 0.5) * dy;
1135 // do interpolation
1136 if (oldInterp)
1137 z = ((TGraphDelaunay*)fDelaunay)->ComputeZ(x, y);
1138 else
1139 z = ((TGraphDelaunay2D*)fDelaunay)->ComputeZ(x, y);
1140
1141 fHistogram->Fill(x, y, z);
1142 }
1143 }
1144
1145 hzmin = GetZminE();
1146 hzmax = GetZmaxE();
1149
1150 if (fMinimum != -1111) fHistogram->SetMinimum(fMinimum);
1151 if (fMaximum != -1111) fHistogram->SetMaximum(fMaximum);
1152
1153 return fHistogram;
1154}
1155
1156
1157////////////////////////////////////////////////////////////////////////////////
1158/// Returns the X maximum
1159
1161{
1162 Double_t v = fX[0];
1163 for (Int_t i = 1; i < fNpoints; i++) if (fX[i] > v) v = fX[i];
1164 return v;
1165}
1166
1167
1168////////////////////////////////////////////////////////////////////////////////
1169/// Returns the X minimum
1170
1172{
1173 Double_t v = fX[0];
1174 for (Int_t i = 1; i < fNpoints; i++) if (fX[i] < v) v = fX[i];
1175 return v;
1176}
1177
1178
1179////////////////////////////////////////////////////////////////////////////////
1180/// Returns the Y maximum
1181
1183{
1184 Double_t v = fY[0];
1185 for (Int_t i = 1; i < fNpoints; i++) if (fY[i] > v) v = fY[i];
1186 return v;
1187}
1188
1189
1190////////////////////////////////////////////////////////////////////////////////
1191/// Returns the Y minimum
1192
1194{
1195 Double_t v = fY[0];
1196 for (Int_t i = 1; i < fNpoints; i++) if (fY[i] < v) v = fY[i];
1197 return v;
1198}
1199
1200
1201////////////////////////////////////////////////////////////////////////////////
1202/// Returns the Z maximum
1203
1205{
1206 Double_t v = fZ[0];
1207 for (Int_t i = 1; i < fNpoints; i++) if (fZ[i] > v) v = fZ[i];
1208 return v;
1209}
1210
1211
1212////////////////////////////////////////////////////////////////////////////////
1213/// Returns the Z minimum
1214
1216{
1217 Double_t v = fZ[0];
1218 for (Int_t i = 1; i < fNpoints; i++) if (fZ[i] < v) v = fZ[i];
1219 return v;
1220}
1221
1222////////////////////////////////////////////////////////////////////////////////
1223/// Get x, y and z values for point number i.
1224/// The function returns -1 in case of an invalid request or the point number otherwise
1225
1227{
1228 if (i < 0 || i >= fNpoints) return -1;
1229 if (!fX || !fY || !fZ) return -1;
1230 x = fX[i];
1231 y = fY[i];
1232 z = fZ[i];
1233 return i;
1234}
1235
1236////////////////////////////////////////////////////////////////////////////////
1237/// Finds the z value at the position (x,y) thanks to
1238/// the Delaunay interpolation.
1239
1241{
1242 if (fNpoints <= 0) {
1243 Error("Interpolate", "Empty TGraph2D");
1244 return 0;
1245 }
1246
1247 if (!fHistogram) GetHistogram("empty");
1248 if (!fDelaunay) {
1250 if (!TestBit(kOldInterpolation) ) {
1251 fDelaunay = hl->FindObject("TGraphDelaunay2D");
1252 if (!fDelaunay) fDelaunay = hl->FindObject("TGraphDelaunay");
1253 }
1254 else {
1255 // if using old implementation
1256 fDelaunay = hl->FindObject("TGraphDelaunay");
1257 if (!fDelaunay) fDelaunay = hl->FindObject("TGraphDelaunay2D");
1258 }
1259 }
1260
1261 if (!fDelaunay) return TMath::QuietNaN();
1262
1264 return ((TGraphDelaunay2D*)fDelaunay)->ComputeZ(x, y);
1265 else if (fDelaunay->IsA() == TGraphDelaunay::Class() )
1266 return ((TGraphDelaunay*)fDelaunay)->ComputeZ(x, y);
1267
1268 // cannot be here
1269 assert(false);
1270 return TMath::QuietNaN();
1271}
1272
1273
1274////////////////////////////////////////////////////////////////////////////////
1275/// Paints this 2D graph with its current attributes
1276
1278{
1279 if (fNpoints <= 0) {
1280 Error("Paint", "Empty TGraph2D");
1281 return;
1282 }
1283
1284 TString opt = option;
1285 opt.ToLower();
1286 if (opt.Contains("p") && !opt.Contains("tri")) {
1287 if (!opt.Contains("pol") &&
1288 !opt.Contains("sph") &&
1289 !opt.Contains("psr")) opt.Append("tri0");
1290 }
1291
1292 if (opt.Contains("line") && !opt.Contains("tri")) opt.Append("tri0");
1293
1294 if (opt.Contains("err") && !opt.Contains("tri")) opt.Append("tri0");
1295
1296 if (opt.Contains("tri0")) {
1297 GetHistogram("empty");
1298 } else if (opt.Contains("old")) {
1299 GetHistogram("old");
1300 } else {
1301 GetHistogram();
1302 }
1303
1312 fHistogram->Paint(opt.Data());
1313}
1314
1315
1316////////////////////////////////////////////////////////////////////////////////
1317/// Print 2D graph values.
1318
1320{
1321 for (Int_t i = 0; i < fNpoints; i++) {
1322 printf("x[%d]=%g, y[%d]=%g, z[%d]=%g\n", i, fX[i], i, fY[i], i, fZ[i]);
1323 }
1324}
1325
1326
1327////////////////////////////////////////////////////////////////////////////////
1328/// Projects a 2-d graph into 1 or 2-d histograms depending on the option parameter.
1329/// option may contain a combination of the characters x,y,z:
1330///
1331/// - option = "x" return the x projection into a TH1D histogram
1332/// - option = "y" return the y projection into a TH1D histogram
1333/// - option = "xy" return the x versus y projection into a TH2D histogram
1334/// - option = "yx" return the y versus x projection into a TH2D histogram
1335
1337{
1338 if (fNpoints <= 0) {
1339 Error("Project", "Empty TGraph2D");
1340 return nullptr;
1341 }
1342
1343 TString opt = option;
1344 opt.ToLower();
1345
1346 Int_t pcase = 0;
1347 if (opt.Contains("x")) pcase = 1;
1348 if (opt.Contains("y")) pcase = 2;
1349 if (opt.Contains("xy")) pcase = 3;
1350 if (opt.Contains("yx")) pcase = 4;
1351
1352 // Create the projection histogram
1353 TH1D *h1 = nullptr;
1354 TH2D *h2 = nullptr;
1355 Int_t nch = strlen(GetName()) + opt.Length() + 2;
1356 char *name = new char[nch];
1357 snprintf(name, nch, "%s_%s", GetName(), option);
1358 nch = strlen(GetTitle()) + opt.Length() + 2;
1359 char *title = new char[nch];
1360 snprintf(title, nch, "%s_%s", GetTitle(), option);
1361
1366
1367 switch (pcase) {
1368 case 1:
1369 // "x"
1370 h1 = new TH1D(name, title, fNpx, hxmin, hxmax);
1371 break;
1372 case 2:
1373 // "y"
1374 h1 = new TH1D(name, title, fNpy, hymin, hymax);
1375 break;
1376 case 3:
1377 // "xy"
1378 h2 = new TH2D(name, title, fNpx, hxmin, hxmax, fNpy, hymin, hymax);
1379 break;
1380 case 4:
1381 // "yx"
1382 h2 = new TH2D(name, title, fNpy, hymin, hymax, fNpx, hxmin, hxmax);
1383 break;
1384 }
1385
1386 delete [] name;
1387 delete [] title;
1388 TH1 *h = h1;
1389 if (h2) h = h2;
1390 if (h == nullptr) return nullptr;
1391
1392 // Fill the projected histogram
1393 Double_t entries = 0;
1394 for (Int_t n = 0; n < fNpoints; n++) {
1395 switch (pcase) {
1396 case 1:
1397 // "x"
1398 h1->Fill(fX[n], fZ[n]);
1399 break;
1400 case 2:
1401 // "y"
1402 h1->Fill(fY[n], fZ[n]);
1403 break;
1404 case 3:
1405 // "xy"
1406 h2->Fill(fX[n], fY[n], fZ[n]);
1407 break;
1408 case 4:
1409 // "yx"
1410 h2->Fill(fY[n], fX[n], fZ[n]);
1411 break;
1412 }
1413 entries += fZ[n];
1414 }
1415 h->SetEntries(entries);
1416 return h;
1417}
1418
1419
1420////////////////////////////////////////////////////////////////////////////////
1421/// Deletes duplicated points.
1422///
1423/// The Delaunay triangulation algorithm assumes that each (x, y) coordinate corresponds to a unique z value,
1424/// meaning duplicate (x, y) points are not allowed. Consequently, when using drawing options that rely on this
1425/// algorithm (e.g., TRI, SURF, etc.), a warning may appear instructing you to remove duplicates.
1426/// This function provides a way to handle such duplicates.
1427///
1428/// Example:
1429/// ~~~ {.cpp}
1430/// g->RemoveDuplicates();
1431/// g->Draw("TRI1");
1432/// ~~~
1433
1435{
1436 for (int i = 0; i < fNpoints; i++) {
1437 double x = fX[i];
1438 double y = fY[i];
1439 for (int j = i + 1; j < fNpoints; j++) {
1440 if (x == fX[j] && y == fY[j]) {
1441 RemovePoint(j);
1442 j--;
1443 }
1444 }
1445 }
1446
1447 return fNpoints;
1448}
1449
1450
1451////////////////////////////////////////////////////////////////////////////////
1452/// Recursively remove object from the list of functions
1453
1455{
1456 if (fFunctions) {
1459 }
1460 if (fHistogram == obj)
1461 fHistogram = nullptr;
1462}
1463
1464
1465////////////////////////////////////////////////////////////////////////////////
1466/// Deletes point number ipoint
1467
1469{
1470 if (ipoint < 0) return -1;
1471 if (ipoint >= fNpoints) return -1;
1472 for (Int_t i = ipoint; i < fNpoints - 1; i++) {
1473 fX[i] = fX[i+1];
1474 fY[i] = fY[i+1];
1475 fZ[i] = fZ[i+1];
1476 }
1477 fNpoints--;
1478 if (fHistogram) {
1479 delete fHistogram;
1480 fHistogram = nullptr;
1481 fDelaunay = nullptr;
1482 }
1483 return ipoint;
1484}
1485
1486
1487////////////////////////////////////////////////////////////////////////////////
1488/// Saves primitive as a C++ statement(s) on output stream out
1489
1490void TGraph2D::SavePrimitive(std::ostream &out, Option_t *option)
1491{
1492 TString arrx = SavePrimitiveVector(out, "graph2d_x", fNpoints, fX, kTRUE);
1493 TString arry = SavePrimitiveVector(out, "graph2d_y", fNpoints, fY);
1494 TString arrz = SavePrimitiveVector(out, "graph2d_z", fNpoints, fZ);
1495
1496 SavePrimitiveConstructor(out, Class(), "graph2d",
1497 TString::Format("%d, %s.data(), %s.data(), %s.data()", fNpoints, arrx.Data(), arry.Data(), arrz.Data()), kFALSE);
1498
1499 if (strcmp(GetName(), "Graph2D"))
1500 out << " graph2d->SetName(\"" << TString(GetName()).ReplaceSpecialCppChars() << "\");\n";
1501
1502 TString title = GetTitle();
1503 if (fHistogram)
1504 title = TString(fHistogram->GetTitle()) + ";" + fHistogram->GetXaxis()->GetTitle() + ";" +
1506
1507 out << " graph2d->SetTitle(\"" << title.ReplaceSpecialCppChars() << "\");\n";
1508
1509 if (!fDirectory)
1510 out << " graph2d->SetDirectory(nullptr);\n";
1511
1512 SaveFillAttributes(out, "graph2d", 0, 1001);
1513 SaveLineAttributes(out, "graph2d", 1, 1, 1);
1514 SaveMarkerAttributes(out, "graph2d", 1, 1, 1);
1515
1516 TH1::SavePrimitiveFunctions(out, "graph2d", fFunctions);
1517
1518 SavePrimitiveDraw(out, "graph2d", option);
1519}
1520
1521////////////////////////////////////////////////////////////////////////////////
1522/// Multiply the values of a TGraph2D by a constant c1.
1523///
1524/// If option contains "x" the x values are scaled
1525/// If option contains "y" the y values are scaled
1526/// If option contains "z" the z values are scaled
1527/// If option contains "xyz" all three x, y and z values are scaled
1528
1530{
1531 TString opt = option; opt.ToLower();
1532 if (opt.Contains("x")) {
1533 for (Int_t i=0; i<GetN(); i++)
1534 GetX()[i] *= c1;
1535 }
1536 if (opt.Contains("y")) {
1537 for (Int_t i=0; i<GetN(); i++)
1538 GetY()[i] *= c1;
1539 }
1540 if (opt.Contains("z")) {
1541 for (Int_t i=0; i<GetN(); i++)
1542 GetZ()[i] *= c1;
1543 }
1544}
1545
1546////////////////////////////////////////////////////////////////////////////////
1547/// Set number of points in the 2D graph.
1548/// Existing coordinates are preserved.
1549/// New coordinates above fNpoints are preset to 0.
1550
1552{
1553 if (n < 0) n = 0;
1554 if (n == fNpoints) return;
1555 if (n > fNpoints) SetPoint(n, 0, 0, 0);
1556 fNpoints = n;
1557}
1558
1559
1560////////////////////////////////////////////////////////////////////////////////
1561/// By default when an 2D graph is created, it is added to the list
1562/// of 2D graph objects in the current directory in memory.
1563/// This method removes reference to this 2D graph from current directory and add
1564/// reference to new directory dir. dir can be 0 in which case the
1565/// 2D graph does not belong to any directory.
1566
1568{
1569 if (fDirectory == dir) return;
1570 if (fDirectory) fDirectory->Remove(this);
1571 fDirectory = dir;
1572 if (fDirectory) fDirectory->Append(this);
1573}
1574
1575
1576////////////////////////////////////////////////////////////////////////////////
1577/// Sets the histogram to be filled.
1578/// If the 2D graph needs to be save in a TFile the following set should be
1579/// followed to read it back:
1580/// 1. Create TGraph2D
1581/// 2. Call g->SetHistogram(h), and do whatever you need to do
1582/// 3. Save g and h to the TFile, exit
1583/// 4. Open the TFile, retrieve g and h
1584/// 5. Call h->SetDirectory(0)
1585/// 6. Call g->SetHistogram(h) again
1586/// 7. Carry on as normal
1587///
1588/// By default use the new interpolation routine based on Triangles
1589/// If the option "old" the old interpolation is used
1590
1592{
1593 TString opt = option;
1594 opt.ToLower();
1595 Bool_t oldInterp = opt.Contains("old");
1596
1597 fUserHisto = kTRUE;
1598 fHistogram = (TH2D*)h;
1599 fNpx = h->GetNbinsX();
1600 fNpy = h->GetNbinsY();
1602}
1603
1604
1605////////////////////////////////////////////////////////////////////////////////
1606/// Sets the extra space (in %) around interpolated area for the 2D histogram
1607
1609{
1610 if (m < 0 || m > 1) {
1611 Warning("SetMargin", "The margin must be >= 0 && <= 1, fMargin set to 0.1");
1612 fMargin = 0.1;
1613 } else {
1614 fMargin = m;
1615 }
1616 if (fHistogram) {
1617 delete fHistogram;
1618 fHistogram = nullptr;
1619 fDelaunay = nullptr;
1620 }
1621}
1622
1623
1624////////////////////////////////////////////////////////////////////////////////
1625/// Sets the histogram bin height for points lying outside the TGraphDelaunay
1626/// convex hull ie: the bins in the margin.
1627
1629{
1630 fZout = z;
1631 if (fHistogram) {
1632 delete fHistogram;
1633 fHistogram = nullptr;
1634 fDelaunay = nullptr;
1635 }
1636}
1637
1638
1639////////////////////////////////////////////////////////////////////////////////
1640/// Set maximum.
1641
1643{
1644 fMaximum = maximum;
1645 TH1 * h = GetHistogram();
1646 if (h) h->SetMaximum(maximum);
1647}
1648
1649
1650////////////////////////////////////////////////////////////////////////////////
1651/// Set minimum.
1652
1654{
1655 fMinimum = minimum;
1656 TH1 * h = GetHistogram();
1657 if (h) h->SetMinimum(minimum);
1658}
1659
1660
1661////////////////////////////////////////////////////////////////////////////////
1662/// Changes the name of this 2D graph
1663
1664void TGraph2D::SetName(const char *name)
1665{
1666 // 2D graphs are named objects in a THashList.
1667 // We must update the hashlist if we change the name
1668 if (fDirectory) fDirectory->Remove(this);
1669 fName = name;
1670 if (fDirectory) fDirectory->Append(this);
1671}
1672
1673
1674////////////////////////////////////////////////////////////////////////////////
1675/// Change the name and title of this 2D graph
1676///
1677
1678void TGraph2D::SetNameTitle(const char *name, const char *title)
1679{
1680 // 2D graphs are named objects in a THashList.
1681 // We must update the hashlist if we change the name
1682 if (fDirectory) fDirectory->Remove(this);
1683 fName = name;
1684 SetTitle(title);
1685 if (fDirectory) fDirectory->Append(this);
1686}
1687
1688
1689////////////////////////////////////////////////////////////////////////////////
1690/// Sets the number of bins along X used to draw the function
1691
1693{
1694 if (npx < 4) {
1695 Warning("SetNpx", "Number of points must be >4 && < 500, fNpx set to 4");
1696 fNpx = 4;
1697 } else if (npx > 500) {
1698 Warning("SetNpx", "Number of points must be >4 && < 500, fNpx set to 500");
1699 fNpx = 500;
1700 } else {
1701 fNpx = npx;
1702 }
1703 if (fHistogram) {
1704 delete fHistogram;
1705 fHistogram = nullptr;
1706 fDelaunay = nullptr;
1707 }
1708}
1709
1710
1711////////////////////////////////////////////////////////////////////////////////
1712/// Sets the number of bins along Y used to draw the function
1713
1715{
1716 if (npy < 4) {
1717 Warning("SetNpy", "Number of points must be >4 && < 500, fNpy set to 4");
1718 fNpy = 4;
1719 } else if (npy > 500) {
1720 Warning("SetNpy", "Number of points must be >4 && < 500, fNpy set to 500");
1721 fNpy = 500;
1722 } else {
1723 fNpy = npy;
1724 }
1725 if (fHistogram) {
1726 delete fHistogram;
1727 fHistogram = nullptr;
1728 fDelaunay = nullptr;
1729 }
1730}
1731
1732
1733////////////////////////////////////////////////////////////////////////////////
1734/// Sets point number n.
1735/// If n is greater than the current size, the arrays are automatically
1736/// extended.
1737
1739{
1740 if (n < 0) return;
1741
1742 if (!fX || !fY || !fZ || n >= fSize) {
1743 // re-allocate the object
1744 Int_t newN = TMath::Max(2 * fSize, n + 1);
1745 Double_t *savex = new Double_t [newN];
1746 Double_t *savey = new Double_t [newN];
1747 Double_t *savez = new Double_t [newN];
1748 if (fX && fSize) {
1749 memcpy(savex, fX, fSize * sizeof(Double_t));
1750 memset(&savex[fSize], 0, (newN - fSize)*sizeof(Double_t));
1751 delete [] fX;
1752 }
1753 if (fY && fSize) {
1754 memcpy(savey, fY, fSize * sizeof(Double_t));
1755 memset(&savey[fSize], 0, (newN - fSize)*sizeof(Double_t));
1756 delete [] fY;
1757 }
1758 if (fZ && fSize) {
1759 memcpy(savez, fZ, fSize * sizeof(Double_t));
1760 memset(&savez[fSize], 0, (newN - fSize)*sizeof(Double_t));
1761 delete [] fZ;
1762 }
1763 fX = savex;
1764 fY = savey;
1765 fZ = savez;
1766 fSize = newN;
1767 }
1768 fX[n] = x;
1769 fY[n] = y;
1770 fZ[n] = z;
1771 fNpoints = TMath::Max(fNpoints, n + 1);
1772}
1773
1774
1775////////////////////////////////////////////////////////////////////////////////
1776/// Sets the 2D graph title.
1777///
1778/// This method allows to change the global title and the axis' titles of a 2D
1779/// graph. If `g` is the 2D graph one can do:
1780///
1781/// ~~~ {.cpp}
1782/// g->SetTitle("Graph title; X axis title; Y axis title; Z axis title");
1783/// ~~~
1784
1785void TGraph2D::SetTitle(const char* title)
1786{
1787 fTitle = title;
1788 if (fHistogram) fHistogram->SetTitle(title);
1789}
1790
1791
1792////////////////////////////////////////////////////////////////////////////////
1793/// Stream a class object
1794
1796{
1797 if (b.IsReading()) {
1798 UInt_t R__s, R__c;
1799 Version_t R__v = b.ReadVersion(&R__s, &R__c);
1800 b.ReadClassBuffer(TGraph2D::Class(), this, R__v, R__s, R__c);
1801
1803 } else {
1804 b.WriteClassBuffer(TGraph2D::Class(), this);
1805 }
1806}
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define g(i)
Definition RSha256.hxx:105
#define h(i)
Definition RSha256.hxx:106
#define e(i)
Definition RSha256.hxx:103
short Version_t
Class version identifier (short)
Definition RtypesCore.h:79
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
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.
#define gDirectory
Definition TDirectory.h:385
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t np
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t format
char name[80]
Definition TGX11.cxx:110
float xmin
float ymin
float xmax
float ymax
#define gROOT
Definition TROOT.h:411
void Printf(const char *fmt,...)
Formats a string in a circular formatting buffer and prints the string.
Definition TString.cxx:2509
R__EXTERN TSystem * gSystem
Definition TSystem.h:572
#define gPad
#define snprintf
Definition civetweb.c:1579
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition DataRange.h:35
Fill Area Attributes class.
Definition TAttFill.h:20
virtual Color_t GetFillColor() const
Return the fill area color.
Definition TAttFill.h:31
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition TAttFill.h:32
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:38
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition TAttFill.h:40
virtual void SaveFillAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1001)
Save fill attributes as C++ statement(s) on output stream out.
Definition TAttFill.cxx:238
Line Attributes class.
Definition TAttLine.h:20
virtual Color_t GetLineColor() const
Return the line color.
Definition TAttLine.h:35
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition TAttLine.h:44
virtual Width_t GetLineWidth() const
Return the line width.
Definition TAttLine.h:37
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition TAttLine.h:45
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:42
virtual Style_t GetLineStyle() const
Return the line style.
Definition TAttLine.h:36
virtual void SaveLineAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t widdef=1)
Save line attributes as C++ statement(s) on output stream out.
Definition TAttLine.cxx:274
Marker Attributes class.
Definition TAttMarker.h:20
virtual void SaveMarkerAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t sizdef=1)
Save line attributes as C++ statement(s) on output stream out.
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition TAttMarker.h:33
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition TAttMarker.h:39
virtual Color_t GetMarkerColor() const
Return the marker color.
Definition TAttMarker.h:32
virtual Size_t GetMarkerSize() const
Return the marker size.
Definition TAttMarker.h:34
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition TAttMarker.h:41
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition TAttMarker.h:46
Class to manage histogram axis.
Definition TAxis.h:32
const char * GetTitle() const override
Returns title of object.
Definition TAxis.h:137
Double_t GetXmax() const
Definition TAxis.h:142
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition TAxis.h:166
Double_t GetXmin() const
Definition TAxis.h:141
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
Buffer base class used for serializing objects.
Definition TBuffer.h:43
TDirectory::TContext keeps track and restore the current directory.
Definition TDirectory.h:89
Describe directory structure in memory.
Definition TDirectory.h:45
virtual void Append(TObject *obj, Bool_t replace=kFALSE)
Append object to this directory.
virtual TObject * Remove(TObject *)
Remove an object from the in-memory list.
A 2-Dim function with parameters.
Definition TF2.h:29
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Graphics object made of three arrays X, Y and Z with the same number of points each.
Definition TGraph2D.h:41
static TClass * Class()
Int_t fMaxIter
Maximum number of iterations to find Delaunay triangles.
Definition TGraph2D.h:48
TF2 * GetFunction(const char *name) const
Return pointer to function with name.
Definition TGraph2D.cxx:990
TObject * FindObject(const char *name) const override
search object named name in the list of functions
Definition TGraph2D.cxx:788
void SetNameTitle(const char *name, const char *title) override
Change the name and title of this 2D graph.
virtual Double_t GetYminE() const
Definition TGraph2D.h:143
TGraph2D()
Graph2D default constructor.
Definition TGraph2D.cxx:242
void Build(Int_t n)
Creates the 2D graph basic data structure.
Definition TGraph2D.cxx:599
Double_t Interpolate(Double_t x, Double_t y)
Finds the z value at the position (x,y) thanks to the Delaunay interpolation.
Int_t fNpoints
Number of points in the data set.
Definition TGraph2D.h:45
virtual Double_t GetZminE() const
Definition TGraph2D.h:145
virtual Double_t GetErrorZ(Int_t bin) const
This function is called by Graph2DFitChisquare.
Definition TGraph2D.cxx:954
Double_t GetMinimum() const
Definition TGraph2D.h:117
Double_t * GetY() const
Definition TGraph2D.h:123
virtual void FitPanel()
Display a GUI panel with all graph fit options.
Definition TGraph2D.cxx:856
virtual void Apply(TF2 *f)
Apply function f to all the data points f may be a 2-D function TF2 or 3-d function TF3 The Z values ...
Definition TGraph2D.cxx:655
void SetMarginBinsContent(Double_t z=0.)
Sets the histogram bin height for points lying outside the TGraphDelaunay convex hull ie: the bins in...
Int_t fNpx
Number of bins along X in fHistogram.
Definition TGraph2D.h:46
void SavePrimitive(std::ostream &out, Option_t *option="") override
Saves primitive as a C++ statement(s) on output stream out.
virtual Double_t GetYmaxE() const
Definition TGraph2D.h:142
Double_t GetYmin() const
Returns the Y minimum.
Bool_t fUserHisto
Definition TGraph2D.h:67
Double_t GetZmin() const
Returns the Z minimum.
void Browse(TBrowser *) override
Browse.
Definition TGraph2D.cxx:668
Double_t * GetX() const
Definition TGraph2D.h:122
virtual Double_t GetZmaxE() const
Definition TGraph2D.h:144
void RecursiveRemove(TObject *obj) override
Recursively remove object from the list of functions.
Double_t * fZ
[fNpoints]
Definition TGraph2D.h:52
void Streamer(TBuffer &) override
Stream a class object.
virtual Double_t GetErrorX(Int_t bin) const
This function is called by Graph2DFitChisquare.
Definition TGraph2D.cxx:934
Double_t fMargin
Extra space (in %) around interpolated area for fHistogram.
Definition TGraph2D.h:55
Double_t fMinimum
Minimum value for plotting along z.
Definition TGraph2D.h:53
Double_t GetXmin() const
Returns the X minimum.
void Print(Option_t *chopt="") const override
Print 2D graph values.
TH2D * GetHistogram(Option_t *option="")
By default returns a pointer to the Delaunay histogram.
TVirtualHistPainter * fPainter
!Pointer to histogram painter
Definition TGraph2D.h:61
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="")
Fits this graph with function with name fname Predefined functions such as gaus, expo and poln are au...
Definition TGraph2D.cxx:810
Double_t GetZmax() const
Returns the Z maximum.
Int_t RemovePoint(Int_t ipoint)
Deletes point number ipoint.
Double_t fMaximum
Maximum value for plotting along z.
Definition TGraph2D.h:54
~TGraph2D() override
TGraph2D destructor.
Definition TGraph2D.cxx:547
TAxis * GetZaxis() const
Get z axis of the graph.
Definition TGraph2D.cxx:902
void SetMargin(Double_t m=0.1)
Sets the extra space (in %) around interpolated area for the 2D histogram.
void SetName(const char *name) override
Changes the name of this 2D graph.
void SetNpy(Int_t npx=40)
Sets the number of bins along Y used to draw the function.
Double_t fZout
fHistogram bin height for points lying outside the interpolated area
Definition TGraph2D.h:56
void SetTitle(const char *title="") override
Sets the 2D graph title.
TH2D * fHistogram
!2D histogram of z values linearly interpolated on the triangles
Definition TGraph2D.h:58
virtual Double_t GetXminE() const
Definition TGraph2D.h:141
TList * GetContourList(Double_t contour)
Returns the X and Y graphs building a contour.
Definition TGraph2D.cxx:915
Double_t GetXmax() const
Returns the X maximum.
virtual Int_t GetPoint(Int_t i, Double_t &x, Double_t &y, Double_t &z) const
Get x, y and z values for point number i.
TAxis * GetYaxis() const
Get y axis of the graph.
Definition TGraph2D.cxx:891
virtual Double_t GetErrorY(Int_t bin) const
This function is called by Graph2DFitChisquare.
Definition TGraph2D.cxx:944
TObject * fDelaunay
! Pointer to Delaunay interpolator object
Definition TGraph2D.h:59
TH1 * Project(Option_t *option="x") const
Projects a 2-d graph into 1 or 2-d histograms depending on the option parameter.
Int_t GetN() const
Definition TGraph2D.h:121
virtual Double_t GetXmaxE() const
Definition TGraph2D.h:140
virtual void Set(Int_t n)
Set number of points in the 2D graph.
void SetMinimum(Double_t minimum=-1111)
Set minimum.
TGraph2D & operator=(const TGraph2D &)
Graph2D operator "=".
Definition TGraph2D.cxx:556
Double_t * fX
[fNpoints]
Definition TGraph2D.h:50
Double_t * fY
[fNpoints] Data set to be plotted
Definition TGraph2D.h:51
virtual void SetHistogram(TH2 *h, Option_t *option="")
Sets the histogram to be filled.
void SetMaximum(Double_t maximum=-1111)
Set maximum.
virtual void Add(TF2 *f, Double_t c1=1)
Performs the operation: z = z + c1*f(x,y,z) Errors are not recalculated.
Definition TGraph2D.cxx:639
Double_t GetYmax() const
Returns the Y maximum.
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Computes distance from point px,py to a graph.
Definition TGraph2D.cxx:727
virtual void SetDirectory(TDirectory *dir)
By default when an 2D graph is created, it is added to the list of 2D graph objects in the current di...
TDirectory * fDirectory
!Pointer to directory holding this 2D graph
Definition TGraph2D.h:60
void Clear(Option_t *option="") override
Free all memory allocated by this object.
Definition TGraph2D.cxx:678
virtual void Scale(Double_t c1=1., Option_t *option="z")
Multiply the values of a TGraph2D by a constant c1.
virtual void DirectoryAutoAdd(TDirectory *)
Perform the automatic addition of the graph to the given directory.
Definition TGraph2D.cxx:712
void CreateInterpolator(Bool_t oldInterp)
Add a TGraphDelaunay in the list of the fHistogram's functions.
Definition TGraph2D.cxx:963
Int_t fNpy
Number of bins along Y in fHistogram.
Definition TGraph2D.h:47
TList * fFunctions
Pointer to list of functions (fits and user)
Definition TGraph2D.h:57
virtual void SetPoint(Int_t point, Double_t x, Double_t y, Double_t z)
Sets point number n.
void Paint(Option_t *option="") override
Paints this 2D graph with its current attributes.
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Executes action corresponding to one event.
Definition TGraph2D.cxx:779
void Draw(Option_t *option="P0") override
Specific drawing options can be used to paint a TGraph2D:
Definition TGraph2D.cxx:759
Int_t fSize
!Real size of fX, fY and fZ
Definition TGraph2D.h:49
@ kOldInterpolation
Definition TGraph2D.h:70
Int_t RemoveDuplicates()
Deletes duplicated points.
TAxis * GetXaxis() const
Get x axis of the graph.
Definition TGraph2D.cxx:880
Double_t * GetZ() const
Definition TGraph2D.h:124
void SetNpx(Int_t npx=40)
Sets the number of bins along X used to draw the function.
TGraphDelaunay2D generates a Delaunay triangulation of a TGraph2D.
static TClass * Class()
TGraphDelaunay generates a Delaunay triangulation of a TGraph2D.
static TClass * Class()
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:927
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:109
TAxis * GetZaxis()
Definition TH1.h:574
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a line.
Definition TH1.cxx:2794
void SetTitle(const char *title) override
Change/set the title.
Definition TH1.cxx:6753
@ kNoStats
Don't draw stats box.
Definition TH1.h:404
TAxis * GetXaxis()
Definition TH1.h:572
TVirtualHistPainter * GetPainter(Option_t *option="")
Return pointer to painter.
Definition TH1.cxx:4488
virtual Double_t GetMaximum(Double_t maxval=FLT_MAX) const
Return maximum value smaller than maxval of bins in the range, unless the value has been overridden b...
Definition TH1.cxx:8573
virtual void SetMaximum(Double_t maximum=-1111)
Definition TH1.h:653
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3315
TAxis * GetYaxis()
Definition TH1.h:573
virtual void SetMinimum(Double_t minimum=-1111)
Definition TH1.h:654
virtual Double_t GetEntries() const
Return the current number of entries.
Definition TH1.cxx:4400
TList * GetListOfFunctions() const
Definition TH1.h:489
void Paint(Option_t *option="") override
Control routine to paint any kind of histograms.
Definition TH1.cxx:6241
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute action corresponding to one event.
Definition TH1.cxx:3211
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition TH1.cxx:9048
static Bool_t AddDirectoryStatus()
Static function: cannot be inlined on Windows/NT.
Definition TH1.cxx:741
static void SavePrimitiveFunctions(std::ostream &out, const char *varname, TList *lst)
Save list of functions Also can be used by TGraph classes.
Definition TH1.cxx:7458
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:356
Service class for 2-D histogram classes.
Definition TH2.h:30
Int_t Fill(Double_t) override
Invalid Fill method.
Definition TH2.cxx:356
A doubly linked list.
Definition TList.h:38
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition TList.cxx:575
void RecursiveRemove(TObject *obj) override
Remove object from this collection and recursively remove the object from all other objects (and coll...
Definition TList.cxx:761
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
Definition TList.cxx:467
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
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
TString fTitle
Definition TNamed.h:33
TString fName
Definition TNamed.h:32
Mother of all ROOT objects.
Definition TObject.h:41
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
static TString SavePrimitiveVector(std::ostream &out, const char *prefix, Int_t len, Double_t *arr, Bool_t empty_line=kFALSE)
Save array in the output stream "out" as vector.
Definition TObject.cxx:788
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition TObject.cxx:203
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual TClass * IsA() const
Definition TObject.h:246
void MakeZombie()
Definition TObject.h:53
static void SavePrimitiveDraw(std::ostream &out, const char *variable_name, Option_t *option=nullptr)
Save invocation of primitive Draw() method Skipped if option contains "nodraw" string.
Definition TObject.cxx:822
static void SavePrimitiveConstructor(std::ostream &out, TClass *cl, const char *variable_name, const char *constructor_agrs="", Bool_t empty_line=kTRUE)
Save object constructor in the output stream "out".
Definition TObject.cxx:771
void ResetBit(UInt_t f)
Definition TObject.h:201
@ kCanDelete
if object in a list can be deleted
Definition TObject.h:68
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition TObject.h:78
@ kMustCleanup
if object destructor must call RecursiveRemove()
Definition TObject.h:70
Longptr_t ExecPlugin(int nargs)
Int_t LoadPlugin()
Load the plugin library for this handler.
Basic string class.
Definition TString.h:138
Ssiz_t Length() const
Definition TString.h:425
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
const char * Data() const
Definition TString.h:384
TString & Append(const char *cs)
Definition TString.h:580
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
virtual Bool_t ExpandPathName(TString &path)
Expand a pathname getting rid of special shell characters like ~.
Definition TSystem.cxx:1285
virtual TList * GetContourList(Double_t contour) const =0
TLine * line
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
TFitResultPtr FitObject(TH1 *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
fitting function for a TH1 (called from TH1::Fit)
Definition HFitImpl.cxx:977
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition HFitImpl.cxx:685
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754.
Definition TMath.h:913
Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec)
Comparing floating points.
Definition TMath.h:429
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124
th1 Draw()
TMarker m
Definition textangle.C:8