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