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