ROOT  6.06/09
Reference Guide
TMultiGraph.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Author: Rene Brun 12/10/2000
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #include "TROOT.h"
13 #include "TEnv.h"
14 #include "TBrowser.h"
15 #include "TMultiGraph.h"
16 #include "TGraph.h"
17 #include "TH1.h"
18 #include "TH2.h"
19 #include "TPolyLine3D.h"
20 #include "TVirtualPad.h"
21 #include "Riostream.h"
22 #include "TVirtualFitter.h"
23 #include "TPluginManager.h"
24 #include "TClass.h"
25 #include "TMath.h"
26 #include "TSystem.h"
27 #include <stdlib.h>
28 
29 #include "HFitInterface.h"
30 #include "Fit/DataRange.h"
31 #include "Math/MinimizerOptions.h"
32 
33 #include <ctype.h>
34 
35 extern void H1LeastSquareSeqnd(Int_t n, Double_t *a, Int_t idim, Int_t &ifail, Int_t k, Double_t *b);
36 
38 
39 
40 ////////////////////////////////////////////////////////////////////////////////
41 
42 /** \class TMultiGraph
43  \ingroup Hist
44 A TMultiGraph is a collection of TGraph (or derived) objects. It allows to
45 manipulate a set of graphs as a single entity. In particular, when drawn,
46 the X and Y axis ranges are automatically computed such as all the graphs
47 will be visible.
48 <p>
49 <tt>TMultiGraph::Add</tt> should be used to add a new graph to the list.
50 <p>
51 The TMultiGraph owns the objects in the list.
52 <p>
53 The drawing options are the same as for TGraph.
54 Like for TGraph, the painting is performed thanks to the TGraphPainter
55 class. All details about the various painting options are given in this class.
56 
57 Example:
58 <pre>
59  TGraph *gr1 = new TGraph(...
60  TGraphErrors *gr2 = new TGraphErrors(...
61  TMultiGraph *mg = new TMultiGraph();
62  mg->Add(gr1,"lp");
63  mg->Add(gr2,"cp");
64  mg->Draw("a");
65 </pre>
66 A special option <tt>3D</tt> allows to draw the graphs in a 3D space. See the
67 following example:
68 
69 Begin_Macro(source)
70 {
71  c0 = new TCanvas("c1","multigraph L3",200,10,700,500);
72  c0->SetFrameFillColor(30);
73 
74  TMultiGraph *mg = new TMultiGraph();
75 
76  TGraph *gr1 = new TGraph(); gr1->SetLineColor(kBlue);
77  TGraph *gr2 = new TGraph(); gr2->SetLineColor(kRed);
78  TGraph *gr3 = new TGraph(); gr3->SetLineColor(kGreen);
79  TGraph *gr4 = new TGraph(); gr4->SetLineColor(kOrange);
80 
81  Double_t dx = 6.28/100;
82  Double_t x = -3.14;
83 
84  for (int i=0; i<=100; i++) {
85  x = x+dx;
86  gr1->SetPoint(i,x,2.*TMath::Sin(x));
87  gr2->SetPoint(i,x,TMath::Cos(x));
88  gr3->SetPoint(i,x,TMath::Cos(x*x));
89  gr4->SetPoint(i,x,TMath::Cos(x*x*x));
90  }
91 
92  mg->Add(gr4); gr4->SetTitle("Cos(x*x*x)"); gr4->SetLineWidth(3);
93  mg->Add(gr3); gr3->SetTitle("Cos(x*x)") ; gr3->SetLineWidth(3);
94  mg->Add(gr2); gr2->SetTitle("Cos(x)") ; gr2->SetLineWidth(3);
95  mg->Add(gr1); gr1->SetTitle("2*Sin(x)") ; gr1->SetLineWidth(3);
96 
97  mg->Draw("a fb l3d");
98  return c0;
99 }
100 End_Macro
101 
102 <p>
103 The number of graphs in a multigraph can be retrieve with:
104 <pre>
105 mg->GetListOfGraphs()->GetSize();
106 </pre>
107 <p>
108 The drawing option for each TGraph may be specified as an optional
109 second argument of the <tt>Add</tt> function.
110 <p>
111 If a draw option is specified, it will be used to draw the graph,
112 otherwise the graph will be drawn with the option specified in
113 <tt>TMultiGraph::Draw</tt>.
114 <p>
115 The following example shows how to fit a TMultiGraph.
116 
117 Begin_Macro(source)
118 {
119  TCanvas *c1 = new TCanvas("c1","c1",600,400);
120 
121  Double_t px1[2] = {2.,4.};
122  Double_t dx1[2] = {0.1,0.1};
123  Double_t py1[2] = {2.1,4.0};
124  Double_t dy1[2] = {0.3,0.2};
125 
126  Double_t px2[2] = {3.,5.};
127  Double_t dx2[2] = {0.1,0.1};
128  Double_t py2[2] = {3.2,4.8};
129  Double_t dy2[2] = {0.3,0.2};
130 
131  gStyle->SetOptFit(0001);
132 
133  TGraphErrors *g1 = new TGraphErrors(2,px1,py1,dx1,dy1);
134  g1->SetMarkerStyle(21);
135  g1->SetMarkerColor(2);
136 
137  TGraphErrors *g2 = new TGraphErrors(2,px2,py2,dx2,dy2);
138  g2->SetMarkerStyle(22);
139  g2->SetMarkerColor(3);
140 
141  TMultiGraph *g = new TMultiGraph();
142  g->Add(g1);
143  g->Add(g2);
144 
145  g->Draw("AP");
146 
147  g->Fit("pol1","FQ");
148  return c1;
149 }
150 End_Macro
151 
152 <p>
153 The axis titles can be modified the following way:
154 <p>
155 <pre>
156  [...]
157  TMultiGraph *mg = new TMultiGraph;
158  mg->SetTitle("title;xaxis title; yaxis title");
159  mg->Add(g1);
160  mg->Add(g2);
161  mg->Draw("apl");
162 </pre>
163 <p>
164 When the graphs in a TMultiGraph are fitted, the fit parameters boxes
165 overlap. The following example shows how to make them all visible.
166 
167 
168 Begin_Macro(source)
169 ../../../tutorials/graphs/multigraph.C
170 End_Macro
171 
172 <p>
173 The axis limits can be changed the like for TGraph. The same methods apply on
174 the multigraph.
175 Note the two differents ways to change limits on X and Y axis.
176 
177 Begin_Macro(source)
178 {
179  TCanvas *c2 = new TCanvas("c2","c2",600,400);
180 
181  TGraph *g[3];
182  Double_t x[10] = {0,1,2,3,4,5,6,7,8,9};
183  Double_t y[10] = {1,2,3,4,5,5,4,3,2,1};
184  TMultiGraph *mg = new TMultiGraph();
185  for (int i=0; i<3; i++) {
186  g[i] = new TGraph(10, x, y);
187  g[i]->SetMarkerStyle(20);
188  g[i]->SetMarkerColor(i+2);
189  for (int j=0; j<10; j++) y[j] = y[j]-1;
190  mg->Add(g[i]);
191  }
192  mg->Draw("APL");
193  mg->GetXaxis()->SetTitle("E_{#gamma} (GeV)");
194  mg->GetYaxis()->SetTitle("Coefficients");
195 
196  // Change the axis limits
197  gPad->Modified();
198  mg->GetXaxis()->SetLimits(1.5,7.5);
199  mg->SetMinimum(0.);
200  mg->SetMaximum(10.);
201 
202  return c2;
203 }
204 End_Macro
205 
206 <p>
207 The method TPad::BuildLegend is able to extract the graphs inside a
208 multigraph. The following example demonstrate this.
209 
210 Begin_Macro(source)
211 {
212  TCanvas *c3 = new TCanvas("c3","c3",600, 400);
213 
214  TMultiGraph * mg = new TMultiGraph("mg","mg");
215 
216  const Int_t size = 10;
217 
218  double px[size];
219  double py1[size];
220  double py2[size];
221  double py3[size];
222 
223  for ( int i = 0; i < size ; ++i ) {
224  px[i] = i;
225  py1[i] = size - i;
226  py2[i] = size - 0.5 * i;
227  py3[i] = size - 0.6 * i;
228  }
229 
230  TGraph * gr1 = new TGraph( size, px, py1 );
231  gr1->SetName("gr1");
232  gr1->SetTitle("graph 1");
233  gr1->SetMarkerStyle(21);
234  gr1->SetDrawOption("AP");
235  gr1->SetLineColor(2);
236  gr1->SetLineWidth(4);
237  gr1->SetFillStyle(0);
238 
239  TGraph * gr2 = new TGraph( size, px, py2 );
240  gr2->SetName("gr2");
241  gr2->SetTitle("graph 2");
242  gr2->SetMarkerStyle(22);
243  gr2->SetMarkerColor(2);
244  gr2->SetDrawOption("P");
245  gr2->SetLineColor(3);
246  gr2->SetLineWidth(4);
247  gr2->SetFillStyle(0);
248 
249  TGraph * gr3 = new TGraph( size, px, py3 );
250  gr3->SetName("gr3");
251  gr3->SetTitle("graph 3");
252  gr3->SetMarkerStyle(23);
253  gr3->SetLineColor(4);
254  gr3->SetLineWidth(4);
255  gr3->SetFillStyle(0);
256 
257  mg->Add( gr1 );
258  mg->Add( gr2 );
259 
260  gr3->Draw("ALP");
261  mg->Draw("LP");
262  c3->BuildLegend();
263 
264  return c3;
265 }
266 End_Macro
267 */
268 
269 
270 ////////////////////////////////////////////////////////////////////////////////
271 /// TMultiGraph default constructor.
272 
274 {
275  fGraphs = 0;
276  fFunctions = 0;
277  fHistogram = 0;
278  fMaximum = -1111;
279  fMinimum = -1111;
280 }
281 
282 
283 ////////////////////////////////////////////////////////////////////////////////
284 /// Constructor with name and title.
285 
286 TMultiGraph::TMultiGraph(const char *name, const char *title)
287  : TNamed(name,title)
288 {
289  fGraphs = 0;
290  fFunctions = 0;
291  fHistogram = 0;
292  fMaximum = -1111;
293  fMinimum = -1111;
294 }
295 
296 
297 ////////////////////////////////////////////////////////////////////////////////
298 /// Copy constructor.
299 
301  TNamed (mg),
302  fGraphs(mg.fGraphs),
303  fFunctions(mg.fFunctions),
304  fHistogram(mg.fHistogram),
305  fMaximum(mg.fMaximum),
306  fMinimum(mg.fMinimum)
307 {
308 }
309 
310 
311 ////////////////////////////////////////////////////////////////////////////////
312 /// Assignement operator.
313 
315 {
316  if (this!=&mg) {
317  TNamed::operator=(mg);
318  fGraphs=mg.fGraphs;
321  fMaximum=mg.fMaximum;
322  fMinimum=mg.fMinimum;
323  }
324  return *this;
325 }
326 
327 
328 ////////////////////////////////////////////////////////////////////////////////
329 /// TMultiGraph destructor.
330 
332 {
333  if (!fGraphs) return;
334  TGraph *g;
335  TIter next(fGraphs);
336  while ((g = (TGraph*) next())) {
338  }
339  fGraphs->Delete();
340  delete fGraphs;
341  fGraphs = 0;
342  delete fHistogram;
343  fHistogram = 0;
344  if (fFunctions) {
346  //special logic to support the case where the same object is
347  //added multiple times in fFunctions.
348  //This case happens when the same object is added with different
349  //drawing modes
350  TObject *obj;
351  while ((obj = fFunctions->First())) {
352  while (fFunctions->Remove(obj)) { }
353  delete obj;
354  }
355  delete fFunctions;
356  }
357 }
358 
359 
360 ////////////////////////////////////////////////////////////////////////////////
361 /// Add a new graph to the list of graphs.
362 /// Note that the graph is now owned by the TMultigraph.
363 /// Deleting the TMultiGraph object will automatically delete the graphs.
364 /// You should not delete the graphs when the TMultigraph is still active.
365 
366 void TMultiGraph::Add(TGraph *graph, Option_t *chopt)
367 {
368  if (!fGraphs) fGraphs = new TList();
369  graph->SetBit(kMustCleanup);
370  fGraphs->Add(graph,chopt);
371 }
372 
373 
374 ////////////////////////////////////////////////////////////////////////////////
375 /// Add all the graphs in "multigraph" to the list of graphs.
376 /// If "chopt" is defined all the graphs in "multigraph" will be added with
377 /// the "chopt" option.
378 /// If "chopt" is undefined each graph will be added with the option it had
379 /// in "multigraph".
380 
381 void TMultiGraph::Add(TMultiGraph *multigraph, Option_t *chopt)
382 {
383  TList *graphlist = multigraph->GetListOfGraphs();
384  if (!graphlist) return;
385 
386  if (!fGraphs) fGraphs = new TList();
387 
388  TObjOptLink *lnk = (TObjOptLink*)graphlist->FirstLink();
389  TObject *obj = 0;
390 
391  while (lnk) {
392  obj = lnk->GetObject();
393  if (!strlen(chopt)) fGraphs->Add(obj,lnk->GetOption());
394  else fGraphs->Add(obj,chopt);
395  lnk = (TObjOptLink*)lnk->Next();
396  }
397 }
398 
399 
400 ////////////////////////////////////////////////////////////////////////////////
401 /// Browse multigraph.
402 
404 {
405  TString opt = gEnv->GetValue("TGraph.BrowseOption", "");
406  if (opt.IsNull()) {
407  opt = b ? b->GetDrawOption() : "alp";
408  opt = (opt == "") ? "alp" : opt.Data();
409  }
410  Draw(opt.Data());
411  gPad->Update();
412 }
413 
414 
415 ////////////////////////////////////////////////////////////////////////////////
416 /// Compute distance from point px,py to each graph.
417 
419 {
420  // Are we on the axis?
421  const Int_t kMaxDiff = 10;
422  Int_t distance = 9999;
423  if (fHistogram) {
424  distance = fHistogram->DistancetoPrimitive(px,py);
425  if (distance <= 0) return distance;
426  }
427 
428  // Loop on the list of graphs
429  if (!fGraphs) return distance;
430  TGraph *g;
431  TIter next(fGraphs);
432  while ((g = (TGraph*) next())) {
433  Int_t dist = g->DistancetoPrimitive(px,py);
434  if (dist <= 0) return 0;
435  if (dist < kMaxDiff) {gPad->SetSelected(g); return dist;}
436  }
437  return distance;
438 }
439 
440 
441 ////////////////////////////////////////////////////////////////////////////////
442 /// Draw this multigraph with its current attributes.
443 ///
444 /// Options to draw a graph are described in TGraphPainter.
445 ///
446 /// The drawing option for each TGraph may be specified as an optional
447 /// second argument of the Add function. You can use GetGraphDrawOption
448 /// to return this option.
449 /// If a draw option is specified, it will be used to draw the graph,
450 /// otherwise the graph will be drawn with the option specified in
451 /// TMultiGraph::Draw. Use GetDrawOption to return the option specified
452 /// when drawing the TMultiGraph.
453 
454 void TMultiGraph::Draw(Option_t *option)
455 {
456  TString opt = option;
457  opt.ToLower();
458 
459  if (gPad) {
460  if (!gPad->IsEditable()) gROOT->MakeDefCanvas();
461  if (opt.Contains("a")) gPad->Clear();
462  }
463  AppendPad(option);
464 }
465 
466 
467 ////////////////////////////////////////////////////////////////////////////////
468 /// Fit this graph with function with name fname.
469 ///
470 /// interface to TF1::Fit(TF1 *f1...
471 
472 TFitResultPtr TMultiGraph::Fit(const char *fname, Option_t *option, Option_t *, Axis_t xmin, Axis_t xmax)
473 {
474  char *linear;
475  linear= (char*)strstr(fname, "++");
476  TF1 *f1=0;
477  if (linear)
478  f1=new TF1(fname, fname, xmin, xmax);
479  else {
480  f1 = (TF1*)gROOT->GetFunction(fname);
481  if (!f1) { Printf("Unknown function: %s",fname); return -1; }
482  }
483 
484  return Fit(f1,option,"",xmin,xmax);
485 }
486 
487 
488 ////////////////////////////////////////////////////////////////////////////////
489 /// Fit this multigraph with function f1.
490 ///
491 /// In this function all graphs of the multigraph are fitted simultaneously
492 ///
493 /// f1 is an already predefined function created by TF1.
494 /// Predefined functions such as gaus, expo and poln are automatically
495 /// created by ROOT.
496 ///
497 /// The list of fit options is given in parameter option.
498 /// option = "W" Set all errors to 1
499 /// = "U" Use a User specified fitting algorithm (via SetFCN)
500 /// = "Q" Quiet mode (minimum printing)
501 /// = "V" Verbose mode (default is between Q and V)
502 /// = "B" Use this option when you want to fix one or more parameters
503 /// and the fitting function is like "gaus","expo","poln","landau".
504 /// = "R" Use the Range specified in the function range
505 /// = "N" Do not store the graphics function, do not draw
506 /// = "0" Do not plot the result of the fit. By default the fitted function
507 /// is drawn unless the option"N" above is specified.
508 /// = "+" Add this new fitted function to the list of fitted functions
509 /// (by default, any previous function is deleted)
510 /// = "C" In case of linear fitting, not calculate the chisquare
511 /// (saves time)
512 /// = "F" If fitting a polN, switch to minuit fitter
513 /// = "ROB" In case of linear fitting, compute the LTS regression
514 /// coefficients (robust(resistant) regression), using
515 /// the default fraction of good points
516 /// "ROB=0.x" - compute the LTS regression coefficients, using
517 /// 0.x as a fraction of good points
518 ///
519 /// When the fit is drawn (by default), the parameter goption may be used
520 /// to specify a list of graphics options. See TGraph::Paint for a complete
521 /// list of these options.
522 ///
523 /// In order to use the Range option, one must first create a function
524 /// with the expression to be fitted. For example, if your graph
525 /// has a defined range between -4 and 4 and you want to fit a gaussian
526 /// only in the interval 1 to 3, you can do:
527 /// TF1 *f1 = new TF1("f1","gaus",1,3);
528 /// graph->Fit("f1","R");
529 ///
530 ///
531 /// who is calling this function
532 /// ============================
533 /// Note that this function is called when calling TGraphErrors::Fit
534 /// or TGraphAsymmErrors::Fit ot TGraphBentErrors::Fit
535 /// see the discussion below on the errors calulation.
536 ///
537 /// Setting initial conditions
538 /// ==========================
539 /// Parameters must be initialized before invoking the Fit function.
540 /// The setting of the parameter initial values is automatic for the
541 /// predefined functions : poln, expo, gaus, landau. One can however disable
542 /// this automatic computation by specifying the option "B".
543 /// You can specify boundary limits for some or all parameters via
544 /// f1->SetParLimits(p_number, parmin, parmax);
545 /// if parmin>=parmax, the parameter is fixed
546 /// Note that you are not forced to fix the limits for all parameters.
547 /// For example, if you fit a function with 6 parameters, you can do:
548 /// func->SetParameters(0,3.1,1.e-6,0.1,-8,100);
549 /// func->SetParLimits(4,-10,-4);
550 /// func->SetParLimits(5, 1,1);
551 /// With this setup, parameters 0->3 can vary freely
552 /// Parameter 4 has boundaries [-10,-4] with initial value -8
553 /// Parameter 5 is fixed to 100.
554 ///
555 /// Fit range
556 /// =========
557 /// The fit range can be specified in two ways:
558 /// - specify rxmax > rxmin (default is rxmin=rxmax=0)
559 /// - specify the option "R". In this case, the function will be taken
560 /// instead of the full graph range.
561 ///
562 /// Changing the fitting function
563 /// =============================
564 /// By default a chi2 fitting function is used for fitting the TGraphs's.
565 /// The function is implemented in FitUtil::EvaluateChi2.
566 /// In case of TGraphErrors an effective chi2 is used
567 /// (see TGraphErrors fit in TGraph::Fit) and is implemented in
568 /// FitUtil::EvaluateChi2Effective
569 /// To specify a User defined fitting function, specify option "U" and
570 /// call the following functions:
571 /// TVirtualFitter::Fitter(mygraph)->SetFCN(MyFittingFunction)
572 /// where MyFittingFunction is of type:
573 /// extern void MyFittingFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
574 ///
575 /// Access to the fit result
576 /// ========================
577 /// The function returns a TFitResultPtr which can hold a pointer to a TFitResult object.
578 /// By default the TFitResultPtr contains only the status of the fit and it converts
579 /// automatically to an integer. If the option "S" is instead used, TFitResultPtr contains
580 /// the TFitResult and behaves as a smart pointer to it. For example one can do:
581 /// TFitResultPtr r = graph->Fit("myFunc","S");
582 /// TMatrixDSym cov = r->GetCovarianceMatrix(); // to access the covariance matrix
583 /// Double_t par0 = r->Parameter(0); // retrieve the value for the parameter 0
584 /// Double_t err0 = r->ParError(0); // retrieve the error for the parameter 0
585 /// r->Print("V"); // print full information of fit including covariance matrix
586 /// r->Write(); // store the result in a file
587 ///
588 /// The fit parameters, error and chi2 (but not covariance matrix) can be retrieved also
589 /// from the fitted function.
590 ///
591 ///
592 /// Associated functions
593 /// ====================
594 /// One or more object (typically a TF1*) can be added to the list
595 /// of functions (fFunctions) associated to each graph.
596 /// When TGraph::Fit is invoked, the fitted function is added to this list.
597 /// Given a graph gr, one can retrieve an associated function
598 /// with: TF1 *myfunc = gr->GetFunction("myfunc");
599 ///
600 /// If the graph is made persistent, the list of
601 /// associated functions is also persistent. Given a pointer (see above)
602 /// to an associated function myfunc, one can retrieve the function/fit
603 /// parameters with calls such as:
604 /// Double_t chi2 = myfunc->GetChisquare();
605 /// Double_t par0 = myfunc->GetParameter(0); //value of 1st parameter
606 /// Double_t err0 = myfunc->GetParError(0); //error on first parameter
607 ///
608 /// Fit Statistics
609 /// ==============
610 /// You can change the statistics box to display the fit parameters with
611 /// the TStyle::SetOptFit(mode) method. This mode has four digits.
612 /// mode = pcev (default = 0111)
613 /// v = 1; print name/values of parameters
614 /// e = 1; print errors (if e=1, v must be 1)
615 /// c = 1; print Chisquare/Number of degress of freedom
616 /// p = 1; print Probability
617 ///
618 /// For example: gStyle->SetOptFit(1011);
619 /// prints the fit probability, parameter names/values, and errors.
620 /// You can change the position of the statistics box with these lines
621 /// (where g is a pointer to the TGraph):
622 ///
623 /// Root > TPaveStats *st = (TPaveStats*)g->GetListOfFunctions()->FindObject("stats")
624 /// Root > st->SetX1NDC(newx1); //new x start position
625 /// Root > st->SetX2NDC(newx2); //new x end position
626 
627 TFitResultPtr TMultiGraph::Fit(TF1 *f1, Option_t *option, Option_t *goption, Axis_t rxmin, Axis_t rxmax)
628 {
629  // internal multigraph fitting methods
630  Foption_t fitOption;
632 
633  // create range and minimizer options with default values
634  ROOT::Fit::DataRange range(rxmin,rxmax);
636  return ROOT::Fit::FitObject(this, f1 , fitOption , minOption, goption, range);
637 
638 }
639 
640 ////////////////////////////////////////////////////////////////////////////////
641 /// Display a panel with all histogram fit options.
642 /// See class TFitPanel for example
643 
645 {
646  if (!gPad)
647  gROOT->MakeDefCanvas();
648 
649  if (!gPad) {
650  Error("FitPanel", "Unable to create a default canvas");
651  return;
652  }
653 
654  // use plugin manager to create instance of TFitEditor
655  TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
656  if (handler && handler->LoadPlugin() != -1) {
657  if (handler->ExecPlugin(2, gPad, this) == 0)
658  Error("FitPanel", "Unable to crate the FitPanel");
659  }
660  else
661  Error("FitPanel", "Unable to find the FitPanel plug-in");
662 }
663 
664 ////////////////////////////////////////////////////////////////////////////////
665 /// Return the draw option for the TGraph gr in this TMultiGraph.
666 /// The return option is the one specified when calling TMultiGraph::Add(gr,option).
667 
669 {
670  if (!fGraphs || !gr) return "";
672  TObject *obj;
673  while ((obj = next())) {
674  if (obj == (TObject*)gr) return next.GetOption();
675  }
676  return "";
677 }
678 
679 
680 ////////////////////////////////////////////////////////////////////////////////
681 /// Compute Initial values of parameters for a gaussian.
682 
684 {
685  Double_t allcha, sumx, sumx2, x, val, rms, mean;
686  Int_t bin;
687  const Double_t sqrtpi = 2.506628;
688 
689  // Compute mean value and RMS of the graph in the given range
690  Int_t np = 0;
691  allcha = sumx = sumx2 = 0;
692  TGraph *g;
693  TIter next(fGraphs);
694  Double_t *px, *py;
695  Int_t npp; //number of points in each graph
696  while ((g = (TGraph*) next())) {
697  px=g->GetX();
698  py=g->GetY();
699  npp=g->GetN();
700  for (bin=0; bin<npp; bin++) {
701  x=px[bin];
702  if (x<xmin || x>xmax) continue;
703  np++;
704  val=py[bin];
705  sumx+=val*x;
706  sumx2+=val*x*x;
707  allcha+=val;
708  }
709  }
710  if (np == 0 || allcha == 0) return;
711  mean = sumx/allcha;
712  rms = TMath::Sqrt(sumx2/allcha - mean*mean);
713 
714  Double_t binwidx = TMath::Abs((xmax-xmin)/np);
715  if (rms == 0) rms = 1;
717  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
718  f1->SetParameter(0,binwidx*allcha/(sqrtpi*rms));
719  f1->SetParameter(1,mean);
720  f1->SetParameter(2,rms);
721  f1->SetParLimits(2,0,10*rms);
722 }
723 
724 
725 ////////////////////////////////////////////////////////////////////////////////
726 /// Compute Initial values of parameters for an exponential.
727 
729 {
730  Double_t constant, slope;
731  Int_t ifail;
732 
733  LeastSquareLinearFit(-1, constant, slope, ifail, xmin, xmax);
734 
736  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
737  f1->SetParameter(0,constant);
738  f1->SetParameter(1,slope);
739 }
740 
741 
742 ////////////////////////////////////////////////////////////////////////////////
743 /// Compute Initial values of parameters for a polynom.
744 
746 {
747  Double_t fitpar[25];
748 
750  TF1 *f1 = (TF1*)grFitter->GetUserFunc();
751  Int_t npar = f1->GetNpar();
752 
753  LeastSquareFit(npar, fitpar, xmin, xmax);
754 
755  for (Int_t i=0;i<npar;i++) f1->SetParameter(i, fitpar[i]);
756 }
757 
758 
759 ////////////////////////////////////////////////////////////////////////////////
760 /// Least squares lpolynomial fitting without weights.
761 ///
762 /// m number of parameters
763 /// a array of parameters
764 /// first 1st point number to fit (default =0)
765 /// last last point number to fit (default=fNpoints-1)
766 ///
767 /// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
768 
770 {
771  const Double_t zero = 0.;
772  const Double_t one = 1.;
773  const Int_t idim = 20;
774 
775  Double_t b[400] /* was [20][20] */;
776  Int_t i, k, l, ifail, bin;
777  Double_t power;
778  Double_t da[20], xk, yk;
779 
780 
781  //count the total number of points to fit
782  TGraph *g;
783  TIter next(fGraphs);
784  Double_t *px, *py;
785  Int_t n=0;
786  Int_t npp;
787  while ((g = (TGraph*) next())) {
788  px=g->GetX();
789  npp=g->GetN();
790  for (bin=0; bin<npp; bin++) {
791  xk=px[bin];
792  if (xk < xmin || xk > xmax) continue;
793  n++;
794  }
795  }
796  if (m <= 2) {
797  LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
798  return;
799  }
800  if (m > idim || m > n) return;
801  da[0] = zero;
802  for (l = 2; l <= m; ++l) {
803  b[l-1] = zero;
804  b[m + l*20 - 21] = zero;
805  da[l-1] = zero;
806  }
807  Int_t np = 0;
808 
809  next.Reset();
810  while ((g = (TGraph*) next())) {
811  px=g->GetX();
812  py=g->GetY();
813  npp=g->GetN();
814 
815  for (k = 0; k <= npp; ++k) {
816  xk = px[k];
817  if (xk < xmin || xk > xmax) continue;
818  np++;
819  yk = py[k];
820  power = one;
821  da[0] += yk;
822  for (l = 2; l <= m; ++l) {
823  power *= xk;
824  b[l-1] += power;
825  da[l-1] += power*yk;
826  }
827  for (l = 2; l <= m; ++l) {
828  power *= xk;
829  b[m + l*20 - 21] += power;
830  }
831  }
832  }
833  b[0] = Double_t(np);
834  for (i = 3; i <= m; ++i) {
835  for (k = i; k <= m; ++k) {
836  b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
837  }
838  }
839  H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
840 
841  if (ifail < 0) {
842  //a[0] = fY[0];
843  py=((TGraph *)fGraphs->First())->GetY();
844  a[0]=py[0];
845  for (i=1; i<m; ++i) a[i] = 0;
846  return;
847  }
848  for (i=0; i<m; ++i) a[i] = da[i];
849 }
850 
851 
852 ////////////////////////////////////////////////////////////////////////////////
853 /// Least square linear fit without weights.
854 ///
855 /// Fit a straight line (a0 + a1*x) to the data in this graph.
856 /// ndata: number of points to fit
857 /// first: first point number to fit
858 /// last: last point to fit O(ndata should be last-first
859 /// ifail: return parameter indicating the status of the fit (ifail=0, fit is OK)
860 ///
861 /// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
862 
865 {
866  Double_t xbar, ybar, x2bar;
867  Int_t i;
868  Double_t xybar;
869  Double_t fn, xk, yk;
870  Double_t det;
871 
872  ifail = -2;
873  xbar = ybar = x2bar = xybar = 0;
874  Int_t np = 0;
875  TGraph *g;
876  TIter next(fGraphs);
877  Double_t *px, *py;
878  Int_t npp;
879  while ((g = (TGraph*) next())) {
880  px=g->GetX();
881  py=g->GetY();
882  npp=g->GetN();
883  for (i = 0; i < npp; ++i) {
884  xk = px[i];
885  if (xk < xmin || xk > xmax) continue;
886  np++;
887  yk = py[i];
888  if (ndata < 0) {
889  if (yk <= 0) yk = 1e-9;
890  yk = TMath::Log(yk);
891  }
892  xbar += xk;
893  ybar += yk;
894  x2bar += xk*xk;
895  xybar += xk*yk;
896  }
897  }
898  fn = Double_t(np);
899  det = fn*x2bar - xbar*xbar;
900  ifail = -1;
901  if (det <= 0) {
902  if (fn > 0) a0 = ybar/fn;
903  else a0 = 0;
904  a1 = 0;
905  return;
906  }
907  ifail = 0;
908  a0 = (x2bar*ybar - xbar*xybar) / det;
909  a1 = (fn*xybar - xbar*ybar) / det;
910 }
911 
912 
913 ////////////////////////////////////////////////////////////////////////////////
914 /// Return 1 if the point (x,y) is inside one of the graphs 0 otherwise.
915 
917 {
918  Int_t in = 0;
919  if (!fGraphs) return in;
920  TGraph *g;
921  TIter next(fGraphs);
922  while ((g = (TGraph*) next())) {
923  in = g->IsInside(x, y);
924  if (in) return in;
925  }
926  return in;
927 }
928 
929 
930 ////////////////////////////////////////////////////////////////////////////////
931 /// Returns a pointer to the histogram used to draw the axis.
932 /// Takes into account the two following cases.
933 /// 1- option 'A' was specified in TMultiGraph::Draw. Return fHistogram
934 /// 2- user had called TPad::DrawFrame. return pointer to hframe histogram
935 
937 {
938  if (fHistogram) return fHistogram;
939  if (!gPad) return 0;
940  gPad->Modified();
941  gPad->Update();
942  if (fHistogram) return fHistogram;
943  TH1F *h1 = (TH1F*)gPad->FindObject("hframe");
944  return h1;
945 }
946 
947 
948 ////////////////////////////////////////////////////////////////////////////////
949 /// Return pointer to function with name.
950 ///
951 /// Functions such as TGraph::Fit store the fitted function in the list of
952 /// functions of this graph.
953 
954 TF1 *TMultiGraph::GetFunction(const char *name) const
955 {
956  if (!fFunctions) return 0;
957  return (TF1*)fFunctions->FindObject(name);
958 }
959 
960 ////////////////////////////////////////////////////////////////////////////////
961 /// Return pointer to list of functions.
962 /// If pointer is null create the list
963 
965 {
966  if (!fFunctions) fFunctions = new TList();
967  return fFunctions;
968 }
969 
970 
971 ////////////////////////////////////////////////////////////////////////////////
972 /// Get x axis of the graph.
973 /// This method returns a valid axis only after the TMultigraph has been drawn.
974 
976 {
977  if (!gPad) return 0;
978  TH1 *h = GetHistogram();
979  if (!h) return 0;
980  return h->GetXaxis();
981 }
982 
983 
984 ////////////////////////////////////////////////////////////////////////////////
985 /// Get y axis of the graph.
986 /// This method returns a valid axis only after the TMultigraph has been drawn.
987 
989 {
990  if (!gPad) return 0;
991  TH1 *h = GetHistogram();
992  if (!h) return 0;
993  return h->GetYaxis();
994 }
995 
996 
997 ////////////////////////////////////////////////////////////////////////////////
998 /// Paint all the graphs of this multigraph.
999 
1000 void TMultiGraph::Paint(Option_t *option)
1002  const TPickerStackGuard pushGuard(this);
1003 
1004  if (!fGraphs) return;
1005  if (fGraphs->GetSize() == 0) return;
1006 
1007  char *l;
1008  Int_t nch = strlen(option);
1009 
1010  TString chopt = option;
1011  chopt.ToUpper();
1012 
1013  l = (char*)strstr(chopt.Data(),"3D");
1014  if (l) {
1015  l = (char*)strstr(chopt.Data(),"L");
1016  if (l) PaintPolyLine3D(chopt.Data());
1017  return;
1018  }
1019 
1020  l = (char*)strstr(chopt.Data(),"PADS");
1021  if (l) {
1022  chopt.ReplaceAll("PADS","");
1023  PaintPads(chopt.Data());
1024  return;
1025  }
1026 
1027  TGraph *g;
1028 
1029  l = (char*)strstr(chopt.Data(),"A");
1030  if (l) {
1031  *l = ' ';
1032  TIter next(fGraphs);
1033  Int_t npt = 100;
1034  Double_t maximum, minimum, rwxmin, rwxmax, rwymin, rwymax, uxmin, uxmax, dx, dy;
1035  rwxmin = gPad->GetUxmin();
1036  rwxmax = gPad->GetUxmax();
1037  rwymin = gPad->GetUymin();
1038  rwymax = gPad->GetUymax();
1039  char *xtitle = 0;
1040  char *ytitle = 0;
1041  Int_t firstx = 0;
1042  Int_t lastx = 0;
1043  Bool_t timedisplay = kFALSE;
1044  char *timeformat = 0;
1045 
1046  if (fHistogram) {
1047  //cleanup in case of a previous unzoom
1048  if (fHistogram->GetMinimum() >= fHistogram->GetMaximum()) {
1049  nch = strlen(fHistogram->GetXaxis()->GetTitle());
1050  firstx = fHistogram->GetXaxis()->GetFirst();
1051  lastx = fHistogram->GetXaxis()->GetLast();
1052  timedisplay = fHistogram->GetXaxis()->GetTimeDisplay();
1053  if (nch) {
1054  xtitle = new char[nch+1];
1055  strlcpy(xtitle,fHistogram->GetXaxis()->GetTitle(),nch+1);
1056  }
1057  nch = strlen(fHistogram->GetYaxis()->GetTitle());
1058  if (nch) {
1059  ytitle = new char[nch+1];
1060  strlcpy(ytitle,fHistogram->GetYaxis()->GetTitle(),nch+1);
1061  }
1062  nch = strlen(fHistogram->GetXaxis()->GetTimeFormat());
1063  if (nch) {
1064  timeformat = new char[nch+1];
1065  strlcpy(timeformat,fHistogram->GetXaxis()->GetTimeFormat(),nch+1);
1066  }
1067  delete fHistogram;
1068  fHistogram = 0;
1069  }
1070  }
1071  if (fHistogram) {
1072  minimum = fHistogram->GetYaxis()->GetXmin();
1073  maximum = fHistogram->GetYaxis()->GetXmax();
1074  uxmin = gPad->PadtoX(rwxmin);
1075  uxmax = gPad->PadtoX(rwxmax);
1076  } else {
1077  Bool_t initialrangeset = kFALSE;
1078  while ((g = (TGraph*) next())) {
1079  if (g->GetN() <= 0) continue;
1080  if (initialrangeset) {
1081  Double_t rx1,ry1,rx2,ry2;
1082  g->ComputeRange(rx1, ry1, rx2, ry2);
1083  if (rx1 < rwxmin) rwxmin = rx1;
1084  if (ry1 < rwymin) rwymin = ry1;
1085  if (rx2 > rwxmax) rwxmax = rx2;
1086  if (ry2 > rwymax) rwymax = ry2;
1087  } else {
1088  g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1089  initialrangeset = kTRUE;
1090  }
1091  if (g->GetN() > npt) npt = g->GetN();
1092  }
1093  if (rwxmin == rwxmax) rwxmax += 1.;
1094  if (rwymin == rwymax) rwymax += 1.;
1095  dx = 0.05*(rwxmax-rwxmin);
1096  dy = 0.05*(rwymax-rwymin);
1097  uxmin = rwxmin - dx;
1098  uxmax = rwxmax + dx;
1099  if (gPad->GetLogy()) {
1100  if (rwymin <= 0) rwymin = 0.001*rwymax;
1101  minimum = rwymin/(1+0.5*TMath::Log10(rwymax/rwymin));
1102  maximum = rwymax*(1+0.2*TMath::Log10(rwymax/rwymin));
1103  } else {
1104  minimum = rwymin - dy;
1105  maximum = rwymax + dy;
1106  }
1107  if (minimum < 0 && rwymin >= 0) minimum = 0;
1108  if (maximum > 0 && rwymax <= 0) maximum = 0;
1109  }
1110 
1111  if (fMinimum != -1111) rwymin = minimum = fMinimum;
1112  if (fMaximum != -1111) rwymax = maximum = fMaximum;
1113  if (uxmin < 0 && rwxmin >= 0) {
1114  if (gPad->GetLogx()) uxmin = 0.9*rwxmin;
1115  //else uxmin = 0;
1116  }
1117  if (uxmax > 0 && rwxmax <= 0) {
1118  if (gPad->GetLogx()) uxmax = 1.1*rwxmax;
1119  //else uxmax = 0;
1120  }
1121  if (minimum < 0 && rwymin >= 0) {
1122  if (gPad->GetLogy()) minimum = 0.9*rwymin;
1123  //else minimum = 0;
1124  }
1125  if (maximum > 0 && rwymax <= 0) {
1126  if (gPad->GetLogy()) maximum = 1.1*rwymax;
1127  //else maximum = 0;
1128  }
1129  if (minimum <= 0 && gPad->GetLogy()) minimum = 0.001*maximum;
1130  if (uxmin <= 0 && gPad->GetLogx()) {
1131  if (uxmax > 1000) uxmin = 1;
1132  else uxmin = 0.001*uxmax;
1133  }
1134  rwymin = minimum;
1135  rwymax = maximum;
1136  if (fHistogram) {
1137  fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1138  }
1139 
1140  // Create a temporary histogram to draw the axis
1141  if (!fHistogram) {
1142  // the graph is created with at least as many channels as there are points
1143  // to permit zooming on the full range
1144  rwxmin = uxmin;
1145  rwxmax = uxmax;
1146  fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1147  if (!fHistogram) return;
1148  fHistogram->SetMinimum(rwymin);
1150  fHistogram->SetMaximum(rwymax);
1151  fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1153  if (xtitle) {fHistogram->GetXaxis()->SetTitle(xtitle); delete [] xtitle;}
1154  if (ytitle) {fHistogram->GetYaxis()->SetTitle(ytitle); delete [] ytitle;}
1155  if (firstx != lastx) fHistogram->GetXaxis()->SetRange(firstx,lastx);
1156  if (timedisplay) {fHistogram->GetXaxis()->SetTimeDisplay(timedisplay);}
1157  if (timeformat) {fHistogram->GetXaxis()->SetTimeFormat(timeformat); delete [] timeformat;}
1158  }
1159  fHistogram->Paint("0");
1160  }
1161 
1162  TGraph *gfit = 0;
1163  if (fGraphs) {
1165  TObject *obj = 0;
1166 
1167  chopt.ReplaceAll("A","");
1168 
1169  while (lnk) {
1170 
1171  obj = lnk->GetObject();
1172 
1173  gPad->PushSelectableObject(obj);
1174 
1175  if (!gPad->PadInHighlightMode() || (gPad->PadInHighlightMode() && obj == gPad->GetSelected())) {
1176  TString opt = lnk->GetOption();
1177  if (!opt.IsWhitespace())
1178  obj->Paint(opt.ReplaceAll("A","").Data());
1179  else {
1180  if (!chopt.IsWhitespace()) obj->Paint(chopt.Data());
1181  else obj->Paint("L");
1182  }
1183  }
1184 
1185  lnk = (TObjOptLink*)lnk->Next();
1186  }
1187 
1188  gfit = (TGraph*)obj; // pick one TGraph in the list to paint the fit parameters.
1189  }
1190 
1191  TObject *f;
1192  TF1 *fit = 0;
1193  if (fFunctions) {
1195  while ((f = (TObject*) next())) {
1196  if (f->InheritsFrom(TF1::Class())) {
1197  if (f->TestBit(TF1::kNotDraw) == 0) f->Paint("lsame");
1198  fit = (TF1*)f;
1199  } else {
1200  f->Paint();
1201  }
1202  }
1203  }
1204 
1205  if (fit) gfit->PaintStats(fit);
1206 }
1207 
1208 
1209 ////////////////////////////////////////////////////////////////////////////////
1210 /// Divides the active pad and draws all Graphs in the Multigraph separately.
1211 
1212 void TMultiGraph::PaintPads(Option_t *option)
1214  TIter next(fGraphs);
1215  Int_t neededPads = fGraphs->GetSize();
1216  Int_t existingPads = 0;
1217  TString opt = (TString)option;
1218 
1219  TVirtualPad *curPad = gPad;
1220  TObject *obj;
1221  TIter nextPad(curPad->GetListOfPrimitives());
1222 
1223  while ((obj = nextPad())) {
1224  if (obj->InheritsFrom(TVirtualPad::Class())) existingPads++;
1225  }
1226  if (existingPads < neededPads) {
1227  curPad->Clear();
1228  Int_t nx = (Int_t)TMath::Sqrt((Double_t)neededPads);
1229  if (nx*nx < neededPads) nx++;
1230  Int_t ny = nx;
1231  if (((nx*ny)-nx) >= neededPads) ny--;
1232  curPad->Divide(nx,ny);
1233  }
1234  Int_t i = 0;
1235  TGraph *g;
1236 
1238  obj = 0;
1239 
1240  while (lnk) {
1241  g = (TGraph*)lnk->GetObject();
1242  i++;
1243  curPad->cd(i);
1244  TString apopt = lnk->GetOption();
1245  if (strlen(apopt)) {
1246  g->Draw((apopt.Append("A")).Data());
1247  } else {
1248  if (strlen(opt)) g->Draw(opt.Append("A"));
1249  else g->Draw("LA");
1250  }
1251  lnk = (TObjOptLink*)lnk->Next();
1252  }
1253 
1254  curPad->cd();
1255 }
1256 
1257 
1258 ////////////////////////////////////////////////////////////////////////////////
1259 /// Paint all the graphs of this multigraph as 3D lines.
1260 
1263  Int_t i, npt=0;
1264  char *l;
1265  Double_t rwxmin=0., rwxmax=0., rwymin=0., rwymax=0.;
1266  TIter next(fGraphs);
1267  TGraph *g;
1268 
1269  g = (TGraph*) next();
1270  if (g) {
1271  g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1272  npt = g->GetN();
1273  }
1274 
1275  while ((g = (TGraph*) next())) {
1276  Double_t rx1,ry1,rx2,ry2;
1277  g->ComputeRange(rx1, ry1, rx2, ry2);
1278  if (rx1 < rwxmin) rwxmin = rx1;
1279  if (ry1 < rwymin) rwymin = ry1;
1280  if (rx2 > rwxmax) rwxmax = rx2;
1281  if (ry2 > rwymax) rwymax = ry2;
1282  if (g->GetN() > npt) npt = g->GetN();
1283  }
1284 
1285  Int_t ndiv = fGraphs->GetSize();
1286  TH2F* frame = new TH2F("frame","", ndiv, 0., (Double_t)(ndiv),
1287  10, rwxmin, rwxmax);
1288 
1289  TAxis *Xaxis = frame->GetXaxis();
1290  Xaxis->SetNdivisions(-ndiv);
1291  next.Reset();
1292  for (i=ndiv; i>=1; i--) {
1293  g = (TGraph*) next();
1294  Xaxis->SetBinLabel(i, g->GetTitle());
1295  }
1296 
1297  frame->SetStats(kFALSE);
1298  frame->SetMinimum(rwymin);
1299  frame->SetMaximum(rwymax);
1300 
1301  l = (char*)strstr(option,"A");
1302  if (l) frame->Paint("lego0,fb,bb");
1303  l = (char*)strstr(option,"BB");
1304  if (!l) frame->Paint("lego0,fb,a,same");
1305 
1306  Double_t *x, *y;
1307  Double_t xyz1[3], xyz2[3];
1308 
1309  next.Reset();
1310  Int_t j = ndiv;
1311  while ((g = (TGraph*) next())) {
1312  npt = g->GetN();
1313  x = g->GetX();
1314  y = g->GetY();
1315  gPad->SetLineColor(g->GetLineColor());
1316  gPad->SetLineWidth(g->GetLineWidth());
1317  gPad->SetLineStyle(g->GetLineStyle());
1318  gPad->TAttLine::Modify();
1319  for (i=0; i<npt-1; i++) {
1320  xyz1[0] = j-0.5;
1321  xyz1[1] = x[i];
1322  xyz1[2] = y[i];
1323  xyz2[0] = j-0.5;
1324  xyz2[1] = x[i+1];
1325  xyz2[2] = y[i+1];
1326  gPad->PaintLine3D(xyz1, xyz2);
1327  }
1328  j--;
1329  }
1330 
1331  l = (char*)strstr(option,"FB");
1332  if (!l) frame->Paint("lego0,bb,a,same");
1333  delete frame;
1334 }
1335 
1336 
1337 ////////////////////////////////////////////////////////////////////////////////
1338 /// Print the list of graphs.
1339 
1340 void TMultiGraph::Print(Option_t *option) const
1342  TGraph *g;
1343  if (fGraphs) {
1344  TIter next(fGraphs);
1345  while ((g = (TGraph*) next())) {
1346  g->Print(option);
1347  }
1348  }
1349 }
1350 
1351 
1352 ////////////////////////////////////////////////////////////////////////////////
1353 /// Recursively remove this object from a list. Typically implemented
1354 /// by classes that can contain multiple references to a same object.
1355 
1358  if (!fGraphs) return;
1359  TObject *objr = fGraphs->Remove(obj);
1360  if (!objr) return;
1361  delete fHistogram; fHistogram = 0;
1362  if (gPad) gPad->Modified();
1363 }
1364 
1365 
1366 ////////////////////////////////////////////////////////////////////////////////
1367 /// Save primitive as a C++ statement(s) on output stream out.
1368 
1369 void TMultiGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
1371  char quote = '"';
1372  out<<" "<<std::endl;
1373  if (gROOT->ClassSaved(TMultiGraph::Class())) {
1374  out<<" ";
1375  } else {
1376  out<<" TMultiGraph *";
1377  }
1378  out<<"multigraph = new TMultiGraph();"<<std::endl;
1379  out<<" multigraph->SetName("<<quote<<GetName()<<quote<<");"<<std::endl;
1380  out<<" multigraph->SetTitle("<<quote<<GetTitle()<<quote<<");"<<std::endl;
1381 
1382  if (fGraphs) {
1384  TObject *g;
1385 
1386  while (lnk) {
1387  g = lnk->GetObject();
1388  g->SavePrimitive(out, Form("multigraph%s",lnk->GetOption()));
1389  lnk = (TObjOptLink*)lnk->Next();
1390  }
1391  }
1392  const char *l = strstr(option,"th2poly");
1393  if (l) {
1394  out<<" "<<l+7<<"->AddBin(multigraph);"<<std::endl;
1395  } else {
1396  out<<" multigraph->Draw(" <<quote<<option<<quote<<");"<<std::endl;
1397  }
1398  TAxis *xaxis = GetXaxis();
1399  TAxis *yaxis = GetYaxis();
1400 
1401  if (xaxis) xaxis->SaveAttributes(out, "multigraph","->GetXaxis()");
1402  if (yaxis) yaxis->SaveAttributes(out, "multigraph","->GetYaxis()");
1403 }
1404 
1405 
1406 ////////////////////////////////////////////////////////////////////////////////
1407 /// Set multigraph maximum.
1408 
1409 void TMultiGraph::SetMaximum(Double_t maximum)
1411  fMaximum = maximum;
1412  if (fHistogram) fHistogram->SetMaximum(maximum);
1413 }
1414 
1415 
1416 ////////////////////////////////////////////////////////////////////////////////
1417 /// Set multigraph minimum.
1418 
1419 void TMultiGraph::SetMinimum(Double_t minimum)
1421  fMinimum = minimum;
1422  if (fHistogram) fHistogram->SetMinimum(minimum);
1423 }
const int nx
Definition: kalman.C:16
const int ndata
void PaintPads(Option_t *chopt="")
Divides the active pad and draws all Graphs in the Multigraph separately.
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:429
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:52
virtual Style_t GetLineStyle() const
Definition: TAttLine.h:48
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:4450
virtual void SetMinimum(Double_t minimum=-1111)
Set multigraph minimum.
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.
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
virtual void SaveAttributes(std::ostream &out, const char *name, const char *subname)
Save axis attributes as C++ statement(s) on output stream out.
Definition: TAxis.cxx:631
virtual void Paint(Option_t *option="")
Control routine to paint any kind of histograms.
Definition: TH1.cxx:5798
float xmin
Definition: THbookFile.cxx:93
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:404
virtual void FitPanel()
Display a panel with all histogram fit options.
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:394
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:154
Double_t Log(Double_t x)
Definition: TMath.h:526
Option_t * GetDrawOption() const
Get option used by the graphics system to draw this object.
Definition: TBrowser.h:108
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:487
virtual void SetTimeFormat(const char *format="")
Change the format used for time plotting.
Definition: TAxis.cxx:943
ClassImp(TSeqCollection) Int_t TSeqCollection TIter next(this)
Return index of object in collection.
TList * GetListOfGraphs() const
Definition: TMultiGraph.h:71
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:8266
THist< 2, float > TH2F
Definition: THist.h:321
const char Option_t
Definition: RtypesCore.h:62
virtual void LeastSquareLinearFit(Int_t ndata, Double_t &a0, Double_t &a1, Int_t &ifail, Double_t xmin, Double_t xmax)
Least square linear fit without weights.
double Axis_t
Definition: RtypesCore.h:72
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:635
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:1771
Double_t distance(const TPoint2 &p1, const TPoint2 &p2)
Definition: CsgOps.cxx:467
TH1 * h
Definition: legend2.C:5
virtual void InitPolynom(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for a polynom.
virtual void SetNdivisions(Int_t n=510, Bool_t optim=kTRUE)
Set the number of divisions for this axis.
Definition: TAttAxis.cxx:211
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:37
virtual Bool_t GetTimeDisplay() const
Definition: TAxis.h:130
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1101
THist< 1, float > TH1F
Definition: THist.h:315
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:395
#define gROOT
Definition: TROOT.h:340
Int_t LoadPlugin()
Load the plugin library for this handler.
Basic string class.
Definition: TString.h:137
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:570
virtual TObject * GetUserFunc() const
Option_t * GetOption() const
Returns the object option stored in the list.
Definition: TList.cxx:959
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1088
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual void Paint(Option_t *chopt="")
Paint all the graphs of this multigraph.
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:496
virtual void RecursiveRemove(TObject *obj)
Recursively remove this object from a list.
Int_t GetN() const
Definition: TGraph.h:132
Long_t ExecPlugin(int nargs, const T &...params)
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:740
virtual TVirtualPad * cd(Int_t subpadnumber=0)=0
void Reset()
Definition: TCollection.h:161
Iterator of linked list.
Definition: TList.h:187
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:732
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:164
const char * Data() const
Definition: TString.h:349
Double_t * GetY() const
Definition: TGraph.h:140
Double_t x[n]
Definition: legend1.C:17
void Class()
Definition: Class.C:29
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
const int ny
Definition: kalman.C:17
Double_t Log10(Double_t x)
Definition: TMath.h:529
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:640
TString & Append(const char *cs)
Definition: TString.h:492
std::vector< std::vector< double > > Data
virtual ~TMultiGraph()
TMultiGraph destructor.
TList * GetListOfFunctions()
Return pointer to list of functions.
virtual void Draw(Option_t *chopt="")
Draw this multigraph with its current attributes.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
TH1F * h1
Definition: legend1.C:5
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition: TVirtualPad.h:59
TMultiGraph & operator=(const TMultiGraph &)
Assignement operator.
TAxis * GetXaxis() const
Get x axis of the graph.
Double_t GetXmin() const
Definition: TAxis.h:137
virtual void LeastSquareFit(Int_t m, Double_t *a, Double_t xmin, Double_t xmax)
Least squares lpolynomial fitting without weights.
char * out
Definition: TBase64.cxx:29
virtual void SetTimeDisplay(Int_t value)
Definition: TAxis.h:161
A doubly linked list.
Definition: TList.h:47
virtual const char * GetTimeFormat() const
Definition: TAxis.h:131
virtual void Print(Option_t *chopt="") const
Print the list of graphs.
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a line.
Definition: TH1.cxx:2612
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:41
virtual void SetRange(Int_t first=0, Int_t last=0)
Set the viewing range for the axis from bin first to last.
Definition: TAxis.cxx:831
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set limits for parameter ipar.
Definition: TF1.cxx:3206
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:40
Double_t * GetX() const
Definition: TGraph.h:139
Class to manage histogram axis.
Definition: TAxis.h:36
virtual Option_t * GetGraphDrawOption(const TGraph *gr) const
Return the draw option for the TGraph gr in this TMultiGraph.
virtual void Browse(TBrowser *b)
Browse multigraph.
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:256
virtual Int_t GetValue(const char *name, Int_t dflt)
Returns the integer value for a resource.
Definition: TEnv.cxx:494
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:674
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:33
Bool_t TestBit(UInt_t f) const
Definition: TObject.h:173
TMarker * m
Definition: textangle.C:8
char * Form(const char *fmt,...)
TList * fFunctions
Definition: TMultiGraph.h:41
virtual void Clear(Option_t *option="")=0
static TVirtualFitter * GetFitter()
static: return the current Fitter
TLine * l
Definition: textangle.C:4
Bool_t IsWhitespace() const
Definition: TString.h:388
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
TAxis * GetYaxis()
Definition: TH1.h:320
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:133
virtual TList * GetListOfPrimitives() const =0
float xmax
Definition: THbookFile.cxx:93
virtual void InitExpo(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for an exponential.
Bool_t IsNull() const
Definition: TString.h:387
virtual TObjLink * FirstLink() const
Definition: TList.h:101
virtual Color_t GetLineColor() const
Definition: TAttLine.h:47
TGraphErrors * gr
Definition: legend1.C:25
#define Printf
Definition: TGeoToOCC.h:18
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to each graph.
virtual void Print(Option_t *chopt="") const
Print graph values.
Definition: TGraph.cxx:1943
virtual void Paint(Option_t *option="")
This method must be overridden if a class wants to paint itself.
Definition: TObject.cxx:563
virtual Int_t GetSize() const
Definition: TCollection.h:95
#define ClassImp(name)
Definition: Rtypes.h:279
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:34
double f(double x)
double Double_t
Definition: RtypesCore.h:55
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:1934
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a graph.
Definition: TGraph.cxx:781
R__EXTERN TEnv * gEnv
Definition: TEnv.h:174
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
Double_t GetXmax() const
Definition: TAxis.h:138
Double_t y[n]
Definition: legend1.C:17
The TH1 histogram class.
Definition: TH1.h:80
Double_t fMaximum
Definition: TMultiGraph.h:43
virtual void SetBinLabel(Int_t bin, const char *label)
Set label for bin.
Definition: TAxis.cxx:793
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis.
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition: HFitImpl.cxx:673
#define name(a, b)
Definition: linkTestLib0.cpp:5
virtual void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0)=0
Abstract Base Class for Fitting.
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:440
void PaintPolyLine3D(Option_t *chopt="")
Paint all the graphs of this multigraph as 3D lines.
Mother of all ROOT objects.
Definition: TObject.h:58
virtual TObject * First() const
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:556
TList * fGraphs
Definition: TMultiGraph.h:40
virtual void SetMaximum(Double_t maximum=-1111)
Set multigraph maximum.
virtual void Add(TObject *obj)
Definition: TList.h:81
TH1F * fHistogram
Definition: TMultiGraph.h:42
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:567
1-Dim function class
Definition: TF1.h:149
A Graph is a graphics object made of two arrays X and Y with npoints each.
Definition: TGraph.h:53
TF1 * f1
Definition: legend1.C:11
#define gPad
Definition: TVirtualPad.h:288
void ResetBit(UInt_t f)
Definition: TObject.h:172
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:431
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:935
virtual void Add(TGraph *graph, Option_t *chopt="")
Add a new graph to the list of graphs.
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
const Bool_t kTRUE
Definition: Rtypes.h:91
virtual Double_t GetMaximum(Double_t maxval=FLT_MAX) const
Return maximum value smaller than maxval of bins in the range, unless the value has been overridden b...
Definition: TH1.cxx:7921
virtual Width_t GetLineWidth() const
Definition: TAttLine.h:49
virtual void SetTitle(const char *title="")
Change (i.e. set) the title of the TNamed.
Definition: TNamed.cxx:152
TObject * obj
TF1 * GetFunction(const char *name) const
Return pointer to function with name.
virtual void InitGaus(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for a gaussian.
const Int_t n
Definition: legend1.C:16
TAxis * GetYaxis() const
Get y axis of the graph.
Double_t fMinimum
Definition: TMultiGraph.h:44
virtual Int_t IsInside(Double_t x, Double_t y) const
Return 1 if the point (x,y) is inside one of the graphs 0 otherwise.
virtual Int_t GetNpar() const
Definition: TF1.h:349
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8320
virtual Double_t GetMinimum(Double_t minval=-FLT_MAX) const
Return minimum value larger than minval of bins in the range, unless the value has been overridden by...
Definition: TH1.cxx:8006
TAxis * GetXaxis()
Definition: TH1.h:319
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TObject.cxx:702