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 (!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 char *linear;
1231 linear = (char*) strstr(fname, "++");
1232 if (linear) {
1233 TF1 f1(fname, fname, xmin, xmax);
1234 return Fit(&f1, option, "", xmin, xmax);
1235 }
1236 TF1 * f1 = (TF1*)gROOT->GetFunction(fname);
1237 if (!f1) {
1238 Printf("Unknown function: %s", fname);
1239 return -1;
1240 }
1241 return Fit(f1, option, "", xmin, xmax);
1242}
1243
1244////////////////////////////////////////////////////////////////////////////////
1245/// Display a GUI panel with all graph fit options.
1246///
1247/// See class TFitEditor for example
1248
1250{
1251 if (!gPad)
1252 gROOT->MakeDefCanvas();
1253
1254 if (!gPad) {
1255 Error("FitPanel", "Unable to create a default canvas");
1256 return;
1257 }
1258
1259 // use plugin manager to create instance of TFitEditor
1260 TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
1261 if (handler && handler->LoadPlugin() != -1) {
1262 if (handler->ExecPlugin(2, gPad, this) == 0)
1263 Error("FitPanel", "Unable to crate the FitPanel");
1264 } else
1265 Error("FitPanel", "Unable to find the FitPanel plug-in");
1266}
1267
1268////////////////////////////////////////////////////////////////////////////////
1269/// Return graph correlation factor
1270
1272{
1273 Double_t rms1 = GetRMS(1);
1274 if (rms1 == 0) return 0;
1275 Double_t rms2 = GetRMS(2);
1276 if (rms2 == 0) return 0;
1277 return GetCovariance() / rms1 / rms2;
1278}
1279
1280////////////////////////////////////////////////////////////////////////////////
1281/// Return covariance of vectors x,y
1282
1284{
1285 if (fNpoints <= 0) return 0;
1286 Double_t sum = fNpoints, sumx = 0, sumy = 0, sumxy = 0;
1287
1288 for (Int_t i = 0; i < fNpoints; i++) {
1289 sumx += fX[i];
1290 sumy += fY[i];
1291 sumxy += fX[i] * fY[i];
1292 }
1293 return sumxy / sum - sumx / sum * sumy / sum;
1294}
1295
1296////////////////////////////////////////////////////////////////////////////////
1297/// Return mean value of X (axis=1) or Y (axis=2)
1298
1300{
1301 if (axis < 1 || axis > 2) return 0;
1302 if (fNpoints <= 0) return 0;
1303 Double_t sumx = 0;
1304 for (Int_t i = 0; i < fNpoints; i++) {
1305 if (axis == 1) sumx += fX[i];
1306 else sumx += fY[i];
1307 }
1308 return sumx / fNpoints;
1309}
1310
1311////////////////////////////////////////////////////////////////////////////////
1312/// Return RMS of X (axis=1) or Y (axis=2)
1313
1315{
1316 if (axis < 1 || axis > 2) return 0;
1317 if (fNpoints <= 0) return 0;
1318 Double_t sumx = 0, sumx2 = 0;
1319 for (Int_t i = 0; i < fNpoints; i++) {
1320 if (axis == 1) {
1321 sumx += fX[i];
1322 sumx2 += fX[i] * fX[i];
1323 } else {
1324 sumx += fY[i];
1325 sumx2 += fY[i] * fY[i];
1326 }
1327 }
1328 Double_t x = sumx / fNpoints;
1329 Double_t rms2 = TMath::Abs(sumx2 / fNpoints - x * x);
1330 return TMath::Sqrt(rms2);
1331}
1332
1333////////////////////////////////////////////////////////////////////////////////
1334/// It always returns a negative value. Real implementation in TGraphErrors
1335
1337{
1338 return -1;
1339}
1340
1341////////////////////////////////////////////////////////////////////////////////
1342/// It always returns a negative value. Real implementation in TGraphErrors
1343
1345{
1346 return -1;
1347}
1348
1349////////////////////////////////////////////////////////////////////////////////
1350/// It always returns a negative value. Real implementation in TGraphErrors
1351/// and TGraphAsymmErrors
1352
1354{
1355 return -1;
1356}
1357
1358////////////////////////////////////////////////////////////////////////////////
1359/// It always returns a negative value. Real implementation in TGraphErrors
1360/// and TGraphAsymmErrors
1361
1363{
1364 return -1;
1365}
1366
1367////////////////////////////////////////////////////////////////////////////////
1368/// It always returns a negative value. Real implementation in TGraphErrors
1369/// and TGraphAsymmErrors
1370
1372{
1373 return -1;
1374}
1375
1376////////////////////////////////////////////////////////////////////////////////
1377/// It always returns a negative value. Real implementation in TGraphErrors
1378/// and TGraphAsymmErrors
1379
1381{
1382 return -1;
1383}
1384
1385////////////////////////////////////////////////////////////////////////////////
1386/// Return pointer to function with name.
1387///
1388/// Functions such as TGraph::Fit store the fitted function in the list of
1389/// functions of this graph.
1390
1391TF1 *TGraph::GetFunction(const char *name) const
1392{
1393 if (!fFunctions) return nullptr;
1394 return (TF1*)fFunctions->FindObject(name);
1395}
1396
1397////////////////////////////////////////////////////////////////////////////////
1398/// Returns a pointer to the histogram used to draw the axis
1399/// Takes into account the two following cases.
1400/// 1. option 'A' was specified in TGraph::Draw. Return fHistogram
1401/// 2. user had called TPad::DrawFrame. return pointer to hframe histogram
1402
1404{
1405 Double_t rwxmin, rwxmax, rwymin, rwymax, maximum, minimum, dx, dy;
1406 Double_t uxmin, uxmax;
1407
1408 ComputeRange(rwxmin, rwymin, rwxmax, rwymax); //this is redefined in TGraphErrors
1409
1410 // (if fHistogram exist) && (if the log scale is on) &&
1411 // (if the computed range minimum is > 0) && (if the fHistogram minimum is zero)
1412 // then it means fHistogram limits have been computed in linear scale
1413 // therefore they might be too strict and cut some points. In that case the
1414 // fHistogram limits should be recomputed ie: the existing fHistogram
1415 // should not be returned.
1416 TH1F *historg = nullptr;
1417 if (fHistogram) {
1418 if (!TestBit(kResetHisto)) {
1419 if (gPad && gPad->GetLogx()) {
1420 if (rwxmin <= 0 || fHistogram->GetXaxis()->GetXmin() != 0) return fHistogram;
1421 } else if (gPad && gPad->GetLogy()) {
1422 if (rwymin <= 0 || fHistogram->GetMinimum() != 0) return fHistogram;
1423 } else {
1424 return fHistogram;
1425 }
1426 } else {
1427 const_cast <TGraph*>(this)->ResetBit(kResetHisto);
1428 }
1429 historg = fHistogram;
1430 }
1431
1432 if (rwxmin == rwxmax) rwxmax += 1.;
1433 if (rwymin == rwymax) rwymax += 1.;
1434 dx = 0.1 * (rwxmax - rwxmin);
1435 dy = 0.1 * (rwymax - rwymin);
1436 uxmin = rwxmin - dx;
1437 uxmax = rwxmax + dx;
1438 minimum = rwymin - dy;
1439 maximum = rwymax + dy;
1440
1441 if (fMinimum != -1111) minimum = fMinimum;
1442 if (fMaximum != -1111) maximum = fMaximum;
1443
1444 // the graph is created with at least as many channels as there are points
1445 // to permit zooming on the full range
1446 if (uxmin < 0 && rwxmin >= 0) {
1447 if (gPad && gPad->GetLogx()) uxmin = 0.9 * rwxmin;
1448 else uxmin = 0;
1449 }
1450 if (uxmax > 0 && rwxmax <= 0) {
1451 if (gPad && gPad->GetLogx()) uxmax = 1.1 * rwxmax;
1452 else uxmax = 0;
1453 }
1454
1455 if (minimum < 0 && rwymin >= 0) minimum = 0.9 * rwymin;
1456
1457 if (minimum <= 0 && gPad && gPad->GetLogy()) minimum = 0.001 * maximum;
1458 if (uxmin <= 0 && gPad && gPad->GetLogx()) {
1459 if (uxmax > 1000) uxmin = 1;
1460 else uxmin = 0.001 * uxmax;
1461 }
1462
1463 rwxmin = uxmin;
1464 rwxmax = uxmax;
1465 Int_t npt = 100;
1466 if (fNpoints > npt) npt = fNpoints;
1467 const char *gname = GetName();
1468 if (!gname[0]) gname = "Graph";
1469 // do not add the histogram to gDirectory
1470 // use local TDirectory::TContect that will set temporarly gDirectory to a nullptr and
1471 // will avoid that histogram is added in the global directory
1472 {
1473 TDirectory::TContext ctx(nullptr);
1474 ((TGraph*)this)->fHistogram = new TH1F(gname, GetTitle(), npt, rwxmin, rwxmax);
1475 }
1476 if (!fHistogram) return nullptr;
1477 fHistogram->SetMinimum(minimum);
1479 fHistogram->SetMaximum(maximum);
1480 fHistogram->GetYaxis()->SetLimits(minimum, maximum);
1481 // Restore the axis attributes if needed
1482 if (historg) {
1483 fHistogram->GetXaxis()->SetTitle(historg->GetXaxis()->GetTitle());
1496
1497 fHistogram->GetYaxis()->SetTitle(historg->GetYaxis()->GetTitle());
1510 delete historg;
1511 }
1512 return fHistogram;
1513}
1514
1515////////////////////////////////////////////////////////////////////////////////
1516/// Get x and y values for point number i.
1517/// The function returns -1 in case of an invalid request or the point number otherwise
1518
1520{
1521 if (i < 0 || i >= fNpoints || !fX || !fY) return -1;
1522 x = fX[i];
1523 y = fY[i];
1524 return i;
1525}
1526
1527////////////////////////////////////////////////////////////////////////////////
1528/// Get x value for point i.
1529
1531{
1532 if (i < 0 || i >= fNpoints || !fX)
1533 return -1.;
1534
1535 return fX[i];
1536}
1537
1538////////////////////////////////////////////////////////////////////////////////
1539/// Get y value for point i.
1540
1542{
1543 if (i < 0 || i >= fNpoints || !fY)
1544 return -1.;
1545
1546 return fY[i];
1547}
1548
1549////////////////////////////////////////////////////////////////////////////////
1550/// Get x axis of the graph.
1551
1553{
1554 TH1 *h = GetHistogram();
1555 if (!h) return 0;
1556 return h->GetXaxis();
1557}
1558
1559////////////////////////////////////////////////////////////////////////////////
1560/// Get y axis of the graph.
1561
1563{
1564 TH1 *h = GetHistogram();
1565 if (!h) return 0;
1566 return h->GetYaxis();
1567}
1568
1569////////////////////////////////////////////////////////////////////////////////
1570/// Implementation to get information on point of graph at cursor position
1571/// Adapted from class TH1
1572
1574{
1575 // localize point
1576 Int_t ipoint = -2;
1577 Int_t i;
1578 // start with a small window (in case the mouse is very close to one point)
1579 for (i = 0; i < fNpoints; i++) {
1580 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
1581 Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
1582
1583 if (dpx * dpx + dpy * dpy < 25) {
1584 ipoint = i;
1585 break;
1586 }
1587 }
1588
1589 Double_t x = gPad->PadtoX(gPad->AbsPixeltoX(px));
1590 Double_t y = gPad->PadtoY(gPad->AbsPixeltoY(py));
1591
1592 if (ipoint == -2)
1593 return Form("x=%g, y=%g", x, y);
1594
1595 Double_t xval = fX[ipoint];
1596 Double_t yval = fY[ipoint];
1597
1598 return Form("x=%g, y=%g, point=%d, xval=%g, yval=%g", x, y, ipoint, xval, yval);
1599}
1600
1601////////////////////////////////////////////////////////////////////////////////
1602/// Compute Initial values of parameters for a gaussian.
1603
1605{
1606 Double_t allcha, sumx, sumx2, x, val, rms, mean;
1607 Int_t bin;
1608 const Double_t sqrtpi = 2.506628;
1609
1610 // Compute mean value and RMS of the graph in the given range
1611 if (xmax <= xmin) {
1612 xmin = fX[0];
1613 xmax = fX[fNpoints-1];
1614 }
1615 Int_t np = 0;
1616 allcha = sumx = sumx2 = 0;
1617 for (bin = 0; bin < fNpoints; bin++) {
1618 x = fX[bin];
1619 if (x < xmin || x > xmax) continue;
1620 np++;
1621 val = fY[bin];
1622 sumx += val * x;
1623 sumx2 += val * x * x;
1624 allcha += val;
1625 }
1626 if (np == 0 || allcha == 0) return;
1627 mean = sumx / allcha;
1628 rms = TMath::Sqrt(sumx2 / allcha - mean * mean);
1629 Double_t binwidx = TMath::Abs((xmax - xmin) / np);
1630 if (rms == 0) rms = 1;
1632 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1633 f1->SetParameter(0, binwidx * allcha / (sqrtpi * rms));
1634 f1->SetParameter(1, mean);
1635 f1->SetParameter(2, rms);
1636 f1->SetParLimits(2, 0, 10 * rms);
1637}
1638
1639////////////////////////////////////////////////////////////////////////////////
1640/// Compute Initial values of parameters for an exponential.
1641
1643{
1644 Double_t constant, slope;
1645 Int_t ifail;
1646 if (xmax <= xmin) {
1647 xmin = fX[0];
1648 xmax = fX[fNpoints-1];
1649 }
1650 Int_t nchanx = fNpoints;
1651
1652 LeastSquareLinearFit(-nchanx, constant, slope, ifail, xmin, xmax);
1653
1655 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1656 f1->SetParameter(0, constant);
1657 f1->SetParameter(1, slope);
1658}
1659
1660////////////////////////////////////////////////////////////////////////////////
1661/// Compute Initial values of parameters for a polynom.
1662
1664{
1665 Double_t fitpar[25];
1666
1668 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1669 Int_t npar = f1->GetNpar();
1670 if (xmax <= xmin) {
1671 xmin = fX[0];
1672 xmax = fX[fNpoints-1];
1673 }
1674
1675 LeastSquareFit(npar, fitpar, xmin, xmax);
1676
1677 for (Int_t i = 0; i < npar; i++) f1->SetParameter(i, fitpar[i]);
1678}
1679
1680////////////////////////////////////////////////////////////////////////////////
1681/// Insert a new point at the mouse position
1682
1684{
1685 Int_t px = gPad->GetEventX();
1686 Int_t py = gPad->GetEventY();
1687
1688 //localize point where to insert
1689 Int_t ipoint = -2;
1690 Int_t i, d = 0;
1691 // start with a small window (in case the mouse is very close to one point)
1692 for (i = 0; i < fNpoints - 1; i++) {
1693 d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1694 if (d < 5) {
1695 ipoint = i + 1;
1696 break;
1697 }
1698 }
1699 if (ipoint == -2) {
1700 //may be we are far from one point, try again with a larger window
1701 for (i = 0; i < fNpoints - 1; i++) {
1702 d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1703 if (d < 10) {
1704 ipoint = i + 1;
1705 break;
1706 }
1707 }
1708 }
1709 if (ipoint == -2) {
1710 //distinguish between first and last point
1711 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[0]));
1712 Int_t dpy = py - gPad->YtoAbsPixel(gPad->XtoPad(fY[0]));
1713 if (dpx * dpx + dpy * dpy < 25) ipoint = 0;
1714 else ipoint = fNpoints;
1715 }
1716
1717
1718 InsertPointBefore(ipoint, gPad->AbsPixeltoX(px), gPad->AbsPixeltoY(py));
1719
1720 gPad->Modified();
1721 return ipoint;
1722}
1723
1724
1725////////////////////////////////////////////////////////////////////////////////
1726/// Insert a new point with coordinates (x,y) before the point number `ipoint`.
1727
1729{
1730 if (ipoint < 0) {
1731 Error("TGraph", "Inserted point index should be >= 0");
1732 return;
1733 }
1734
1735 if (ipoint > fNpoints) {
1736 Error("TGraph", "Inserted point index should be <= %d", fNpoints);
1737 return;
1738 }
1739
1740 if (ipoint == fNpoints) {
1741 SetPoint(ipoint, x, y);
1742 return;
1743 }
1744
1745 Double_t **ps = ExpandAndCopy(fNpoints + 1, ipoint);
1746 CopyAndRelease(ps, ipoint, fNpoints++, ipoint + 1);
1747
1748 // To avoid redefinitions in descendant classes
1749 FillZero(ipoint, ipoint + 1);
1750
1751 fX[ipoint] = x;
1752 fY[ipoint] = y;
1753}
1754
1755
1756////////////////////////////////////////////////////////////////////////////////
1757/// Integrate the TGraph data within a given (index) range.
1758/// Note that this function computes the area of the polygon enclosed by the points of the TGraph.
1759/// The polygon segments, which are defined by the points of the TGraph, do not need to form a closed polygon,
1760/// since the last polygon segment, which closes the polygon, is taken as the line connecting the last TGraph point
1761/// with the first one. It is clear that the order of the point is essential in defining the polygon.
1762/// Also note that the segments should not intersect.
1763///
1764/// NB:
1765/// - if last=-1 (default) last is set to the last point.
1766/// - if (first <0) the first point (0) is taken.
1767///
1768/// ### Method:
1769///
1770/// There are many ways to calculate the surface of a polygon. It all depends on what kind of data
1771/// you have to deal with. The most evident solution would be to divide the polygon in triangles and
1772/// calculate the surface of them. But this can quickly become complicated as you will have to test
1773/// every segments of every triangles and check if they are intersecting with a current polygon's
1774/// segment or if it goes outside the polygon. Many calculations that would lead to many problems...
1775///
1776/// ### The solution (implemented by R.Brun)
1777/// Fortunately for us, there is a simple way to solve this problem, as long as the polygon's
1778/// segments don't intersect.
1779/// It takes the x coordinate of the current vertex and multiply it by the y coordinate of the next
1780/// vertex. Then it subtracts from it the result of the y coordinate of the current vertex multiplied
1781/// by the x coordinate of the next vertex. Then divide the result by 2 to get the surface/area.
1782///
1783/// ### Sources
1784/// - http://forums.wolfram.com/mathgroup/archive/1998/Mar/msg00462.html
1785/// - http://stackoverflow.com/questions/451426/how-do-i-calculate-the-surface-area-of-a-2d-polygon
1786
1788{
1789 if (first < 0) first = 0;
1790 if (last < 0) last = fNpoints - 1;
1791 if (last >= fNpoints) last = fNpoints - 1;
1792 if (first >= last) return 0;
1793 Int_t np = last - first + 1;
1794 Double_t sum = 0.0;
1795 //for(Int_t i=first;i<=last;i++) {
1796 // Int_t j = first + (i-first+1)%np;
1797 // sum += TMath::Abs(fX[i]*fY[j]);
1798 // sum -= TMath::Abs(fY[i]*fX[j]);
1799 //}
1800 for (Int_t i = first; i <= last; i++) {
1801 Int_t j = first + (i - first + 1) % np;
1802 sum += (fY[i] + fY[j]) * (fX[j] - fX[i]);
1803 }
1804 return 0.5 * TMath::Abs(sum);
1805}
1806
1807////////////////////////////////////////////////////////////////////////////////
1808/// Return 1 if the point (x,y) is inside the polygon defined by
1809/// the graph vertices 0 otherwise.
1810///
1811/// Algorithm:
1812///
1813/// The loop is executed with the end-point coordinates of a line segment
1814/// (X1,Y1)-(X2,Y2) and the Y-coordinate of a horizontal line.
1815/// The counter inter is incremented if the line (X1,Y1)-(X2,Y2) intersects
1816/// the horizontal line. In this case XINT is set to the X-coordinate of the
1817/// intersection point. If inter is an odd number, then the point x,y is within
1818/// the polygon.
1819
1821{
1822 return (Int_t)TMath::IsInside(x, y, fNpoints, fX, fY);
1823}
1824
1825////////////////////////////////////////////////////////////////////////////////
1826/// Least squares polynomial fitting without weights.
1827///
1828/// \param [in] m number of parameters
1829/// \param [in] a array of parameters
1830/// \param [in] xmin 1st point number to fit (default =0)
1831/// \param [in] xmax last point number to fit (default=fNpoints-1)
1832///
1833/// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
1834
1836{
1837 const Double_t zero = 0.;
1838 const Double_t one = 1.;
1839 const Int_t idim = 20;
1840
1841 Double_t b[400] /* was [20][20] */;
1842 Int_t i, k, l, ifail;
1843 Double_t power;
1844 Double_t da[20], xk, yk;
1845 Int_t n = fNpoints;
1846 if (xmax <= xmin) {
1847 xmin = fX[0];
1848 xmax = fX[fNpoints-1];
1849 }
1850
1851 if (m <= 2) {
1852 LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
1853 return;
1854 }
1855 if (m > idim || m > n) return;
1856 da[0] = zero;
1857 for (l = 2; l <= m; ++l) {
1858 b[l-1] = zero;
1859 b[m + l*20 - 21] = zero;
1860 da[l-1] = zero;
1861 }
1862 Int_t np = 0;
1863 for (k = 0; k < fNpoints; ++k) {
1864 xk = fX[k];
1865 if (xk < xmin || xk > xmax) continue;
1866 np++;
1867 yk = fY[k];
1868 power = one;
1869 da[0] += yk;
1870 for (l = 2; l <= m; ++l) {
1871 power *= xk;
1872 b[l-1] += power;
1873 da[l-1] += power * yk;
1874 }
1875 for (l = 2; l <= m; ++l) {
1876 power *= xk;
1877 b[m + l*20 - 21] += power;
1878 }
1879 }
1880 b[0] = Double_t(np);
1881 for (i = 3; i <= m; ++i) {
1882 for (k = i; k <= m; ++k) {
1883 b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
1884 }
1885 }
1886 H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
1887
1888 if (ifail < 0) {
1889 a[0] = fY[0];
1890 for (i = 1; i < m; ++i) a[i] = 0;
1891 return;
1892 }
1893 for (i = 0; i < m; ++i) a[i] = da[i];
1894}
1895
1896////////////////////////////////////////////////////////////////////////////////
1897/// Least square linear fit without weights.
1898///
1899/// Fit a straight line (a0 + a1*x) to the data in this graph.
1900///
1901/// \param [in] ndata if ndata<0, fits the logarithm of the graph (used in InitExpo() to set
1902/// the initial parameter values for a fit with exponential function.
1903/// \param [in] a0 constant
1904/// \param [in] a1 slope
1905/// \param [in] ifail return parameter indicating the status of the fit (ifail=0, fit is OK)
1906/// \param [in] xmin, xmax fitting range
1907///
1908/// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
1909
1911{
1912 Double_t xbar, ybar, x2bar;
1913 Int_t i;
1914 Double_t xybar;
1915 Double_t fn, xk, yk;
1916 Double_t det;
1917 if (xmax <= xmin) {
1918 xmin = fX[0];
1919 xmax = fX[fNpoints-1];
1920 }
1921
1922 ifail = -2;
1923 xbar = ybar = x2bar = xybar = 0;
1924 Int_t np = 0;
1925 for (i = 0; i < fNpoints; ++i) {
1926 xk = fX[i];
1927 if (xk < xmin || xk > xmax) continue;
1928 np++;
1929 yk = fY[i];
1930 if (ndata < 0) {
1931 if (yk <= 0) yk = 1e-9;
1932 yk = TMath::Log(yk);
1933 }
1934 xbar += xk;
1935 ybar += yk;
1936 x2bar += xk * xk;
1937 xybar += xk * yk;
1938 }
1939 fn = Double_t(np);
1940 det = fn * x2bar - xbar * xbar;
1941 ifail = -1;
1942 if (det <= 0) {
1943 if (fn > 0) a0 = ybar / fn;
1944 else a0 = 0;
1945 a1 = 0;
1946 return;
1947 }
1948 ifail = 0;
1949 a0 = (x2bar * ybar - xbar * xybar) / det;
1950 a1 = (fn * xybar - xbar * ybar) / det;
1951}
1952
1953////////////////////////////////////////////////////////////////////////////////
1954/// Draw this graph with its current attributes.
1955
1957{
1959 if (painter) painter->PaintHelper(this, option);
1960}
1961
1962////////////////////////////////////////////////////////////////////////////////
1963/// Draw the (x,y) as a graph.
1964
1965void TGraph::PaintGraph(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
1966{
1968 if (painter) painter->PaintGraph(this, npoints, x, y, chopt);
1969}
1970
1971////////////////////////////////////////////////////////////////////////////////
1972/// Draw the (x,y) as a histogram.
1973
1974void TGraph::PaintGrapHist(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
1975{
1977 if (painter) painter->PaintGrapHist(this, npoints, x, y, chopt);
1978}
1979
1980////////////////////////////////////////////////////////////////////////////////
1981/// Draw the stats
1982
1984{
1986 if (painter) painter->PaintStats(this, fit);
1987}
1988
1989////////////////////////////////////////////////////////////////////////////////
1990/// Print graph values.
1991
1993{
1994 for (Int_t i = 0; i < fNpoints; i++) {
1995 printf("x[%d]=%g, y[%d]=%g\n", i, fX[i], i, fY[i]);
1996 }
1997}
1998
1999////////////////////////////////////////////////////////////////////////////////
2000/// Recursively remove object from the list of functions
2001
2003{
2004 if (fFunctions) {
2006 }
2007 if (fHistogram == obj)
2008 fHistogram = nullptr;
2009}
2010
2011////////////////////////////////////////////////////////////////////////////////
2012/// Delete point close to the mouse position
2013
2015{
2016 Int_t px = gPad->GetEventX();
2017 Int_t py = gPad->GetEventY();
2018
2019 //localize point to be deleted
2020 Int_t ipoint = -2;
2021 Int_t i;
2022 // start with a small window (in case the mouse is very close to one point)
2023 for (i = 0; i < fNpoints; i++) {
2024 Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
2025 Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
2026 if (dpx * dpx + dpy * dpy < 100) {
2027 ipoint = i;
2028 break;
2029 }
2030 }
2031 return RemovePoint(ipoint);
2032}
2033
2034////////////////////////////////////////////////////////////////////////////////
2035/// Delete point number ipoint
2036
2038{
2039 if (ipoint < 0) return -1;
2040 if (ipoint >= fNpoints) return -1;
2041
2042 Double_t **ps = ShrinkAndCopy(fNpoints - 1, ipoint);
2043 CopyAndRelease(ps, ipoint + 1, fNpoints--, ipoint);
2044 if (gPad) gPad->Modified();
2045 return ipoint;
2046}
2047
2048////////////////////////////////////////////////////////////////////////////////
2049/// Save the graph as .csv, .tsv or .txt. In case of any other extension, fall
2050/// back to TObject::SaveAs
2051///
2052/// The result can be immediately imported into Excel, gnuplot, Python or whatever,
2053/// without the needing to install pyroot, etc.
2054///
2055/// \param filename the name of the file where to store the graph
2056/// \param option some tuning options
2057///
2058/// The file extension defines the delimiter used:
2059/// - `.csv` : comma
2060/// - `.tsv` : tab
2061/// - `.txt` : space
2062///
2063/// If option = "title" a title line is generated with the axis titles.
2064
2065void TGraph::SaveAs(const char *filename, Option_t *option) const
2066{
2067 char del = '\0';
2068 TString ext = "";
2069 TString fname = filename;
2070 TString opt = option;
2071
2072 if (filename) {
2073 if (fname.EndsWith(".csv")) {del = ','; ext = "csv";}
2074 else if (fname.EndsWith(".tsv")) {del = '\t'; ext = "tsv";}
2075 else if (fname.EndsWith(".txt")) {del = ' '; ext = "txt";}
2076 }
2077 if (del) {
2078 std::ofstream out;
2079 out.open(filename, std::ios::out);
2080 if (!out.good ()) {
2081 Error("SaveAs", "cannot open file: %s", filename);
2082 return;
2083 }
2085 if(opt.Contains("title"))
2086 out << "# " << GetXaxis()->GetTitle() << "\tex\t" << GetYaxis()->GetTitle() << "\tey" << std::endl;
2087 double *ex = this->GetEX();
2088 double *ey = this->GetEY();
2089 for(int i=0 ; i<fNpoints ; i++)
2090 out << fX[i] << del << (ex?ex[i]:0) << del << fY[i] << del << (ey?ey[i]:0) << std::endl;
2092 if(opt.Contains("title"))
2093 out << "# " << GetXaxis()->GetTitle() << "\texl\t" << "\texh\t" << GetYaxis()->GetTitle() << "\teyl" << "\teyh" << std::endl;
2094 double *exl = this->GetEXlow();
2095 double *exh = this->GetEXhigh();
2096 double *eyl = this->GetEYlow();
2097 double *eyh = this->GetEYhigh();
2098 for(int i=0 ; i<fNpoints ; i++)
2099 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;
2100 } else {
2101 if(opt.Contains("title"))
2102 out << "# " << GetXaxis()->GetTitle() << "\t" << GetYaxis()->GetTitle() << std::endl;
2103 for (int i=0 ; i<fNpoints ; i++)
2104 out << fX[i] << del << fY[i] << std::endl;
2105 }
2106 out.close();
2107 Info("SaveAs", "%s file: %s has been generated", ext.Data(), filename);
2108 } else {
2110 }
2111}
2112
2113
2114////////////////////////////////////////////////////////////////////////////////
2115/// Save primitive as a C++ statement(s) on output stream out
2116
2117void TGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
2118{
2119 char quote = '"';
2120 out << " " << std::endl;
2121 static Int_t frameNumber = 0;
2122 frameNumber++;
2123
2124 if (fNpoints >= 1) {
2125 Int_t i;
2126 TString fXName = TString(GetName()) + Form("_fx%d",frameNumber);
2127 TString fYName = TString(GetName()) + Form("_fy%d",frameNumber);
2128 out << " Double_t " << fXName << "[" << fNpoints << "] = {" << std::endl;
2129 for (i = 0; i < fNpoints-1; i++) out << " " << fX[i] << "," << std::endl;
2130 out << " " << fX[fNpoints-1] << "};" << std::endl;
2131 out << " Double_t " << fYName << "[" << fNpoints << "] = {" << std::endl;
2132 for (i = 0; i < fNpoints-1; i++) out << " " << fY[i] << "," << std::endl;
2133 out << " " << fY[fNpoints-1] << "};" << std::endl;
2134 if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2135 else out << " TGraph *";
2136 out << "graph = new TGraph(" << fNpoints << "," << fXName << "," << fYName << ");" << std::endl;
2137 } else {
2138 if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2139 else out << " TGraph *";
2140 out << "graph = new TGraph();" << std::endl;
2141 }
2142
2143 out << " graph->SetName(" << quote << GetName() << quote << ");" << std::endl;
2144 out << " graph->SetTitle(" << quote << GetTitle() << quote << ");" << std::endl;
2145
2146 SaveFillAttributes(out, "graph", 0, 1001);
2147 SaveLineAttributes(out, "graph", 1, 1, 1);
2148 SaveMarkerAttributes(out, "graph", 1, 1, 1);
2149
2150 if (fHistogram) {
2151 TString hname = fHistogram->GetName();
2152 hname += frameNumber;
2153 fHistogram->SetName(Form("Graph_%s", hname.Data()));
2154 fHistogram->SavePrimitive(out, "nodraw");
2155 out << " graph->SetHistogram(" << fHistogram->GetName() << ");" << std::endl;
2156 out << " " << std::endl;
2157 }
2158
2159 // save list of functions
2160 TIter next(fFunctions);
2161 TObject *obj;
2162 while ((obj = next())) {
2163 obj->SavePrimitive(out, Form("nodraw #%d\n",++frameNumber));
2164 if (obj->InheritsFrom("TPaveStats")) {
2165 out << " graph->GetListOfFunctions()->Add(ptstats);" << std::endl;
2166 out << " ptstats->SetParent(graph->GetListOfFunctions());" << std::endl;
2167 } else {
2168 TString objname;
2169 objname.Form("%s%d",obj->GetName(),frameNumber);
2170 if (obj->InheritsFrom("TF1")) {
2171 out << " " << objname << "->SetParent(graph);\n";
2172 }
2173 out << " graph->GetListOfFunctions()->Add("
2174 << objname << ");" << std::endl;
2175 }
2176 }
2177
2178 const char *l;
2179 l = strstr(option, "multigraph");
2180 if (l) {
2181 out << " multigraph->Add(graph," << quote << l + 10 << quote << ");" << std::endl;
2182 return;
2183 }
2184 l = strstr(option, "th2poly");
2185 if (l) {
2186 out << " " << l + 7 << "->AddBin(graph);" << std::endl;
2187 return;
2188 }
2189 out << " graph->Draw(" << quote << option << quote << ");" << std::endl;
2190}
2191
2192////////////////////////////////////////////////////////////////////////////////
2193/// Multiply the values of a TGraph by a constant c1.
2194///
2195/// If option contains "x" the x values are scaled
2196/// If option contains "y" the y values are scaled
2197/// If option contains "xy" both x and y values are scaled
2198
2200{
2201 TString opt = option; opt.ToLower();
2202 if (opt.Contains("x")) {
2203 for (Int_t i=0; i<GetN(); i++)
2204 GetX()[i] *= c1;
2205 }
2206 if (opt.Contains("y")) {
2207 for (Int_t i=0; i<GetN(); i++)
2208 GetY()[i] *= c1;
2209 }
2210}
2211
2212////////////////////////////////////////////////////////////////////////////////
2213/// Set number of points in the graph
2214/// Existing coordinates are preserved
2215/// New coordinates above fNpoints are preset to 0.
2216
2218{
2219 if (n < 0) n = 0;
2220 if (n == fNpoints) return;
2221 Double_t **ps = Allocate(n);
2223 if (n > fNpoints) {
2225 }
2226 fNpoints = n;
2227}
2228
2229////////////////////////////////////////////////////////////////////////////////
2230/// Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
2231
2233{
2234 return TestBit(kNotEditable) ? kFALSE : kTRUE;
2235}
2236
2237////////////////////////////////////////////////////////////////////////////////
2238/// if editable=kFALSE, the graph cannot be modified with the mouse
2239/// by default a TGraph is editable
2240
2242{
2243 if (editable) ResetBit(kNotEditable);
2244 else SetBit(kNotEditable);
2245}
2246
2247////////////////////////////////////////////////////////////////////////////////
2248/// Set highlight (enable/disble) mode for the graph
2249/// by default highlight mode is disable
2250
2252{
2253 if (IsHighlight() == set) return;
2254
2256 if (!painter) return;
2257 SetBit(kIsHighlight, set);
2258 painter->SetHighlight(this);
2259}
2260
2261////////////////////////////////////////////////////////////////////////////////
2262/// Set the maximum of the graph.
2263
2265{
2266 fMaximum = maximum;
2267 GetHistogram()->SetMaximum(maximum);
2268}
2269
2270////////////////////////////////////////////////////////////////////////////////
2271/// Set the minimum of the graph.
2272
2274{
2275 fMinimum = minimum;
2276 GetHistogram()->SetMinimum(minimum);
2277}
2278
2279////////////////////////////////////////////////////////////////////////////////
2280/// Set x and y values for point number i.
2281
2283{
2284 if (i < 0) return;
2286
2287 if (i >= fMaxSize) {
2288 Double_t **ps = ExpandAndCopy(i + 1, fNpoints);
2289 CopyAndRelease(ps, 0, 0, 0);
2290 }
2291 if (i >= fNpoints) {
2292 // points above i can be not initialized
2293 // set zero up to i-th point to avoid redefinition
2294 // of this method in descendant classes
2295 FillZero(fNpoints, i + 1);
2296 fNpoints = i + 1;
2297 }
2298 fX[i] = x;
2299 fY[i] = y;
2300 if (gPad) gPad->Modified();
2301}
2302
2303////////////////////////////////////////////////////////////////////////////////
2304/// Set x value for point i.
2305
2307{
2308 SetPoint(i, x, GetPointY(i));
2309}
2310
2311////////////////////////////////////////////////////////////////////////////////
2312/// Set y value for point i.
2313
2315{
2316 SetPoint(i, GetPointX(i), y);
2317}
2318
2319////////////////////////////////////////////////////////////////////////////////
2320/// Set graph name.
2321void TGraph::SetName(const char *name)
2322{
2323 fName = name;
2325}
2326
2327////////////////////////////////////////////////////////////////////////////////
2328/// Change (i.e. set) the title
2329///
2330/// if title is in the form `stringt;stringx;stringy;stringz`
2331/// the graph title is set to `stringt`, the x axis title to `stringx`,
2332/// the y axis title to `stringy`, and the z axis title to `stringz`.
2333///
2334/// To insert the character `;` in one of the titles, one should use `#;`
2335/// or `#semicolon`.
2336
2337void TGraph::SetTitle(const char* title)
2338{
2339 fTitle = title;
2340 fTitle.ReplaceAll("#;",2,"#semicolon",10);
2341 Int_t p = fTitle.Index(";");
2342
2343 if (p>0) {
2344 if (!fHistogram) GetHistogram();
2345 fHistogram->SetTitle(title);
2346 Int_t n = fTitle.Length()-p;
2347 if (p>0) fTitle.Remove(p,n);
2348 fTitle.ReplaceAll("#semicolon",10,"#;",2);
2349 } else {
2350 if (fHistogram) fHistogram->SetTitle(title);
2351 }
2352}
2353
2354////////////////////////////////////////////////////////////////////////////////
2355/// Set graph name and title
2356
2357void TGraph::SetNameTitle(const char *name, const char *title)
2358{
2359 SetName(name);
2360 SetTitle(title);
2361}
2362
2363////////////////////////////////////////////////////////////////////////////////
2364/// Set statistics option on/off.
2365///
2366/// By default, the statistics box is drawn.
2367/// The paint options can be selected via gStyle->SetOptStats.
2368/// This function sets/resets the kNoStats bit in the graph object.
2369/// It has priority over the Style option.
2370
2372{
2374 if (!stats) {
2376 //remove the "stats" object from the list of functions
2377 if (fFunctions) {
2378 TObject *obj = fFunctions->FindObject("stats");
2379 if (obj) {
2380 fFunctions->Remove(obj);
2381 delete obj;
2382 }
2383 }
2384 }
2385}
2386
2387////////////////////////////////////////////////////////////////////////////////
2388/// if size*2 <= fMaxSize allocate new arrays of size points,
2389/// copy points [0,oend).
2390/// Return newarray (passed or new instance if it was zero
2391/// and allocations are needed)
2392
2394{
2395 if (size * 2 > fMaxSize || !fMaxSize) {
2396 return 0;
2397 }
2398 Double_t **newarrays = Allocate(size);
2399 CopyPoints(newarrays, 0, oend, 0);
2400 return newarrays;
2401}
2402
2403////////////////////////////////////////////////////////////////////////////////
2404/// Sorts the points of this TGraph using in-place quicksort (see e.g. older glibc).
2405/// To compare two points the function parameter greaterfunc is used (see TGraph::CompareX for an
2406/// example of such a method, which is also the default comparison function for Sort). After
2407/// the sort, greaterfunc(this, i, j) will return kTRUE for all i>j if ascending == kTRUE, and
2408/// kFALSE otherwise.
2409///
2410/// The last two parameters are used for the recursive quick sort, stating the range to be sorted
2411///
2412/// Examples:
2413/// ~~~ {.cpp}
2414/// // sort points along x axis
2415/// graph->Sort();
2416/// // sort points along their distance to origin
2417/// graph->Sort(&TGraph::CompareRadius);
2418///
2419/// Bool_t CompareErrors(const TGraph* gr, Int_t i, Int_t j) {
2420/// const TGraphErrors* ge=(const TGraphErrors*)gr;
2421/// return (ge->GetEY()[i]>ge->GetEY()[j]); }
2422/// // sort using the above comparison function, largest errors first
2423/// graph->Sort(&CompareErrors, kFALSE);
2424/// ~~~
2425
2426void TGraph::Sort(Bool_t (*greaterfunc)(const TGraph*, Int_t, Int_t) /*=TGraph::CompareX()*/,
2427 Bool_t ascending /*=kTRUE*/, Int_t low /* =0 */, Int_t high /* =-1111 */)
2428{
2429
2430 // set the bit in case of an ascending =sort in X
2431 if (greaterfunc == TGraph::CompareX && ascending && low == 0 && high == -1111)
2433
2434 if (high == -1111) high = GetN() - 1;
2435 // Termination condition
2436 if (high <= low) return;
2437
2438 int left, right;
2439 left = low; // low is the pivot element
2440 right = high;
2441 while (left < right) {
2442 // move left while item < pivot
2443 while (left <= high && greaterfunc(this, left, low) != ascending)
2444 left++;
2445 // move right while item > pivot
2446 while (right > low && greaterfunc(this, right, low) == ascending)
2447 right--;
2448 if (left < right && left < high && right > low)
2449 SwapPoints(left, right);
2450 }
2451 // right is final position for the pivot
2452 if (right > low)
2453 SwapPoints(low, right);
2454 Sort(greaterfunc, ascending, low, right - 1);
2455 Sort(greaterfunc, ascending, right + 1, high);
2456}
2457
2458////////////////////////////////////////////////////////////////////////////////
2459/// Stream an object of class TGraph.
2460
2462{
2463 if (b.IsReading()) {
2464 UInt_t R__s, R__c;
2465 Version_t R__v = b.ReadVersion(&R__s, &R__c);
2466 if (R__v > 2) {
2467 b.ReadClassBuffer(TGraph::Class(), this, R__v, R__s, R__c);
2468 if (fHistogram) fHistogram->SetDirectory(nullptr);
2469 TIter next(fFunctions);
2470 TObject *obj;
2471 while ((obj = next())) {
2472 if (obj->InheritsFrom(TF1::Class())) {
2473 TF1 *f1 = (TF1*)obj;
2474 f1->SetParent(this);
2475 }
2476 }
2478 return;
2479 }
2480 //====process old versions before automatic schema evolution
2485 b >> fNpoints;
2487 fX = new Double_t[fNpoints];
2488 fY = new Double_t[fNpoints];
2489 if (R__v < 2) {
2490 Float_t *x = new Float_t[fNpoints];
2491 Float_t *y = new Float_t[fNpoints];
2492 b.ReadFastArray(x, fNpoints);
2493 b.ReadFastArray(y, fNpoints);
2494 for (Int_t i = 0; i < fNpoints; i++) {
2495 fX[i] = x[i];
2496 fY[i] = y[i];
2497 }
2498 delete [] y;
2499 delete [] x;
2500 } else {
2501 b.ReadFastArray(fX, fNpoints);
2502 b.ReadFastArray(fY, fNpoints);
2503 }
2504 b >> fFunctions;
2505 b >> fHistogram;
2506 if (fHistogram) fHistogram->SetDirectory(nullptr);
2507 if (R__v < 2) {
2508 Float_t mi, ma;
2509 b >> mi;
2510 b >> ma;
2511 fMinimum = mi;
2512 fMaximum = ma;
2513 } else {
2514 b >> fMinimum;
2515 b >> fMaximum;
2516 }
2517 b.CheckByteCount(R__s, R__c, TGraph::IsA());
2518 //====end of old versions
2519
2520 } else {
2521 b.WriteClassBuffer(TGraph::Class(), this);
2522 }
2523}
2524
2525////////////////////////////////////////////////////////////////////////////////
2526/// Swap points.
2527
2529{
2530 SwapValues(fX, pos1, pos2);
2531 SwapValues(fY, pos1, pos2);
2532}
2533
2534////////////////////////////////////////////////////////////////////////////////
2535/// Swap values.
2536
2538{
2539 Double_t tmp = arr[pos1];
2540 arr[pos1] = arr[pos2];
2541 arr[pos2] = tmp;
2542}
2543
2544////////////////////////////////////////////////////////////////////////////////
2545/// Set current style settings in this graph
2546/// This function is called when either TCanvas::UseCurrentStyle
2547/// or TROOT::ForceStyle have been invoked.
2548
2550{
2551 if (gStyle->IsReading()) {
2560 } else {
2569 }
2571
2572 TIter next(GetListOfFunctions());
2573 TObject *obj;
2574
2575 while ((obj = next())) {
2576 obj->UseCurrentStyle();
2577 }
2578}
2579
2580////////////////////////////////////////////////////////////////////////////////
2581/// Adds all graphs from the collection to this graph.
2582/// Returns the total number of points in the result or -1 in case of an error.
2583
2585{
2586 TIter next(li);
2587 while (TObject* o = next()) {
2588 TGraph *g = dynamic_cast<TGraph*>(o);
2589 if (!g) {
2590 Error("Merge",
2591 "Cannot merge - an object which doesn't inherit from TGraph found in the list");
2592 return -1;
2593 }
2594 DoMerge(g);
2595 }
2596 return GetN();
2597}
2598
2599////////////////////////////////////////////////////////////////////////////////
2600/// protected function to perform the merge operation of a graph
2601
2603{
2604 Double_t x = 0, y = 0;
2605 for (Int_t i = 0 ; i < g->GetN(); i++) {
2606 g->GetPoint(i, x, y);
2607 SetPoint(GetN(), x, y);
2608 }
2609 return kTRUE;
2610}
2611
2612////////////////////////////////////////////////////////////////////////////////
2613/// Move all graph points on specified values dx,dy
2614/// If log argument specified, calculation done in logarithmic scale like:
2615/// new_value = exp( log(old_value) + delta );
2616
2618{
2619 Double_t x = 0, y = 0;
2620 for (Int_t i = 0 ; i < GetN(); i++) {
2621 GetPoint(i, x, y);
2622 if (!logx) {
2623 x += dx;
2624 } else if (x > 0) {
2625 x = TMath::Exp(TMath::Log(x) + dx);
2626 }
2627 if (!logy) {
2628 y += dy;
2629 } else if (y > 0) {
2630 y = TMath::Exp(TMath::Log(y) + dy);
2631 }
2632 SetPoint(i, x, y);
2633 }
2634}
2635
2636
2637////////////////////////////////////////////////////////////////////////////////
2638/// Find zero of a continuous function.
2639/// This function finds a real zero of the continuous real
2640/// function Y(X) in a given interval (A,B). See accompanying
2641/// notes for details of the argument list and calling sequence
2642
2644 , Int_t maxiterations)
2645{
2646 static Double_t a, b, ya, ytest, y1, x1, h;
2647 static Int_t j1, it, j3, j2;
2648 Double_t yb, x2;
2649 yb = 0;
2650
2651 // Calculate Y(X) at X=AZ.
2652 if (k <= 0) {
2653 a = AZ;
2654 b = BZ;
2655 X = a;
2656 j1 = 1;
2657 it = 1;
2658 k = j1;
2659 return;
2660 }
2661
2662 // Test whether Y(X) is sufficiently small.
2663
2664 if (TMath::Abs(Y) <= E2) {
2665 k = 2;
2666 return;
2667 }
2668
2669 // Calculate Y(X) at X=BZ.
2670
2671 if (j1 == 1) {
2672 ya = Y;
2673 X = b;
2674 j1 = 2;
2675 return;
2676 }
2677 // Test whether the signs of Y(AZ) and Y(BZ) are different.
2678 // if not, begin the binary subdivision.
2679
2680 if (j1 != 2) goto L100;
2681 if (ya * Y < 0) goto L120;
2682 x1 = a;
2683 y1 = ya;
2684 j1 = 3;
2685 h = b - a;
2686 j2 = 1;
2687 x2 = a + 0.5 * h;
2688 j3 = 1;
2689 it++; //*-*- Check whether (maxiterations) function values have been calculated.
2690 if (it >= maxiterations) k = j1;
2691 else X = x2;
2692 return;
2693
2694 // Test whether a bracket has been found .
2695 // If not,continue the search
2696
2697L100:
2698 if (j1 > 3) goto L170;
2699 if (ya*Y >= 0) {
2700 if (j3 >= j2) {
2701 h = 0.5 * h;
2702 j2 = 2 * j2;
2703 a = x1;
2704 ya = y1;
2705 x2 = a + 0.5 * h;
2706 j3 = 1;
2707 } else {
2708 a = X;
2709 ya = Y;
2710 x2 = X + h;
2711 j3++;
2712 }
2713 it++;
2714 if (it >= maxiterations) k = j1;
2715 else X = x2;
2716 return;
2717 }
2718
2719 // The first bracket has been found.calculate the next X by the
2720 // secant method based on the bracket.
2721
2722L120:
2723 b = X;
2724 yb = Y;
2725 j1 = 4;
2726L130:
2727 if (TMath::Abs(ya) > TMath::Abs(yb)) {
2728 x1 = a;
2729 y1 = ya;
2730 X = b;
2731 Y = yb;
2732 } else {
2733 x1 = b;
2734 y1 = yb;
2735 X = a;
2736 Y = ya;
2737 }
2738
2739 // Use the secant method based on the function values y1 and Y.
2740 // check that x2 is inside the interval (a,b).
2741
2742L150:
2743 x2 = X - Y * (X - x1) / (Y - y1);
2744 x1 = X;
2745 y1 = Y;
2746 ytest = 0.5 * TMath::Min(TMath::Abs(ya), TMath::Abs(yb));
2747 if ((x2 - a)*(x2 - b) < 0) {
2748 it++;
2749 if (it >= maxiterations) k = j1;
2750 else X = x2;
2751 return;
2752 }
2753
2754 // Calculate the next value of X by bisection . Check whether
2755 // the maximum accuracy has been achieved.
2756
2757L160:
2758 x2 = 0.5 * (a + b);
2759 ytest = 0;
2760 if ((x2 - a)*(x2 - b) >= 0) {
2761 k = 2;
2762 return;
2763 }
2764 it++;
2765 if (it >= maxiterations) k = j1;
2766 else X = x2;
2767 return;
2768
2769
2770 // Revise the bracket (a,b).
2771
2772L170:
2773 if (j1 != 4) return;
2774 if (ya * Y < 0) {
2775 b = X;
2776 yb = Y;
2777 } else {
2778 a = X;
2779 ya = Y;
2780 }
2781
2782 // Use ytest to decide the method for the next value of X.
2783
2784 if (ytest <= 0) goto L130;
2785 if (TMath::Abs(Y) - ytest <= 0) goto L150;
2786 goto L160;
2787}
#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:4846
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:2447
void Printf(const char *fmt,...)
Formats a string in a circular formatting buffer and prints the string.
Definition: TString.cxx:2461
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:479
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:1053
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 void SetParent(TObject *p=0)
Definition: TF1.h:670
virtual Int_t GetNpar() const
Definition: TF1.h:486
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set lower and upper limits for parameter ipar.
Definition: TF1.cxx:3541
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:1530
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:1787
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:1820
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:1835
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:2549
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2282
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:2584
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:2393
virtual Double_t GetRMS(Int_t axis=1) const
Return RMS of X (axis=1) or Y (axis=2)
Definition: TGraph.cxx:1314
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:1956
@ 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:1362
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:2617
virtual Double_t GetErrorYlow(Int_t bin) const
It always returns a negative value.
Definition: TGraph.cxx:1380
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:1992
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:2264
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:1371
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:2065
Int_t GetN() const
Definition: TGraph.h:126
TF1 * GetFunction(const char *name) const
Return pointer to function with name.
Definition: TGraph.cxx:1391
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:1910
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:2426
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:1573
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:2199
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:1283
static void SwapValues(Double_t *arr, Int_t pos1, Int_t pos2)
Swap values.
Definition: TGraph.cxx:2537
void Streamer(TBuffer &) override
Stream an object of class TGraph.
Definition: TGraph.cxx:2461
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:2643
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:1249
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:2602
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:1728
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:1642
virtual Int_t RemovePoint()
Delete point close to the mouse position.
Definition: TGraph.cxx:2014
virtual void InitGaus(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for a gaussian.
Definition: TGraph.cxx:1604
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:2321
virtual void SetHighlight(Bool_t set=kTRUE)
Set highlight (enable/disble) mode for the graph by default highlight mode is disable.
Definition: TGraph.cxx:2251
virtual void SwapPoints(Int_t pos1, Int_t pos2)
Swap points.
Definition: TGraph.cxx:2528
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:1552
Bool_t GetEditable() const
Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
Definition: TGraph.cxx:2232
virtual Double_t * GetEXhigh() const
Definition: TGraph.h:137
virtual Double_t GetCorrelationFactor() const
Return graph correlation factor.
Definition: TGraph.cxx:1271
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:1299
Double_t * fX
[fNpoints] array of X points
Definition: TGraph.h:47
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:1983
TAxis * GetYaxis() const
Get y axis of the graph.
Definition: TGraph.cxx:1562
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:2371
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1403
virtual Double_t GetErrorY(Int_t bin) const
It always returns a negative value. Real implementation in TGraphErrors.
Definition: TGraph.cxx:1344
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:2117
virtual Double_t GetPointY(Int_t i) const
Get y value for point i.
Definition: TGraph.cxx:1541
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:1965
void SetTitle(const char *title="") override
Change (i.e.
Definition: TGraph.cxx:2337
virtual Int_t InsertPoint()
Insert a new point at the mouse position.
Definition: TGraph.cxx:1683
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:2002
virtual void SetPointY(Int_t i, Double_t y)
Set y value for point i.
Definition: TGraph.cxx:2314
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:1974
void SetNameTitle(const char *name="", const char *title="") override
Set graph name and title.
Definition: TGraph.cxx:2357
virtual void SetPointX(Int_t i, Double_t x)
Set x value for point i.
Definition: TGraph.cxx:2306
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:2217
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:1519
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:2241
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum of the graph.
Definition: TGraph.cxx:2273
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:1353
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:1663
virtual Double_t GetErrorX(Int_t bin) const
It always returns a negative value. Real implementation in TGraphErrors.
Definition: TGraph.cxx:1336
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:8812
void SetTitle(const char *title) override
See GetStatOverflows for more information.
Definition: TH1.cxx:6708
void UseCurrentStyle() override
Copy current attributes from/to current style.
Definition: TH1.cxx:7388
@ kNoStats
Don't draw stats box.
Definition: TH1.h:163
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
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:7155
void SetName(const char *name) override
Change the name of this histogram.
Definition: TH1.cxx:8835
TObject * Clone(const char *newname="") const override
Make a complete copy of the underlying object.
Definition: TH1.cxx:2737
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:37
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:359
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
virtual void UseCurrentStyle()
Set current style settings in this object This function is called when either TCanvas::UseCurrentStyl...
Definition: TObject.cxx:717
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:107
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TObject.cxx:666
virtual void SaveAs(const char *filename="", Option_t *option="") const
Save this object in the file specified by filename.
Definition: TObject.cxx:601
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
void MakeZombie()
Definition: TObject.h:49
void ResetBit(UInt_t f)
Definition: TObject.h:186
@ kCanDelete
if object in a list can be deleted
Definition: TObject.h:58
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition: TObject.h:68
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:867
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:1150
Int_t Atoi() const
Return integer value of string.
Definition: TString.cxx:1946
Bool_t EndsWith(const char *pat, ECaseCompare cmp=kExact) const
Return true if string ends with the specified string.
Definition: TString.cxx:2202
Double_t Atof() const
Return floating-point value contained in string.
Definition: TString.cxx:2012
Bool_t IsFloat() const
Returns kTRUE if string contains a floating point or integer number.
Definition: TString.cxx:1816
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:1788
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1163
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:2336
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2314
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:971
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:1025
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:1159
Double_t Exp(Double_t x)
Definition: TMath.h:677
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:629
Double_t Log(Double_t x)
Definition: TMath.h:710
Double_t Sqrt(Double_t x)
Definition: TMath.h:641
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:176
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMathBase.h:358
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMathBase.h:274
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Definition: first.py:1
auto * m
Definition: textangle.C:8
auto * l
Definition: textangle.C:4
auto * a
Definition: textangle.C:12
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345