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