Logo ROOT  
Reference Guide
TGraph.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id$
2// Author: Rene Brun, Olivier Couet 12/12/94
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
13#include "TROOT.h"
14#include "TBuffer.h"
15#include "TEnv.h"
16#include "TGraph.h"
17#include "TGraphErrors.h"
18#include "TGraphAsymmErrors.h"
19#include "TGraphBentErrors.h"
20#include "TH1.h"
21#include "TF1.h"
22#include "TStyle.h"
23#include "TMath.h"
24#include "TVectorD.h"
25#include "Foption.h"
26#include "TRandom.h"
27#include "TSpline.h"
28#include "TVirtualFitter.h"
29#include "TVirtualPad.h"
31#include "TBrowser.h"
32#include "TSystem.h"
33#include "TPluginManager.h"
34#include "strtok.h"
35
36#include <cstdlib>
37#include <string>
38#include <cassert>
39#include <iostream>
40#include <fstream>
41#include <cstring>
42
43#include "HFitInterface.h"
44#include "Fit/DataRange.h"
46
47extern void H1LeastSquareSeqnd(Int_t n, Double_t *a, Int_t idim, Int_t &ifail, Int_t k, Double_t *b);
48
50
51////////////////////////////////////////////////////////////////////////////////
52
53/** \class TGraph
54 \ingroup Graphs
55A TGraph is an object made of two arrays X and Y with npoints each.
56The TGraph painting is performed thanks to the TGraphPainter
57class. All details about the various painting options are given in this class.
58
59#### Notes
60
61 - Unlike histogram or tree (or even TGraph2D), TGraph objects
62 are not automatically attached to the current TFile, in order to keep the
63 management and size of the TGraph as small as possible.
64 - The TGraph constructors do not have the TGraph title and name as parameters.
65 A TGraph has the default title and name "Graph". To change the default title
66 and name `SetTitle` and `SetName` should be called on the TGraph after its creation.
67 TGraph was a light weight object to start with, like TPolyline or TPolyMarker.
68 That’s why it did not have any title and name parameters in the constructors.
69
70#### Example
71
72The picture below gives an example:
73
74Begin_Macro(source)
75{
76 double x[100], y[100];
77 int n = 20;
78 for (int i=0;i<n;i++) {
79 x[i] = i*0.1;
80 y[i] = 10*sin(x[i]+0.2);
81 }
82 auto g = new TGraph(n,x,y);
83 g->Draw("AC*");
84}
85End_Macro
86
87#### Default X-Points
88
89If one doesn't specify the points in the x-axis, they will get the default values 0, 1, 2, 3, (etc. depending
90on the length of the y-points):
91
92Begin_Macro(source)
93{
94 double y[6] = {3, 8, 1, 10, 5, 7};
95 auto g = new TGraph(6,y);
96 g->Draw();
97}
98End_Macro
99
100*/
101
102////////////////////////////////////////////////////////////////////////////////
103/// Graph default constructor.
104
106{
107 fNpoints = -1; //will be reset to 0 in CtorAllocate
108 if (!CtorAllocate()) return;
109}
110
111////////////////////////////////////////////////////////////////////////////////
112/// Constructor with only the number of points set
113/// the arrays x and y will be set later
114
116 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
117{
118 fNpoints = n;
119 if (!CtorAllocate()) return;
120 FillZero(0, fNpoints);
121}
122
123////////////////////////////////////////////////////////////////////////////////
124/// Graph normal constructor with ints.
125
127 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
128{
129 if (!x || !y) {
130 fNpoints = 0;
131 } else {
132 fNpoints = n;
133 }
134 if (!CtorAllocate()) return;
135 for (Int_t i = 0; i < n; i++) {
136 fX[i] = (Double_t)x[i];
137 fY[i] = (Double_t)y[i];
138 }
139}
140
141////////////////////////////////////////////////////////////////////////////////
142/// Graph normal constructor with floats.
143
145 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
146{
147 if (!x || !y) {
148 fNpoints = 0;
149 } else {
150 fNpoints = n;
151 }
152 if (!CtorAllocate()) return;
153 for (Int_t i = 0; i < n; i++) {
154 fX[i] = x[i];
155 fY[i] = y[i];
156 }
157}
158
159////////////////////////////////////////////////////////////////////////////////
160/// Default X-Points constructor. The points along the x-axis get the default
161/// values `start`, `start+step`, `start+2*step`, `start+3*step`, etc ...
162
164 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
165{
166 if (!y) {
167 fNpoints = 0;
168 } else {
169 fNpoints = n;
170 }
171 if (!CtorAllocate()) return;
172 for (Int_t i = 0; i < n; i++) {
173 fX[i] = start+i*step;
174 fY[i] = y[i];
175 }
176}
177
178////////////////////////////////////////////////////////////////////////////////
179/// Graph normal constructor with doubles.
180
182 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
183{
184 if (!x || !y) {
185 fNpoints = 0;
186 } else {
187 fNpoints = n;
188 }
189 if (!CtorAllocate()) return;
190 n = fNpoints * sizeof(Double_t);
191 memcpy(fX, x, n);
192 memcpy(fY, y, n);
193}
194
195////////////////////////////////////////////////////////////////////////////////
196/// Copy constructor for this graph
197
200{
204 else fFunctions = new TList;
205 if (gr.fHistogram) {
207 fHistogram->SetDirectory(nullptr);
208 } else {
209 fHistogram = nullptr;
210 }
213 if (!fMaxSize) {
214 fX = fY = nullptr;
215 return;
216 } else {
217 fX = new Double_t[fMaxSize];
218 fY = new Double_t[fMaxSize];
219 }
220
221 Int_t n = gr.GetN() * sizeof(Double_t);
222 memcpy(fX, gr.fX, n);
223 memcpy(fY, gr.fY, n);
224}
225
226////////////////////////////////////////////////////////////////////////////////
227/// Equal operator for this graph
228
230{
231 if (this != &gr) {
236
239
240 // delete list of functions and their contents before copying it
241 if (fFunctions) {
242 // delete previous lists of functions
243 if (!fFunctions->IsEmpty()) {
245 // use TList::Remove to take into account the case the same object is
246 // added multiple times in the list
247 TObject *obj;
248 while ((obj = fFunctions->First())) {
249 while (fFunctions->Remove(obj)) { }
250 delete obj;
251 }
252 }
253 delete fFunctions;
254 }
255
257 else fFunctions = new TList;
258
259 if (fHistogram) delete fHistogram;
260 if (gr.fHistogram) {
261 fHistogram = new TH1F(*(gr.fHistogram));
262 fHistogram->SetDirectory(nullptr);
263 } else {
264 fHistogram = nullptr;
265 }
266
269 if (fX) delete [] fX;
270 if (fY) delete [] fY;
271 if (!fMaxSize) {
272 fX = fY = nullptr;
273 return *this;
274 } else {
275 fX = new Double_t[fMaxSize];
276 fY = new Double_t[fMaxSize];
277 }
278
279 Int_t n = gr.GetN() * sizeof(Double_t);
280 if (n > 0) {
281 memcpy(fX, gr.fX, n);
282 memcpy(fY, gr.fY, n);
283 }
284 }
285 return *this;
286}
287
288////////////////////////////////////////////////////////////////////////////////
289/// Graph constructor with two vectors of floats in input
290/// A graph is build with the X coordinates taken from vx and Y coord from vy
291/// The number of points in the graph is the minimum of number of points
292/// in vx and vy.
293
294TGraph::TGraph(const TVectorF &vx, const TVectorF &vy)
295 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
296{
297 fNpoints = TMath::Min(vx.GetNrows(), vy.GetNrows());
298 if (!CtorAllocate()) return;
299 Int_t ivxlow = vx.GetLwb();
300 Int_t ivylow = vy.GetLwb();
301 for (Int_t i = 0; i < fNpoints; i++) {
302 fX[i] = vx(i + ivxlow);
303 fY[i] = vy(i + ivylow);
304 }
305}
306
307////////////////////////////////////////////////////////////////////////////////
308/// Graph constructor with two vectors of doubles in input
309/// A graph is build with the X coordinates taken from vx and Y coord from vy
310/// The number of points in the graph is the minimum of number of points
311/// in vx and vy.
312
313TGraph::TGraph(const TVectorD &vx, const TVectorD &vy)
314 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
315{
316 fNpoints = TMath::Min(vx.GetNrows(), vy.GetNrows());
317 if (!CtorAllocate()) return;
318 Int_t ivxlow = vx.GetLwb();
319 Int_t ivylow = vy.GetLwb();
320 for (Int_t i = 0; i < fNpoints; i++) {
321 fX[i] = vx(i + ivxlow);
322 fY[i] = vy(i + ivylow);
323 }
324}
325
326////////////////////////////////////////////////////////////////////////////////
327/// Graph constructor importing its parameters from the TH1 object passed as argument
328
330 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
331{
332 if (!h) {
333 Error("TGraph", "Pointer to histogram is null");
334 fNpoints = 0;
335 return;
336 }
337 if (h->GetDimension() != 1) {
338 Error("TGraph", "Histogram must be 1-D; h %s is %d-D", h->GetName(), h->GetDimension());
339 fNpoints = 0;
340 } else {
341 fNpoints = ((TH1*)h)->GetXaxis()->GetNbins();
342 }
343
344 if (!CtorAllocate()) return;
345
346 TAxis *xaxis = ((TH1*)h)->GetXaxis();
347 for (Int_t i = 0; i < fNpoints; i++) {
348 fX[i] = xaxis->GetBinCenter(i + 1);
349 fY[i] = h->GetBinContent(i + 1);
350 }
351 h->TAttLine::Copy(*this);
352 h->TAttFill::Copy(*this);
353 h->TAttMarker::Copy(*this);
354
355 std::string gname = "Graph_from_" + std::string(h->GetName());
356 SetName(gname.c_str());
357 SetTitle(h->GetTitle());
358}
359
360////////////////////////////////////////////////////////////////////////////////
361/// Graph constructor importing its parameters from the TF1 object passed as argument
362/// - if option =="" (default), a TGraph is created with points computed
363/// at the fNpx points of f.
364/// - if option =="d", a TGraph is created with points computed with the derivatives
365/// at the fNpx points of f.
366/// - if option =="i", a TGraph is created with points computed with the integral
367/// at the fNpx points of f.
368/// - if option =="I", a TGraph is created with points computed with the integral
369/// at the fNpx+1 points of f and the integral is normalized to 1.
370
372 : TNamed("Graph", "Graph"), TAttLine(), TAttFill(0, 1000), TAttMarker()
373{
374 char coption = ' ';
375 if (!f) {
376 Error("TGraph", "Pointer to function is null");
377 fNpoints = 0;
378 } else {
379 fNpoints = f->GetNpx();
380 if (option) coption = *option;
381 if (coption == 'i' || coption == 'I') fNpoints++;
382 }
383 if (!CtorAllocate()) return;
384
385 Double_t xmin = f->GetXmin();
386 Double_t xmax = f->GetXmax();
387 Double_t dx = (xmax - xmin) / fNpoints;
388 Double_t integ = 0;
389 Int_t i;
390 for (i = 0; i < fNpoints; i++) {
391 if (coption == 'i' || coption == 'I') {
392 fX[i] = xmin + i * dx;
393 if (i == 0) fY[i] = 0;
394 else fY[i] = integ + ((TF1*)f)->Integral(fX[i] - dx, fX[i]);
395 integ = fY[i];
396 } else if (coption == 'd' || coption == 'D') {
397 fX[i] = xmin + (i + 0.5) * dx;
398 fY[i] = ((TF1*)f)->Derivative(fX[i]);
399 } else {
400 fX[i] = xmin + (i + 0.5) * dx;
401 fY[i] = ((TF1*)f)->Eval(fX[i]);
402 }
403 }
404 if (integ != 0 && coption == 'I') {
405 for (i = 1; i < fNpoints; i++) fY[i] /= integ;
406 }
407
408 f->TAttLine::Copy(*this);
409 f->TAttFill::Copy(*this);
410 f->TAttMarker::Copy(*this);
411
412 SetName(f->GetName());
413 SetTitle(f->GetTitle());
414}
415
416////////////////////////////////////////////////////////////////////////////////
417/// Graph constructor reading input from filename.
418///
419/// `filename` is assumed to contain at least two columns of numbers.
420/// The string format is by default `"%lg %lg"`.
421/// This is a standard c formatting for `scanf()`.
422///
423/// If columns of numbers should be skipped, a `"%*lg"` or `"%*s"` for each column
424/// can be added, e.g. `"%lg %*lg %lg"` would read x-values from the first and
425/// y-values from the third column.
426///
427/// For files separated by a specific delimiter different from ' ' and '\\t' (e.g.
428/// ';' in csv files) you can avoid using `%*s` to bypass this delimiter by explicitly
429/// specify the `option` argument,
430/// e.g. option=`" \\t,;"` for columns of figures separated by any of these characters
431/// (' ', '\\t', ',', ';')
432/// used once (e.g. `"1;1"`) or in a combined way (`" 1;,;; 1"`).
433/// Note in that case, the instantiation is about two times slower.
434
435TGraph::TGraph(const char *filename, const char *format, Option_t *option)
436 : TNamed("Graph", filename), TAttLine(), TAttFill(0, 1000), TAttMarker()
437{
438 Double_t x, y;
439 TString fname = filename;
440 gSystem->ExpandPathName(fname);
441
442 std::ifstream infile(fname.Data());
443 if (!infile.good()) {
444 MakeZombie();
445 Error("TGraph", "Cannot open file: %s, TGraph is Zombie", filename);
446 fNpoints = 0;
447 return;
448 } else {
449 fNpoints = 100; //initial number of points
450 }
451 if (!CtorAllocate()) return;
452 std::string line;
453 Int_t np = 0;
454
455 // No delimiters specified (standard constructor).
456 if (strcmp(option, "") == 0) {
457
458 while (std::getline(infile, line, '\n')) {
459 if (2 != sscanf(line.c_str(), format, &x, &y)) {
460 continue; //skip empty and ill-formed lines
461 }
462 SetPoint(np, x, y);
463 np++;
464 }
465 Set(np);
466
467 // A delimiter has been specified in "option"
468 } else {
469
470 // Checking format and creating its boolean counterpart
471 TString format_ = TString(format) ;
472 format_.ReplaceAll(" ", "") ;
473 format_.ReplaceAll("\t", "") ;
474 format_.ReplaceAll("lg", "") ;
475 format_.ReplaceAll("s", "") ;
476 format_.ReplaceAll("%*", "0") ;
477 format_.ReplaceAll("%", "1") ;
478 if (!format_.IsDigit()) {
479 Error("TGraph", "Incorrect input format! Allowed formats are {\"%%lg\",\"%%*lg\" or \"%%*s\"}");
480 return;
481 }
482 Int_t ntokens = format_.Length() ;
483 if (ntokens < 2) {
484 Error("TGraph", "Incorrect input format! Only %d tag(s) in format whereas 2 \"%%lg\" tags are expected!", ntokens);
485 return;
486 }
487 Int_t ntokensToBeSaved = 0 ;
488 Bool_t * isTokenToBeSaved = new Bool_t [ntokens] ;
489 for (Int_t idx = 0; idx < ntokens; idx++) {
490 isTokenToBeSaved[idx] = TString::Format("%c", format_[idx]).Atoi() ; //atoi(&format_[idx]) does not work for some reason...
491 if (isTokenToBeSaved[idx] == 1) {
492 ntokensToBeSaved++ ;
493 }
494 }
495 if (ntokens >= 2 && ntokensToBeSaved != 2) { //first condition not to repeat the previous error message
496 Error("TGraph", "Incorrect input format! There are %d \"%%lg\" tag(s) in format whereas 2 and only 2 are expected!", ntokensToBeSaved);
497 delete [] isTokenToBeSaved ;
498 return;
499 }
500
501 // Initializing loop variables
502 Bool_t isLineToBeSkipped = kFALSE ; //empty and ill-formed lines
503 char * token = NULL ;
504 TString token_str = "" ;
505 Int_t token_idx = 0 ;
506 Double_t * value = new Double_t [2] ; //x,y buffers
507 Int_t value_idx = 0 ;
508
509 // Looping
510 char *rest;
511 while (std::getline(infile, line, '\n')) {
512 if (line != "") {
513 if (line[line.size() - 1] == char(13)) { // removing DOS CR character
514 line.erase(line.end() - 1, line.end()) ;
515 }
516 //token = R__STRTOK_R(const_cast<char *>(line.c_str()), option, rest);
517 token = R__STRTOK_R(const_cast<char *>(line.c_str()), option, &rest);
518 while (token != NULL && value_idx < 2) {
519 if (isTokenToBeSaved[token_idx]) {
520 token_str = TString(token) ;
521 token_str.ReplaceAll("\t", "") ;
522 if (!token_str.IsFloat()) {
523 isLineToBeSkipped = kTRUE ;
524 break ;
525 } else {
526 value[value_idx] = token_str.Atof() ;
527 value_idx++ ;
528 }
529 }
530 token = R__STRTOK_R(NULL, option, &rest); // next token
531 token_idx++ ;
532 }
533 if (!isLineToBeSkipped && value_idx == 2) {
534 x = value[0] ;
535 y = value[1] ;
536 SetPoint(np, x, y) ;
537 np++ ;
538 }
539 }
540 isLineToBeSkipped = kFALSE ;
541 token = NULL ;
542 token_idx = 0 ;
543 value_idx = 0 ;
544 }
545 Set(np) ;
546
547 // Cleaning
548 delete [] isTokenToBeSaved ;
549 delete [] value ;
550 delete token ;
551 }
552 infile.close();
553}
554
555////////////////////////////////////////////////////////////////////////////////
556/// Graph default destructor.
557
559{
560 delete [] fX;
561 delete [] fY;
562 if (fFunctions) {
564 //special logic to support the case where the same object is
565 //added multiple times in fFunctions.
566 //This case happens when the same object is added with different
567 //drawing modes
568 TObject *obj;
569 while ((obj = fFunctions->First())) {
570 while (fFunctions->Remove(obj)) { }
571 delete obj;
572 }
573 delete fFunctions;
574 fFunctions = nullptr; //to avoid accessing a deleted object in RecursiveRemove
575 }
576 delete fHistogram;
577}
578
579////////////////////////////////////////////////////////////////////////////////
580/// Allocate internal data structures for `newsize` points.
581
583{
584 return AllocateArrays(2, newsize);
585}
586
587////////////////////////////////////////////////////////////////////////////////
588/// Allocate arrays.
589
591{
592 if (arraySize < 0) {
593 arraySize = 0;
594 }
595 Double_t **newarrays = new Double_t*[Narrays];
596 if (!arraySize) {
597 for (Int_t i = 0; i < Narrays; ++i)
598 newarrays[i] = 0;
599 } else {
600 for (Int_t i = 0; i < Narrays; ++i)
601 newarrays[i] = new Double_t[arraySize];
602 }
603 fMaxSize = arraySize;
604 return newarrays;
605}
606
607////////////////////////////////////////////////////////////////////////////////
608/// Apply function f to all the data points
609/// f may be a 1-D function TF1 or 2-d function TF2
610/// The Y values of the graph are replaced by the new values computed
611/// using the function
612
614{
616
617 for (Int_t i = 0; i < fNpoints; i++) {
618 fY[i] = f->Eval(fX[i], fY[i]);
619 }
620 if (gPad) gPad->Modified();
621}
622
623////////////////////////////////////////////////////////////////////////////////
624/// Browse
625
627{
628 TString opt = gEnv->GetValue("TGraph.BrowseOption", "");
629 if (opt.IsNull()) {
630 opt = b ? b->GetDrawOption() : "alp";
631 opt = (opt == "") ? "alp" : opt.Data();
632 }
633 Draw(opt.Data());
634 gPad->Update();
635}
636
637////////////////////////////////////////////////////////////////////////////////
638/// Return the chisquare of this graph with respect to f1.
639/// The chisquare is computed as the sum of the quantity below at each point:
640/// \f[
641/// \frac{(y-f1(x))^{2}}{ey^{2}+(\frac{1}{2}(exl+exh)f1'(x))^{2}}
642/// \f]
643/// where x and y are the graph point coordinates and f1'(x) is the derivative of function f1(x).
644/// This method to approximate the uncertainty in y because of the errors in x, is called
645/// "effective variance" method.
646/// In case of a pure TGraph, the denominator is 1.
647/// In case of a TGraphErrors or TGraphAsymmErrors the errors are taken
648/// into account.
649/// By default the range of the graph is used whatever function range.
650/// Use option "R" to use the function range
651
653{
654 if (!func) {
655 Error("Chisquare","Function pointer is Null - return -1");
656 return -1;
657 }
658
659 TString opt(option); opt.ToUpper();
660 bool useRange = opt.Contains("R");
661
662 return ROOT::Fit::Chisquare(*this, *func,useRange);
663}
664
665////////////////////////////////////////////////////////////////////////////////
666/// Return kTRUE if point number "left"'s argument (angle with respect to positive
667/// x-axis) is bigger than that of point number "right". Can be used by Sort.
668
670{
671 Double_t xl = 0, yl = 0, xr = 0, yr = 0;
672 gr->GetPoint(left, xl, yl);
673 gr->GetPoint(right, xr, yr);
674 return (TMath::ATan2(yl, xl) > TMath::ATan2(yr, xr));
675}
676
677////////////////////////////////////////////////////////////////////////////////
678/// Return kTRUE if fX[left] > fX[right]. Can be used by Sort.
679
681{
682 return gr->fX[left] > gr->fX[right];
683}
684
685////////////////////////////////////////////////////////////////////////////////
686/// Return kTRUE if fY[left] > fY[right]. Can be used by Sort.
687
689{
690 return gr->fY[left] > gr->fY[right];
691}
692
693////////////////////////////////////////////////////////////////////////////////
694/// Return kTRUE if point number "left"'s distance to origin is bigger than
695/// that of point number "right". Can be used by Sort.
696
698{
699 return gr->fX[left] * gr->fX[left] + gr->fY[left] * gr->fY[left]
700 > gr->fX[right] * gr->fX[right] + gr->fY[right] * gr->fY[right];
701}
702
703////////////////////////////////////////////////////////////////////////////////
704/// Compute the x/y range of the points in this graph
705
707{
708 if (fNpoints <= 0) {
709 xmin = xmax = ymin = ymax = 0;
710 return;
711 }
712 xmin = xmax = fX[0];
713 ymin = ymax = fY[0];
714
715 Double_t xminl = 0; // Positive minimum. Used in case of log scale along X axis.
716 Double_t yminl = 0; // Positive minimum. Used in case of log scale along Y axis.
717
718 for (Int_t i = 1; i < fNpoints; i++) {
719 if (fX[i] < xmin) xmin = fX[i];
720 if (fX[i] > xmax) xmax = fX[i];
721 if (fY[i] < ymin) ymin = fY[i];
722 if (fY[i] > ymax) ymax = fY[i];
723 if (ymin>0 && (yminl==0 || ymin<yminl)) yminl = ymin;
724 if (xmin>0 && (xminl==0 || xmin<xminl)) xminl = xmin;
725 }
726
727 if (gPad && gPad->GetLogy() && yminl>0) ymin = yminl;
728 if (gPad && gPad->GetLogx() && xminl>0) xmin = xminl;
729}
730
731////////////////////////////////////////////////////////////////////////////////
732/// Copy points from fX and fY to arrays[0] and arrays[1]
733/// or to fX and fY if arrays == 0 and ibegin != iend.
734/// If newarrays is non null, replace fX, fY with pointers from newarrays[0,1].
735/// Delete newarrays, old fX and fY
736
737void TGraph::CopyAndRelease(Double_t **newarrays, Int_t ibegin, Int_t iend,
738 Int_t obegin)
739{
740 CopyPoints(newarrays, ibegin, iend, obegin);
741 if (newarrays) {
742 delete[] fX;
743 fX = newarrays[0];
744 delete[] fY;
745 fY = newarrays[1];
746 delete[] newarrays;
747 }
748}
749
750////////////////////////////////////////////////////////////////////////////////
751/// Copy points from fX and fY to arrays[0] and arrays[1]
752/// or to fX and fY if arrays == 0 and ibegin != iend.
753
755 Int_t obegin)
756{
757 if (ibegin < 0 || iend <= ibegin || obegin < 0) { // Error;
758 return kFALSE;
759 }
760 if (!arrays && ibegin == obegin) { // No copying is needed
761 return kFALSE;
762 }
763 Int_t n = (iend - ibegin) * sizeof(Double_t);
764 if (arrays) {
765 memmove(&arrays[0][obegin], &fX[ibegin], n);
766 memmove(&arrays[1][obegin], &fY[ibegin], n);
767 } else {
768 memmove(&fX[obegin], &fX[ibegin], n);
769 memmove(&fY[obegin], &fY[ibegin], n);
770 }
771 return kTRUE;
772}
773
774////////////////////////////////////////////////////////////////////////////////
775/// In constructors set fNpoints than call this method.
776/// Return kFALSE if the graph will contain no points.
777///Note: This function should be called only from the constructor
778/// since it does not delete previously existing arrays
779
781{
782 fHistogram = nullptr;
783 fMaximum = -1111;
784 fMinimum = -1111;
786 fFunctions = new TList;
787 if (fNpoints <= 0) {
788 fNpoints = 0;
789 fMaxSize = 0;
790 fX = nullptr;
791 fY = nullptr;
792 return kFALSE;
793 } else {
795 fX = new Double_t[fMaxSize];
796 fY = new Double_t[fMaxSize];
797 }
798 return kTRUE;
799}
800
801////////////////////////////////////////////////////////////////////////////////
802/// Draw this graph with its current attributes.
803///
804/// The options to draw a graph are described in TGraphPainter class.
805
807{
808 TString opt = option;
809 opt.ToLower();
810
811 if (opt.Contains("same")) {
812 opt.ReplaceAll("same", "");
813 }
814
815 // in case of option *, set marker style to 3 (star) and replace
816 // * option by option P.
817 Ssiz_t pos;
818 if ((pos = opt.Index("*")) != kNPOS) {
820 opt.Replace(pos, 1, "p");
821 }
822
823 // If no option is specified, it is defined as "alp" in case there
824 // no current pad or if the current pad as no axis defined.
825 if (!option || !strlen(option)) {
826 if (gPad) {
827 if (!gPad->GetListOfPrimitives()->FindObject("TFrame")) opt = "alp";
828 } else {
829 opt = "alp";
830 }
831 }
832
833 if (gPad) {
834 if (!gPad->IsEditable()) gROOT->MakeDefCanvas();
835 if (opt.Contains("a")) gPad->Clear();
836 }
837
838 AppendPad(opt);
839
840 gPad->IncrementPaletteColor(1, opt);
841
842}
843
844////////////////////////////////////////////////////////////////////////////////
845/// Compute distance from point px,py to a graph.
846///
847/// Compute the closest distance of approach from point px,py to this line.
848/// The distance is computed in pixels units.
849
851{
853 if (painter) return painter->DistancetoPrimitiveHelper(this, px, py);
854 else return 0;
855}
856
857////////////////////////////////////////////////////////////////////////////////
858/// Draw this graph with new attributes.
859
861{
862 TGraph *newgraph = new TGraph(n, x, y);
863 TAttLine::Copy(*newgraph);
864 TAttFill::Copy(*newgraph);
865 TAttMarker::Copy(*newgraph);
866 newgraph->SetBit(kCanDelete);
867 newgraph->AppendPad(option);
868}
869
870////////////////////////////////////////////////////////////////////////////////
871/// Draw this graph with new attributes.
872
874{
875 TGraph *newgraph = new TGraph(n, x, y);
876 TAttLine::Copy(*newgraph);
877 TAttFill::Copy(*newgraph);
878 TAttMarker::Copy(*newgraph);
879 newgraph->SetBit(kCanDelete);
880 newgraph->AppendPad(option);
881}
882
883////////////////////////////////////////////////////////////////////////////////
884/// Draw this graph with new attributes.
885
887{
888 const Double_t *xx = x;
889 const Double_t *yy = y;
890 if (!xx) xx = fX;
891 if (!yy) yy = fY;
892 TGraph *newgraph = new TGraph(n, xx, yy);
893 TAttLine::Copy(*newgraph);
894 TAttFill::Copy(*newgraph);
895 TAttMarker::Copy(*newgraph);
896 newgraph->SetBit(kCanDelete);
897 newgraph->AppendPad(option);
898}
899
900////////////////////////////////////////////////////////////////////////////////
901/// Display a panel with all graph drawing options.
902
904{
906 if (painter) painter->DrawPanelHelper(this);
907}
908
909////////////////////////////////////////////////////////////////////////////////
910/// Interpolate points in this graph at x using a TSpline.
911///
912/// - if spline==0 and option="" a linear interpolation between the two points
913/// close to x is computed. If x is outside the graph range, a linear
914/// extrapolation is computed.
915/// - if spline==0 and option="S" a TSpline3 object is created using this graph
916/// and the interpolated value from the spline is returned.
917/// the internally created spline is deleted on return.
918/// - if spline is specified, it is used to return the interpolated value.
919///
920/// If the points are sorted in X a binary search is used (significantly faster)
921/// One needs to set the bit TGraph::SetBit(TGraph::kIsSortedX) before calling
922/// TGraph::Eval to indicate that the graph is sorted in X.
923
925{
926
927 if (spline) {
928 //spline interpolation using the input spline
929 return spline->Eval(x);
930 }
931
932 if (fNpoints == 0) return 0;
933 if (fNpoints == 1) return fY[0];
934
935 if (option && *option) {
936 TString opt = option;
937 opt.ToLower();
938 // create a TSpline every time when using option "s" and no spline pointer is given
939 if (opt.Contains("s")) {
940
941 // points must be sorted before using a TSpline
942 std::vector<Double_t> xsort(fNpoints);
943 std::vector<Double_t> ysort(fNpoints);
944 std::vector<Int_t> indxsort(fNpoints);
945 TMath::Sort(fNpoints, fX, &indxsort[0], false);
946 for (Int_t i = 0; i < fNpoints; ++i) {
947 xsort[i] = fX[ indxsort[i] ];
948 ysort[i] = fY[ indxsort[i] ];
949 }
950
951 // spline interpolation creating a new spline
952 TSpline3 s("", &xsort[0], &ysort[0], fNpoints);
953 Double_t result = s.Eval(x);
954 return result;
955 }
956 }
957 //linear interpolation
958 //In case x is < fX[0] or > fX[fNpoints-1] return the extrapolated point
959
960 //find points in graph around x assuming points are not sorted
961 // (if point are sorted use a binary search)
962 Int_t low = -1;
963 Int_t up = -1;
966 if (low == -1) {
967 // use first two points for doing an extrapolation
968 low = 0;
969 }
970 if (fX[low] == x) return fY[low];
971 if (low == fNpoints-1) low--; // for extrapolating
972 up = low+1;
973 }
974 else {
975 // case TGraph is not sorted
976
977 // find neighbours simply looping all points
978 // and find also the 2 adjacent points: (low2 < low < x < up < up2 )
979 // needed in case x is outside the graph ascissa interval
980 Int_t low2 = -1;
981 Int_t up2 = -1;
982
983 for (Int_t i = 0; i < fNpoints; ++i) {
984 if (fX[i] < x) {
985 if (low == -1 || fX[i] > fX[low]) {
986 low2 = low;
987 low = i;
988 } else if (low2 == -1) low2 = i;
989 } else if (fX[i] > x) {
990 if (up == -1 || fX[i] < fX[up]) {
991 up2 = up;
992 up = i;
993 } else if (up2 == -1) up2 = i;
994 } else // case x == fX[i]
995 return fY[i]; // no interpolation needed
996 }
997
998 // treat cases when x is outside graph min max abscissa
999 if (up == -1) {
1000 up = low;
1001 low = low2;
1002 }
1003 if (low == -1) {
1004 low = up;
1005 up = up2;
1006 }
1007 }
1008 // do now the linear interpolation
1009 assert(low != -1 && up != -1);
1010
1011 if (fX[low] == fX[up]) return fY[low];
1012 Double_t yn = fY[up] + (x - fX[up]) * (fY[low] - fY[up]) / (fX[low] - fX[up]);
1013 return yn;
1014}
1015
1016////////////////////////////////////////////////////////////////////////////////
1017/// Execute action corresponding to one event.
1018///
1019/// This member function is called when a graph is clicked with the locator
1020///
1021/// If Left button clicked on one of the line end points, this point
1022/// follows the cursor until button is released.
1023///
1024/// if Middle button clicked, the line is moved parallel to itself
1025/// until the button is released.
1026
1028{
1030 if (painter) painter->ExecuteEventHelper(this, event, px, py);
1031}
1032
1033////////////////////////////////////////////////////////////////////////////////
1034/// If array sizes <= newsize, expand storage to 2*newsize.
1035
1037{
1038 Double_t **ps = ExpandAndCopy(newsize, fNpoints);
1039 CopyAndRelease(ps, 0, 0, 0);
1040}
1041
1042////////////////////////////////////////////////////////////////////////////////
1043/// If graph capacity is less than newsize points then make array sizes
1044/// equal to least multiple of step to contain newsize points.
1045
1046void TGraph::Expand(Int_t newsize, Int_t step)
1047{
1048 if (newsize <= fMaxSize) {
1049 return;
1050 }
1051 Double_t **ps = Allocate(step * (newsize / step + (newsize % step ? 1 : 0)));
1052 CopyAndRelease(ps, 0, fNpoints, 0);
1053}
1054
1055////////////////////////////////////////////////////////////////////////////////
1056/// if size > fMaxSize allocate new arrays of 2*size points and copy iend first
1057/// points.
1058/// Return pointer to new arrays.
1059
1061{
1062 if (size <= fMaxSize) {
1063 return 0;
1064 }
1065 Double_t **newarrays = Allocate(2 * size);
1066 CopyPoints(newarrays, 0, iend, 0);
1067 return newarrays;
1068}
1069
1070////////////////////////////////////////////////////////////////////////////////
1071/// Set zero values for point arrays in the range [begin, end)
1072/// Should be redefined in descendant classes
1073
1075{
1076 memset(fX + begin, 0, (end - begin)*sizeof(Double_t));
1077 memset(fY + begin, 0, (end - begin)*sizeof(Double_t));
1078}
1079
1080////////////////////////////////////////////////////////////////////////////////
1081/// Search object named name in the list of functions
1082
1084{
1085 if (fFunctions) return fFunctions->FindObject(name);
1086 return 0;
1087}
1088
1089////////////////////////////////////////////////////////////////////////////////
1090/// Search object obj in the list of functions
1091
1093{
1094 if (fFunctions) return fFunctions->FindObject(obj);
1095 return 0;
1096}
1097
1098////////////////////////////////////////////////////////////////////////////////
1099/// Fit this graph with function f1.
1100///
1101/// \param[in] f1 pointer to the function object
1102/// \param[in] option string defining the fit options (see table below).
1103/// \param[in] goption specify a list of graphics options. See TGraph::Draw and TGraphPainter for a complete list of these possible options.
1104/// \param[in] rxmin lower fitting range
1105/// \param[in] rxmax upper fitting range
1106///
1107/// \anchor GFitOpt
1108/// ### Graph Fitting Options
1109/// The list of fit options is given in parameter option.
1110///
1111/// option | description
1112/// -------|------------
1113/// "S" | The full result of the fit is returned in the `TFitResultPtr`. This is needed to get the covariance matrix of the fit. See `TFitResult` and the base class `ROOT::Math::FitResult`.
1114/// "W" | Ignore all point errors when fitting a TGraphErrors or TGraphAsymmErrors
1115/// "F" | Uses the default minimizer (e.g. Minuit) when fitting a linear function (e.g. polN) instead of the linear fitter.
1116/// "U" | Uses a user specified objective function (e.g. user providedlikelihood function) defined using `TVirtualFitter::SetFCN`
1117/// "E" | Performs a better parameter errors estimation using the Minos technique for all fit parameters.
1118/// "M" | Uses the IMPROVE algorithm (available only in TMinuit). This algorithm attempts improve the found local minimum by searching for a better one.
1119/// "Q" | Quiet mode (minimum printing)
1120/// "V" | Verbose mode (default is between Q and V)
1121/// "+" | Adds this new fitted function to the list of fitted functions. By default, the previous function is deleted and only the last one is kept.
1122/// "N" | Does not store the graphics function, does not draw the histogram with the function after fitting.
1123/// "0" | Does not draw the histogram and the fitted function after fitting, but in contrast to option "N", it stores the fitted function in the histogram list of functions.
1124/// "R" | Fit using a fitting range specified in the function range with `TF1::SetRange`.
1125/// "B" | Use this option when you want to fix one or more parameters and the fitting function is a predefined one (e.g gaus, expo,..), otherwise in case of pre-defined functions, some default initial values and limits are set.
1126/// "C" | In case of linear fitting, do no calculate the chisquare (saves CPU time).
1127/// "G" | Uses the gradient implemented in `TF1::GradientPar` for the minimization. This allows to use Automatic Differentiation when it is supported by the provided TF1 function.
1128/// "EX0" | When fitting a TGraphErrors or TGraphAsymErrors do not consider errors in the X coordinates
1129/// "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
1130///
1131///
1132/// This function is used for fitting also the derived TGraph classes such as TGraphErrors or TGraphAsymmErrors.
1133/// See the note below on how the errors are used when fitting a TGraphErrors or TGraphAsymmErrors.
1134///
1135/// The fitting of the TGraph, i.e simple data points without any error associated, is performed using the
1136/// un-weighted least-square (chi-square) method.
1137///
1138///
1139///\anchor GFitErrors
1140/// ### TGraphErrors fit:
1141///
1142/// In case of a TGraphErrors or TGraphAsymmErrors object, when `x` errors are present, the error along x,
1143/// is projected along the y-direction by calculating the function at the points `x-ex_low` and
1144/// `x+ex_high`, where `ex_low` and `ex_high` are the corresponding lower and upper error in x.
1145/// The chi-square is then computed as the sum of the quantity below at each data point:
1146///
1147/// \f[
1148/// \frac{(y-f(x))^{2}}{ey^{2}+(\frac{1}{2}(exl+exh)f'(x))^{2}}
1149/// \f]
1150///
1151/// where `x` and `y` are the point coordinates, and `f'(x)` is the derivative of the
1152/// function `f(x)`.
1153///
1154/// In case of asymmetric errors, if the function lies below (above) the data point, `ey` is `ey_low` (`ey_high`).
1155///
1156/// The approach used to approximate the uncertainty in y because of the
1157/// errors in x is to make it equal the error in x times the slope of the line.
1158/// This approach is called "effective variance method" and
1159/// the implementation is provided in the function FitUtil::EvaluateChi2Effective
1160///
1161/// \anchor GFitLinear
1162/// ### Linear fitting:
1163/// When the fitting function is linear (contains the `++` sign) or the fitting
1164/// function is a polynomial, a linear fitter is initialised.
1165/// To create a linear function, use the following syntax: linear parts
1166/// separated by `++` sign.
1167/// Example: to fit the parameters of the function `p0*x + p1*sin(x)`, you can create a
1168/// TF1 object as
1169///
1170/// TF1 *f1 = new TF1("f1", "x++sin(x)", xmin, xmax);
1171///
1172/// For such a TF1 you don't have to set the initial conditions and the linear fitter is used.
1173/// Going via the linear fitter for functions, linear in parameters, gives a
1174/// considerable advantage in speed.
1175/// When using the linear fitting it is also possible to perform a robust fitting with the
1176/// Least Trimmed Square (LTS) regression algorithm, by using the fit option `ROB`.
1177/// See the tutorial `fitLinearRobust.C`.
1178///
1179/// ### Notes on TGraph/TGraphErrors Fitting:
1180///
1181/// 1. By using the "effective variance" method a simple linear regression
1182/// becomes a non-linear case, which takes several iterations
1183/// instead of 0 as in the linear case.
1184/// 2. The effective variance technique assumes that there is no correlation
1185/// between the x and y coordinate.
1186/// 3. The standard chi2 (least square) method without error in the coordinates (x) can
1187/// be forced by using option "EX0"
1188/// 4. The linear fitter doesn't take into account the errors in x. When fitting a
1189/// TGraphErrors with a linear functions the errors in x will not be considered.
1190/// If errors in x are important, use option "F" for linear function fitting.
1191/// 5. When fitting a TGraph (i.e. no errors associated with each point),
1192/// a correction is applied to the errors on the parameters with the following
1193/// formula:
1194/// `parameter_error *= sqrt(chisquare/(ndf-1))`
1195///
1196/// ### General Fitting documentation
1197///
1198/// See in TH1::Fit for the documentation of
1199/// - [Fit Result](\ref HFitRes)
1200/// - [Fit Status](\ref HFitStatus)
1201/// - [Fit Statistics Box](\ref HFitStatBox)
1202/// - [Fitting in a Range](\ref HFitRange)
1203/// - [Setting Initial Conditions](\ref HFitInitial)
1204
1206{
1207 Foption_t fitOption;
1209 // create range and minimizer options with default values
1210 ROOT::Fit::DataRange range(rxmin, rxmax);
1212 return ROOT::Fit::FitObject(this, f1 , fitOption , minOption, goption, range);
1213}
1214
1215////////////////////////////////////////////////////////////////////////////////
1216/// Fit this graph with function with name `fname`.
1217///
1218/// This is a different interface to TGraph fitting using TGraph::Fit(TF1 *f1,Option_t *, Option_t *, Axis_t, Axis_t)
1219/// See there for the details about fitting a TGraph.
1220///
1221/// The parameter `fname` is the name of an already predefined function created by TF1 or TF2
1222/// Predefined functions such as gaus, expo and poln are automatically
1223/// created by ROOT.
1224///
1225/// The parameter `fname` can also be a formula, accepted by the linear fitter (linear parts divided
1226/// by "++" sign), for example "x++sin(x)" for fitting "[0]*x+[1]*sin(x)"
1227
1229{
1230 const char *linear = fname ? strstr(fname, "++") : nullptr;
1231 if (linear) {
1232 TF1 f1(fname, fname, xmin, xmax);
1233 return Fit(&f1, option, "", xmin, xmax);
1234 }
1235 TF1 * f1 = (TF1*)gROOT->GetFunction(fname);
1236 if (!f1) {
1237 Printf("Unknown function: %s", fname);
1238 return -1;
1239 }
1240 return Fit(f1, option, "", xmin, xmax);
1241}
1242
1243////////////////////////////////////////////////////////////////////////////////
1244/// Display a GUI panel with all graph fit options.
1245///
1246/// See class TFitEditor for example
1247
1249{
1250 if (!gPad)
1251 gROOT->MakeDefCanvas();
1252
1253 if (!gPad) {
1254 Error("FitPanel", "Unable to create a default canvas");
1255 return;
1256 }
1257
1258 // use plugin manager to create instance of TFitEditor
1259 TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
1260 if (handler && handler->LoadPlugin() != -1) {
1261 if (handler->ExecPlugin(2, gPad, this) == 0)
1262 Error("FitPanel", "Unable to crate the FitPanel");
1263 } else
1264 Error("FitPanel", "Unable to find the FitPanel plug-in");
1265}
1266
1267////////////////////////////////////////////////////////////////////////////////
1268/// Return graph correlation factor
1269
1271{
1272 Double_t rms1 = GetRMS(1);
1273 if (rms1 == 0) return 0;
1274 Double_t rms2 = GetRMS(2);
1275 if (rms2 == 0) return 0;
1276 return GetCovariance() / rms1 / rms2;
1277}
1278
1279////////////////////////////////////////////////////////////////////////////////
1280/// Return covariance of vectors x,y
1281
1283{
1284 if (fNpoints <= 0) return 0;
1285 Double_t sum = fNpoints, sumx = 0, sumy = 0, sumxy = 0;
1286
1287 for (Int_t i = 0; i < fNpoints; i++) {
1288 sumx += fX[i];
1289 sumy += fY[i];
1290 sumxy += fX[i] * fY[i];
1291 }
1292 return sumxy / sum - sumx / sum * sumy / sum;
1293}
1294
1295////////////////////////////////////////////////////////////////////////////////
1296/// Return mean value of X (axis=1) or Y (axis=2)
1297
1299{
1300 if (axis < 1 || axis > 2) return 0;
1301 if (fNpoints <= 0) return 0;
1302 Double_t sumx = 0;
1303 for (Int_t i = 0; i < fNpoints; i++) {
1304 if (axis == 1) sumx += fX[i];
1305 else sumx += fY[i];
1306 }
1307 return sumx / fNpoints;
1308}
1309
1310////////////////////////////////////////////////////////////////////////////////
1311/// Return RMS of X (axis=1) or Y (axis=2)
1312
1314{
1315 if (axis < 1 || axis > 2) return 0;
1316 if (fNpoints <= 0) return 0;
1317 Double_t sumx = 0, sumx2 = 0;
1318 for (Int_t i = 0; i < fNpoints; i++) {
1319 if (axis == 1) {
1320 sumx += fX[i];
1321 sumx2 += fX[i] * fX[i];
1322 } else {
1323 sumx += fY[i];
1324 sumx2 += fY[i] * fY[i];
1325 }
1326 }
1327 Double_t x = sumx / fNpoints;
1328 Double_t rms2 = TMath::Abs(sumx2 / fNpoints - x * x);
1329 return TMath::Sqrt(rms2);
1330}
1331
1332////////////////////////////////////////////////////////////////////////////////
1333/// It always returns a negative value. Real implementation in TGraphErrors
1334
1336{
1337 return -1;
1338}
1339
1340////////////////////////////////////////////////////////////////////////////////
1341/// It always returns a negative value. Real implementation in TGraphErrors
1342
1344{
1345 return -1;
1346}
1347
1348////////////////////////////////////////////////////////////////////////////////
1349/// It always returns a negative value. Real implementation in TGraphErrors
1350/// and TGraphAsymmErrors
1351
1353{
1354 return -1;
1355}
1356
1357////////////////////////////////////////////////////////////////////////////////
1358/// It always returns a negative value. Real implementation in TGraphErrors
1359/// and TGraphAsymmErrors
1360
1362{
1363 return -1;
1364}
1365
1366////////////////////////////////////////////////////////////////////////////////
1367/// It always returns a negative value. Real implementation in TGraphErrors
1368/// and TGraphAsymmErrors
1369
1371{
1372 return -1;
1373}
1374
1375////////////////////////////////////////////////////////////////////////////////
1376/// It always returns a negative value. Real implementation in TGraphErrors
1377/// and TGraphAsymmErrors
1378
1380{
1381 return -1;
1382}
1383
1384////////////////////////////////////////////////////////////////////////////////
1385/// Return pointer to function with name.
1386///
1387/// Functions such as TGraph::Fit store the fitted function in the list of
1388/// functions of this graph.
1389
1390TF1 *TGraph::GetFunction(const char *name) const
1391{
1392 if (!fFunctions) return nullptr;
1393 return (TF1*)fFunctions->FindObject(name);
1394}
1395
1396////////////////////////////////////////////////////////////////////////////////
1397/// Returns a pointer to the histogram used to draw the axis
1398/// Takes into account the two following cases.
1399/// 1. option 'A' was specified in TGraph::Draw. Return fHistogram
1400/// 2. user had called TPad::DrawFrame. return pointer to hframe histogram
1401
1403{
1404 Double_t rwxmin, rwxmax, rwymin, rwymax, maximum, minimum, dx, dy;
1405 Double_t uxmin, uxmax;
1406
1407 ComputeRange(rwxmin, rwymin, rwxmax, rwymax); //this is redefined in TGraphErrors
1408
1409 // (if fHistogram exist) && (if the log scale is on) &&
1410 // (if the computed range minimum is > 0) && (if the fHistogram minimum is zero)
1411 // then it means fHistogram limits have been computed in linear scale
1412 // therefore they might be too strict and cut some points. In that case the
1413 // fHistogram limits should be recomputed ie: the existing fHistogram
1414 // should not be returned.
1415 TH1F *historg = nullptr;
1416 if (fHistogram) {
1417 if (!TestBit(kResetHisto)) {
1418 if (gPad && gPad->GetLogx()) {
1419 if (rwxmin <= 0 || fHistogram->GetXaxis()->GetXmin() != 0) return fHistogram;
1420 } else if (gPad && gPad->GetLogy()) {
1421 if (rwymin <= 0 || fHistogram->GetMinimum() != 0) return fHistogram;
1422 } else {
1423 return fHistogram;
1424 }
1425 } else {
1426 const_cast <TGraph*>(this)->ResetBit(kResetHisto);
1427 }
1428 historg = fHistogram;
1429 }
1430
1431 if (rwxmin == rwxmax) rwxmax += 1.;
1432 if (rwymin == rwymax) rwymax += 1.;
1433 dx = 0.1 * (rwxmax - rwxmin);
1434 dy = 0.1 * (rwymax - rwymin);
1435 uxmin = rwxmin - dx;
1436 uxmax = rwxmax + dx;
1437 minimum = rwymin - dy;
1438 maximum = rwymax + dy;
1439
1440 if (fMinimum != -1111) minimum = fMinimum;
1441 if (fMaximum != -1111) maximum = fMaximum;
1442
1443 // the graph is created with at least as many channels as there are points
1444 // to permit zooming on the full range
1445 if (uxmin < 0 && rwxmin >= 0) {
1446 if (gPad && gPad->GetLogx()) uxmin = 0.9 * rwxmin;
1447 else uxmin = 0;
1448 }
1449 if (uxmax > 0 && rwxmax <= 0) {
1450 if (gPad && gPad->GetLogx()) uxmax = 1.1 * rwxmax;
1451 else uxmax = 0;
1452 }
1453
1454 if (minimum < 0 && rwymin >= 0) minimum = 0.9 * rwymin;
1455
1456 if (minimum <= 0 && gPad && gPad->GetLogy()) minimum = 0.001 * maximum;
1457 if (uxmin <= 0 && gPad && gPad->GetLogx()) {
1458 if (uxmax > 1000) uxmin = 1;
1459 else uxmin = 0.001 * uxmax;
1460 }
1461
1462 rwxmin = uxmin;
1463 rwxmax = uxmax;
1464 Int_t npt = 100;
1465 if (fNpoints > npt) npt = fNpoints;
1466 const char *gname = GetName();
1467 if (!gname[0]) gname = "Graph";
1468 // do not add the histogram to gDirectory
1469 // use local TDirectory::TContect that will set temporarly gDirectory to a nullptr and
1470 // will avoid that histogram is added in the global directory
1471 {
1472 TDirectory::TContext ctx(nullptr);
1473 ((TGraph*)this)->fHistogram = new TH1F(gname, GetTitle(), npt, rwxmin, rwxmax);
1474 }
1475 if (!fHistogram) return nullptr;
1476 fHistogram->SetMinimum(minimum);
1478 fHistogram->SetMaximum(maximum);
1479 fHistogram->GetYaxis()->SetLimits(minimum, maximum);
1480 // Restore the axis attributes if needed
1481 if (historg) {
1482 fHistogram->GetXaxis()->SetTitle(historg->GetXaxis()->GetTitle());
1495
1496 fHistogram->GetYaxis()->SetTitle(historg->GetYaxis()->GetTitle());
1509 delete historg;
1510 }
1511 return fHistogram;
1512}
1513
1514////////////////////////////////////////////////////////////////////////////////
1515/// Get x and y values for point number i.
1516/// The function returns -1 in case of an invalid request or the point number otherwise
1517
1519{
1520 if (i < 0 || i >= fNpoints || !fX || !fY) return -1;
1521 x = fX[i];
1522 y = fY[i];
1523 return i;
1524}
1525
1526////////////////////////////////////////////////////////////////////////////////
1527/// Get x value for point i.
1528
1530{
1531 if (i < 0 || i >= fNpoints || !fX)
1532 return -1.;
1533
1534 return fX[i];
1535}
1536
1537////////////////////////////////////////////////////////////////////////////////
1538/// Get y value for point i.
1539
1541{
1542 if (i < 0 || i >= fNpoints || !fY)
1543 return -1.;
1544
1545 return fY[i];
1546}
1547
1548////////////////////////////////////////////////////////////////////////////////
1549/// Get x axis of the graph.
1550
1552{
1553 TH1 *h = GetHistogram();
1554 if (!h) return 0;
1555 return h->GetXaxis();
1556}
1557
1558////////////////////////////////////////////////////////////////////////////////
1559/// Get y axis of the graph.
1560
1562{
1563 TH1 *h = GetHistogram();
1564 if (!h) return 0;
1565 return h->GetYaxis();
1566}
1567
1568////////////////////////////////////////////////////////////////////////////////
1569/// Implementation to get information on point of graph at cursor position
1570/// Adapted from class TH1
1571
1573{
1574 // localize point
1575 Int_t ipoint = -2;
1576 Int_t i;
1577 // start with a small window (in case the mouse is very close to one point)
1578 for (i = 0; i < fNpoints; i++) {
1579 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
1580 Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
1581
1582 if (dpx * dpx + dpy * dpy < 25) {
1583 ipoint = i;
1584 break;
1585 }
1586 }
1587
1588 Double_t x = gPad->PadtoX(gPad->AbsPixeltoX(px));
1589 Double_t y = gPad->PadtoY(gPad->AbsPixeltoY(py));
1590
1591 if (ipoint == -2)
1592 return Form("x=%g, y=%g", x, y);
1593
1594 Double_t xval = fX[ipoint];
1595 Double_t yval = fY[ipoint];
1596
1597 return Form("x=%g, y=%g, point=%d, xval=%g, yval=%g", x, y, ipoint, xval, yval);
1598}
1599
1600////////////////////////////////////////////////////////////////////////////////
1601/// Compute Initial values of parameters for a gaussian.
1602
1604{
1605 Double_t allcha, sumx, sumx2, x, val, rms, mean;
1606 Int_t bin;
1607 const Double_t sqrtpi = 2.506628;
1608
1609 // Compute mean value and RMS of the graph in the given range
1610 if (xmax <= xmin) {
1611 xmin = fX[0];
1612 xmax = fX[fNpoints-1];
1613 }
1614 Int_t np = 0;
1615 allcha = sumx = sumx2 = 0;
1616 for (bin = 0; bin < fNpoints; bin++) {
1617 x = fX[bin];
1618 if (x < xmin || x > xmax) continue;
1619 np++;
1620 val = fY[bin];
1621 sumx += val * x;
1622 sumx2 += val * x * x;
1623 allcha += val;
1624 }
1625 if (np == 0 || allcha == 0) return;
1626 mean = sumx / allcha;
1627 rms = TMath::Sqrt(sumx2 / allcha - mean * mean);
1628 Double_t binwidx = TMath::Abs((xmax - xmin) / np);
1629 if (rms == 0) rms = 1;
1631 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1632 f1->SetParameter(0, binwidx * allcha / (sqrtpi * rms));
1633 f1->SetParameter(1, mean);
1634 f1->SetParameter(2, rms);
1635 f1->SetParLimits(2, 0, 10 * rms);
1636}
1637
1638////////////////////////////////////////////////////////////////////////////////
1639/// Compute Initial values of parameters for an exponential.
1640
1642{
1643 Double_t constant, slope;
1644 Int_t ifail;
1645 if (xmax <= xmin) {
1646 xmin = fX[0];
1647 xmax = fX[fNpoints-1];
1648 }
1649 Int_t nchanx = fNpoints;
1650
1651 LeastSquareLinearFit(-nchanx, constant, slope, ifail, xmin, xmax);
1652
1654 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1655 f1->SetParameter(0, constant);
1656 f1->SetParameter(1, slope);
1657}
1658
1659////////////////////////////////////////////////////////////////////////////////
1660/// Compute Initial values of parameters for a polynom.
1661
1663{
1664 Double_t fitpar[25];
1665
1667 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1668 Int_t npar = f1->GetNpar();
1669 if (xmax <= xmin) {
1670 xmin = fX[0];
1671 xmax = fX[fNpoints-1];
1672 }
1673
1674 LeastSquareFit(npar, fitpar, xmin, xmax);
1675
1676 for (Int_t i = 0; i < npar; i++) f1->SetParameter(i, fitpar[i]);
1677}
1678
1679////////////////////////////////////////////////////////////////////////////////
1680/// Insert a new point at the mouse position
1681
1683{
1684 Int_t px = gPad->GetEventX();
1685 Int_t py = gPad->GetEventY();
1686
1687 //localize point where to insert
1688 Int_t ipoint = -2;
1689 Int_t i, d = 0;
1690 // start with a small window (in case the mouse is very close to one point)
1691 for (i = 0; i < fNpoints - 1; i++) {
1692 d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1693 if (d < 5) {
1694 ipoint = i + 1;
1695 break;
1696 }
1697 }
1698 if (ipoint == -2) {
1699 //may be we are far from one point, try again with a larger window
1700 for (i = 0; i < fNpoints - 1; i++) {
1701 d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1702 if (d < 10) {
1703 ipoint = i + 1;
1704 break;
1705 }
1706 }
1707 }
1708 if (ipoint == -2) {
1709 //distinguish between first and last point
1710 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[0]));
1711 Int_t dpy = py - gPad->YtoAbsPixel(gPad->XtoPad(fY[0]));
1712 if (dpx * dpx + dpy * dpy < 25) ipoint = 0;
1713 else ipoint = fNpoints;
1714 }
1715
1716
1717 InsertPointBefore(ipoint, gPad->AbsPixeltoX(px), gPad->AbsPixeltoY(py));
1718
1719 gPad->Modified();
1720 return ipoint;
1721}
1722
1723
1724////////////////////////////////////////////////////////////////////////////////
1725/// Insert a new point with coordinates (x,y) before the point number `ipoint`.
1726
1728{
1729 if (ipoint < 0) {
1730 Error("TGraph", "Inserted point index should be >= 0");
1731 return;
1732 }
1733
1734 if (ipoint > fNpoints) {
1735 Error("TGraph", "Inserted point index should be <= %d", fNpoints);
1736 return;
1737 }
1738
1739 if (ipoint == fNpoints) {
1740 SetPoint(ipoint, x, y);
1741 return;
1742 }
1743
1744 Double_t **ps = ExpandAndCopy(fNpoints + 1, ipoint);
1745 CopyAndRelease(ps, ipoint, fNpoints++, ipoint + 1);
1746
1747 // To avoid redefinitions in descendant classes
1748 FillZero(ipoint, ipoint + 1);
1749
1750 fX[ipoint] = x;
1751 fY[ipoint] = y;
1752}
1753
1754
1755////////////////////////////////////////////////////////////////////////////////
1756/// Integrate the TGraph data within a given (index) range.
1757/// Note that this function computes the area of the polygon enclosed by the points of the TGraph.
1758/// The polygon segments, which are defined by the points of the TGraph, do not need to form a closed polygon,
1759/// since the last polygon segment, which closes the polygon, is taken as the line connecting the last TGraph point
1760/// with the first one. It is clear that the order of the point is essential in defining the polygon.
1761/// Also note that the segments should not intersect.
1762///
1763/// NB:
1764/// - if last=-1 (default) last is set to the last point.
1765/// - if (first <0) the first point (0) is taken.
1766///
1767/// ### Method:
1768///
1769/// There are many ways to calculate the surface of a polygon. It all depends on what kind of data
1770/// you have to deal with. The most evident solution would be to divide the polygon in triangles and
1771/// calculate the surface of them. But this can quickly become complicated as you will have to test
1772/// every segments of every triangles and check if they are intersecting with a current polygon's
1773/// segment or if it goes outside the polygon. Many calculations that would lead to many problems...
1774///
1775/// ### The solution (implemented by R.Brun)
1776/// Fortunately for us, there is a simple way to solve this problem, as long as the polygon's
1777/// segments don't intersect.
1778/// It takes the x coordinate of the current vertex and multiply it by the y coordinate of the next
1779/// vertex. Then it subtracts from it the result of the y coordinate of the current vertex multiplied
1780/// by the x coordinate of the next vertex. Then divide the result by 2 to get the surface/area.
1781///
1782/// ### Sources
1783/// - http://forums.wolfram.com/mathgroup/archive/1998/Mar/msg00462.html
1784/// - http://stackoverflow.com/questions/451426/how-do-i-calculate-the-surface-area-of-a-2d-polygon
1785
1787{
1788 if (first < 0) first = 0;
1789 if (last < 0) last = fNpoints - 1;
1790 if (last >= fNpoints) last = fNpoints - 1;
1791 if (first >= last) return 0;
1792 Int_t np = last - first + 1;
1793 Double_t sum = 0.0;
1794 //for(Int_t i=first;i<=last;i++) {
1795 // Int_t j = first + (i-first+1)%np;
1796 // sum += TMath::Abs(fX[i]*fY[j]);
1797 // sum -= TMath::Abs(fY[i]*fX[j]);
1798 //}
1799 for (Int_t i = first; i <= last; i++) {
1800 Int_t j = first + (i - first + 1) % np;
1801 sum += (fY[i] + fY[j]) * (fX[j] - fX[i]);
1802 }
1803 return 0.5 * TMath::Abs(sum);
1804}
1805
1806////////////////////////////////////////////////////////////////////////////////
1807/// Return 1 if the point (x,y) is inside the polygon defined by
1808/// the graph vertices 0 otherwise.
1809///
1810/// Algorithm:
1811///
1812/// The loop is executed with the end-point coordinates of a line segment
1813/// (X1,Y1)-(X2,Y2) and the Y-coordinate of a horizontal line.
1814/// The counter inter is incremented if the line (X1,Y1)-(X2,Y2) intersects
1815/// the horizontal line. In this case XINT is set to the X-coordinate of the
1816/// intersection point. If inter is an odd number, then the point x,y is within
1817/// the polygon.
1818
1820{
1821 return (Int_t)TMath::IsInside(x, y, fNpoints, fX, fY);
1822}
1823
1824////////////////////////////////////////////////////////////////////////////////
1825/// Least squares polynomial fitting without weights.
1826///
1827/// \param [in] m number of parameters
1828/// \param [in] a array of parameters
1829/// \param [in] xmin 1st point number to fit (default =0)
1830/// \param [in] xmax last point number to fit (default=fNpoints-1)
1831///
1832/// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
1833
1835{
1836 const Double_t zero = 0.;
1837 const Double_t one = 1.;
1838 const Int_t idim = 20;
1839
1840 Double_t b[400] /* was [20][20] */;
1841 Int_t i, k, l, ifail;
1842 Double_t power;
1843 Double_t da[20], xk, yk;
1844 Int_t n = fNpoints;
1845 if (xmax <= xmin) {
1846 xmin = fX[0];
1847 xmax = fX[fNpoints-1];
1848 }
1849
1850 if (m <= 2) {
1851 LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
1852 return;
1853 }
1854 if (m > idim || m > n) return;
1855 da[0] = zero;
1856 for (l = 2; l <= m; ++l) {
1857 b[l-1] = zero;
1858 b[m + l*20 - 21] = zero;
1859 da[l-1] = zero;
1860 }
1861 Int_t np = 0;
1862 for (k = 0; k < fNpoints; ++k) {
1863 xk = fX[k];
1864 if (xk < xmin || xk > xmax) continue;
1865 np++;
1866 yk = fY[k];
1867 power = one;
1868 da[0] += yk;
1869 for (l = 2; l <= m; ++l) {
1870 power *= xk;
1871 b[l-1] += power;
1872 da[l-1] += power * yk;
1873 }
1874 for (l = 2; l <= m; ++l) {
1875 power *= xk;
1876 b[m + l*20 - 21] += power;
1877 }
1878 }
1879 b[0] = Double_t(np);
1880 for (i = 3; i <= m; ++i) {
1881 for (k = i; k <= m; ++k) {
1882 b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
1883 }
1884 }
1885 H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
1886
1887 if (ifail < 0) {
1888 a[0] = fY[0];
1889 for (i = 1; i < m; ++i) a[i] = 0;
1890 return;
1891 }
1892 for (i = 0; i < m; ++i) a[i] = da[i];
1893}
1894
1895////////////////////////////////////////////////////////////////////////////////
1896/// Least square linear fit without weights.
1897///
1898/// Fit a straight line (a0 + a1*x) to the data in this graph.
1899///
1900/// \param [in] ndata if ndata<0, fits the logarithm of the graph (used in InitExpo() to set
1901/// the initial parameter values for a fit with exponential function.
1902/// \param [in] a0 constant
1903/// \param [in] a1 slope
1904/// \param [in] ifail return parameter indicating the status of the fit (ifail=0, fit is OK)
1905/// \param [in] xmin, xmax fitting range
1906///
1907/// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
1908
1910{
1911 Double_t xbar, ybar, x2bar;
1912 Int_t i;
1913 Double_t xybar;
1914 Double_t fn, xk, yk;
1915 Double_t det;
1916 if (xmax <= xmin) {
1917 xmin = fX[0];
1918 xmax = fX[fNpoints-1];
1919 }
1920
1921 ifail = -2;
1922 xbar = ybar = x2bar = xybar = 0;
1923 Int_t np = 0;
1924 for (i = 0; i < fNpoints; ++i) {
1925 xk = fX[i];
1926 if (xk < xmin || xk > xmax) continue;
1927 np++;
1928 yk = fY[i];
1929 if (ndata < 0) {
1930 if (yk <= 0) yk = 1e-9;
1931 yk = TMath::Log(yk);
1932 }
1933 xbar += xk;
1934 ybar += yk;
1935 x2bar += xk * xk;
1936 xybar += xk * yk;
1937 }
1938 fn = Double_t(np);
1939 det = fn * x2bar - xbar * xbar;
1940 ifail = -1;
1941 if (det <= 0) {
1942 if (fn > 0) a0 = ybar / fn;
1943 else a0 = 0;
1944 a1 = 0;
1945 return;
1946 }
1947 ifail = 0;
1948 a0 = (x2bar * ybar - xbar * xybar) / det;
1949 a1 = (fn * xybar - xbar * ybar) / det;
1950}
1951
1952////////////////////////////////////////////////////////////////////////////////
1953/// Draw this graph with its current attributes.
1954
1956{
1958 if (painter) painter->PaintHelper(this, option);
1959}
1960
1961////////////////////////////////////////////////////////////////////////////////
1962/// Draw the (x,y) as a graph.
1963
1964void TGraph::PaintGraph(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
1965{
1967 if (painter) painter->PaintGraph(this, npoints, x, y, chopt);
1968}
1969
1970////////////////////////////////////////////////////////////////////////////////
1971/// Draw the (x,y) as a histogram.
1972
1973void TGraph::PaintGrapHist(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
1974{
1976 if (painter) painter->PaintGrapHist(this, npoints, x, y, chopt);
1977}
1978
1979////////////////////////////////////////////////////////////////////////////////
1980/// Draw the stats
1981
1983{
1985 if (painter) painter->PaintStats(this, fit);
1986}
1987
1988////////////////////////////////////////////////////////////////////////////////
1989/// Print graph values.
1990
1992{
1993 for (Int_t i = 0; i < fNpoints; i++) {
1994 printf("x[%d]=%g, y[%d]=%g\n", i, fX[i], i, fY[i]);
1995 }
1996}
1997
1998////////////////////////////////////////////////////////////////////////////////
1999/// Recursively remove object from the list of functions
2000
2002{
2003 if (fFunctions) {
2005 }
2006 if (fHistogram == obj)
2007 fHistogram = nullptr;
2008}
2009
2010////////////////////////////////////////////////////////////////////////////////
2011/// Delete point close to the mouse position
2012
2014{
2015 Int_t px = gPad->GetEventX();
2016 Int_t py = gPad->GetEventY();
2017
2018 //localize point to be deleted
2019 Int_t ipoint = -2;
2020 Int_t i;
2021 // start with a small window (in case the mouse is very close to one point)
2022 for (i = 0; i < fNpoints; i++) {
2023 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
2024 Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
2025 if (dpx * dpx + dpy * dpy < 100) {
2026 ipoint = i;
2027 break;
2028 }
2029 }
2030 return RemovePoint(ipoint);
2031}
2032
2033////////////////////////////////////////////////////////////////////////////////
2034/// Delete point number ipoint
2035
2037{
2038 if (ipoint < 0) return -1;
2039 if (ipoint >= fNpoints) return -1;
2040
2041 Double_t **ps = ShrinkAndCopy(fNpoints - 1, ipoint);
2042 CopyAndRelease(ps, ipoint + 1, fNpoints--, ipoint);
2043 if (gPad) gPad->Modified();
2044 return ipoint;
2045}
2046
2047////////////////////////////////////////////////////////////////////////////////
2048/// Save the graph as .csv, .tsv or .txt. In case of any other extension, fall
2049/// back to TObject::SaveAs
2050///
2051/// The result can be immediately imported into Excel, gnuplot, Python or whatever,
2052/// without the needing to install pyroot, etc.
2053///
2054/// \param filename the name of the file where to store the graph
2055/// \param option some tuning options
2056///
2057/// The file extension defines the delimiter used:
2058/// - `.csv` : comma
2059/// - `.tsv` : tab
2060/// - `.txt` : space
2061///
2062/// If option = "title" a title line is generated with the axis titles.
2063
2064void TGraph::SaveAs(const char *filename, Option_t *option) const
2065{
2066 char del = '\0';
2067 TString ext = "";
2068 TString fname = filename;
2069 TString opt = option;
2070
2071 if (filename) {
2072 if (fname.EndsWith(".csv")) {del = ','; ext = "csv";}
2073 else if (fname.EndsWith(".tsv")) {del = '\t'; ext = "tsv";}
2074 else if (fname.EndsWith(".txt")) {del = ' '; ext = "txt";}
2075 }
2076 if (del) {
2077 std::ofstream out;
2078 out.open(filename, std::ios::out);
2079 if (!out.good ()) {
2080 Error("SaveAs", "cannot open file: %s", filename);
2081 return;
2082 }
2084 if(opt.Contains("title"))
2085 out << "# " << GetXaxis()->GetTitle() << "\tex\t" << GetYaxis()->GetTitle() << "\tey" << std::endl;
2086 double *ex = this->GetEX();
2087 double *ey = this->GetEY();
2088 for(int i=0 ; i<fNpoints ; i++)
2089 out << fX[i] << del << (ex?ex[i]:0) << del << fY[i] << del << (ey?ey[i]:0) << std::endl;
2091 if(opt.Contains("title"))
2092 out << "# " << GetXaxis()->GetTitle() << "\texl\t" << "\texh\t" << GetYaxis()->GetTitle() << "\teyl" << "\teyh" << std::endl;
2093 double *exl = this->GetEXlow();
2094 double *exh = this->GetEXhigh();
2095 double *eyl = this->GetEYlow();
2096 double *eyh = this->GetEYhigh();
2097 for(int i=0 ; i<fNpoints ; i++)
2098 out << fX[i] << del << (exl?exl[i]:0) << del << (exh?exh[i]:0) << del << fY[i] << del << (eyl?eyl[i]:0) << del << (eyh?eyh[i]:0) << std::endl;
2099 } else {
2100 if(opt.Contains("title"))
2101 out << "# " << GetXaxis()->GetTitle() << "\t" << GetYaxis()->GetTitle() << std::endl;
2102 for (int i=0 ; i<fNpoints ; i++)
2103 out << fX[i] << del << fY[i] << std::endl;
2104 }
2105 out.close();
2106 Info("SaveAs", "%s file: %s has been generated", ext.Data(), filename);
2107 } else {
2109 }
2110}
2111
2112
2113////////////////////////////////////////////////////////////////////////////////
2114/// Save primitive as a C++ statement(s) on output stream out
2115
2116void TGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
2117{
2118 char quote = '"';
2119 out << " " << std::endl;
2120 static Int_t frameNumber = 0;
2121 frameNumber++;
2122
2123 if (fNpoints >= 1) {
2124 Int_t i;
2125 TString fXName = TString(GetName()) + Form("_fx%d",frameNumber);
2126 TString fYName = TString(GetName()) + Form("_fy%d",frameNumber);
2127 out << " Double_t " << fXName << "[" << fNpoints << "] = {" << std::endl;
2128 for (i = 0; i < fNpoints-1; i++) out << " " << fX[i] << "," << std::endl;
2129 out << " " << fX[fNpoints-1] << "};" << std::endl;
2130 out << " Double_t " << fYName << "[" << fNpoints << "] = {" << std::endl;
2131 for (i = 0; i < fNpoints-1; i++) out << " " << fY[i] << "," << std::endl;
2132 out << " " << fY[fNpoints-1] << "};" << std::endl;
2133 if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2134 else out << " TGraph *";
2135 out << "graph = new TGraph(" << fNpoints << "," << fXName << "," << fYName << ");" << std::endl;
2136 } else {
2137 if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2138 else out << " TGraph *";
2139 out << "graph = new TGraph();" << std::endl;
2140 }
2141
2142 out << " graph->SetName(" << quote << GetName() << quote << ");" << std::endl;
2143 out << " graph->SetTitle(" << quote << GetTitle() << quote << ");" << std::endl;
2144
2145 SaveFillAttributes(out, "graph", 0, 1001);
2146 SaveLineAttributes(out, "graph", 1, 1, 1);
2147 SaveMarkerAttributes(out, "graph", 1, 1, 1);
2148
2149 if (fHistogram) {
2150 TString hname = fHistogram->GetName();
2151 hname += frameNumber;
2152 fHistogram->SetName(Form("Graph_%s", hname.Data()));
2153 fHistogram->SavePrimitive(out, "nodraw");
2154 out << " graph->SetHistogram(" << fHistogram->GetName() << ");" << std::endl;
2155 out << " " << std::endl;
2156 }
2157
2158 // save list of functions
2159 TIter next(fFunctions);
2160 TObject *obj;
2161 while ((obj = next())) {
2162 obj->SavePrimitive(out, Form("nodraw #%d\n",++frameNumber));
2163 if (obj->InheritsFrom("TPaveStats")) {
2164 out << " graph->GetListOfFunctions()->Add(ptstats);" << std::endl;
2165 out << " ptstats->SetParent(graph->GetListOfFunctions());" << std::endl;
2166 } else {
2167 TString objname;
2168 objname.Form("%s%d",obj->GetName(),frameNumber);
2169 if (obj->InheritsFrom("TF1")) {
2170 out << " " << objname << "->SetParent(graph);\n";
2171 }
2172 out << " graph->GetListOfFunctions()->Add("
2173 << objname << ");" << std::endl;
2174 }
2175 }
2176
2177 const char *soption = option ? option : "";
2178 const char *l = strstr(soption, "multigraph");
2179 if (l) {
2180 out << " multigraph->Add(graph," << quote << l + 10 << quote << ");" << std::endl;
2181 return;
2182 }
2183 l = strstr(soption, "th2poly");
2184 if (l) {
2185 out << " " << l + 7 << "->AddBin(graph);" << std::endl;
2186 return;
2187 }
2188 out << " graph->Draw(" << quote << soption << quote << ");" << std::endl;
2189}
2190
2191////////////////////////////////////////////////////////////////////////////////
2192/// Multiply the values of a TGraph by a constant c1.
2193///
2194/// If option contains "x" the x values are scaled
2195/// If option contains "y" the y values are scaled
2196/// If option contains "xy" both x and y values are scaled
2197
2199{
2200 TString opt = option; opt.ToLower();
2201 if (opt.Contains("x")) {
2202 for (Int_t i=0; i<GetN(); i++)
2203 GetX()[i] *= c1;
2204 }
2205 if (opt.Contains("y")) {
2206 for (Int_t i=0; i<GetN(); i++)
2207 GetY()[i] *= c1;
2208 }
2209}
2210
2211////////////////////////////////////////////////////////////////////////////////
2212/// Set number of points in the graph
2213/// Existing coordinates are preserved
2214/// New coordinates above fNpoints are preset to 0.
2215
2217{
2218 if (n < 0) n = 0;
2219 if (n == fNpoints) return;
2220 Double_t **ps = Allocate(n);
2222 if (n > fNpoints) {
2224 }
2225 fNpoints = n;
2226}
2227
2228////////////////////////////////////////////////////////////////////////////////
2229/// Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
2230
2232{
2233 return TestBit(kNotEditable) ? kFALSE : kTRUE;
2234}
2235
2236////////////////////////////////////////////////////////////////////////////////
2237/// if editable=kFALSE, the graph cannot be modified with the mouse
2238/// by default a TGraph is editable
2239
2241{
2242 if (editable) ResetBit(kNotEditable);
2243 else SetBit(kNotEditable);
2244}
2245
2246////////////////////////////////////////////////////////////////////////////////
2247/// Set highlight (enable/disble) mode for the graph
2248/// by default highlight mode is disable
2249
2251{
2252 if (IsHighlight() == set) return;
2253
2255 if (!painter) return;
2256 SetBit(kIsHighlight, set);
2257 painter->SetHighlight(this);
2258}
2259
2260////////////////////////////////////////////////////////////////////////////////
2261/// Set the maximum of the graph.
2262
2264{
2265 fMaximum = maximum;
2266 GetHistogram()->SetMaximum(maximum);
2267}
2268
2269////////////////////////////////////////////////////////////////////////////////
2270/// Set the minimum of the graph.
2271
2273{
2274 fMinimum = minimum;
2275 GetHistogram()->SetMinimum(minimum);
2276}
2277
2278////////////////////////////////////////////////////////////////////////////////
2279/// Set x and y values for point number i.
2280
2282{
2283 if (i < 0) return;
2285
2286 if (i >= fMaxSize) {
2287 Double_t **ps = ExpandAndCopy(i + 1, fNpoints);
2288 CopyAndRelease(ps, 0, 0, 0);
2289 }
2290 if (i >= fNpoints) {
2291 // points above i can be not initialized
2292 // set zero up to i-th point to avoid redefinition
2293 // of this method in descendant classes
2294 FillZero(fNpoints, i + 1);
2295 fNpoints = i + 1;
2296 }
2297 fX[i] = x;
2298 fY[i] = y;
2299 if (gPad) gPad->Modified();
2300}
2301
2302////////////////////////////////////////////////////////////////////////////////
2303/// Set x value for point i.
2304
2306{
2307 SetPoint(i, x, GetPointY(i));
2308}
2309
2310////////////////////////////////////////////////////////////////////////////////
2311/// Set y value for point i.
2312
2314{
2315 SetPoint(i, GetPointX(i), y);
2316}
2317
2318////////////////////////////////////////////////////////////////////////////////
2319/// Set graph name.
2320void TGraph::SetName(const char *name)
2321{
2322 fName = name;
2324}
2325
2326////////////////////////////////////////////////////////////////////////////////
2327/// Change (i.e. set) the title
2328///
2329/// if title is in the form `stringt;stringx;stringy;stringz`
2330/// the graph title is set to `stringt`, the x axis title to `stringx`,
2331/// the y axis title to `stringy`, and the z axis title to `stringz`.
2332///
2333/// To insert the character `;` in one of the titles, one should use `#;`
2334/// or `#semicolon`.
2335
2336void TGraph::SetTitle(const char* title)
2337{
2338 fTitle = title;
2339 fTitle.ReplaceAll("#;",2,"#semicolon",10);
2340 Int_t p = fTitle.Index(";");
2341
2342 if (p>0) {
2343 if (!fHistogram) GetHistogram();
2344 fHistogram->SetTitle(title);
2345 Int_t n = fTitle.Length()-p;
2346 if (p>0) fTitle.Remove(p,n);
2347 fTitle.ReplaceAll("#semicolon",10,"#;",2);
2348 } else {
2349 if (fHistogram) fHistogram->SetTitle(title);
2350 }
2351}
2352
2353////////////////////////////////////////////////////////////////////////////////
2354/// Set graph name and title
2355
2356void TGraph::SetNameTitle(const char *name, const char *title)
2357{
2358 SetName(name);
2359 SetTitle(title);
2360}
2361
2362////////////////////////////////////////////////////////////////////////////////
2363/// Set statistics option on/off.
2364///
2365/// By default, the statistics box is drawn.
2366/// The paint options can be selected via gStyle->SetOptStats.
2367/// This function sets/resets the kNoStats bit in the graph object.
2368/// It has priority over the Style option.
2369
2371{
2373 if (!stats) {
2375 //remove the "stats" object from the list of functions
2376 if (fFunctions) {
2377 TObject *obj = fFunctions->FindObject("stats");
2378 if (obj) {
2379 fFunctions->Remove(obj);
2380 delete obj;
2381 }
2382 }
2383 }
2384}
2385
2386////////////////////////////////////////////////////////////////////////////////
2387/// if size*2 <= fMaxSize allocate new arrays of size points,
2388/// copy points [0,oend).
2389/// Return newarray (passed or new instance if it was zero
2390/// and allocations are needed)
2391
2393{
2394 if (size * 2 > fMaxSize || !fMaxSize) {
2395 return 0;
2396 }
2397 Double_t **newarrays = Allocate(size);
2398 CopyPoints(newarrays, 0, oend, 0);
2399 return newarrays;
2400}
2401
2402////////////////////////////////////////////////////////////////////////////////
2403/// Sorts the points of this TGraph using in-place quicksort (see e.g. older glibc).
2404/// To compare two points the function parameter greaterfunc is used (see TGraph::CompareX for an
2405/// example of such a method, which is also the default comparison function for Sort). After
2406/// the sort, greaterfunc(this, i, j) will return kTRUE for all i>j if ascending == kTRUE, and
2407/// kFALSE otherwise.
2408///
2409/// The last two parameters are used for the recursive quick sort, stating the range to be sorted
2410///
2411/// Examples:
2412/// ~~~ {.cpp}
2413/// // sort points along x axis
2414/// graph->Sort();
2415/// // sort points along their distance to origin
2416/// graph->Sort(&TGraph::CompareRadius);
2417///
2418/// Bool_t CompareErrors(const TGraph* gr, Int_t i, Int_t j) {
2419/// const TGraphErrors* ge=(const TGraphErrors*)gr;
2420/// return (ge->GetEY()[i]>ge->GetEY()[j]); }
2421/// // sort using the above comparison function, largest errors first
2422/// graph->Sort(&CompareErrors, kFALSE);
2423/// ~~~
2424
2425void TGraph::Sort(Bool_t (*greaterfunc)(const TGraph*, Int_t, Int_t) /*=TGraph::CompareX()*/,
2426 Bool_t ascending /*=kTRUE*/, Int_t low /* =0 */, Int_t high /* =-1111 */)
2427{
2428
2429 // set the bit in case of an ascending =sort in X
2430 if (greaterfunc == TGraph::CompareX && ascending && low == 0 && high == -1111)
2432
2433 if (high == -1111) high = GetN() - 1;
2434 // Termination condition
2435 if (high <= low) return;
2436
2437 int left, right;
2438 left = low; // low is the pivot element
2439 right = high;
2440 while (left < right) {
2441 // move left while item < pivot
2442 while (left <= high && greaterfunc(this, left, low) != ascending)
2443 left++;
2444 // move right while item > pivot
2445 while (right > low && greaterfunc(this, right, low) == ascending)
2446 right--;
2447 if (left < right && left < high && right > low)
2448 SwapPoints(left, right);
2449 }
2450 // right is final position for the pivot
2451 if (right > low)
2452 SwapPoints(low, right);
2453 Sort(greaterfunc, ascending, low, right - 1);
2454 Sort(greaterfunc, ascending, right + 1, high);
2455}
2456
2457////////////////////////////////////////////////////////////////////////////////
2458/// Stream an object of class TGraph.
2459
2461{
2462 if (b.IsReading()) {
2463 UInt_t R__s, R__c;
2464 Version_t R__v = b.ReadVersion(&R__s, &R__c);
2465 if (R__v > 2) {
2466 b.ReadClassBuffer(TGraph::Class(), this, R__v, R__s, R__c);
2467 if (fHistogram) fHistogram->SetDirectory(nullptr);
2468 TIter next(fFunctions);
2469 TObject *obj;
2470 while ((obj = next())) {
2471 if (obj->InheritsFrom(TF1::Class())) {
2472 TF1 *f1 = (TF1*)obj;
2473 f1->SetParent(this);
2474 }
2475 }
2477 return;
2478 }
2479 //====process old versions before automatic schema evolution
2484 b >> fNpoints;
2486 fX = new Double_t[fNpoints];
2487 fY = new Double_t[fNpoints];
2488 if (R__v < 2) {
2489 Float_t *x = new Float_t[fNpoints];
2490 Float_t *y = new Float_t[fNpoints];
2491 b.ReadFastArray(x, fNpoints);
2492 b.ReadFastArray(y, fNpoints);
2493 for (Int_t i = 0; i < fNpoints; i++) {
2494 fX[i] = x[i];
2495 fY[i] = y[i];
2496 }
2497 delete [] y;
2498 delete [] x;
2499 } else {
2500 b.ReadFastArray(fX, fNpoints);
2501 b.ReadFastArray(fY, fNpoints);
2502 }
2503 b >> fFunctions;
2504 b >> fHistogram;
2505 if (fHistogram) fHistogram->SetDirectory(nullptr);
2506 if (R__v < 2) {
2507 Float_t mi, ma;
2508 b >> mi;
2509 b >> ma;
2510 fMinimum = mi;
2511 fMaximum = ma;
2512 } else {
2513 b >> fMinimum;
2514 b >> fMaximum;
2515 }
2516 b.CheckByteCount(R__s, R__c, TGraph::IsA());
2517 //====end of old versions
2518
2519 } else {
2520 b.WriteClassBuffer(TGraph::Class(), this);
2521 }
2522}
2523
2524////////////////////////////////////////////////////////////////////////////////
2525/// Swap points.
2526
2528{
2529 SwapValues(fX, pos1, pos2);
2530 SwapValues(fY, pos1, pos2);
2531}
2532
2533////////////////////////////////////////////////////////////////////////////////
2534/// Swap values.
2535
2537{
2538 Double_t tmp = arr[pos1];
2539 arr[pos1] = arr[pos2];
2540 arr[pos2] = tmp;
2541}
2542
2543////////////////////////////////////////////////////////////////////////////////
2544/// Set current style settings in this graph
2545/// This function is called when either TCanvas::UseCurrentStyle
2546/// or TROOT::ForceStyle have been invoked.
2547
2549{
2550 if (gStyle->IsReading()) {
2559 } else {
2568 }
2570
2571 TIter next(GetListOfFunctions());
2572 TObject *obj;
2573
2574 while ((obj = next())) {
2575 obj->UseCurrentStyle();
2576 }
2577}
2578
2579////////////////////////////////////////////////////////////////////////////////
2580/// Adds all graphs from the collection to this graph.
2581/// Returns the total number of points in the result or -1 in case of an error.
2582
2584{
2585 TIter next(li);
2586 while (TObject* o = next()) {
2587 TGraph *g = dynamic_cast<TGraph*>(o);
2588 if (!g) {
2589 Error("Merge",
2590 "Cannot merge - an object which doesn't inherit from TGraph found in the list");
2591 return -1;
2592 }
2593 DoMerge(g);
2594 }
2595 return GetN();
2596}
2597
2598////////////////////////////////////////////////////////////////////////////////
2599/// protected function to perform the merge operation of a graph
2600
2602{
2603 Double_t x = 0, y = 0;
2604 for (Int_t i = 0 ; i < g->GetN(); i++) {
2605 g->GetPoint(i, x, y);
2606 SetPoint(GetN(), x, y);
2607 }
2608 return kTRUE;
2609}
2610
2611////////////////////////////////////////////////////////////////////////////////
2612/// Move all graph points on specified values dx,dy
2613/// If log argument specified, calculation done in logarithmic scale like:
2614/// new_value = exp( log(old_value) + delta );
2615
2617{
2618 Double_t x = 0, y = 0;
2619 for (Int_t i = 0 ; i < GetN(); i++) {
2620 GetPoint(i, x, y);
2621 if (!logx) {
2622 x += dx;
2623 } else if (x > 0) {
2624 x = TMath::Exp(TMath::Log(x) + dx);
2625 }
2626 if (!logy) {
2627 y += dy;
2628 } else if (y > 0) {
2629 y = TMath::Exp(TMath::Log(y) + dy);
2630 }
2631 SetPoint(i, x, y);
2632 }
2633}
2634
2635
2636////////////////////////////////////////////////////////////////////////////////
2637/// Find zero of a continuous function.
2638/// This function finds a real zero of the continuous real
2639/// function Y(X) in a given interval (A,B). See accompanying
2640/// notes for details of the argument list and calling sequence
2641
2643 , Int_t maxiterations)
2644{
2645 static Double_t a, b, ya, ytest, y1, x1, h;
2646 static Int_t j1, it, j3, j2;
2647 Double_t yb, x2;
2648 yb = 0;
2649
2650 // Calculate Y(X) at X=AZ.
2651 if (k <= 0) {
2652 a = AZ;
2653 b = BZ;
2654 X = a;
2655 j1 = 1;
2656 it = 1;
2657 k = j1;
2658 return;
2659 }
2660
2661 // Test whether Y(X) is sufficiently small.
2662
2663 if (TMath::Abs(Y) <= E2) {
2664 k = 2;
2665 return;
2666 }
2667
2668 // Calculate Y(X) at X=BZ.
2669
2670 if (j1 == 1) {
2671 ya = Y;
2672 X = b;
2673 j1 = 2;
2674 return;
2675 }
2676 // Test whether the signs of Y(AZ) and Y(BZ) are different.
2677 // if not, begin the binary subdivision.
2678
2679 if (j1 != 2) goto L100;
2680 if (ya * Y < 0) goto L120;
2681 x1 = a;
2682 y1 = ya;
2683 j1 = 3;
2684 h = b - a;
2685 j2 = 1;
2686 x2 = a + 0.5 * h;
2687 j3 = 1;
2688 it++; //*-*- Check whether (maxiterations) function values have been calculated.
2689 if (it >= maxiterations) k = j1;
2690 else X = x2;
2691 return;
2692
2693 // Test whether a bracket has been found .
2694 // If not,continue the search
2695
2696L100:
2697 if (j1 > 3) goto L170;
2698 if (ya*Y >= 0) {
2699 if (j3 >= j2) {
2700 h = 0.5 * h;
2701 j2 = 2 * j2;
2702 a = x1;
2703 ya = y1;
2704 x2 = a + 0.5 * h;
2705 j3 = 1;
2706 } else {
2707 a = X;
2708 ya = Y;
2709 x2 = X + h;
2710 j3++;
2711 }
2712 it++;
2713 if (it >= maxiterations) k = j1;
2714 else X = x2;
2715 return;
2716 }
2717
2718 // The first bracket has been found.calculate the next X by the
2719 // secant method based on the bracket.
2720
2721L120:
2722 b = X;
2723 yb = Y;
2724 j1 = 4;
2725L130:
2726 if (TMath::Abs(ya) > TMath::Abs(yb)) {
2727 x1 = a;
2728 y1 = ya;
2729 X = b;
2730 Y = yb;
2731 } else {
2732 x1 = b;
2733 y1 = yb;
2734 X = a;
2735 Y = ya;
2736 }
2737
2738 // Use the secant method based on the function values y1 and Y.
2739 // check that x2 is inside the interval (a,b).
2740
2741L150:
2742 x2 = X - Y * (X - x1) / (Y - y1);
2743 x1 = X;
2744 y1 = Y;
2745 ytest = 0.5 * TMath::Min(TMath::Abs(ya), TMath::Abs(yb));
2746 if ((x2 - a)*(x2 - b) < 0) {
2747 it++;
2748 if (it >= maxiterations) k = j1;
2749 else X = x2;
2750 return;
2751 }
2752
2753 // Calculate the next value of X by bisection . Check whether
2754 // the maximum accuracy has been achieved.
2755
2756L160:
2757 x2 = 0.5 * (a + b);
2758 ytest = 0;
2759 if ((x2 - a)*(x2 - b) >= 0) {
2760 k = 2;
2761 return;
2762 }
2763 it++;
2764 if (it >= maxiterations) k = j1;
2765 else X = x2;
2766 return;
2767
2768
2769 // Revise the bracket (a,b).
2770
2771L170:
2772 if (j1 != 4) return;
2773 if (ya * Y < 0) {
2774 b = X;
2775 yb = Y;
2776 } else {
2777 a = X;
2778 ya = Y;
2779 }
2780
2781 // Use ytest to decide the method for the next value of X.
2782
2783 if (ytest <= 0) goto L130;
2784 if (TMath::Abs(Y) - ytest <= 0) goto L150;
2785 goto L160;
2786}
#define d(i)
Definition: RSha256.hxx:102
#define f(i)
Definition: RSha256.hxx:104
#define h(i)
Definition: RSha256.hxx:106
#define e(i)
Definition: RSha256.hxx:103
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
const Ssiz_t kNPOS
Definition: RtypesCore.h:124
int Int_t
Definition: RtypesCore.h:45
short Version_t
Definition: RtypesCore.h:65
const Bool_t kFALSE
Definition: RtypesCore.h:101
float Float_t
Definition: RtypesCore.h:57
double Double_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:375
R__EXTERN TEnv * gEnv
Definition: TEnv.h:170
winID h TVirtualViewer3D TVirtualGLPainter p
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 del
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 b
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 char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
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 x2
Option_t Option_t TPoint TPoint const char x1
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
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 g
Option_t Option_t TPoint TPoint const char y1
char name[80]
Definition: TGX11.cxx:110
void H1LeastSquareSeqnd(Int_t n, Double_t *a, Int_t idim, Int_t &ifail, Int_t k, Double_t *b)
Extracted from CERN Program library routine DSEQN.
Definition: TH1.cxx:4836
float xmin
Definition: THbookFile.cxx:95
float ymin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
float ymax
Definition: THbookFile.cxx:95
Binding & operator=(OUT(*fun)(void))
#define gROOT
Definition: TROOT.h:404
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2452
void Printf(const char *fmt,...)
Formats a string in a circular formatting buffer and prints the string.
Definition: TString.cxx:2466
R__EXTERN TStyle * gStyle
Definition: TStyle.h:414
R__EXTERN TSystem * gSystem
Definition: TSystem.h:559
#define gPad
Definition: TVirtualPad.h:288
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:35
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:36
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:301
virtual Style_t GetTitleFont() const
Definition: TAttAxis.h:47
virtual Float_t GetLabelOffset() const
Definition: TAttAxis.h:40
virtual void SetLabelSize(Float_t size=0.04)
Set size of axis labels.
Definition: TAttAxis.cxx:206
virtual Style_t GetLabelFont() const
Definition: TAttAxis.h:39
virtual void SetTitleFont(Style_t font=62)
Set the title font.
Definition: TAttAxis.cxx:330
virtual void SetLabelOffset(Float_t offset=0.005)
Set distance between the axis and the labels.
Definition: TAttAxis.cxx:194
virtual void SetLabelFont(Style_t font=62)
Set labels' font.
Definition: TAttAxis.cxx:183
virtual void SetTitleSize(Float_t size=0.04)
Set size of axis title.
Definition: TAttAxis.cxx:312
virtual Float_t GetTitleSize() const
Definition: TAttAxis.h:44
virtual Float_t GetLabelSize() const
Definition: TAttAxis.h:41
virtual Float_t GetTitleOffset() const
Definition: TAttAxis.h:43
virtual void SetNdivisions(Int_t n=510, Bool_t optim=kTRUE)
Set the number of divisions for this axis.
Definition: TAttAxis.cxx:236
Fill Area Attributes class.
Definition: TAttFill.h:19
virtual void Streamer(TBuffer &)
virtual Color_t GetFillColor() const
Return the fill area color.
Definition: TAttFill.h:30
void Copy(TAttFill &attfill) const
Copy this fill attributes to a new TAttFill.
Definition: TAttFill.cxx:204
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:236
Line Attributes class.
Definition: TAttLine.h:18
virtual void Streamer(TBuffer &)
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
void Copy(TAttLine &attline) const
Copy this line attributes to a new TAttLine.
Definition: TAttLine.cxx:175
Int_t DistancetoLine(Int_t px, Int_t py, Double_t xp1, Double_t yp1, Double_t xp2, Double_t yp2)
Compute distance from point px,py to a line.
Definition: TAttLine.cxx:209
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:273
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.
Definition: TAttMarker.cxx:345
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
void Copy(TAttMarker &attmarker) const
Copy this marker attributes to a new TAttMarker.
Definition: TAttMarker.cxx:241
virtual void Streamer(TBuffer &)
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 Bool_t GetTimeDisplay() const
Definition: TAxis.h:127
Bool_t GetRotateTitle() const
Definition: TAxis.h:125
const char * GetTitle() const override
Returns title of object.
Definition: TAxis.h:130
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition: TAxis.cxx:478
Bool_t GetCenterTitle() const
Definition: TAxis.h:115
Bool_t GetNoExponent() const
Definition: TAxis.h:123
virtual void SetTimeDisplay(Int_t value)
Definition: TAxis.h:162
void RotateTitle(Bool_t rotate=kTRUE)
Rotate title by 180 degrees.
Definition: TAxis.h:194
void CenterTitle(Bool_t center=kTRUE)
Center axis title.
Definition: TAxis.h:185
void SetNoExponent(Bool_t noExponent=kTRUE)
Set the NoExponent flag By default, an exponent of the form 10^N is used when the label value are eit...
Definition: TAxis.h:224
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:155
virtual const char * GetTimeFormat() const
Definition: TAxis.h:128
virtual void SetTimeFormat(const char *format="")
Change the format used for time plotting.
Definition: TAxis.cxx:1052
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:37
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
Collection abstract base class.
Definition: TCollection.h:65
virtual Bool_t IsEmpty() const
Definition: TCollection.h:188
TObject * Clone(const char *newname="") const override
Make a clone of an collection using the Streamer facility.
TDirectory::TContext keeps track and restore the current directory.
Definition: TDirectory.h:89
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition: TEnv.cxx:491
1-Dim function class
Definition: TF1.h:213
static TClass * Class()
virtual Int_t GetNpar() const
Definition: TF1.h:486
virtual void SetParent(TObject *p=nullptr)
Definition: TF1.h:670
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set lower and upper limits for parameter ipar.
Definition: TF1.cxx:3545
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:639
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:32
static TClass * Class()
static TClass * Class()
static TClass * Class()
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
virtual Double_t GetPointX(Int_t i) const
Get x value for point i.
Definition: TGraph.cxx:1529
static TClass * Class()
virtual Double_t Integral(Int_t first=0, Int_t last=-1) const
Integrate the TGraph data within a given (index) range.
Definition: TGraph.cxx:1786
Int_t fNpoints
Number of points <= fMaxSize.
Definition: TGraph.h:46
virtual Int_t IsInside(Double_t x, Double_t y) const
Return 1 if the point (x,y) is inside the polygon defined by the graph vertices 0 otherwise.
Definition: TGraph.cxx:1819
virtual void LeastSquareFit(Int_t m, Double_t *a, Double_t xmin=0, Double_t xmax=0)
Least squares polynomial fitting without weights.
Definition: TGraph.cxx:1834
virtual Double_t Chisquare(TF1 *f1, Option_t *option="") const
Return the chisquare of this graph with respect to f1.
Definition: TGraph.cxx:652
void UseCurrentStyle() override
Set current style settings in this graph This function is called when either TCanvas::UseCurrentStyle...
Definition: TGraph.cxx:2548
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2281
Double_t * GetY() const
Definition: TGraph.h:134
virtual Int_t Merge(TCollection *list)
Adds all graphs from the collection to this graph.
Definition: TGraph.cxx:2583
Int_t fMaxSize
!Current dimension of arrays fX and fY
Definition: TGraph.h:45
~TGraph() override
Graph default destructor.
Definition: TGraph.cxx:558
Double_t ** ShrinkAndCopy(Int_t size, Int_t iend)
if size*2 <= fMaxSize allocate new arrays of size points, copy points [0,oend).
Definition: TGraph.cxx:2392
virtual Double_t GetRMS(Int_t axis=1) const
Return RMS of X (axis=1) or Y (axis=2)
Definition: TGraph.cxx:1313
TH1F * fHistogram
Pointer to histogram used for drawing axis.
Definition: TGraph.h:50
void Paint(Option_t *chopt="") override
Draw this graph with its current attributes.
Definition: TGraph.cxx:1955
@ kNotEditable
Bit set if graph is non editable.
Definition: TGraph.h:73
@ kIsHighlight
Bit set if graph is highlight.
Definition: TGraph.h:75
@ kIsSortedX
Graph is sorted in X points.
Definition: TGraph.h:74
@ kClipFrame
Clip to the frame boundary.
Definition: TGraph.h:71
@ kResetHisto
fHistogram must be reset in GetHistogram
Definition: TGraph.h:72
@ kNoStats
Don't draw stats box.
Definition: TGraph.h:70
virtual Double_t GetErrorXlow(Int_t bin) const
It always returns a negative value.
Definition: TGraph.cxx:1361
virtual void MovePoints(Double_t dx, Double_t dy, Bool_t logx=kFALSE, Bool_t logy=kFALSE)
Move all graph points on specified values dx,dy If log argument specified, calculation done in logari...
Definition: TGraph.cxx:2616
virtual Double_t GetErrorYlow(Int_t bin) const
It always returns a negative value.
Definition: TGraph.cxx:1379
virtual void CopyAndRelease(Double_t **newarrays, Int_t ibegin, Int_t iend, Int_t obegin)
Copy points from fX and fY to arrays[0] and arrays[1] or to fX and fY if arrays == 0 and ibegin !...
Definition: TGraph.cxx:737
Double_t GetMinimum() const
Definition: TGraph.h:146
void Print(Option_t *chopt="") const override
Print graph values.
Definition: TGraph.cxx:1991
virtual Double_t * GetEYlow() const
Definition: TGraph.h:140
virtual void SetMaximum(Double_t maximum=-1111)
Set the maximum of the graph.
Definition: TGraph.cxx:2263
virtual Double_t * GetEX() const
Definition: TGraph.h:135
static Bool_t CompareY(const TGraph *gr, Int_t left, Int_t right)
Return kTRUE if fY[left] > fY[right]. Can be used by Sort.
Definition: TGraph.cxx:688
static Bool_t CompareRadius(const TGraph *gr, Int_t left, Int_t right)
Return kTRUE if point number "left"'s distance to origin is bigger than that of point number "right".
Definition: TGraph.cxx:697
virtual Double_t GetErrorYhigh(Int_t bin) const
It always returns a negative value.
Definition: TGraph.cxx:1370
TClass * IsA() const override
Definition: TGraph.h:196
static Bool_t CompareX(const TGraph *gr, Int_t left, Int_t right)
Return kTRUE if fX[left] > fX[right]. Can be used by Sort.
Definition: TGraph.cxx:680
void SaveAs(const char *filename, Option_t *option="") const override
Save the graph as .csv, .tsv or .txt.
Definition: TGraph.cxx:2064
Int_t GetN() const
Definition: TGraph.h:126
TF1 * GetFunction(const char *name) const
Return pointer to function with name.
Definition: TGraph.cxx:1390
virtual void LeastSquareLinearFit(Int_t n, Double_t &a0, Double_t &a1, Int_t &ifail, Double_t xmin=0, Double_t xmax=0)
Least square linear fit without weights.
Definition: TGraph.cxx:1909
Double_t * fY
[fNpoints] array of Y points
Definition: TGraph.h:48
Bool_t CtorAllocate()
In constructors set fNpoints than call this method.
Definition: TGraph.cxx:780
virtual void DrawGraph(Int_t n, const Int_t *x, const Int_t *y, Option_t *option="")
Draw this graph with new attributes.
Definition: TGraph.cxx:860
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Axis_t xmin=0, Axis_t xmax=0)
Fit this graph with function with name fname.
Definition: TGraph.cxx:1228
virtual void Sort(Bool_t(*greater)(const TGraph *, Int_t, Int_t)=&TGraph::CompareX, Bool_t ascending=kTRUE, Int_t low=0, Int_t high=-1111)
Sorts the points of this TGraph using in-place quicksort (see e.g.
Definition: TGraph.cxx:2425
static Bool_t CompareArg(const TGraph *gr, Int_t left, Int_t right)
Return kTRUE if point number "left"'s argument (angle with respect to positive x-axis) is bigger than...
Definition: TGraph.cxx:669
virtual void ComputeRange(Double_t &xmin, Double_t &ymin, Double_t &xmax, Double_t &ymax) const
Compute the x/y range of the points in this graph.
Definition: TGraph.cxx:706
char * GetObjectInfo(Int_t px, Int_t py) const override
Implementation to get information on point of graph at cursor position Adapted from class TH1.
Definition: TGraph.cxx:1572
Double_t ** AllocateArrays(Int_t Narrays, Int_t arraySize)
Allocate arrays.
Definition: TGraph.cxx:590
virtual void Scale(Double_t c1=1., Option_t *option="y")
Multiply the values of a TGraph by a constant c1.
Definition: TGraph.cxx:2198
TList * fFunctions
Pointer to list of functions (fits and user)
Definition: TGraph.h:49
virtual Double_t GetCovariance() const
Return covariance of vectors x,y.
Definition: TGraph.cxx:1282
static void SwapValues(Double_t *arr, Int_t pos1, Int_t pos2)
Swap values.
Definition: TGraph.cxx:2536
void Streamer(TBuffer &) override
Stream an object of class TGraph.
Definition: TGraph.cxx:2460
void Zero(Int_t &k, Double_t AZ, Double_t BZ, Double_t E2, Double_t &X, Double_t &Y, Int_t maxiterations)
Find zero of a continuous function.
Definition: TGraph.cxx:2642
virtual Double_t ** Allocate(Int_t newsize)
Allocate internal data structures for newsize points.
Definition: TGraph.cxx:582
virtual void FitPanel()
Display a GUI panel with all graph fit options.
Definition: TGraph.cxx:1248
void Browse(TBrowser *b) override
Browse.
Definition: TGraph.cxx:626
virtual Bool_t DoMerge(const TGraph *g)
protected function to perform the merge operation of a graph
Definition: TGraph.cxx:2601
virtual void InsertPointBefore(Int_t ipoint, Double_t x, Double_t y)
Insert a new point with coordinates (x,y) before the point number ipoint.
Definition: TGraph.cxx:1727
TList * GetListOfFunctions() const
Definition: TGraph.h:120
void ExecuteEvent(Int_t event, Int_t px, Int_t py) override
Execute action corresponding to one event.
Definition: TGraph.cxx:1027
Double_t * GetX() const
Definition: TGraph.h:133
virtual Double_t Eval(Double_t x, TSpline *spline=nullptr, Option_t *option="") const
Interpolate points in this graph at x using a TSpline.
Definition: TGraph.cxx:924
virtual void InitExpo(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for an exponential.
Definition: TGraph.cxx:1641
virtual Int_t RemovePoint()
Delete point close to the mouse position.
Definition: TGraph.cxx:2013
virtual void InitGaus(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for a gaussian.
Definition: TGraph.cxx:1603
virtual Bool_t IsHighlight() const
Definition: TGraph.h:161
virtual void Apply(TF1 *f)
Apply function f to all the data points f may be a 1-D function TF1 or 2-d function TF2 The Y values ...
Definition: TGraph.cxx:613
void SetName(const char *name="") override
Set graph name.
Definition: TGraph.cxx:2320
virtual void SetHighlight(Bool_t set=kTRUE)
Set highlight (enable/disble) mode for the graph by default highlight mode is disable.
Definition: TGraph.cxx:2250
virtual void SwapPoints(Int_t pos1, Int_t pos2)
Swap points.
Definition: TGraph.cxx:2527
void Draw(Option_t *chopt="") override
Draw this graph with its current attributes.
Definition: TGraph.cxx:806
TAxis * GetXaxis() const
Get x axis of the graph.
Definition: TGraph.cxx:1551
Bool_t GetEditable() const
Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
Definition: TGraph.cxx:2231
virtual Double_t * GetEXhigh() const
Definition: TGraph.h:137
virtual Double_t GetCorrelationFactor() const
Return graph correlation factor.
Definition: TGraph.cxx:1270
virtual void FillZero(Int_t begin, Int_t end, Bool_t from_ctor=kTRUE)
Set zero values for point arrays in the range [begin, end) Should be redefined in descendant classes.
Definition: TGraph.cxx:1074
virtual Double_t * GetEYhigh() const
Definition: TGraph.h:139
Double_t ** ExpandAndCopy(Int_t size, Int_t iend)
if size > fMaxSize allocate new arrays of 2*size points and copy iend first points.
Definition: TGraph.cxx:1060
virtual void Expand(Int_t newsize)
If array sizes <= newsize, expand storage to 2*newsize.
Definition: TGraph.cxx:1036
virtual Double_t GetMean(Int_t axis=1) const
Return mean value of X (axis=1) or Y (axis=2)
Definition: TGraph.cxx:1298
Double_t * fX
[fNpoints] array of X points
Definition: TGraph.h:47
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:1982
TAxis * GetYaxis() const
Get y axis of the graph.
Definition: TGraph.cxx:1561
TObject * FindObject(const char *name) const override
Search object named name in the list of functions.
Definition: TGraph.cxx:1083
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TGraph.cxx:2370
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1402
virtual Double_t GetErrorY(Int_t bin) const
It always returns a negative value. Real implementation in TGraphErrors.
Definition: TGraph.cxx:1343
virtual Double_t * GetEY() const
Definition: TGraph.h:136
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
Definition: TGraph.cxx:2116
virtual Double_t GetPointY(Int_t i) const
Get y value for point i.
Definition: TGraph.cxx:1540
Double_t fMinimum
Minimum value for plotting along y.
Definition: TGraph.h:51
void PaintGraph(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
Draw the (x,y) as a graph.
Definition: TGraph.cxx:1964
void SetTitle(const char *title="") override
Change (i.e.
Definition: TGraph.cxx:2336
virtual Int_t InsertPoint()
Insert a new point at the mouse position.
Definition: TGraph.cxx:1682
virtual Double_t * GetEXlow() const
Definition: TGraph.h:138
void RecursiveRemove(TObject *obj) override
Recursively remove object from the list of functions.
Definition: TGraph.cxx:2001
virtual void SetPointY(Int_t i, Double_t y)
Set y value for point i.
Definition: TGraph.cxx:2313
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a graph.
Definition: TGraph.cxx:850
virtual void DrawPanel()
Display a panel with all graph drawing options.
Definition: TGraph.cxx:903
void PaintGrapHist(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
Draw the (x,y) as a histogram.
Definition: TGraph.cxx:1973
void SetNameTitle(const char *name="", const char *title="") override
Set graph name and title.
Definition: TGraph.cxx:2356
virtual void SetPointX(Int_t i, Double_t x)
Set x value for point i.
Definition: TGraph.cxx:2305
virtual void Set(Int_t n)
Set number of points in the graph Existing coordinates are preserved New coordinates above fNpoints a...
Definition: TGraph.cxx:2216
virtual Int_t GetPoint(Int_t i, Double_t &x, Double_t &y) const
Get x and y values for point number i.
Definition: TGraph.cxx:1518
virtual void SetEditable(Bool_t editable=kTRUE)
if editable=kFALSE, the graph cannot be modified with the mouse by default a TGraph is editable
Definition: TGraph.cxx:2240
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum of the graph.
Definition: TGraph.cxx:2272
TGraph()
Graph default constructor.
Definition: TGraph.cxx:105
virtual Bool_t CopyPoints(Double_t **newarrays, Int_t ibegin, Int_t iend, Int_t obegin)
Copy points from fX and fY to arrays[0] and arrays[1] or to fX and fY if arrays == 0 and ibegin !...
Definition: TGraph.cxx:754
Double_t fMaximum
Maximum value for plotting along y.
Definition: TGraph.h:52
virtual Double_t GetErrorXhigh(Int_t bin) const
It always returns a negative value.
Definition: TGraph.cxx:1352
TGraph & operator=(const TGraph &)
Equal operator for this graph.
Definition: TGraph.cxx:229
virtual void InitPolynom(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for a polynom.
Definition: TGraph.cxx:1662
virtual Double_t GetErrorX(Int_t bin) const
It always returns a negative value. Real implementation in TGraphErrors.
Definition: TGraph.cxx:1335
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:574
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition: TH1.cxx:8804
void SetTitle(const char *title) override
Change (i.e.
Definition: TH1.cxx:6700
void UseCurrentStyle() override
Copy current attributes from/to current style.
Definition: TH1.cxx:7380
@ kNoStats
Don't draw stats box.
Definition: TH1.h:163
TAxis * GetXaxis()
Definition: TH1.h:319
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:397
TAxis * GetYaxis()
Definition: TH1.h:320
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:398
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
Definition: TH1.cxx:7147
void SetName(const char *name) override
Change the name of this histogram.
Definition: TH1.cxx:8827
TObject * Clone(const char *newname="") const override
Make a complete copy of the underlying object.
Definition: TH1.cxx:2727
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:578
void RecursiveRemove(TObject *obj) override
Remove object from this collection and recursively remove the object from all other objects (and coll...
Definition: TList.cxx:764
TObject * Remove(TObject *obj) override
Remove object from the list.
Definition: TList.cxx:822
TObject * First() const override
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:659
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
void Streamer(TBuffer &) override
Stream an object of class TObject.
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
TString fTitle
Definition: TNamed.h:33
TString fName
Definition: TNamed.h:32
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
Mother of all ROOT objects.
Definition: TObject.h:41
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:433
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:201
virtual void UseCurrentStyle()
Set current style settings in this object This function is called when either TCanvas::UseCurrentStyl...
Definition: TObject.cxx:789
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:177
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TObject.cxx:739
virtual void SaveAs(const char *filename="", Option_t *option="") const
Save this object in the file specified by filename.
Definition: TObject.cxx:674
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:768
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:519
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:963
void MakeZombie()
Definition: TObject.h:53
void ResetBit(UInt_t f)
Definition: TObject.h:200
@ kCanDelete
if object in a list can be deleted
Definition: TObject.h:62
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition: TObject.h:72
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:937
Longptr_t ExecPlugin(int nargs, const T &... params)
Int_t LoadPlugin()
Load the plugin library for this handler.
Class to create third splines to interpolate knots Arbitrary conditions can be introduced for first a...
Definition: TSpline.h:201
Base class for spline implementation containing the Draw/Paint methods.
Definition: TSpline.h:31
virtual Double_t Eval(Double_t x) const =0
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:1155
Int_t Atoi() const
Return integer value of string.
Definition: TString.cxx:1951
Bool_t EndsWith(const char *pat, ECaseCompare cmp=kExact) const
Return true if string ends with the specified string.
Definition: TString.cxx:2207
Double_t Atof() const
Return floating-point value contained in string.
Definition: TString.cxx:2017
Bool_t IsFloat() const
Returns kTRUE if string contains a floating point or integer number.
Definition: TString.cxx:1821
TString & Replace(Ssiz_t pos, Ssiz_t n, const char *s)
Definition: TString.h:682
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:1793
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1168
Bool_t IsNull() const
Definition: TString.h:407
TString & Remove(Ssiz_t pos)
Definition: TString.h:673
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:2341
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2319
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition: TString.h:639
void SetHistFillColor(Color_t color=1)
Definition: TStyle.h:363
Color_t GetHistLineColor() const
Definition: TStyle.h:225
Bool_t IsReading() const
Definition: TStyle.h:283
void SetHistLineStyle(Style_t styl=0)
Definition: TStyle.h:366
Style_t GetHistFillStyle() const
Definition: TStyle.h:226
Color_t GetHistFillColor() const
Definition: TStyle.h:224
void SetHistLineColor(Color_t color=1)
Definition: TStyle.h:364
Style_t GetHistLineStyle() const
Definition: TStyle.h:227
void SetHistFillStyle(Style_t styl=0)
Definition: TStyle.h:365
Width_t GetHistLineWidth() const
Definition: TStyle.h:228
void SetHistLineWidth(Width_t width=1)
Definition: TStyle.h:367
virtual Bool_t ExpandPathName(TString &path)
Expand a pathname getting rid of special shell characters like ~.
Definition: TSystem.cxx:1274
TVectorT.
Definition: TVectorT.h:27
Int_t GetNrows() const
Definition: TVectorT.h:75
Int_t GetLwb() const
Definition: TVectorT.h:73
Abstract Base Class for Fitting.
static TVirtualFitter * GetFitter()
static: return the current Fitter
virtual TObject * GetUserFunc() const
Abstract interface to a histogram painter.
virtual void DrawPanelHelper(TGraph *theGraph)=0
virtual void ExecuteEventHelper(TGraph *theGraph, Int_t event, Int_t px, Int_t py)=0
virtual void SetHighlight(TGraph *theGraph)=0
virtual void PaintGrapHist(TGraph *theGraph, Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)=0
virtual Int_t DistancetoPrimitiveHelper(TGraph *theGraph, Int_t px, Int_t py)=0
virtual void PaintGraph(TGraph *theGraph, Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)=0
virtual void PaintStats(TGraph *theGraph, TF1 *fit)=0
virtual void PaintHelper(TGraph *theGraph, Option_t *option)=0
static TVirtualGraphPainter * GetPainter()
Static function returning a pointer to the current graph painter.
TLine * line
Double_t y[n]
Definition: legend1.C:17
return c1
Definition: legend1.C:41
Double_t x[n]
Definition: legend1.C:17
Double_t ey[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
TGraphErrors * gr
Definition: legend1.C:25
Double_t ex[n]
Definition: legend1.C:17
TF1 * f1
Definition: legend1.C:11
def fit(model, train_loader, val_loader, num_epochs, batch_size, optimizer, criterion, save_best, scheduler)
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:965
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition: HFitImpl.cxx:684
double Chisquare(const TH1 &h1, TF1 &f1, bool useRange, bool usePL=false)
compute the chi2 value for an histogram given a function (see TH1::Chisquare for the documentation)
Definition: HFitImpl.cxx:1019
static constexpr double s
static constexpr double ps
Bool_t IsInside(T xp, T yp, Int_t np, T *x, T *y)
Function which returns kTRUE if point xp,yp lies inside the polygon defined by the np points in array...
Definition: TMath.h:1231
Double_t Exp(Double_t x)
Returns the base-e exponential function of x, which is e raised to the power x.
Definition: TMath.h:707
Double_t ATan2(Double_t y, Double_t x)
Returns the principal value of the arc tangent of y/x, expressed in radians.
Definition: TMath.h:644
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition: TMath.h:754
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition: TMath.h:660
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition: TMathBase.h:198
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Sort the n elements of the array a of generic templated type Element.
Definition: TMathBase.h:431
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Binary search in an array of n values to locate value.
Definition: TMathBase.h:347
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123
Definition: first.py:1
TMarker m
Definition: textangle.C:8
TLine l
Definition: textangle.C:4
TArc a
Definition: textangle.C:12
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345