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