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"
27 #include "TVirtualGraphPainter.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"
42 #include "Math/MinimizerOptions.h"
43 
44 extern 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 Hist
52 A TGraph is an object made of two arrays X and Y with npoints each.
53 The TGraph painting is performed thanks to the TGraphPainter
54 class. 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 
67 The picture below gives an example:
68 
69 Begin_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 }
81 End_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 
108 TGraph::TGraph(Int_t n, const Int_t *x, const Int_t *y)
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 {
164  fNpoints = gr.fNpoints;
165  fMaxSize = gr.fMaxSize;
167  else fFunctions = new TList;
168  if (gr.fHistogram) {
170  fHistogram->SetDirectory(nullptr);
171  } else {
172  fHistogram = nullptr;
173  }
174  fMinimum = gr.fMinimum;
175  fMaximum = gr.fMaximum;
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 
200  fNpoints = gr.fNpoints;
201  fMaxSize = gr.fMaxSize;
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 
230  fMinimum = gr.fMinimum;
231  fMaximum = gr.fMaximum;
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 
257 TGraph::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 
276 TGraph::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 
334 TGraph::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 
398 TGraph::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 
615 Double_t TGraph::Chisquare(TF1 *func, Option_t * option) const
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 
700 void 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 {
757  fMaxSize = fNpoints;
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 
769 void TGraph::Draw(Option_t *option)
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) {
782  SetMarkerStyle(3);
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 
823 void 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 
836 void 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 
849 void 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 
887 Double_t TGraph::Eval(Double_t x, TSpline *spline, Option_t *option) const
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;
927  if (TestBit(TGraph::kIsSortedX) ) {
928  low = TMath::BinarySearch(fNpoints, fX, x);
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 
999 void TGraph::Expand(Int_t newsize)
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 
1009 void 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 
1046 TObject *TGraph::FindObject(const char *name) const
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 
1316 TFitResultPtr 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);
1322  ROOT::Math::MinimizerOptions minOption;
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 /// This function is called by GraphFitChisquare.
1417 /// It always returns a negative value. Real implementation in TGraphErrors
1418 
1420 {
1421  return -1;
1422 }
1423 
1424 ////////////////////////////////////////////////////////////////////////////////
1425 /// This function is called by GraphFitChisquare.
1426 /// It always returns a negative value. Real implementation in TGraphErrors
1427 
1429 {
1430  return -1;
1431 }
1432 
1433 ////////////////////////////////////////////////////////////////////////////////
1434 /// This function is called by GraphFitChisquare.
1435 /// It always returns a negative value. Real implementation in TGraphErrors
1436 /// and TGraphAsymmErrors
1437 
1439 {
1440  return -1;
1441 }
1442 
1443 ////////////////////////////////////////////////////////////////////////////////
1444 /// This function is called by GraphFitChisquare.
1445 /// It always returns a negative value. Real implementation in TGraphErrors
1446 /// and TGraphAsymmErrors
1447 
1449 {
1450  return -1;
1451 }
1452 
1453 ////////////////////////////////////////////////////////////////////////////////
1454 /// This function is called by GraphFitChisquare.
1455 /// It always returns a negative value. Real implementation in TGraphErrors
1456 /// and TGraphAsymmErrors
1457 
1459 {
1460  return -1;
1461 }
1462 
1463 ////////////////////////////////////////////////////////////////////////////////
1464 /// This function is called by GraphFitChisquare.
1465 /// It always returns a negative value. Real implementation in TGraphErrors
1466 /// and TGraphAsymmErrors
1467 
1469 {
1470  return -1;
1471 }
1472 
1473 ////////////////////////////////////////////////////////////////////////////////
1474 /// Return pointer to function with name.
1475 ///
1476 /// Functions such as TGraph::Fit store the fitted function in the list of
1477 /// functions of this graph.
1478 
1479 TF1 *TGraph::GetFunction(const char *name) const
1480 {
1481  if (!fFunctions) return nullptr;
1482  return (TF1*)fFunctions->FindObject(name);
1483 }
1484 
1485 ////////////////////////////////////////////////////////////////////////////////
1486 /// Returns a pointer to the histogram used to draw the axis
1487 /// Takes into account the two following cases.
1488 /// 1. option 'A' was specified in TGraph::Draw. Return fHistogram
1489 /// 2. user had called TPad::DrawFrame. return pointer to hframe histogram
1490 
1492 {
1493  Double_t rwxmin, rwxmax, rwymin, rwymax, maximum, minimum, dx, dy;
1494  Double_t uxmin, uxmax;
1495 
1496  ComputeRange(rwxmin, rwymin, rwxmax, rwymax); //this is redefined in TGraphErrors
1497 
1498  // (if fHistogram exist) && (if the log scale is on) &&
1499  // (if the computed range minimum is > 0) && (if the fHistogram minimum is zero)
1500  // then it means fHistogram limits have been computed in linear scale
1501  // therefore they might be too strict and cut some points. In that case the
1502  // fHistogram limits should be recomputed ie: the existing fHistogram
1503  // should not be returned.
1504  TH1F *historg = nullptr;
1505  if (fHistogram) {
1506  if (!TestBit(kResetHisto)) {
1507  if (gPad && gPad->GetLogx()) {
1508  if (rwxmin <= 0 || fHistogram->GetXaxis()->GetXmin() != 0) return fHistogram;
1509  } else if (gPad && gPad->GetLogy()) {
1510  if (rwymin <= 0 || fHistogram->GetMinimum() != 0) return fHistogram;
1511  } else {
1512  return fHistogram;
1513  }
1514  } else {
1515  const_cast <TGraph*>(this)->ResetBit(kResetHisto);
1516  }
1517  historg = fHistogram;
1518  }
1519 
1520  if (rwxmin == rwxmax) rwxmax += 1.;
1521  if (rwymin == rwymax) rwymax += 1.;
1522  dx = 0.1 * (rwxmax - rwxmin);
1523  dy = 0.1 * (rwymax - rwymin);
1524  uxmin = rwxmin - dx;
1525  uxmax = rwxmax + dx;
1526  minimum = rwymin - dy;
1527  maximum = rwymax + dy;
1528 
1529  if (fMinimum != -1111) minimum = fMinimum;
1530  if (fMaximum != -1111) maximum = fMaximum;
1531 
1532  // the graph is created with at least as many channels as there are points
1533  // to permit zooming on the full range
1534  if (uxmin < 0 && rwxmin >= 0) {
1535  if (gPad && gPad->GetLogx()) uxmin = 0.9 * rwxmin;
1536  else uxmin = 0;
1537  }
1538  if (uxmax > 0 && rwxmax <= 0) {
1539  if (gPad && gPad->GetLogx()) uxmax = 1.1 * rwxmax;
1540  else uxmax = 0;
1541  }
1542 
1543  if (minimum < 0 && rwymin >= 0) minimum = 0.9 * rwymin;
1544 
1545  if (minimum <= 0 && gPad && gPad->GetLogy()) minimum = 0.001 * maximum;
1546  if (uxmin <= 0 && gPad && gPad->GetLogx()) {
1547  if (uxmax > 1000) uxmin = 1;
1548  else uxmin = 0.001 * uxmax;
1549  }
1550 
1551  rwxmin = uxmin;
1552  rwxmax = uxmax;
1553  Int_t npt = 100;
1554  if (fNpoints > npt) npt = fNpoints;
1555  const char *gname = GetName();
1556  if (!gname[0]) gname = "Graph";
1557  // do not add the histogram to gDirectory
1558  // use local TDirectory::TContect that will set temporarly gDirectory to a nullptr and
1559  // will avoid that histogram is added in the global directory
1560  {
1561  TDirectory::TContext ctx(nullptr);
1562  ((TGraph*)this)->fHistogram = new TH1F(gname, GetTitle(), npt, rwxmin, rwxmax);
1563  }
1564  if (!fHistogram) return nullptr;
1565  fHistogram->SetMinimum(minimum);
1567  fHistogram->SetMaximum(maximum);
1568  fHistogram->GetYaxis()->SetLimits(minimum, maximum);
1569  // Restore the axis attributes if needed
1570  if (historg) {
1571  fHistogram->GetXaxis()->SetTitle(historg->GetXaxis()->GetTitle());
1584 
1585  fHistogram->GetYaxis()->SetTitle(historg->GetYaxis()->GetTitle());
1598  delete historg;
1599  }
1600  return fHistogram;
1601 }
1602 
1603 ////////////////////////////////////////////////////////////////////////////////
1604 /// Get x and y values for point number i.
1605 /// The function returns -1 in case of an invalid request or the point number otherwise
1606 
1608 {
1609  if (i < 0 || i >= fNpoints || !fX || !fY) return -1;
1610  x = fX[i];
1611  y = fY[i];
1612  return i;
1613 }
1614 
1615 ////////////////////////////////////////////////////////////////////////////////
1616 /// Get x value for point i.
1617 
1619 {
1620  if (i < 0 || i >= fNpoints || !fX)
1621  return -1.;
1622 
1623  return fX[i];
1624 }
1625 
1626 ////////////////////////////////////////////////////////////////////////////////
1627 /// Get y value for point i.
1628 
1630 {
1631  if (i < 0 || i >= fNpoints || !fY)
1632  return -1.;
1633 
1634  return fY[i];
1635 }
1636 
1637 ////////////////////////////////////////////////////////////////////////////////
1638 /// Get x axis of the graph.
1639 
1641 {
1642  TH1 *h = GetHistogram();
1643  if (!h) return 0;
1644  return h->GetXaxis();
1645 }
1646 
1647 ////////////////////////////////////////////////////////////////////////////////
1648 /// Get y axis of the graph.
1649 
1651 {
1652  TH1 *h = GetHistogram();
1653  if (!h) return 0;
1654  return h->GetYaxis();
1655 }
1656 
1657 ////////////////////////////////////////////////////////////////////////////////
1658 /// Implementation to get information on point of graph at cursor position
1659 /// Adapted from class TH1
1660 
1661 char *TGraph::GetObjectInfo(Int_t px, Int_t py) const
1662 {
1663  // localize point
1664  Int_t ipoint = -2;
1665  Int_t i;
1666  // start with a small window (in case the mouse is very close to one point)
1667  for (i = 0; i < fNpoints; i++) {
1668  Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
1669  Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
1670 
1671  if (dpx * dpx + dpy * dpy < 25) {
1672  ipoint = i;
1673  break;
1674  }
1675  }
1676 
1677  Double_t x = gPad->PadtoX(gPad->AbsPixeltoX(px));
1678  Double_t y = gPad->PadtoY(gPad->AbsPixeltoY(py));
1679 
1680  if (ipoint == -2)
1681  return Form("x=%g, y=%g", x, y);
1682 
1683  Double_t xval = fX[ipoint];
1684  Double_t yval = fY[ipoint];
1685 
1686  return Form("x=%g, y=%g, point=%d, xval=%g, yval=%g", x, y, ipoint, xval, yval);
1687 }
1688 
1689 ////////////////////////////////////////////////////////////////////////////////
1690 /// Compute Initial values of parameters for a gaussian.
1691 
1693 {
1694  Double_t allcha, sumx, sumx2, x, val, rms, mean;
1695  Int_t bin;
1696  const Double_t sqrtpi = 2.506628;
1697 
1698  // Compute mean value and RMS of the graph in the given range
1699  if (xmax <= xmin) {
1700  xmin = fX[0];
1701  xmax = fX[fNpoints-1];
1702  }
1703  Int_t np = 0;
1704  allcha = sumx = sumx2 = 0;
1705  for (bin = 0; bin < fNpoints; bin++) {
1706  x = fX[bin];
1707  if (x < xmin || x > xmax) continue;
1708  np++;
1709  val = fY[bin];
1710  sumx += val * x;
1711  sumx2 += val * x * x;
1712  allcha += val;
1713  }
1714  if (np == 0 || allcha == 0) return;
1715  mean = sumx / allcha;
1716  rms = TMath::Sqrt(sumx2 / allcha - mean * mean);
1717  Double_t binwidx = TMath::Abs((xmax - xmin) / np);
1718  if (rms == 0) rms = 1;
1720  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1721  f1->SetParameter(0, binwidx * allcha / (sqrtpi * rms));
1722  f1->SetParameter(1, mean);
1723  f1->SetParameter(2, rms);
1724  f1->SetParLimits(2, 0, 10 * rms);
1725 }
1726 
1727 ////////////////////////////////////////////////////////////////////////////////
1728 /// Compute Initial values of parameters for an exponential.
1729 
1731 {
1732  Double_t constant, slope;
1733  Int_t ifail;
1734  if (xmax <= xmin) {
1735  xmin = fX[0];
1736  xmax = fX[fNpoints-1];
1737  }
1738  Int_t nchanx = fNpoints;
1739 
1740  LeastSquareLinearFit(-nchanx, constant, slope, ifail, xmin, xmax);
1741 
1743  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1744  f1->SetParameter(0, constant);
1745  f1->SetParameter(1, slope);
1746 }
1747 
1748 ////////////////////////////////////////////////////////////////////////////////
1749 /// Compute Initial values of parameters for a polynom.
1750 
1752 {
1753  Double_t fitpar[25];
1754 
1756  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
1757  Int_t npar = f1->GetNpar();
1758  if (xmax <= xmin) {
1759  xmin = fX[0];
1760  xmax = fX[fNpoints-1];
1761  }
1762 
1763  LeastSquareFit(npar, fitpar, xmin, xmax);
1764 
1765  for (Int_t i = 0; i < npar; i++) f1->SetParameter(i, fitpar[i]);
1766 }
1767 
1768 ////////////////////////////////////////////////////////////////////////////////
1769 /// Insert a new point at the mouse position
1770 
1772 {
1773  Int_t px = gPad->GetEventX();
1774  Int_t py = gPad->GetEventY();
1775 
1776  //localize point where to insert
1777  Int_t ipoint = -2;
1778  Int_t i, d = 0;
1779  // start with a small window (in case the mouse is very close to one point)
1780  for (i = 0; i < fNpoints - 1; i++) {
1781  d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1782  if (d < 5) {
1783  ipoint = i + 1;
1784  break;
1785  }
1786  }
1787  if (ipoint == -2) {
1788  //may be we are far from one point, try again with a larger window
1789  for (i = 0; i < fNpoints - 1; i++) {
1790  d = DistancetoLine(px, py, gPad->XtoPad(fX[i]), gPad->YtoPad(fY[i]), gPad->XtoPad(fX[i+1]), gPad->YtoPad(fY[i+1]));
1791  if (d < 10) {
1792  ipoint = i + 1;
1793  break;
1794  }
1795  }
1796  }
1797  if (ipoint == -2) {
1798  //distinguish between first and last point
1799  Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[0]));
1800  Int_t dpy = py - gPad->YtoAbsPixel(gPad->XtoPad(fY[0]));
1801  if (dpx * dpx + dpy * dpy < 25) ipoint = 0;
1802  else ipoint = fNpoints;
1803  }
1804 
1805 
1806  InsertPointBefore(ipoint, gPad->AbsPixeltoX(px), gPad->AbsPixeltoY(py));
1807 
1808  gPad->Modified();
1809  return ipoint;
1810 }
1811 
1812 
1813 ////////////////////////////////////////////////////////////////////////////////
1814 /// Insert a new point with coordinates (x,y) before the point number `ipoint`.
1815 
1817 {
1818  if (ipoint < 0) {
1819  Error("TGraph", "Inserted point index should be >= 0");
1820  return;
1821  }
1822 
1823  if (ipoint > fNpoints) {
1824  Error("TGraph", "Inserted point index should be <= %d", fNpoints);
1825  return;
1826  }
1827 
1828  if (ipoint == fNpoints) {
1829  SetPoint(ipoint, x, y);
1830  return;
1831  }
1832 
1833  Double_t **ps = ExpandAndCopy(fNpoints + 1, ipoint);
1834  CopyAndRelease(ps, ipoint, fNpoints++, ipoint + 1);
1835 
1836  // To avoid redefinitions in descendant classes
1837  FillZero(ipoint, ipoint + 1);
1838 
1839  fX[ipoint] = x;
1840  fY[ipoint] = y;
1841 }
1842 
1843 
1844 ////////////////////////////////////////////////////////////////////////////////
1845 /// Integrate the TGraph data within a given (index) range.
1846 /// Note that this function computes the area of the polygon enclosed by the points of the TGraph.
1847 /// The polygon segments, which are defined by the points of the TGraph, do not need to form a closed polygon,
1848 /// since the last polygon segment, which closes the polygon, is taken as the line connecting the last TGraph point
1849 /// with the first one. It is clear that the order of the point is essential in defining the polygon.
1850 /// Also note that the segments should not intersect.
1851 ///
1852 /// NB:
1853 /// - if last=-1 (default) last is set to the last point.
1854 /// - if (first <0) the first point (0) is taken.
1855 ///
1856 /// ### Method:
1857 ///
1858 /// There are many ways to calculate the surface of a polygon. It all depends on what kind of data
1859 /// you have to deal with. The most evident solution would be to divide the polygon in triangles and
1860 /// calculate the surface of them. But this can quickly become complicated as you will have to test
1861 /// every segments of every triangles and check if they are intersecting with a current polygon's
1862 /// segment or if it goes outside the polygon. Many calculations that would lead to many problems...
1863 ///
1864 /// ### The solution (implemented by R.Brun)
1865 /// Fortunately for us, there is a simple way to solve this problem, as long as the polygon's
1866 /// segments don't intersect.
1867 /// It takes the x coordinate of the current vertex and multiply it by the y coordinate of the next
1868 /// vertex. Then it subtracts from it the result of the y coordinate of the current vertex multiplied
1869 /// by the x coordinate of the next vertex. Then divide the result by 2 to get the surface/area.
1870 ///
1871 /// ### Sources
1872 /// - http://forums.wolfram.com/mathgroup/archive/1998/Mar/msg00462.html
1873 /// - http://stackoverflow.com/questions/451426/how-do-i-calculate-the-surface-area-of-a-2d-polygon
1874 
1876 {
1877  if (first < 0) first = 0;
1878  if (last < 0) last = fNpoints - 1;
1879  if (last >= fNpoints) last = fNpoints - 1;
1880  if (first >= last) return 0;
1881  Int_t np = last - first + 1;
1882  Double_t sum = 0.0;
1883  //for(Int_t i=first;i<=last;i++) {
1884  // Int_t j = first + (i-first+1)%np;
1885  // sum += TMath::Abs(fX[i]*fY[j]);
1886  // sum -= TMath::Abs(fY[i]*fX[j]);
1887  //}
1888  for (Int_t i = first; i <= last; i++) {
1889  Int_t j = first + (i - first + 1) % np;
1890  sum += (fY[i] + fY[j]) * (fX[j] - fX[i]);
1891  }
1892  return 0.5 * TMath::Abs(sum);
1893 }
1894 
1895 ////////////////////////////////////////////////////////////////////////////////
1896 /// Return 1 if the point (x,y) is inside the polygon defined by
1897 /// the graph vertices 0 otherwise.
1898 ///
1899 /// Algorithm:
1900 ///
1901 /// The loop is executed with the end-point coordinates of a line segment
1902 /// (X1,Y1)-(X2,Y2) and the Y-coordinate of a horizontal line.
1903 /// The counter inter is incremented if the line (X1,Y1)-(X2,Y2) intersects
1904 /// the horizontal line. In this case XINT is set to the X-coordinate of the
1905 /// intersection point. If inter is an odd number, then the point x,y is within
1906 /// the polygon.
1907 
1909 {
1910  return (Int_t)TMath::IsInside(x, y, fNpoints, fX, fY);
1911 }
1912 
1913 ////////////////////////////////////////////////////////////////////////////////
1914 /// Least squares polynomial fitting without weights.
1915 ///
1916 /// \param [in] m number of parameters
1917 /// \param [in] a array of parameters
1918 /// \param [in] xmin 1st point number to fit (default =0)
1919 /// \param [in] xmax last point number to fit (default=fNpoints-1)
1920 ///
1921 /// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
1922 
1924 {
1925  const Double_t zero = 0.;
1926  const Double_t one = 1.;
1927  const Int_t idim = 20;
1928 
1929  Double_t b[400] /* was [20][20] */;
1930  Int_t i, k, l, ifail;
1931  Double_t power;
1932  Double_t da[20], xk, yk;
1933  Int_t n = fNpoints;
1934  if (xmax <= xmin) {
1935  xmin = fX[0];
1936  xmax = fX[fNpoints-1];
1937  }
1938 
1939  if (m <= 2) {
1940  LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
1941  return;
1942  }
1943  if (m > idim || m > n) return;
1944  da[0] = zero;
1945  for (l = 2; l <= m; ++l) {
1946  b[l-1] = zero;
1947  b[m + l*20 - 21] = zero;
1948  da[l-1] = zero;
1949  }
1950  Int_t np = 0;
1951  for (k = 0; k < fNpoints; ++k) {
1952  xk = fX[k];
1953  if (xk < xmin || xk > xmax) continue;
1954  np++;
1955  yk = fY[k];
1956  power = one;
1957  da[0] += yk;
1958  for (l = 2; l <= m; ++l) {
1959  power *= xk;
1960  b[l-1] += power;
1961  da[l-1] += power * yk;
1962  }
1963  for (l = 2; l <= m; ++l) {
1964  power *= xk;
1965  b[m + l*20 - 21] += power;
1966  }
1967  }
1968  b[0] = Double_t(np);
1969  for (i = 3; i <= m; ++i) {
1970  for (k = i; k <= m; ++k) {
1971  b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
1972  }
1973  }
1974  H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
1975 
1976  if (ifail < 0) {
1977  a[0] = fY[0];
1978  for (i = 1; i < m; ++i) a[i] = 0;
1979  return;
1980  }
1981  for (i = 0; i < m; ++i) a[i] = da[i];
1982 }
1983 
1984 ////////////////////////////////////////////////////////////////////////////////
1985 /// Least square linear fit without weights.
1986 ///
1987 /// Fit a straight line (a0 + a1*x) to the data in this graph.
1988 ///
1989 /// \param [in] ndata if ndata<0, fits the logarithm of the graph (used in InitExpo() to set
1990 /// the initial parameter values for a fit with exponential function.
1991 /// \param [in] a0 constant
1992 /// \param [in] a1 slope
1993 /// \param [in] ifail return parameter indicating the status of the fit (ifail=0, fit is OK)
1994 /// \param [in] xmin, xmax fitting range
1995 ///
1996 /// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
1997 
1999 {
2000  Double_t xbar, ybar, x2bar;
2001  Int_t i;
2002  Double_t xybar;
2003  Double_t fn, xk, yk;
2004  Double_t det;
2005  if (xmax <= xmin) {
2006  xmin = fX[0];
2007  xmax = fX[fNpoints-1];
2008  }
2009 
2010  ifail = -2;
2011  xbar = ybar = x2bar = xybar = 0;
2012  Int_t np = 0;
2013  for (i = 0; i < fNpoints; ++i) {
2014  xk = fX[i];
2015  if (xk < xmin || xk > xmax) continue;
2016  np++;
2017  yk = fY[i];
2018  if (ndata < 0) {
2019  if (yk <= 0) yk = 1e-9;
2020  yk = TMath::Log(yk);
2021  }
2022  xbar += xk;
2023  ybar += yk;
2024  x2bar += xk * xk;
2025  xybar += xk * yk;
2026  }
2027  fn = Double_t(np);
2028  det = fn * x2bar - xbar * xbar;
2029  ifail = -1;
2030  if (det <= 0) {
2031  if (fn > 0) a0 = ybar / fn;
2032  else a0 = 0;
2033  a1 = 0;
2034  return;
2035  }
2036  ifail = 0;
2037  a0 = (x2bar * ybar - xbar * xybar) / det;
2038  a1 = (fn * xybar - xbar * ybar) / det;
2039 }
2040 
2041 ////////////////////////////////////////////////////////////////////////////////
2042 /// Draw this graph with its current attributes.
2043 
2045 {
2047  if (painter) painter->PaintHelper(this, option);
2048 }
2049 
2050 ////////////////////////////////////////////////////////////////////////////////
2051 /// Draw the (x,y) as a graph.
2052 
2053 void TGraph::PaintGraph(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
2054 {
2056  if (painter) painter->PaintGraph(this, npoints, x, y, chopt);
2057 }
2058 
2059 ////////////////////////////////////////////////////////////////////////////////
2060 /// Draw the (x,y) as a histogram.
2061 
2062 void TGraph::PaintGrapHist(Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)
2063 {
2065  if (painter) painter->PaintGrapHist(this, npoints, x, y, chopt);
2066 }
2067 
2068 ////////////////////////////////////////////////////////////////////////////////
2069 /// Draw the stats
2070 
2072 {
2074  if (painter) painter->PaintStats(this, fit);
2075 }
2076 
2077 ////////////////////////////////////////////////////////////////////////////////
2078 /// Print graph values.
2079 
2081 {
2082  for (Int_t i = 0; i < fNpoints; i++) {
2083  printf("x[%d]=%g, y[%d]=%g\n", i, fX[i], i, fY[i]);
2084  }
2085 }
2086 
2087 ////////////////////////////////////////////////////////////////////////////////
2088 /// Recursively remove object from the list of functions
2089 
2091 {
2092  if (fFunctions) {
2094  }
2095  if (fHistogram == obj)
2096  fHistogram = nullptr;
2097 }
2098 
2099 ////////////////////////////////////////////////////////////////////////////////
2100 /// Delete point close to the mouse position
2101 
2103 {
2104  Int_t px = gPad->GetEventX();
2105  Int_t py = gPad->GetEventY();
2106 
2107  //localize point to be deleted
2108  Int_t ipoint = -2;
2109  Int_t i;
2110  // start with a small window (in case the mouse is very close to one point)
2111  for (i = 0; i < fNpoints; i++) {
2112  Int_t dpx = px - gPad->XtoAbsPixel(gPad->XtoPad(fX[i]));
2113  Int_t dpy = py - gPad->YtoAbsPixel(gPad->YtoPad(fY[i]));
2114  if (dpx * dpx + dpy * dpy < 100) {
2115  ipoint = i;
2116  break;
2117  }
2118  }
2119  return RemovePoint(ipoint);
2120 }
2121 
2122 ////////////////////////////////////////////////////////////////////////////////
2123 /// Delete point number ipoint
2124 
2126 {
2127  if (ipoint < 0) return -1;
2128  if (ipoint >= fNpoints) return -1;
2129 
2130  Double_t **ps = ShrinkAndCopy(fNpoints - 1, ipoint);
2131  CopyAndRelease(ps, ipoint + 1, fNpoints--, ipoint);
2132  if (gPad) gPad->Modified();
2133  return ipoint;
2134 }
2135 
2136 ////////////////////////////////////////////////////////////////////////////////
2137 /// Save primitive as a C++ statement(s) on output stream out
2138 
2139 void TGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
2140 {
2141  char quote = '"';
2142  out << " " << std::endl;
2143  static Int_t frameNumber = 0;
2144  frameNumber++;
2145 
2146  if (fNpoints >= 1) {
2147  Int_t i;
2148  TString fXName = TString(GetName()) + Form("_fx%d",frameNumber);
2149  TString fYName = TString(GetName()) + Form("_fy%d",frameNumber);
2150  out << " Double_t " << fXName << "[" << fNpoints << "] = {" << std::endl;
2151  for (i = 0; i < fNpoints-1; i++) out << " " << fX[i] << "," << std::endl;
2152  out << " " << fX[fNpoints-1] << "};" << std::endl;
2153  out << " Double_t " << fYName << "[" << fNpoints << "] = {" << std::endl;
2154  for (i = 0; i < fNpoints-1; i++) out << " " << fY[i] << "," << std::endl;
2155  out << " " << fY[fNpoints-1] << "};" << std::endl;
2156  if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2157  else out << " TGraph *";
2158  out << "graph = new TGraph(" << fNpoints << "," << fXName << "," << fYName << ");" << std::endl;
2159  } else {
2160  if (gROOT->ClassSaved(TGraph::Class())) out << " ";
2161  else out << " TGraph *";
2162  out << "graph = new TGraph();" << std::endl;
2163  }
2164 
2165  out << " graph->SetName(" << quote << GetName() << quote << ");" << std::endl;
2166  out << " graph->SetTitle(" << quote << GetTitle() << quote << ");" << std::endl;
2167 
2168  SaveFillAttributes(out, "graph", 0, 1001);
2169  SaveLineAttributes(out, "graph", 1, 1, 1);
2170  SaveMarkerAttributes(out, "graph", 1, 1, 1);
2171 
2172  if (fHistogram) {
2173  TString hname = fHistogram->GetName();
2174  hname += frameNumber;
2175  fHistogram->SetName(Form("Graph_%s", hname.Data()));
2176  fHistogram->SavePrimitive(out, "nodraw");
2177  out << " graph->SetHistogram(" << fHistogram->GetName() << ");" << std::endl;
2178  out << " " << std::endl;
2179  }
2180 
2181  // save list of functions
2182  TIter next(fFunctions);
2183  TObject *obj;
2184  while ((obj = next())) {
2185  obj->SavePrimitive(out, Form("nodraw #%d\n",++frameNumber));
2186  if (obj->InheritsFrom("TPaveStats")) {
2187  out << " graph->GetListOfFunctions()->Add(ptstats);" << std::endl;
2188  out << " ptstats->SetParent(graph->GetListOfFunctions());" << std::endl;
2189  } else {
2190  TString objname;
2191  objname.Form("%s%d",obj->GetName(),frameNumber);
2192  if (obj->InheritsFrom("TF1")) {
2193  out << " " << objname << "->SetParent(graph);\n";
2194  }
2195  out << " graph->GetListOfFunctions()->Add("
2196  << objname << ");" << std::endl;
2197  }
2198  }
2199 
2200  const char *l;
2201  l = strstr(option, "multigraph");
2202  if (l) {
2203  out << " multigraph->Add(graph," << quote << l + 10 << quote << ");" << std::endl;
2204  return;
2205  }
2206  l = strstr(option, "th2poly");
2207  if (l) {
2208  out << " " << l + 7 << "->AddBin(graph);" << std::endl;
2209  return;
2210  }
2211  out << " graph->Draw(" << quote << option << quote << ");" << std::endl;
2212 }
2213 
2214 ////////////////////////////////////////////////////////////////////////////////
2215 /// Set number of points in the graph
2216 /// Existing coordinates are preserved
2217 /// New coordinates above fNpoints are preset to 0.
2218 
2220 {
2221  if (n < 0) n = 0;
2222  if (n == fNpoints) return;
2223  Double_t **ps = Allocate(n);
2225  if (n > fNpoints) {
2227  }
2228  fNpoints = n;
2229 }
2230 
2231 ////////////////////////////////////////////////////////////////////////////////
2232 /// Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
2233 
2235 {
2236  return TestBit(kNotEditable) ? kFALSE : kTRUE;
2237 }
2238 
2239 ////////////////////////////////////////////////////////////////////////////////
2240 /// if editable=kFALSE, the graph cannot be modified with the mouse
2241 /// by default a TGraph is editable
2242 
2244 {
2245  if (editable) ResetBit(kNotEditable);
2246  else SetBit(kNotEditable);
2247 }
2248 
2249 ////////////////////////////////////////////////////////////////////////////////
2250 /// Set highlight (enable/disble) mode for the graph
2251 /// by default highlight mode is disable
2252 
2254 {
2255  if (IsHighlight() == set) return;
2256 
2258  if (!painter) return;
2259  SetBit(kIsHighlight, set);
2260  painter->SetHighlight(this);
2261 }
2262 
2263 ////////////////////////////////////////////////////////////////////////////////
2264 /// Set the maximum of the graph.
2265 
2267 {
2268  fMaximum = maximum;
2269  GetHistogram()->SetMaximum(maximum);
2270 }
2271 
2272 ////////////////////////////////////////////////////////////////////////////////
2273 /// Set the minimum of the graph.
2274 
2276 {
2277  fMinimum = minimum;
2278  GetHistogram()->SetMinimum(minimum);
2279 }
2280 
2281 ////////////////////////////////////////////////////////////////////////////////
2282 /// Set x and y values for point number i.
2283 
2285 {
2286  if (i < 0) return;
2288 
2289  if (i >= fMaxSize) {
2290  Double_t **ps = ExpandAndCopy(i + 1, fNpoints);
2291  CopyAndRelease(ps, 0, 0, 0);
2292  }
2293  if (i >= fNpoints) {
2294  // points above i can be not initialized
2295  // set zero up to i-th point to avoid redefinition
2296  // of this method in descendant classes
2297  FillZero(fNpoints, i + 1);
2298  fNpoints = i + 1;
2299  }
2300  fX[i] = x;
2301  fY[i] = y;
2302  if (gPad) gPad->Modified();
2303 }
2304 
2305 ////////////////////////////////////////////////////////////////////////////////
2306 /// Set x value for point i.
2307 
2309 {
2310  SetPoint(i, x, GetPointY(i));
2311 }
2312 
2313 ////////////////////////////////////////////////////////////////////////////////
2314 /// Set y value for point i.
2315 
2317 {
2318  SetPoint(i, GetPointX(i), y);
2319 }
2320 
2321 ////////////////////////////////////////////////////////////////////////////////
2322 /// Set graph name.
2323 void TGraph::SetName(const char *name)
2324 {
2325  fName = name;
2327 }
2328 
2329 ////////////////////////////////////////////////////////////////////////////////
2330 /// Change (i.e. set) the title
2331 ///
2332 /// if title is in the form `stringt;stringx;stringy;stringz`
2333 /// the graph title is set to `stringt`, the x axis title to `stringx`,
2334 /// the y axis title to `stringy`, and the z axis title to `stringz`.
2335 ///
2336 /// To insert the character `;` in one of the titles, one should use `#;`
2337 /// or `#semicolon`.
2338 
2339 void TGraph::SetTitle(const char* title)
2340 {
2341  fTitle = title;
2342  fTitle.ReplaceAll("#;",2,"#semicolon",10);
2343  Int_t p = fTitle.Index(";");
2344 
2345  if (p>0) {
2346  if (!fHistogram) GetHistogram();
2347  fHistogram->SetTitle(title);
2348  Int_t n = fTitle.Length()-p;
2349  if (p>0) fTitle.Remove(p,n);
2350  fTitle.ReplaceAll("#semicolon",10,"#;",2);
2351  } else {
2352  if (fHistogram) fHistogram->SetTitle(title);
2353  }
2354 }
2355 
2356 ////////////////////////////////////////////////////////////////////////////////
2357 /// Set graph name and title
2358 
2359 void TGraph::SetNameTitle(const char *name, const char *title)
2360 {
2361  SetName(name);
2362  SetTitle(title);
2363 }
2364 
2365 ////////////////////////////////////////////////////////////////////////////////
2366 /// if size*2 <= fMaxSize allocate new arrays of size points,
2367 /// copy points [0,oend).
2368 /// Return newarray (passed or new instance if it was zero
2369 /// and allocations are needed)
2370 
2372 {
2373  if (size * 2 > fMaxSize || !fMaxSize) {
2374  return 0;
2375  }
2376  Double_t **newarrays = Allocate(size);
2377  CopyPoints(newarrays, 0, oend, 0);
2378  return newarrays;
2379 }
2380 
2381 ////////////////////////////////////////////////////////////////////////////////
2382 /// Sorts the points of this TGraph using in-place quicksort (see e.g. older glibc).
2383 /// To compare two points the function parameter greaterfunc is used (see TGraph::CompareX for an
2384 /// example of such a method, which is also the default comparison function for Sort). After
2385 /// the sort, greaterfunc(this, i, j) will return kTRUE for all i>j if ascending == kTRUE, and
2386 /// kFALSE otherwise.
2387 ///
2388 /// The last two parameters are used for the recursive quick sort, stating the range to be sorted
2389 ///
2390 /// Examples:
2391 /// ~~~ {.cpp}
2392 /// // sort points along x axis
2393 /// graph->Sort();
2394 /// // sort points along their distance to origin
2395 /// graph->Sort(&TGraph::CompareRadius);
2396 ///
2397 /// Bool_t CompareErrors(const TGraph* gr, Int_t i, Int_t j) {
2398 /// const TGraphErrors* ge=(const TGraphErrors*)gr;
2399 /// return (ge->GetEY()[i]>ge->GetEY()[j]); }
2400 /// // sort using the above comparison function, largest errors first
2401 /// graph->Sort(&CompareErrors, kFALSE);
2402 /// ~~~
2403 
2404 void TGraph::Sort(Bool_t (*greaterfunc)(const TGraph*, Int_t, Int_t) /*=TGraph::CompareX()*/,
2405  Bool_t ascending /*=kTRUE*/, Int_t low /* =0 */, Int_t high /* =-1111 */)
2406 {
2407 
2408  // set the bit in case of an ascending =sort in X
2409  if (greaterfunc == TGraph::CompareX && ascending && low == 0 && high == -1111)
2411 
2412  if (high == -1111) high = GetN() - 1;
2413  // Termination condition
2414  if (high <= low) return;
2415 
2416  int left, right;
2417  left = low; // low is the pivot element
2418  right = high;
2419  while (left < right) {
2420  // move left while item < pivot
2421  while (left <= high && greaterfunc(this, left, low) != ascending)
2422  left++;
2423  // move right while item > pivot
2424  while (right > low && greaterfunc(this, right, low) == ascending)
2425  right--;
2426  if (left < right && left < high && right > low)
2427  SwapPoints(left, right);
2428  }
2429  // right is final position for the pivot
2430  if (right > low)
2431  SwapPoints(low, right);
2432  Sort(greaterfunc, ascending, low, right - 1);
2433  Sort(greaterfunc, ascending, right + 1, high);
2434 }
2435 
2436 ////////////////////////////////////////////////////////////////////////////////
2437 /// Stream an object of class TGraph.
2438 
2439 void TGraph::Streamer(TBuffer &b)
2440 {
2441  if (b.IsReading()) {
2442  UInt_t R__s, R__c;
2443  Version_t R__v = b.ReadVersion(&R__s, &R__c);
2444  if (R__v > 2) {
2445  b.ReadClassBuffer(TGraph::Class(), this, R__v, R__s, R__c);
2446  if (fHistogram) fHistogram->SetDirectory(nullptr);
2447  TIter next(fFunctions);
2448  TObject *obj;
2449  while ((obj = next())) {
2450  if (obj->InheritsFrom(TF1::Class())) {
2451  TF1 *f1 = (TF1*)obj;
2452  f1->SetParent(this);
2453  }
2454  }
2455  fMaxSize = fNpoints;
2456  return;
2457  }
2458  //====process old versions before automatic schema evolution
2459  TNamed::Streamer(b);
2460  TAttLine::Streamer(b);
2461  TAttFill::Streamer(b);
2462  TAttMarker::Streamer(b);
2463  b >> fNpoints;
2464  fMaxSize = fNpoints;
2465  fX = new Double_t[fNpoints];
2466  fY = new Double_t[fNpoints];
2467  if (R__v < 2) {
2468  Float_t *x = new Float_t[fNpoints];
2469  Float_t *y = new Float_t[fNpoints];
2470  b.ReadFastArray(x, fNpoints);
2471  b.ReadFastArray(y, fNpoints);
2472  for (Int_t i = 0; i < fNpoints; i++) {
2473  fX[i] = x[i];
2474  fY[i] = y[i];
2475  }
2476  delete [] y;
2477  delete [] x;
2478  } else {
2479  b.ReadFastArray(fX, fNpoints);
2480  b.ReadFastArray(fY, fNpoints);
2481  }
2482  b >> fFunctions;
2483  b >> fHistogram;
2484  if (fHistogram) fHistogram->SetDirectory(nullptr);
2485  if (R__v < 2) {
2486  Float_t mi, ma;
2487  b >> mi;
2488  b >> ma;
2489  fMinimum = mi;
2490  fMaximum = ma;
2491  } else {
2492  b >> fMinimum;
2493  b >> fMaximum;
2494  }
2495  b.CheckByteCount(R__s, R__c, TGraph::IsA());
2496  //====end of old versions
2497 
2498  } else {
2499  b.WriteClassBuffer(TGraph::Class(), this);
2500  }
2501 }
2502 
2503 ////////////////////////////////////////////////////////////////////////////////
2504 /// Swap points.
2505 
2507 {
2508  SwapValues(fX, pos1, pos2);
2509  SwapValues(fY, pos1, pos2);
2510 }
2511 
2512 ////////////////////////////////////////////////////////////////////////////////
2513 /// Swap values.
2514 
2515 void TGraph::SwapValues(Double_t* arr, Int_t pos1, Int_t pos2)
2516 {
2517  Double_t tmp = arr[pos1];
2518  arr[pos1] = arr[pos2];
2519  arr[pos2] = tmp;
2520 }
2521 
2522 ////////////////////////////////////////////////////////////////////////////////
2523 /// Set current style settings in this graph
2524 /// This function is called when either TCanvas::UseCurrentStyle
2525 /// or TROOT::ForceStyle have been invoked.
2526 
2528 {
2529  if (gStyle->IsReading()) {
2538  } else {
2547  }
2549 
2550  TIter next(GetListOfFunctions());
2551  TObject *obj;
2552 
2553  while ((obj = next())) {
2554  obj->UseCurrentStyle();
2555  }
2556 }
2557 
2558 ////////////////////////////////////////////////////////////////////////////////
2559 /// Adds all graphs from the collection to this graph.
2560 /// Returns the total number of points in the result or -1 in case of an error.
2561 
2563 {
2564  TIter next(li);
2565  while (TObject* o = next()) {
2566  TGraph *g = dynamic_cast<TGraph*>(o);
2567  if (!g) {
2568  Error("Merge",
2569  "Cannot merge - an object which doesn't inherit from TGraph found in the list");
2570  return -1;
2571  }
2572  DoMerge(g);
2573  }
2574  return GetN();
2575 }
2576 
2577 ////////////////////////////////////////////////////////////////////////////////
2578 /// protected function to perform the merge operation of a graph
2579 
2581 {
2582  Double_t x = 0, y = 0;
2583  for (Int_t i = 0 ; i < g->GetN(); i++) {
2584  g->GetPoint(i, x, y);
2585  SetPoint(GetN(), x, y);
2586  }
2587  return kTRUE;
2588 }
2589 
2590 ////////////////////////////////////////////////////////////////////////////////
2591 /// Move all graph points on specified values dx,dy
2592 /// If log argument specified, calculation done in logarithmic scale like:
2593 /// new_value = exp( log(old_value) + delta );
2594 
2596 {
2597  Double_t x = 0, y = 0;
2598  for (Int_t i = 0 ; i < GetN(); i++) {
2599  GetPoint(i, x, y);
2600  if (!logx) {
2601  x += dx;
2602  } else if (x > 0) {
2603  x = TMath::Exp(TMath::Log(x) + dx);
2604  }
2605  if (!logy) {
2606  y += dy;
2607  } else if (y > 0) {
2608  y = TMath::Exp(TMath::Log(y) + dy);
2609  }
2610  SetPoint(i, x, y);
2611  }
2612 }
2613 
2614 
2615 ////////////////////////////////////////////////////////////////////////////////
2616 /// Find zero of a continuous function.
2617 /// This function finds a real zero of the continuous real
2618 /// function Y(X) in a given interval (A,B). See accompanying
2619 /// notes for details of the argument list and calling sequence
2620 
2622  , Int_t maxiterations)
2623 {
2624  static Double_t a, b, ya, ytest, y1, x1, h;
2625  static Int_t j1, it, j3, j2;
2626  Double_t yb, x2;
2627  yb = 0;
2628 
2629  // Calculate Y(X) at X=AZ.
2630  if (k <= 0) {
2631  a = AZ;
2632  b = BZ;
2633  X = a;
2634  j1 = 1;
2635  it = 1;
2636  k = j1;
2637  return;
2638  }
2639 
2640  // Test whether Y(X) is sufficiently small.
2641 
2642  if (TMath::Abs(Y) <= E2) {
2643  k = 2;
2644  return;
2645  }
2646 
2647  // Calculate Y(X) at X=BZ.
2648 
2649  if (j1 == 1) {
2650  ya = Y;
2651  X = b;
2652  j1 = 2;
2653  return;
2654  }
2655  // Test whether the signs of Y(AZ) and Y(BZ) are different.
2656  // if not, begin the binary subdivision.
2657 
2658  if (j1 != 2) goto L100;
2659  if (ya * Y < 0) goto L120;
2660  x1 = a;
2661  y1 = ya;
2662  j1 = 3;
2663  h = b - a;
2664  j2 = 1;
2665  x2 = a + 0.5 * h;
2666  j3 = 1;
2667  it++; //*-*- Check whether (maxiterations) function values have been calculated.
2668  if (it >= maxiterations) k = j1;
2669  else X = x2;
2670  return;
2671 
2672  // Test whether a bracket has been found .
2673  // If not,continue the search
2674 
2675 L100:
2676  if (j1 > 3) goto L170;
2677  if (ya*Y >= 0) {
2678  if (j3 >= j2) {
2679  h = 0.5 * h;
2680  j2 = 2 * j2;
2681  a = x1;
2682  ya = y1;
2683  x2 = a + 0.5 * h;
2684  j3 = 1;
2685  } else {
2686  a = X;
2687  ya = Y;
2688  x2 = X + h;
2689  j3++;
2690  }
2691  it++;
2692  if (it >= maxiterations) k = j1;
2693  else X = x2;
2694  return;
2695  }
2696 
2697  // The first bracket has been found.calculate the next X by the
2698  // secant method based on the bracket.
2699 
2700 L120:
2701  b = X;
2702  yb = Y;
2703  j1 = 4;
2704 L130:
2705  if (TMath::Abs(ya) > TMath::Abs(yb)) {
2706  x1 = a;
2707  y1 = ya;
2708  X = b;
2709  Y = yb;
2710  } else {
2711  x1 = b;
2712  y1 = yb;
2713  X = a;
2714  Y = ya;
2715  }
2716 
2717  // Use the secant method based on the function values y1 and Y.
2718  // check that x2 is inside the interval (a,b).
2719 
2720 L150:
2721  x2 = X - Y * (X - x1) / (Y - y1);
2722  x1 = X;
2723  y1 = Y;
2724  ytest = 0.5 * TMath::Min(TMath::Abs(ya), TMath::Abs(yb));
2725  if ((x2 - a)*(x2 - b) < 0) {
2726  it++;
2727  if (it >= maxiterations) k = j1;
2728  else X = x2;
2729  return;
2730  }
2731 
2732  // Calculate the next value of X by bisection . Check whether
2733  // the maximum accuracy has been achieved.
2734 
2735 L160:
2736  x2 = 0.5 * (a + b);
2737  ytest = 0;
2738  if ((x2 - a)*(x2 - b) >= 0) {
2739  k = 2;
2740  return;
2741  }
2742  it++;
2743  if (it >= maxiterations) k = j1;
2744  else X = x2;
2745  return;
2746 
2747 
2748  // Revise the bracket (a,b).
2749 
2750 L170:
2751  if (j1 != 4) return;
2752  if (ya * Y < 0) {
2753  b = X;
2754  yb = Y;
2755  } else {
2756  a = X;
2757  ya = Y;
2758  }
2759 
2760  // Use ytest to decide the method for the next value of X.
2761 
2762  if (ytest <= 0) goto L130;
2763  if (TMath::Abs(Y) - ytest <= 0) goto L150;
2764  goto L160;
2765 }
TGraph::DrawPanel
virtual void DrawPanel()
Display a panel with all graph drawing options.
Definition: TGraph.cxx:866
l
auto * l
Definition: textangle.C:4
TGraph::LeastSquareFit
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:1923
TAttMarker::SetMarkerSize
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
TGraph::Eval
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
TGraph::GetListOfFunctions
TList * GetListOfFunctions() const
Definition: TGraph.h:118
m
auto * m
Definition: textangle.C:8
n
const Int_t n
Definition: legend1.C:16
TAxis
Class to manage histogram axis.
Definition: TAxis.h:30
first
Definition: first.py:1
TGraph::FindObject
virtual TObject * FindObject(const char *name) const
Search object named name in the list of functions.
Definition: TGraph.cxx:1046
TObject::SavePrimitive
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
TBrowser
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:37
ymax
float ymax
Definition: THbookFile.cxx:95
TGraph::CompareArg
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
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
TGraph::GetCorrelationFactor
virtual Double_t GetCorrelationFactor() const
Return graph correlation factor.
Definition: TGraph.cxx:1353
e
#define e(i)
Definition: RSha256.hxx:103
TGraph::GetErrorYlow
virtual Double_t GetErrorYlow(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1468
TGraph::IsInside
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:1908
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
TGraph::fFunctions
TList * fFunctions
Pointer to list of functions (fits and user)
Definition: TGraph.h:49
TVectorD.h
TGraph::SetTitle
virtual void SetTitle(const char *title="")
Change (i.e.
Definition: TGraph.cxx:2339
TMath::ATan2
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:679
Version_t
short Version_t
Definition: RtypesCore.h:65
TGraph::Draw
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:769
TGraph::kClipFrame
@ kClipFrame
Clip to the frame boundary.
Definition: TGraph.h:70
TH1::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
Definition: TH1.cxx:7122
f
#define f(i)
Definition: RSha256.hxx:104
TAttFill::SetFillStyle
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
TGraph::fHistogram
TH1F * fHistogram
Pointer to histogram used for drawing axis.
Definition: TGraph.h:50
Option_t
const char Option_t
Definition: RtypesCore.h:66
TGraph::DrawGraph
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
Foption_t
Definition: Foption.h:24
TGraph::Merge
virtual Int_t Merge(TCollection *list)
Adds all graphs from the collection to this graph.
Definition: TGraph.cxx:2562
TString::Atoi
Int_t Atoi() const
Return integer value of string.
Definition: TString.cxx:1943
TGraph::fMaxSize
Int_t fMaxSize
!Current dimension of arrays fX and fY
Definition: TGraph.h:45
kNPOS
const Ssiz_t kNPOS
Definition: RtypesCore.h:124
TGraph::CopyPoints
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
TH1::SetMinimum
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:399
TGraph::Apply
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
TSpline::Eval
virtual Double_t Eval(Double_t x) const =0
gEnv
R__EXTERN TEnv * gEnv
Definition: TEnv.h:171
TList::FindObject
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:578
TAxis::SetLimits
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:154
TGraph::Paint
virtual void Paint(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:2044
TString::Data
const char * Data() const
Definition: TString.h:369
TGraph::SetEditable
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:2243
TNamed::operator=
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
TStyle::SetHistLineStyle
void SetHistLineStyle(Style_t styl=0)
Definition: TStyle.h:365
Axis_t
double Axis_t
Definition: RtypesCore.h:85
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
TNamed::GetTitle
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
TAttMarker::GetMarkerSize
virtual Size_t GetMarkerSize() const
Return the marker size.
Definition: TAttMarker.h:33
ROOT::Fit::FitOptionsMake
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition: HFitImpl.cxx:684
TStyle::GetHistFillColor
Color_t GetHistFillColor() const
Definition: TStyle.h:223
TString::Replace
TString & Replace(Ssiz_t pos, Ssiz_t n, const char *s)
Definition: TString.h:682
TF1::GetNpar
virtual Int_t GetNpar() const
Definition: TF1.h:481
TGraph::ShrinkAndCopy
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:2371
TVirtualGraphPainter::PaintHelper
virtual void PaintHelper(TGraph *theGraph, Option_t *option)=0
TGraph.h
xmax
float xmax
Definition: THbookFile.cxx:95
TAttAxis::SetLabelOffset
virtual void SetLabelOffset(Float_t offset=0.005)
Set distance between the axis and the labels.
Definition: TAttAxis.cxx:192
sum
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345
HFitInterface.h
TMath::Log
Double_t Log(Double_t x)
Definition: TMath.h:760
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
TAttMarker::Copy
void Copy(TAttMarker &attmarker) const
Copy this marker attributes to a new TAttMarker.
Definition: TAttMarker.cxx:239
TGraph::GetPointY
virtual Double_t GetPointY(Int_t i) const
Get y value for point i.
Definition: TGraph.cxx:1629
TString::Atof
Double_t Atof() const
Return floating-point value contained in string.
Definition: TString.cxx:2009
TGraph::operator=
TGraph & operator=(const TGraph &)
Equal operator for this graph.
Definition: TGraph.cxx:192
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
TAttAxis::GetLabelSize
virtual Float_t GetLabelSize() const
Definition: TAttAxis.h:41
DataRange.h
TGraph::PaintGrapHist
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:2062
operator=
Binding & operator=(OUT(*fun)(void))
Definition: TRInterface_Binding.h:15
TRandom.h
Float_t
float Float_t
Definition: RtypesCore.h:57
TCollection::Clone
virtual TObject * Clone(const char *newname="") const
Make a clone of an collection using the Streamer facility.
Definition: TCollection.cxx:263
TMath::Exp
Double_t Exp(Double_t x)
Definition: TMath.h:727
TStyle.h
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:162
TGraph::RemovePoint
virtual Int_t RemovePoint()
Delete point close to the mouse position.
Definition: TGraph.cxx:2102
TGraph::SetNameTitle
virtual void SetNameTitle(const char *name="", const char *title="")
Set graph name and title.
Definition: TGraph.cxx:2359
Int_t
int Int_t
Definition: RtypesCore.h:45
TAxis::GetBinCenter
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition: TAxis.cxx:478
TAxis::SetTimeFormat
virtual void SetTimeFormat(const char *format="")
Change the format used for time plotting.
Definition: TAxis.cxx:1020
TBrowser.h
TGeant4Unit::ps
static constexpr double ps
Definition: TGeant4SystemOfUnits.h:165
ROOT::Fit::FitObject
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
TNamed::fName
TString fName
Definition: TNamed.h:32
TObject::AppendPad
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:107
TGraph::fX
Double_t * fX
[fNpoints] array of X points
Definition: TGraph.h:47
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
TAttLine::SaveLineAttributes
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
x
Double_t x[n]
Definition: legend1.C:17
TString::Length
Ssiz_t Length() const
Definition: TString.h:410
H1LeastSquareSeqnd
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:4804
TDirectory::TContext
TDirectory::TContext keeps track and restore the current directory.
Definition: TDirectory.h:89
TSpline
Base class for spline implementation containing the Draw/Paint methods.
Definition: TSpline.h:31
TAttMarker::GetMarkerStyle
virtual Style_t GetMarkerStyle() const
Return the marker style.
Definition: TAttMarker.h:32
TGraph::GetErrorY
virtual Double_t GetErrorY(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1428
TGraph::CompareRadius
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
TEnv::GetValue
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition: TEnv.cxx:491
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
TAttLine::SetLineColor
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
TAxis::GetTimeDisplay
virtual Bool_t GetTimeDisplay() const
Definition: TAxis.h:126
TAxis::SetTimeDisplay
virtual void SetTimeDisplay(Int_t value)
Definition: TAxis.h:161
TMath::Sort
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMathBase.h:362
TString::IsFloat
Bool_t IsFloat() const
Returns kTRUE if string contains a floating point or integer number.
Definition: TString.cxx:1813
TAxis::CenterTitle
void CenterTitle(Bool_t center=kTRUE)
Center axis title.
Definition: TAxis.h:184
TString::Format
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:2333
TGraph::GetErrorXlow
virtual Double_t GetErrorXlow(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1448
TH1::SetName
virtual void SetName(const char *name)
Change the name of this histogram.
Definition: TH1.cxx:8793
TEnv.h
TGraph::ExecuteEvent
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TGraph.cxx:990
TString
Basic string class.
Definition: TString.h:136
TAttAxis::GetLabelOffset
virtual Float_t GetLabelOffset() const
Definition: TAttAxis.h:40
TVirtualFitter.h
TMath::IsInside
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
TGraph::GetYaxis
TAxis * GetYaxis() const
Get y axis of the graph.
Definition: TGraph.cxx:1650
TGraph::InitPolynom
virtual void InitPolynom(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for a polynom.
Definition: TGraph.cxx:1751
TAttFill::SaveFillAttributes
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
MinimizerOptions.h
Foption.h
TVirtualFitter
Abstract Base Class for Fitting.
Definition: TVirtualFitter.h:29
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
TVirtualGraphPainter.h
b
#define b(i)
Definition: RSha256.hxx:100
TVirtualGraphPainter::DrawPanelHelper
virtual void DrawPanelHelper(TGraph *theGraph)=0
bool
TString::ReplaceAll
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
TGraph::SwapPoints
virtual void SwapPoints(Int_t pos1, Int_t pos2)
Swap points.
Definition: TGraph.cxx:2506
TGraph::GetCovariance
virtual Double_t GetCovariance() const
Return covariance of vectors x,y.
Definition: TGraph.cxx:1365
TGraph::SetPointX
virtual void SetPointX(Int_t i, Double_t x)
Set x value for point i.
Definition: TGraph.cxx:2308
TStyle::SetHistFillStyle
void SetHistFillStyle(Style_t styl=0)
Definition: TStyle.h:364
TAttAxis::SetTitleSize
virtual void SetTitleSize(Float_t size=0.04)
Set size of axis title.
Definition: TAttAxis.cxx:303
TAxis::GetTitle
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:129
TNamed::fTitle
TString fTitle
Definition: TNamed.h:33
TVectorT::GetNrows
Int_t GetNrows() const
Definition: TVectorT.h:75
TH1::UseCurrentStyle
void UseCurrentStyle()
Copy current attributes from/to current style.
Definition: TH1.cxx:7355
x1
static const double x1[5]
Definition: RooGaussKronrodIntegrator1D.cxx:346
TGraph::InsertPoint
virtual Int_t InsertPoint()
Insert a new point at the mouse position.
Definition: TGraph.cxx:1771
TROOT.h
TF1::SetParameter
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:634
TH1::SetTitle
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition: TH1.cxx:6678
TGraph::RecursiveRemove
virtual void RecursiveRemove(TObject *obj)
Recursively remove object from the list of functions.
Definition: TGraph.cxx:2090
TH1::Clone
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
Definition: TH1.cxx:2740
TString::ToUpper
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1160
TGraph::Integral
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:1875
TList::First
virtual TObject * First() const
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:659
TObject::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:359
TGraph::DoMerge
virtual Bool_t DoMerge(const TGraph *g)
protected function to perform the merge operation of a graph
Definition: TGraph.cxx:2580
TPluginHandler
Definition: TPluginManager.h:101
TSpline.h
TString::Form
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2311
TAttLine::Copy
void Copy(TAttLine &attline) const
Copy this line attributes to a new TAttLine.
Definition: TAttLine.cxx:175
TAttMarker::SaveMarkerAttributes
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:343
TAttAxis::SetLabelSize
virtual void SetLabelSize(Float_t size=0.04)
Set size of axis labels.
Definition: TAttAxis.cxx:203
TGraph::GetXaxis
TAxis * GetXaxis() const
Get x axis of the graph.
Definition: TGraph.cxx:1640
TObject::UseCurrentStyle
virtual void UseCurrentStyle()
Set current style settings in this object This function is called when either TCanvas::UseCurrentStyl...
Definition: TObject.cxx:717
PyTorch_Generate_CNN_Model.fit
def fit(model, train_loader, val_loader, num_epochs, batch_size, optimizer, criterion, save_best, scheduler)
Definition: PyTorch_Generate_CNN_Model.py:32
TAttFill::Copy
void Copy(TAttFill &attfill) const
Copy this fill attributes to a new TAttFill.
Definition: TAttFill.cxx:204
TBuffer.h
TAttLine
Line Attributes class.
Definition: TAttLine.h:18
gStyle
R__EXTERN TStyle * gStyle
Definition: TStyle.h:412
TObject::ResetBit
void ResetBit(UInt_t f)
Definition: TObject.h:186
TSystem.h
TFitResultPtr
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:32
TH1::GetYaxis
TAxis * GetYaxis()
Definition: TH1.h:321
xmin
float xmin
Definition: THbookFile.cxx:95
TGraph::IsHighlight
virtual Bool_t IsHighlight() const
Definition: TGraph.h:159
TGraph::GetHistogram
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1491
TGraph::Chisquare
virtual Double_t Chisquare(TF1 *f1, Option_t *option="") const
Return the chisquare of this graph with respect to f1.
Definition: TGraph.cxx:615
h
#define h(i)
Definition: RSha256.hxx:106
TStyle::SetHistFillColor
void SetHistFillColor(Color_t color=1)
Definition: TStyle.h:362
TGraph::SetName
virtual void SetName(const char *name="")
Set graph name.
Definition: TGraph.cxx:2323
TAttLine::GetLineColor
virtual Color_t GetLineColor() const
Return the line color.
Definition: TAttLine.h:33
size
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
TObject::SetBit
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
TGraph::Sort
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:2404
gr
TGraphErrors * gr
Definition: legend1.C:25
TGraph::DistancetoPrimitive
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a graph.
Definition: TGraph.cxx:813
TAttAxis::SetTitleFont
virtual void SetTitleFont(Style_t font=62)
Set the title font.
Definition: TAttAxis.cxx:321
TGraph::InitGaus
virtual void InitGaus(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for a gaussian.
Definition: TGraph.cxx:1692
TAttMarker::SetMarkerColor
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
a
auto * a
Definition: textangle.C:12
TString::Remove
TString & Remove(Ssiz_t pos)
Definition: TString.h:673
TNamed
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
TAttLine::GetLineStyle
virtual Style_t GetLineStyle() const
Return the line style.
Definition: TAttLine.h:34
TGraph::Zero
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:2621
TGraph::GetErrorX
virtual Double_t GetErrorX(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1419
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
TGraph::fMaximum
Double_t fMaximum
Maximum value for plotting along y.
Definition: TGraph.h:52
TStyle::GetHistFillStyle
Style_t GetHistFillStyle() const
Definition: TStyle.h:225
TH1::SetMaximum
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:398
TVirtualGraphPainter::SetHighlight
virtual void SetHighlight(TGraph *theGraph)=0
TAttMarker
Marker Attributes class.
Definition: TAttMarker.h:19
ROOT::Fit::DataRange
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:35
TGraph::Print
virtual void Print(Option_t *chopt="") const
Print graph values.
Definition: TGraph.cxx:2080
ROOT::Math::MinimizerOptions
Minimizer options.
Definition: MinimizerOptions.h:32
TGraph::SetMaximum
virtual void SetMaximum(Double_t maximum=-1111)
Set the maximum of the graph.
Definition: TGraph.cxx:2266
TGraph::fMinimum
Double_t fMinimum
Minimum value for plotting along y.
Definition: TGraph.h:51
TVirtualPad.h
UInt_t
unsigned int UInt_t
Definition: RtypesCore.h:46
TAxis::GetTimeFormat
virtual const char * GetTimeFormat() const
Definition: TAxis.h:127
TGraph::FillZero
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
TH1::SetDirectory
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:8770
TAttFill::GetFillStyle
virtual Style_t GetFillStyle() const
Return the fill area style.
Definition: TAttFill.h:31
TGraph::fY
Double_t * fY
[fNpoints] array of Y points
Definition: TGraph.h:48
y
Double_t y[n]
Definition: legend1.C:17
TAxis::GetRotateTitle
Bool_t GetRotateTitle() const
Definition: TAxis.h:124
TGraph::kIsSortedX
@ kIsSortedX
Graph is sorted in X points.
Definition: TGraph.h:73
TObject::MakeZombie
void MakeZombie()
Definition: TObject.h:49
TStyle::SetHistLineWidth
void SetHistLineWidth(Width_t width=1)
Definition: TStyle.h:366
TSystem::ExpandPathName
virtual Bool_t ExpandPathName(TString &path)
Expand a pathname getting rid of special shell characters like ~.
Definition: TSystem.cxx:1273
TGraph::Expand
virtual void Expand(Int_t newsize)
If array sizes <= newsize, expand storage to 2*newsize.
Definition: TGraph.cxx:999
TGraph::CtorAllocate
Bool_t CtorAllocate()
In constructors set fNpoints than call this method.
Definition: TGraph.cxx:743
TGraph::UseCurrentStyle
virtual void UseCurrentStyle()
Set current style settings in this graph This function is called when either TCanvas::UseCurrentStyle...
Definition: TGraph.cxx:2527
TStyle::IsReading
Bool_t IsReading() const
Definition: TStyle.h:282
TNamed::SetTitle
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
line
TLine * line
Definition: entrylistblock_figure1.C:235
TCollection::IsEmpty
virtual Bool_t IsEmpty() const
Definition: TCollection.h:186
TAttMarker::GetMarkerColor
virtual Color_t GetMarkerColor() const
Return the marker color.
Definition: TAttMarker.h:31
TString::IsDigit
Bool_t IsDigit() const
Returns true if all characters in string are digits (0-9) or white spaces, i.e.
Definition: TString.cxx:1785
TGraph::kNotEditable
@ kNotEditable
Bit set if graph is non editable.
Definition: TGraph.h:72
TObject::kCanDelete
@ kCanDelete
if object in a list can be deleted
Definition: TObject.h:58
ROOT::Fit::Chisquare
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
TGraph::SetPoint
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2284
TMath::Min
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
TObject::kInvalidObject
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition: TObject.h:68
TGraph::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
Definition: TGraph.cxx:2139
TPluginHandler::ExecPlugin
Longptr_t ExecPlugin(int nargs, const T &... params)
Definition: TPluginManager.h:165
TVirtualGraphPainter
Abstract interface to a histogram painter.
Definition: TVirtualGraphPainter.h:27
TGraph::GetPoint
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:1607
ymin
float ymin
Definition: THbookFile.cxx:95
TMath::BinarySearch
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMathBase.h:278
TString::Index
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition: TString.h:639
TGraph::LeastSquareLinearFit
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:1998
Printf
void Printf(const char *fmt,...)
TGraph::MovePoints
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:2595
gSystem
R__EXTERN TSystem * gSystem
Definition: TSystem.h:559
TAxis::GetNoExponent
Bool_t GetNoExponent() const
Definition: TAxis.h:122
TVirtualGraphPainter::PaintGrapHist
virtual void PaintGrapHist(TGraph *theGraph, Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)=0
TGraph::ComputeRange
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
TAttLine::SetLineWidth
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
TGraph::Set
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:2219
TVectorT
TVectorT.
Definition: TVectorT.h:27
TVirtualFitter::GetUserFunc
virtual TObject * GetUserFunc() const
Definition: TVirtualFitter.h:84
TGraph::SetHighlight
virtual void SetHighlight(Bool_t set=kTRUE)
Set highlight (enable/disble) mode for the graph by default highlight mode is disable.
Definition: TGraph.cxx:2253
TGraph::CopyAndRelease
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
TStyle::SetHistLineColor
void SetHistLineColor(Color_t color=1)
Definition: TStyle.h:363
TGraph::SwapValues
static void SwapValues(Double_t *arr, Int_t pos1, Int_t pos2)
Swap values.
Definition: TGraph.cxx:2515
TPluginHandler::LoadPlugin
Int_t LoadPlugin()
Load the plugin library for this handler.
Definition: TPluginManager.cxx:254
TVirtualGraphPainter::DistancetoPrimitiveHelper
virtual Int_t DistancetoPrimitiveHelper(TGraph *theGraph, Int_t px, Int_t py)=0
f1
TF1 * f1
Definition: legend1.C:11
TString::IsNull
Bool_t IsNull() const
Definition: TString.h:407
Double_t
double Double_t
Definition: RtypesCore.h:59
TGraph::AllocateArrays
Double_t ** AllocateArrays(Int_t Narrays, Int_t arraySize)
Allocate arrays.
Definition: TGraph.cxx:553
TGraph
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
TSpline3
Class to create third splines to interpolate knots Arbitrary conditions can be introduced for first a...
Definition: TSpline.h:201
TPluginManager.h
TAxis::SetNoExponent
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
TF1.h
TList::Remove
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:822
TGraph::GetObjectInfo
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:1661
ROOT::Fit::kGraph
@ kGraph
Definition: HFitInterface.h:50
TH1F
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:575
TGraph::GetMinimum
Double_t GetMinimum() const
Definition: TGraph.h:144
TAttFill::SetFillColor
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
TAttAxis::GetTitleFont
virtual Style_t GetTitleFont() const
Definition: TAttAxis.h:47
TObject
Mother of all ROOT objects.
Definition: TObject.h:37
TVirtualGraphPainter::PaintStats
virtual void PaintStats(TGraph *theGraph, TF1 *fit)=0
TGraph::GetErrorXhigh
virtual Double_t GetErrorXhigh(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1438
TGraph::Fit
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
TH1
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
TVirtualGraphPainter::PaintGraph
virtual void PaintGraph(TGraph *theGraph, Int_t npoints, const Double_t *x, const Double_t *y, Option_t *chopt)=0
TGraph::GetN
Int_t GetN() const
Definition: TGraph.h:124
name
char name[80]
Definition: TGX11.cxx:110
TGraph::SetPointY
virtual void SetPointY(Int_t i, Double_t y)
Set y value for point i.
Definition: TGraph.cxx:2316
d
#define d(i)
Definition: RSha256.hxx:102
TGraph::GetEditable
Bool_t GetEditable() const
Return kTRUE if kNotEditable bit is not set, kFALSE otherwise.
Definition: TGraph.cxx:2234
TAttAxis::SetNdivisions
virtual void SetNdivisions(Int_t n=510, Bool_t optim=kTRUE)
Set the number of divisions for this axis.
Definition: TAttAxis.cxx:228
TGraph::CompareY
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
TGraph::InsertPointBefore
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:1816
x2
static const double x2[5]
Definition: RooGaussKronrodIntegrator1D.cxx:364
gPad
#define gPad
Definition: TVirtualPad.h:287
TAttFill::GetFillColor
virtual Color_t GetFillColor() const
Return the fill area color.
Definition: TAttFill.h:30
TGraph::GetRMS
virtual Double_t GetRMS(Int_t axis=1) const
Return RMS of X (axis=1) or Y (axis=2)
Definition: TGraph.cxx:1396
TIter
Definition: TCollection.h:233
TAttLine::DistancetoLine
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
TGraph::fNpoints
Int_t fNpoints
Number of points <= fMaxSize.
Definition: TGraph.h:46
TH1::kNoStats
@ kNoStats
Don't draw stats box.
Definition: TH1.h:164
TVirtualGraphPainter::GetPainter
static TVirtualGraphPainter * GetPainter()
Static function returning a pointer to the current graph painter.
Definition: TVirtualGraphPainter.cxx:29
TStyle::GetHistLineWidth
Width_t GetHistLineWidth() const
Definition: TStyle.h:227
TCollection
Collection abstract base class.
Definition: TCollection.h:63
TAttFill
Fill Area Attributes class.
Definition: TAttFill.h:19
TGraph::GetFunction
TF1 * GetFunction(const char *name) const
Return pointer to function with name.
Definition: TGraph.cxx:1479
TAxis::GetCenterTitle
Bool_t GetCenterTitle() const
Definition: TAxis.h:115
TAttLine::GetLineWidth
virtual Width_t GetLineWidth() const
Return the line width.
Definition: TAttLine.h:35
TF1::SetParLimits
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set limits for parameter ipar.
Definition: TF1.cxx:3518
TAttAxis::GetLabelFont
virtual Style_t GetLabelFont() const
Definition: TAttAxis.h:39
TGraph::GetErrorYhigh
virtual Double_t GetErrorYhigh(Int_t bin) const
This function is called by GraphFitChisquare.
Definition: TGraph.cxx:1458
TAttMarker::SetMarkerStyle
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TStyle::GetHistLineColor
Color_t GetHistLineColor() const
Definition: TStyle.h:224
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
TF1::SetParent
virtual void SetParent(TObject *p=0)
Definition: TF1.h:665
TGraph::Allocate
virtual Double_t ** Allocate(Int_t newsize)
Allocate internal data structures for newsize points.
Definition: TGraph.cxx:545
TAxis::RotateTitle
void RotateTitle(Bool_t rotate=kTRUE)
Rotate title by 180 degrees.
Definition: TAxis.h:193
TVirtualGraphPainter::ExecuteEventHelper
virtual void ExecuteEventHelper(TGraph *theGraph, Int_t event, Int_t px, Int_t py)=0
TVirtualFitter::GetFitter
static TVirtualFitter * GetFitter()
static: return the current Fitter
Definition: TVirtualFitter.cxx:209
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
TF1
1-Dim function class
Definition: TF1.h:213
TGraph::kResetHisto
@ kResetHisto
fHistogram must be reset in GetHistogram
Definition: TGraph.h:71
TGraph::GetMean
virtual Double_t GetMean(Int_t axis=1) const
Return mean value of X (axis=1) or Y (axis=2)
Definition: TGraph.cxx:1381
Class
void Class()
Definition: Class.C:29
TString::ToLower
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1147
TGraph::FitPanel
virtual void FitPanel()
Display a GUI panel with all graph fit options.
Definition: TGraph.cxx:1331
TH1::GetXaxis
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:320
TGraph::kIsHighlight
@ kIsHighlight
Bit set if graph is highlight.
Definition: TGraph.h:74
TList::RecursiveRemove
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
TAttLine::SetLineStyle
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
TGraph::TGraph
TGraph()
Graph default constructor.
Definition: TGraph.cxx:87
TAttAxis::SetLabelFont
virtual void SetLabelFont(Style_t font=62)
Set labels' font.
Definition: TAttAxis.cxx:182
TH1.h
TGraph::PaintStats
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:2071
TGraph::SetMinimum
virtual void SetMinimum(Double_t minimum=-1111)
Set the minimum of the graph.
Definition: TGraph.cxx:2275
TGraph::GetPointX
virtual Double_t GetPointX(Int_t i) const
Get x value for point i.
Definition: TGraph.cxx:1618
TGraph::PaintGraph
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:2053
TAttAxis::GetNdivisions
virtual Int_t GetNdivisions() const
Definition: TAttAxis.h:36
TGraph::Browse
virtual void Browse(TBrowser *b)
Browse.
Definition: TGraph.cxx:589
TAttAxis::GetTitleSize
virtual Float_t GetTitleSize() const
Definition: TAttAxis.h:44
TGraph::~TGraph
virtual ~TGraph()
Graph default destructor.
Definition: TGraph.cxx:521
TList
A doubly linked list.
Definition: TList.h:44
TVectorT::GetLwb
Int_t GetLwb() const
Definition: TVectorT.h:73
TStyle::GetHistLineStyle
Style_t GetHistLineStyle() const
Definition: TStyle.h:226
TMath.h
gROOT
#define gROOT
Definition: TROOT.h:406
TGraph::ExpandAndCopy
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
TGraph::InitExpo
virtual void InitExpo(Double_t xmin=0, Double_t xmax=0)
Compute Initial values of parameters for an exponential.
Definition: TGraph.cxx:1730
int
TGraph::CompareX
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
g
#define g(i)
Definition: RSha256.hxx:105
TAttAxis::GetTitleOffset
virtual Float_t GetTitleOffset() const
Definition: TAttAxis.h:43