Logo ROOT  
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 "TVirtualPad.h"
20#include "TVirtualFitter.h"
21#include "TPluginManager.h"
22#include "TMath.h"
23#include "TF1.h"
24#include "strlcpy.h"
25
26#include "HFitInterface.h"
27#include "Fit/DataRange.h"
29
30#include <iostream>
31#include <cstdlib>
32#include <cctype>
33
34extern void H1LeastSquareSeqnd(Int_t n, Double_t *a, Int_t idim, Int_t &ifail, Int_t k, Double_t *b);
35
37
38
39////////////////////////////////////////////////////////////////////////////////
40
41/** \class TMultiGraph
42 \ingroup Graphs
43 \brief A TMultiGraph is a collection of TGraph (or derived) objects.
44
45- [Introduction](\ref MG00)
46- [MultiGraphs' drawing](\ref MG01)
47 - [Setting drawing options](\ref MG01a)
48 - [Titles setting](\ref MG01b)
49 - [The option \"3D\"](\ref MG01c)
50 - [Legend drawing](\ref MG01d)
51 - [Automatic coloring](\ref MG01e)
52 - [Reverse axis](\ref MG01f)
53- [MultiGraphs' fitting](\ref MG02)
54 - [Fit box position](\ref MG02a)
55- [Axis' limits setting](\ref MG03)
56
57
58\anchor MG00
59### Introduction
60
61A TMultiGraph allows to manipulate a set of graphs as a single entity. In particular,
62when drawn, the X and Y axis ranges are automatically computed such as all the graphs
63will be visible.
64
65`TMultiGraph::Add` should be used to add a new graph to the list.
66
67The TMultiGraph owns the objects in the list.
68
69The number of graphs in a multigraph can be retrieve with:
70~~~ {.cpp}
71mg->GetListOfGraphs()->GetEntries();
72~~~
73
74\anchor MG01
75### MultiGraphs' Drawing
76
77The drawing options are the same as for TGraph.
78Like for TGraph, the painting is performed thanks to the TGraphPainter
79class. All details about the various painting options are given in this class.
80
81Example:
82~~~ {.cpp}
83 TGraph *gr1 = new TGraph(...
84 TGraphErrors *gr2 = new TGraphErrors(...
85 TMultiGraph *mg = new TMultiGraph();
86 mg->Add(gr1,"lp");
87 mg->Add(gr2,"cp");
88 mg->Draw("a");
89~~~
90
91\anchor MG01a
92#### Setting drawing options
93
94The drawing option for each TGraph may be specified as an optional
95second argument of the `Add` function.
96
97If a draw option is specified, it will be used to draw the graph,
98otherwise the graph will be drawn with the option specified in
99`TMultiGraph::Draw`
100
101\anchor MG01b
102#### Titles setting
103
104The global title and the axis titles can be modified the following way:
105
106~~~ {.cpp}
107 [...]
108 auto mg = new TMultiGraph;
109 mg->SetTitle("title;xaxis title; yaxis title");
110 mg->Add(g1);
111 mg->Add(g2);
112 mg->Draw("apl");
113~~~
114
115\anchor MG01c
116#### The option "3D"
117
118A special option `3D` allows to draw the graphs in a 3D space. See the
119following example:
120
121Begin_Macro(source)
122{
123 auto c0 = new TCanvas("c1","multigraph L3",200,10,700,500);
124
125 auto mg = new TMultiGraph();
126
127 auto gr1 = new TGraph(); gr1->SetLineColor(kBlue);
128 auto gr2 = new TGraph(); gr2->SetLineColor(kRed);
129 auto gr3 = new TGraph(); gr3->SetLineColor(kGreen);
130 auto gr4 = new TGraph(); gr4->SetLineColor(kOrange);
131
132 Double_t dx = 6.28/1000;
133 Double_t x = -3.14;
134
135 for (int i=0; i<=1000; i++) {
136 x = x+dx;
137 gr1->SetPoint(i,x,2.*TMath::Sin(x));
138 gr2->SetPoint(i,x,TMath::Cos(x));
139 gr3->SetPoint(i,x,TMath::Cos(x*x));
140 gr4->SetPoint(i,x,TMath::Cos(x*x*x));
141 }
142
143 mg->Add(gr4); gr4->SetTitle("Cos(x*x*x)"); gr4->SetLineWidth(3);
144 mg->Add(gr3); gr3->SetTitle("Cos(x*x)") ; gr3->SetLineWidth(3);
145 mg->Add(gr2); gr2->SetTitle("Cos(x)") ; gr2->SetLineWidth(3);
146 mg->Add(gr1); gr1->SetTitle("2*Sin(x)") ; gr1->SetLineWidth(3);
147
148 mg->SetTitle("Multi-graph Title; X-axis Title; Y-axis Title");
149
150 mg->Draw("a fb l3d");
151
152 mg->GetHistogram()->GetXaxis()->SetRangeUser(0.,2.5);
153 gPad->Modified();
154 gPad->Update();
155}
156End_Macro
157
158\anchor MG01d
159#### Legend drawing
160
161The method TPad::BuildLegend is able to extract the graphs inside a
162multigraph. The following example demonstrate this.
163
164Begin_Macro(source)
165{
166 auto c3 = new TCanvas("c3","c3",600, 400);
167
168 auto mg = new TMultiGraph("mg","mg");
169
170 const Int_t size = 10;
171
172 double px[size];
173 double py1[size];
174 double py2[size];
175 double py3[size];
176
177 for ( int i = 0; i < size ; ++i ) {
178 px[i] = i;
179 py1[i] = size - i;
180 py2[i] = size - 0.5 * i;
181 py3[i] = size - 0.6 * i;
182 }
183
184 auto gr1 = new TGraph( size, px, py1 );
185 gr1->SetName("gr1");
186 gr1->SetTitle("graph 1");
187 gr1->SetMarkerStyle(21);
188 gr1->SetDrawOption("AP");
189 gr1->SetLineColor(2);
190 gr1->SetLineWidth(4);
191 gr1->SetFillStyle(0);
192
193 auto gr2 = new TGraph( size, px, py2 );
194 gr2->SetName("gr2");
195 gr2->SetTitle("graph 2");
196 gr2->SetMarkerStyle(22);
197 gr2->SetMarkerColor(2);
198 gr2->SetDrawOption("P");
199 gr2->SetLineColor(3);
200 gr2->SetLineWidth(4);
201 gr2->SetFillStyle(0);
202
203 auto gr3 = new TGraph( size, px, py3 );
204 gr3->SetName("gr3");
205 gr3->SetTitle("graph 3");
206 gr3->SetMarkerStyle(23);
207 gr3->SetLineColor(4);
208 gr3->SetLineWidth(4);
209 gr3->SetFillStyle(0);
210
211 mg->Add( gr1 );
212 mg->Add( gr2 );
213
214 gr3->Draw("ALP");
215 mg->Draw("LP");
216 c3->BuildLegend();
217}
218End_Macro
219
220\anchor MG01e
221#### Automatic coloring
222
223Automatic coloring according to the current palette is available as shown in the
224following example:
225
226Begin_Macro(source)
227../../../tutorials/graphs/multigraphpalettecolor.C
228End_Macro
229
230\anchor MG01f
231#### Reverse axis
232
233\since **ROOT version 6.19/02**
234
235When a TMultiGraph is drawn, the X-axis is drawn with increasing values from left to
236right and the Y-axis from bottom to top. The two options RX and RY allow to change
237this order. The option RX allows to draw the X-axis with increasing values from
238right to left and the RY option allows to draw the Y-axis with increasing values
239from top to bottom. The following example illustrate how to use these options.
240
241Begin_Macro(source)
242{
243 auto *c = new TCanvas();
244 c->Divide(2,1);
245
246 auto *g1 = new TGraphErrors();
247 g1->SetPoint(0,-4,-3);
248 g1->SetPoint(1,1,1);
249 g1->SetPoint(2,2,1);
250 g1->SetPoint(3,3,4);
251 g1->SetPoint(4,5,5);
252 g1->SetPointError(0,1.,2.);
253 g1->SetPointError(1,2,1);
254 g1->SetPointError(2,2,3);
255 g1->SetPointError(3,3,2);
256 g1->SetPointError(4,4,5);
257 g1->SetMarkerStyle(21);
258
259 auto *g2 = new TGraph();
260 g2->SetPoint(0,4,8);
261 g2->SetPoint(1,5,9);
262 g2->SetPoint(2,6,10);
263 g2->SetPoint(3,10,11);
264 g2->SetPoint(4,15,12);
265 g2->SetLineColor(kRed);
266 g2->SetLineWidth(5);
267
268 auto mg = new TMultiGraph();
269 mg->Add(g1,"P");
270 mg->Add(g2,"L");
271
272 c->cd(1); gPad->SetGrid(1,1);
273 mg->Draw("A");
274
275 c->cd(2); gPad->SetGrid(1,1);
276 mg->Draw("A RX RY");
277}
278End_Macro
279
280\anchor MG02
281### MultiGraphs' fitting
282
283The following example shows how to fit a TMultiGraph.
284
285Begin_Macro(source)
286{
287 auto c1 = new TCanvas("c1","c1",600,400);
288
289 Double_t px1[2] = {2.,4.};
290 Double_t dx1[2] = {0.1,0.1};
291 Double_t py1[2] = {2.1,4.0};
292 Double_t dy1[2] = {0.3,0.2};
293
294 Double_t px2[2] = {3.,5.};
295 Double_t dx2[2] = {0.1,0.1};
296 Double_t py2[2] = {3.2,4.8};
297 Double_t dy2[2] = {0.3,0.2};
298
299 gStyle->SetOptFit(0001);
300
301 auto g1 = new TGraphErrors(2,px1,py1,dx1,dy1);
302 g1->SetMarkerStyle(21);
303 g1->SetMarkerColor(2);
304
305 auto g2 = new TGraphErrors(2,px2,py2,dx2,dy2);
306 g2->SetMarkerStyle(22);
307 g2->SetMarkerColor(3);
308
309 auto g = new TMultiGraph();
310 g->Add(g1);
311 g->Add(g2);
312
313 g->Draw("AP");
314
315 g->Fit("pol1","FQ");
316}
317End_Macro
318
319\anchor MG02a
320#### Fit box position
321
322When the graphs in a TMultiGraph are fitted, the fit parameters boxes
323overlap. The following example shows how to make them all visible.
324
325
326Begin_Macro(source)
327../../../tutorials/graphs/multigraph.C
328End_Macro
329
330\anchor MG03
331### Axis' limits setting
332
333The axis limits can be changed the like for TGraph. The same methods apply on
334the multigraph.
335Note the two differents ways to change limits on X and Y axis.
336
337Begin_Macro(source)
338{
339 auto c2 = new TCanvas("c2","c2",600,400);
340
341 TGraph *g[3];
342 Double_t x[10] = {0,1,2,3,4,5,6,7,8,9};
343 Double_t y[10] = {1,2,3,4,5,5,4,3,2,1};
344 auto mg = new TMultiGraph();
345 for (int i=0; i<3; i++) {
346 g[i] = new TGraph(10, x, y);
347 g[i]->SetMarkerStyle(20);
348 g[i]->SetMarkerColor(i+2);
349 for (int j=0; j<10; j++) y[j] = y[j]-1;
350 mg->Add(g[i]);
351 }
352 mg->Draw("APL");
353 mg->GetXaxis()->SetTitle("E_{#gamma} (GeV)");
354 mg->GetYaxis()->SetTitle("Coefficients");
355
356 // Change the axis limits
357 gPad->Modified();
358 mg->GetXaxis()->SetLimits(1.5,7.5);
359 mg->SetMinimum(0.);
360 mg->SetMaximum(10.);
361}
362End_Macro
363*/
364
365
366////////////////////////////////////////////////////////////////////////////////
367/// TMultiGraph default constructor.
368
370{
371 fGraphs = 0;
372 fFunctions = 0;
373 fHistogram = 0;
374 fMaximum = -1111;
375 fMinimum = -1111;
376}
377
378
379////////////////////////////////////////////////////////////////////////////////
380/// Constructor with name and title.
381
382TMultiGraph::TMultiGraph(const char *name, const char *title)
383 : TNamed(name,title)
384{
385 fGraphs = 0;
386 fFunctions = 0;
387 fHistogram = 0;
388 fMaximum = -1111;
389 fMinimum = -1111;
390}
391
392
393////////////////////////////////////////////////////////////////////////////////
394/// Copy constructor.
395
397 TNamed (mg),
398 fGraphs(mg.fGraphs),
399 fFunctions(mg.fFunctions),
400 fHistogram(mg.fHistogram),
401 fMaximum(mg.fMaximum),
402 fMinimum(mg.fMinimum)
403{
404}
405
406
407////////////////////////////////////////////////////////////////////////////////
408/// Assignment operator.
409
411{
412 if (this!=&mg) {
414 fGraphs=mg.fGraphs;
415 fFunctions=mg.fFunctions;
416 fHistogram=mg.fHistogram;
417 fMaximum=mg.fMaximum;
418 fMinimum=mg.fMinimum;
419 }
420 return *this;
421}
422
423
424////////////////////////////////////////////////////////////////////////////////
425/// TMultiGraph destructor.
426
428{
429 if (!fGraphs) return;
430 TGraph *g;
431 TIter next(fGraphs);
432 while ((g = (TGraph*) next())) {
433 g->ResetBit(kMustCleanup);
434 }
435 fGraphs->Delete();
436 delete fGraphs;
437 fGraphs = 0;
438 delete fHistogram;
439 fHistogram = 0;
440 if (fFunctions) {
442 //special logic to support the case where the same object is
443 //added multiple times in fFunctions.
444 //This case happens when the same object is added with different
445 //drawing modes
446 TObject *obj;
447 while ((obj = fFunctions->First())) {
448 while (fFunctions->Remove(obj)) { }
449 delete obj;
450 }
451 delete fFunctions;
452 }
453}
454
455
456////////////////////////////////////////////////////////////////////////////////
457/// Add a new graph to the list of graphs.
458/// Note that the graph is now owned by the TMultigraph.
459/// Deleting the TMultiGraph object will automatically delete the graphs.
460/// You should not delete the graphs when the TMultigraph is still active.
461
463{
464 if (!fGraphs) fGraphs = new TList();
465 graph->SetBit(kMustCleanup);
466 fGraphs->Add(graph,chopt);
467}
468
469
470////////////////////////////////////////////////////////////////////////////////
471/// Add all the graphs in "multigraph" to the list of graphs.
472///
473/// - If "chopt" is defined all the graphs in "multigraph" will be added with
474/// the "chopt" option.
475/// - If "chopt" is undefined each graph will be added with the option it had
476/// in "multigraph".
477
478void TMultiGraph::Add(TMultiGraph *multigraph, Option_t *chopt)
479{
480 TList *graphlist = multigraph->GetListOfGraphs();
481 if (!graphlist) return;
482
483 if (!fGraphs) fGraphs = new TList();
484
485 TObjOptLink *lnk = (TObjOptLink*)graphlist->FirstLink();
486 TObject *obj = 0;
487
488 while (lnk) {
489 obj = lnk->GetObject();
490 if (!strlen(chopt)) fGraphs->Add(obj,lnk->GetOption());
491 else fGraphs->Add(obj,chopt);
492 lnk = (TObjOptLink*)lnk->Next();
493 }
494}
495
496
497////////////////////////////////////////////////////////////////////////////////
498/// Browse multigraph.
499
501{
502 TString opt = gEnv->GetValue("TGraph.BrowseOption", "");
503 if (opt.IsNull()) {
504 opt = b ? b->GetDrawOption() : "alp";
505 opt = (opt == "") ? "alp" : opt.Data();
506 }
507 Draw(opt.Data());
508 gPad->Update();
509}
510
511
512////////////////////////////////////////////////////////////////////////////////
513/// Compute distance from point px,py to each graph.
514
516{
517 // Are we on the axis?
518 const Int_t kMaxDiff = 10;
519 Int_t distance = 9999;
520 if (fHistogram) {
521 distance = fHistogram->DistancetoPrimitive(px,py);
522 if (distance <= 0) return distance;
523 }
524
525 // Loop on the list of graphs
526 if (!fGraphs) return distance;
527 TGraph *g;
528 TIter next(fGraphs);
529 while ((g = (TGraph*) next())) {
530 Int_t dist = g->DistancetoPrimitive(px,py);
531 if (dist <= 0) return 0;
532 if (dist < kMaxDiff) {gPad->SetSelected(g); return dist;}
533 }
534 return distance;
535}
536
537
538////////////////////////////////////////////////////////////////////////////////
539/// Draw this multigraph with its current attributes.
540///
541/// Options to draw a graph are described in TGraphPainter.
542///
543/// The drawing option for each TGraph may be specified as an optional
544/// second argument of the Add function. You can use GetGraphDrawOption
545/// to return this option.
546///
547/// If a draw option is specified, it will be used to draw the graph,
548/// otherwise the graph will be drawn with the option specified in
549/// TMultiGraph::Draw. Use GetDrawOption to return the option specified
550/// when drawing the TMultiGraph.
551
553{
554 TString opt = option;
555 opt.ToLower();
556
557 if (gPad) {
558 if (!gPad->IsEditable()) gROOT->MakeDefCanvas();
559 if (opt.Contains("a")) gPad->Clear();
560 }
561 AppendPad(option);
562}
563
564
565////////////////////////////////////////////////////////////////////////////////
566/// Fit this graph with function with name fname.
567///
568/// interface to TF1::Fit(TF1 *f1...
569
571{
572 char *linear;
573 linear= (char*)strstr(fname, "++");
574 if (linear) {
575 TF1 f1(fname, fname, xmin, xmax);
576 return Fit(&f1,option,"",xmin,xmax);
577 }
578 TF1 * f1 = (TF1*)gROOT->GetFunction(fname);
579 if (!f1) { Printf("Unknown function: %s",fname); return -1; }
580
581 return Fit(f1,option,"",xmin,xmax);
582}
583
584
585////////////////////////////////////////////////////////////////////////////////
586/// Fit this multigraph with function f1.
587///
588/// In this function all graphs of the multigraph are fitted simultaneously
589///
590/// f1 is an already predefined function created by TF1.
591/// Predefined functions such as gaus, expo and poln are automatically
592/// created by ROOT.
593///
594/// The list of fit options is given in parameter `option`which may takes the
595/// following values:
596///
597/// - "W" Ignore all the point errors
598/// - "U" Use a User specified fitting algorithm (via SetFCN)
599/// - "Q" Quiet mode (minimum printing)
600/// - "V" Verbose mode (default is between Q and V)
601/// - "B" Use this option when you want to fix one or more parameters
602/// and the fitting function is like "gaus","expo","poln","landau".
603/// - "R" Use the Range specified in the function range
604/// - "N" Do not store the graphics function, do not draw
605/// - "0" Do not plot the result of the fit. By default the fitted function
606/// is drawn unless the option"N" above is specified.
607/// - "+" Add this new fitted function to the list of fitted functions
608/// (by default, any previous function is deleted)
609/// - "C" In case of linear fitting, not calculate the chisquare (saves time)
610/// - "F" If fitting a polN, switch to minuit fitter
611/// - "ROB" In case of linear fitting, compute the LTS regression
612/// coefficients (robust(resistant) regression), using
613/// the default fraction of good points
614/// - "ROB=0.x" - compute the LTS regression coefficients, using
615/// 0.x as a fraction of good points
616///
617/// When the fit is drawn (by default), the parameter goption may be used
618/// to specify a list of graphics options. See TGraph::Paint for a complete
619/// list of these options.
620///
621/// In order to use the Range option, one must first create a function
622/// with the expression to be fitted. For example, if your graph
623/// has a defined range between -4 and 4 and you want to fit a gaussian
624/// only in the interval 1 to 3, you can do:
625/// ~~~ {.cpp}
626/// TF1 *f1 = new TF1("f1","gaus",1,3);
627/// graph->Fit("f1","R");
628/// ~~~
629///
630/// ### Who is calling this function ?
631///
632/// Note that this function is called when calling TGraphErrors::Fit
633/// or TGraphAsymmErrors::Fit ot TGraphBentErrors::Fit
634/// see the discussion below on the errors calculation.
635///
636/// ### Setting initial conditions
637///
638/// Parameters must be initialized before invoking the Fit function.
639/// The setting of the parameter initial values is automatic for the
640/// predefined functions : poln, expo, gaus, landau. One can however disable
641/// this automatic computation by specifying the option "B".
642/// You can specify boundary limits for some or all parameters via
643/// ~~~ {.cpp}
644/// f1->SetParLimits(p_number, parmin, parmax);
645/// ~~~
646/// if `parmin>=parmax`, the parameter is fixed
647/// Note that you are not forced to fix the limits for all parameters.
648/// For example, if you fit a function with 6 parameters, you can do:
649/// ~~~ {.cpp}
650/// func->SetParameters(0,3.1,1.e-6,0.1,-8,100);
651/// func->SetParLimits(4,-10,-4);
652/// func->SetParLimits(5, 1,1);
653/// ~~~
654/// With this setup, parameters 0->3 can vary freely
655/// Parameter 4 has boundaries [-10,-4] with initial value -8
656/// Parameter 5 is fixed to 100.
657///
658/// ### Fit range
659///
660/// The fit range can be specified in two ways:
661///
662/// - specify rxmax > rxmin (default is rxmin=rxmax=0)
663/// - specify the option "R". In this case, the function will be taken
664/// instead of the full graph range.
665///
666/// ### Changing the fitting function
667///
668/// By default a chi2 fitting function is used for fitting the TGraphs's.
669/// The function is implemented in `FitUtil::EvaluateChi2`.
670/// In case of TGraphErrors an effective chi2 is used
671/// (see TGraphErrors fit in TGraph::Fit) and is implemented in
672/// `FitUtil::EvaluateChi2Effective`
673/// To specify a User defined fitting function, specify option "U" and
674/// call the following function:
675/// ~~~ {.cpp}
676/// TVirtualFitter::Fitter(mygraph)->SetFCN(MyFittingFunction)
677/// ~~~
678/// where MyFittingFunction is of type:
679/// ~~~ {.cpp}
680/// extern void MyFittingFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
681/// ~~~
682///
683/// ### Access to the fit result
684///
685/// The function returns a TFitResultPtr which can hold a pointer to a TFitResult object.
686/// By default the TFitResultPtr contains only the status of the fit and it converts
687/// automatically to an integer. If the option "S" is instead used, TFitResultPtr contains
688/// the TFitResult and behaves as a smart pointer to it. For example one can do:
689/// ~~~ {.cpp}
690/// TFitResultPtr r = graph->Fit("myFunc","S");
691/// TMatrixDSym cov = r->GetCovarianceMatrix(); // to access the covariance matrix
692/// Double_t par0 = r->Parameter(0); // retrieve the value for the parameter 0
693/// Double_t err0 = r->ParError(0); // retrieve the error for the parameter 0
694/// r->Print("V"); // print full information of fit including covariance matrix
695/// r->Write(); // store the result in a file
696/// ~~~
697///
698/// The fit parameters, error and chi2 (but not covariance matrix) can be retrieved also
699/// from the fitted function.
700///
701/// ### Associated functions
702///
703/// One or more object (typically a TF1*) can be added to the list
704/// of functions (fFunctions) associated to each graph.
705/// When TGraph::Fit is invoked, the fitted function is added to this list.
706/// Given a graph gr, one can retrieve an associated function
707/// with:
708/// ~~~ {.cpp}
709/// TF1 *myfunc = gr->GetFunction("myfunc");
710/// ~~~
711///
712/// If the graph is made persistent, the list of
713/// associated functions is also persistent. Given a pointer (see above)
714/// to an associated function myfunc, one can retrieve the function/fit
715/// parameters with calls such as:
716/// ~~~ {.cpp}
717/// Double_t chi2 = myfunc->GetChisquare();
718/// Double_t par0 = myfunc->GetParameter(0); //value of 1st parameter
719/// Double_t err0 = myfunc->GetParError(0); //error on first parameter
720/// ~~~
721///
722/// ### Fit Statistics
723///
724/// You can change the statistics box to display the fit parameters with
725/// the TStyle::SetOptFit(mode) method. This mode has four digits.
726/// mode = pcev (default = 0111)
727///
728/// - v = 1; print name/values of parameters
729/// - e = 1; print errors (if e=1, v must be 1)
730/// - c = 1; print Chisquare/Number of degrees of freedom
731/// - p = 1; print Probability
732///
733/// For example: `gStyle->SetOptFit(1011);`
734/// prints the fit probability, parameter names/values, and errors.
735/// You can change the position of the statistics box with these lines
736/// (where g is a pointer to the TGraph):
737///
738/// ~~~ {.cpp}
739/// Root > TPaveStats *st = (TPaveStats*)g->GetListOfFunctions()->FindObject("stats")
740/// Root > st->SetX1NDC(newx1); //new x start position
741/// Root > st->SetX2NDC(newx2); //new x end position
742/// ~~~
743
745{
746 // internal multigraph fitting methods
747 Foption_t fitOption;
749
750 // create range and minimizer options with default values
751 ROOT::Fit::DataRange range(rxmin,rxmax);
753 return ROOT::Fit::FitObject(this, f1 , fitOption , minOption, goption, range);
754
755}
756
757////////////////////////////////////////////////////////////////////////////////
758/// Display a panel with all histogram fit options.
759/// See class TFitPanel for example
760
762{
763 if (!gPad)
764 gROOT->MakeDefCanvas();
765
766 if (!gPad) {
767 Error("FitPanel", "Unable to create a default canvas");
768 return;
769 }
770
771 // use plugin manager to create instance of TFitEditor
772 TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
773 if (handler && handler->LoadPlugin() != -1) {
774 if (handler->ExecPlugin(2, gPad, this) == 0)
775 Error("FitPanel", "Unable to crate the FitPanel");
776 }
777 else
778 Error("FitPanel", "Unable to find the FitPanel plug-in");
779}
780
781////////////////////////////////////////////////////////////////////////////////
782/// Return the draw option for the TGraph `gr` in this TMultiGraph.
783/// The return option is the one specified when calling TMultiGraph::Add(gr,option).
784
786{
787 if (!fGraphs || !gr) return "";
788 TListIter next(fGraphs);
789 TObject *obj;
790 while ((obj = next())) {
791 if (obj == (TObject*)gr) return next.GetOption();
792 }
793 return "";
794}
795
796
797////////////////////////////////////////////////////////////////////////////////
798/// Compute Initial values of parameters for a gaussian.
799
801{
802 Double_t allcha, sumx, sumx2, x, val, rms, mean;
803 Int_t bin;
804 const Double_t sqrtpi = 2.506628;
805
806 // Compute mean value and RMS of the graph in the given range
807 Int_t np = 0;
808 allcha = sumx = sumx2 = 0;
809 TGraph *g;
810 TIter next(fGraphs);
811 Double_t *px, *py;
812 Int_t npp; //number of points in each graph
813 while ((g = (TGraph*) next())) {
814 px=g->GetX();
815 py=g->GetY();
816 npp=g->GetN();
817 for (bin=0; bin<npp; bin++) {
818 x=px[bin];
819 if (x<xmin || x>xmax) continue;
820 np++;
821 val=py[bin];
822 sumx+=val*x;
823 sumx2+=val*x*x;
824 allcha+=val;
825 }
826 }
827 if (np == 0 || allcha == 0) return;
828 mean = sumx/allcha;
829 rms = TMath::Sqrt(sumx2/allcha - mean*mean);
830
831 Double_t binwidx = TMath::Abs((xmax-xmin)/np);
832 if (rms == 0) rms = 1;
834 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
835 f1->SetParameter(0,binwidx*allcha/(sqrtpi*rms));
836 f1->SetParameter(1,mean);
837 f1->SetParameter(2,rms);
838 f1->SetParLimits(2,0,10*rms);
839}
840
841
842////////////////////////////////////////////////////////////////////////////////
843/// Compute Initial values of parameters for an exponential.
844
846{
847 Double_t constant, slope;
848 Int_t ifail;
849
850 LeastSquareLinearFit(-1, constant, slope, ifail, xmin, xmax);
851
853 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
854 f1->SetParameter(0,constant);
855 f1->SetParameter(1,slope);
856}
857
858
859////////////////////////////////////////////////////////////////////////////////
860/// Compute Initial values of parameters for a polynom.
861
863{
864 Double_t fitpar[25];
865
867 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
868 Int_t npar = f1->GetNpar();
869
870 LeastSquareFit(npar, fitpar, xmin, xmax);
871
872 for (Int_t i=0;i<npar;i++) f1->SetParameter(i, fitpar[i]);
873}
874
875
876////////////////////////////////////////////////////////////////////////////////
877/// Least squares lpolynomial fitting without weights.
878///
879/// - m number of parameters
880/// - a array of parameters
881/// - first 1st point number to fit (default =0)
882/// - last last point number to fit (default=fNpoints-1)
883///
884/// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
885
887{
888 const Double_t zero = 0.;
889 const Double_t one = 1.;
890 const Int_t idim = 20;
891
892 Double_t b[400] /* was [20][20] */;
893 Int_t i, k, l, ifail, bin;
894 Double_t power;
895 Double_t da[20], xk, yk;
896
897
898 //count the total number of points to fit
899 TGraph *g;
900 TIter next(fGraphs);
901 Double_t *px, *py;
902 Int_t n=0;
903 Int_t npp;
904 while ((g = (TGraph*) next())) {
905 px=g->GetX();
906 npp=g->GetN();
907 for (bin=0; bin<npp; bin++) {
908 xk=px[bin];
909 if (xk < xmin || xk > xmax) continue;
910 n++;
911 }
912 }
913 if (m <= 2) {
914 LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
915 return;
916 }
917 if (m > idim || m > n) return;
918 da[0] = zero;
919 for (l = 2; l <= m; ++l) {
920 b[l-1] = zero;
921 b[m + l*20 - 21] = zero;
922 da[l-1] = zero;
923 }
924 Int_t np = 0;
925
926 next.Reset();
927 while ((g = (TGraph*) next())) {
928 px=g->GetX();
929 py=g->GetY();
930 npp=g->GetN();
931
932 for (k = 0; k <= npp; ++k) {
933 xk = px[k];
934 if (xk < xmin || xk > xmax) continue;
935 np++;
936 yk = py[k];
937 power = one;
938 da[0] += yk;
939 for (l = 2; l <= m; ++l) {
940 power *= xk;
941 b[l-1] += power;
942 da[l-1] += power*yk;
943 }
944 for (l = 2; l <= m; ++l) {
945 power *= xk;
946 b[m + l*20 - 21] += power;
947 }
948 }
949 }
950 b[0] = Double_t(np);
951 for (i = 3; i <= m; ++i) {
952 for (k = i; k <= m; ++k) {
953 b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
954 }
955 }
956 H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
957
958 if (ifail < 0) {
959 //a[0] = fY[0];
960 py=((TGraph *)fGraphs->First())->GetY();
961 a[0]=py[0];
962 for (i=1; i<m; ++i) a[i] = 0;
963 return;
964 }
965 for (i=0; i<m; ++i) a[i] = da[i];
966}
967
968
969////////////////////////////////////////////////////////////////////////////////
970/// Least square linear fit without weights.
971///
972/// Fit a straight line (a0 + a1*x) to the data in this graph.
973///
974/// - ndata: number of points to fit
975/// - first: first point number to fit
976/// - last: last point to fit O(ndata should be last-first
977/// - ifail: return parameter indicating the status of the fit (ifail=0, fit is OK)
978///
979/// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
980
983{
984 Double_t xbar, ybar, x2bar;
985 Int_t i;
986 Double_t xybar;
987 Double_t fn, xk, yk;
988 Double_t det;
989
990 ifail = -2;
991 xbar = ybar = x2bar = xybar = 0;
992 Int_t np = 0;
993 TGraph *g;
994 TIter next(fGraphs);
995 Double_t *px, *py;
996 Int_t npp;
997 while ((g = (TGraph*) next())) {
998 px=g->GetX();
999 py=g->GetY();
1000 npp=g->GetN();
1001 for (i = 0; i < npp; ++i) {
1002 xk = px[i];
1003 if (xk < xmin || xk > xmax) continue;
1004 np++;
1005 yk = py[i];
1006 if (ndata < 0) {
1007 if (yk <= 0) yk = 1e-9;
1008 yk = TMath::Log(yk);
1009 }
1010 xbar += xk;
1011 ybar += yk;
1012 x2bar += xk*xk;
1013 xybar += xk*yk;
1014 }
1015 }
1016 fn = Double_t(np);
1017 det = fn*x2bar - xbar*xbar;
1018 ifail = -1;
1019 if (det <= 0) {
1020 if (fn > 0) a0 = ybar/fn;
1021 else a0 = 0;
1022 a1 = 0;
1023 return;
1024 }
1025 ifail = 0;
1026 a0 = (x2bar*ybar - xbar*xybar) / det;
1027 a1 = (fn*xybar - xbar*ybar) / det;
1028}
1029
1030
1031////////////////////////////////////////////////////////////////////////////////
1032/// Return 1 if the point (x,y) is inside one of the graphs 0 otherwise.
1033
1035{
1036 Int_t in = 0;
1037 if (!fGraphs) return in;
1038 TGraph *g;
1039 TIter next(fGraphs);
1040 while ((g = (TGraph*) next())) {
1041 in = g->IsInside(x, y);
1042 if (in) return in;
1043 }
1044 return in;
1045}
1046
1047
1048////////////////////////////////////////////////////////////////////////////////
1049/// Returns a pointer to the histogram used to draw the axis.
1050/// Takes into account following cases.
1051///
1052/// 1. if `fHistogram` exists it is returned
1053/// 2. if `fHistogram` doesn't exists and `gPad` exists `gPad` is updated. That
1054/// may trigger the creation of `fHistogram`. If `fHistogram` still does not
1055/// exit but `hframe` does (if user called `TPad::DrawFrame`) the pointer to
1056/// `hframe` histogram is returned
1057/// 3. after the two previous steps, if `fHistogram` still doesn't exist, then
1058/// it is created.
1059
1061{
1062 if (fHistogram) return fHistogram;
1063
1064 if (gPad) {
1065 gPad->Modified();
1066 gPad->Update();
1067 if (fHistogram) return fHistogram;
1068 TH1F *h1 = (TH1F*)gPad->FindObject("hframe");
1069 if (h1) return h1;
1070 }
1071
1072 Bool_t initialrangeset = kFALSE;
1073 Double_t rwxmin = 0.,rwxmax = 0.,rwymin = 0.,rwymax = 0.;
1074 TGraph *g;
1075 Int_t npt = 100 ;
1076 TIter next(fGraphs);
1077 while ((g = (TGraph*) next())) {
1078 if (g->GetN() <= 0) continue;
1079 if (initialrangeset) {
1080 Double_t rx1,ry1,rx2,ry2;
1081 g->ComputeRange(rx1, ry1, rx2, ry2);
1082 if (rx1 < rwxmin) rwxmin = rx1;
1083 if (ry1 < rwymin) rwymin = ry1;
1084 if (rx2 > rwxmax) rwxmax = rx2;
1085 if (ry2 > rwymax) rwymax = ry2;
1086 } else {
1087 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1088 initialrangeset = kTRUE;
1089 }
1090 if (g->GetN() > npt) npt = g->GetN();
1091 }
1092 if (rwxmin == rwxmax) rwxmax += 1.;
1093 if (rwymin == rwymax) rwymax += 1.;
1094 double dx = 0.05*(rwxmax-rwxmin);
1095 double dy = 0.05*(rwymax-rwymin);
1096 rwxmin = rwxmin - dx;
1097 rwxmax = rwxmax + dx;
1098 rwymin = rwymin - dy;
1099 rwymax = rwymax + dy;
1100 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1101 if (!fHistogram) return 0;
1102 fHistogram->SetMinimum(rwymin);
1104 fHistogram->SetMaximum(rwymax);
1105 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1107 return fHistogram;
1108}
1109
1110
1111////////////////////////////////////////////////////////////////////////////////
1112/// Return pointer to function with name.
1113///
1114/// Functions such as TGraph::Fit store the fitted function in the list of
1115/// functions of this graph.
1116
1118{
1119 if (!fFunctions) return 0;
1120 return (TF1*)fFunctions->FindObject(name);
1121}
1122
1123////////////////////////////////////////////////////////////////////////////////
1124/// Return pointer to list of functions.
1125/// If pointer is null create the list
1126
1128{
1129 if (!fFunctions) fFunctions = new TList();
1130 return fFunctions;
1131}
1132
1133
1134////////////////////////////////////////////////////////////////////////////////
1135/// Get x axis of the graph.
1136/// This method returns a valid axis only after the TMultigraph has been drawn.
1137
1139{
1140 TH1 *h = GetHistogram();
1141 if (!h) return 0;
1142 return h->GetXaxis();
1143}
1144
1145
1146////////////////////////////////////////////////////////////////////////////////
1147/// Get y axis of the graph.
1148/// This method returns a valid axis only after the TMultigraph has been drawn.
1149
1151{
1152 TH1 *h = GetHistogram();
1153 if (!h) return 0;
1154 return h->GetYaxis();
1155}
1156
1157
1158////////////////////////////////////////////////////////////////////////////////
1159/// Paint all the graphs of this multigraph.
1160
1162{
1163 const TPickerStackGuard pushGuard(this);
1164
1165 if (!fGraphs) return;
1166 if (fGraphs->GetSize() == 0) return;
1167
1168 char option[128];
1169 strlcpy(option,choptin,128);
1170 Int_t nch = strlen(choptin);
1171 for (Int_t i=0;i<nch;i++) option[i] = toupper(option[i]);
1172
1173 // Automatic color
1174 char *l1 = strstr(option,"PFC"); // Automatic Fill Color
1175 char *l2 = strstr(option,"PLC"); // Automatic Line Color
1176 char *l3 = strstr(option,"PMC"); // Automatic Marker Color
1177 if (l1 || l2 || l3) {
1178 TString opt1 = option; opt1.ToLower();
1179 if (l1) memcpy(l1," ",3);
1180 if (l2) memcpy(l2," ",3);
1181 if (l3) memcpy(l3," ",3);
1183 TGraph* gAti;
1184 Int_t ngraphs = fGraphs->GetSize();
1185 Int_t ic;
1186 gPad->IncrementPaletteColor(ngraphs, opt1);
1187 for (Int_t i=0;i<ngraphs;i++) {
1188 ic = gPad->NextPaletteColor();
1189 gAti = (TGraph*)(fGraphs->At(i));
1190 if (l1) gAti->SetFillColor(ic);
1191 if (l2) gAti->SetLineColor(ic);
1192 if (l3) gAti->SetMarkerColor(ic);
1193 lnk = (TObjOptLink*)lnk->Next();
1194 }
1195 }
1196
1197 char *l;
1198
1199 TString chopt = option;
1200
1201 l = (char*)strstr(chopt.Data(),"3D");
1202 if (l) {
1203 l = (char*)strstr(chopt.Data(),"L");
1204 if (l) PaintPolyLine3D(chopt.Data());
1205 return;
1206 }
1207
1208 l = (char*)strstr(chopt.Data(),"PADS");
1209 if (l) {
1210 chopt.ReplaceAll("PADS","");
1211 PaintPads(chopt.Data());
1212 return;
1213 }
1214
1215 char *lrx = (char *)strstr(chopt.Data(), "RX"); // Reverse graphs along X axis
1216 char *lry = (char *)strstr(chopt.Data(), "RY"); // Reverse graphs along Y axis
1217 if (lrx || lry) {
1218 PaintReverse(chopt.Data());
1219 return;
1220 }
1221
1222 TGraph *g;
1223
1224 l = (char*)strstr(chopt.Data(),"A");
1225 if (l) {
1226 *l = ' ';
1227 TIter next(fGraphs);
1228 Int_t npt = 100;
1229 Double_t maximum, minimum, rwxmin, rwxmax, rwymin, rwymax, uxmin, uxmax, dx, dy;
1230 rwxmin = gPad->GetUxmin();
1231 rwxmax = gPad->GetUxmax();
1232 rwymin = gPad->GetUymin();
1233 rwymax = gPad->GetUymax();
1234 char *xtitle = 0;
1235 char *ytitle = 0;
1236 Int_t firstx = 0;
1237 Int_t lastx = 0;
1238 Bool_t timedisplay = kFALSE;
1239 char *timeformat = 0;
1240
1241 if (fHistogram) {
1242 //cleanup in case of a previous unzoom and in case one of the TGraph has changed
1244 TGraph* gAti;
1245 Int_t ngraphs = fGraphs->GetSize();
1246 Bool_t reset_hist = kFALSE;
1247 for (Int_t i=0;i<ngraphs;i++) {
1248 gAti = (TGraph*)(fGraphs->At(i));
1249 if(gAti->TestBit(TGraph::kResetHisto)) {reset_hist = kTRUE; break;}
1250 lnk = (TObjOptLink*)lnk->Next();
1251 }
1252 if (fHistogram->GetMinimum() >= fHistogram->GetMaximum() || reset_hist) {
1253 nch = strlen(fHistogram->GetXaxis()->GetTitle());
1254 firstx = fHistogram->GetXaxis()->GetFirst();
1255 lastx = fHistogram->GetXaxis()->GetLast();
1256 timedisplay = fHistogram->GetXaxis()->GetTimeDisplay();
1257 if (nch) {
1258 xtitle = new char[nch+1];
1259 strlcpy(xtitle,fHistogram->GetXaxis()->GetTitle(),nch+1);
1260 }
1261 nch = strlen(fHistogram->GetYaxis()->GetTitle());
1262 if (nch) {
1263 ytitle = new char[nch+1];
1264 strlcpy(ytitle,fHistogram->GetYaxis()->GetTitle(),nch+1);
1265 }
1266 nch = strlen(fHistogram->GetXaxis()->GetTimeFormat());
1267 if (nch) {
1268 timeformat = new char[nch+1];
1269 strlcpy(timeformat,fHistogram->GetXaxis()->GetTimeFormat(),nch+1);
1270 }
1271 delete fHistogram;
1272 fHistogram = 0;
1273 }
1274 }
1275 if (fHistogram) {
1276 minimum = fHistogram->GetYaxis()->GetXmin();
1277 maximum = fHistogram->GetYaxis()->GetXmax();
1278 uxmin = gPad->PadtoX(rwxmin);
1279 uxmax = gPad->PadtoX(rwxmax);
1280 } else {
1281 Bool_t initialrangeset = kFALSE;
1282 while ((g = (TGraph*) next())) {
1283 if (g->GetN() <= 0) continue;
1284 if (initialrangeset) {
1285 Double_t rx1,ry1,rx2,ry2;
1286 g->ComputeRange(rx1, ry1, rx2, ry2);
1287 if (rx1 < rwxmin) rwxmin = rx1;
1288 if (ry1 < rwymin) rwymin = ry1;
1289 if (rx2 > rwxmax) rwxmax = rx2;
1290 if (ry2 > rwymax) rwymax = ry2;
1291 } else {
1292 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1293 initialrangeset = kTRUE;
1294 }
1295 if (g->GetN() > npt) npt = g->GetN();
1296 }
1297 if (rwxmin == rwxmax) rwxmax += 1.;
1298 if (rwymin == rwymax) rwymax += 1.;
1299 dx = 0.05*(rwxmax-rwxmin);
1300 dy = 0.05*(rwymax-rwymin);
1301 uxmin = rwxmin - dx;
1302 uxmax = rwxmax + dx;
1303 if (gPad->GetLogy()) {
1304 if (rwymin <= 0) rwymin = 0.001*rwymax;
1305 minimum = rwymin/(1+0.5*TMath::Log10(rwymax/rwymin));
1306 maximum = rwymax*(1+0.2*TMath::Log10(rwymax/rwymin));
1307 } else {
1308 minimum = rwymin - dy;
1309 maximum = rwymax + dy;
1310 }
1311 if (minimum < 0 && rwymin >= 0) minimum = 0;
1312 if (maximum > 0 && rwymax <= 0) maximum = 0;
1313 }
1314
1315 if (fMinimum != -1111) rwymin = minimum = fMinimum;
1316 if (fMaximum != -1111) rwymax = maximum = fMaximum;
1317 if (uxmin < 0 && rwxmin >= 0) {
1318 if (gPad->GetLogx()) uxmin = 0.9*rwxmin;
1319 //else uxmin = 0;
1320 }
1321 if (uxmax > 0 && rwxmax <= 0) {
1322 if (gPad->GetLogx()) uxmax = 1.1*rwxmax;
1323 //else uxmax = 0;
1324 }
1325 if (minimum < 0 && rwymin >= 0) {
1326 if (gPad->GetLogy()) minimum = 0.9*rwymin;
1327 //else minimum = 0;
1328 }
1329 if (maximum > 0 && rwymax <= 0) {
1330 if (gPad->GetLogy()) maximum = 1.1*rwymax;
1331 //else maximum = 0;
1332 }
1333 if (minimum <= 0 && gPad->GetLogy()) minimum = 0.001*maximum;
1334 if (uxmin <= 0 && gPad->GetLogx()) {
1335 if (uxmax > 1000) uxmin = 1;
1336 else uxmin = 0.001*uxmax;
1337 }
1338 rwymin = minimum;
1339 rwymax = maximum;
1340 if (fHistogram) {
1341 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1342 }
1343
1344 // Create a temporary histogram to draw the axis
1345 if (!fHistogram) {
1346 // the graph is created with at least as many channels as there are points
1347 // to permit zooming on the full range
1348 rwxmin = uxmin;
1349 rwxmax = uxmax;
1350 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1351 if (!fHistogram) return;
1352 fHistogram->SetMinimum(rwymin);
1354 fHistogram->SetMaximum(rwymax);
1355 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1357 if (xtitle) {fHistogram->GetXaxis()->SetTitle(xtitle); delete [] xtitle;}
1358 if (ytitle) {fHistogram->GetYaxis()->SetTitle(ytitle); delete [] ytitle;}
1359 if (firstx != lastx) fHistogram->GetXaxis()->SetRange(firstx,lastx);
1360 if (timedisplay) {fHistogram->GetXaxis()->SetTimeDisplay(timedisplay);}
1361 if (timeformat) {fHistogram->GetXaxis()->SetTimeFormat(timeformat); delete [] timeformat;}
1362 }
1363 fHistogram->Paint("0");
1364 }
1365
1366 TGraph *gfit = nullptr;
1367 if (fGraphs) {
1369 TObject *obj = 0;
1370
1371 chopt.ReplaceAll("A","");
1372
1373 while (lnk) {
1374
1375 obj = lnk->GetObject();
1376
1377 gPad->PushSelectableObject(obj);
1378
1379 if (!gPad->PadInHighlightMode() || (gPad->PadInHighlightMode() && obj == gPad->GetSelected())) {
1380 TString opt = lnk->GetOption();
1381 if (!opt.IsWhitespace())
1382 obj->Paint(opt.ReplaceAll("A","").Data());
1383 else {
1384 if (!chopt.IsWhitespace()) obj->Paint(chopt.Data());
1385 else obj->Paint("L");
1386 }
1387 }
1388
1389 lnk = (TObjOptLink*)lnk->Next();
1390 }
1391
1392 gfit = (TGraph*)obj; // pick one TGraph in the list to paint the fit parameters.
1393 }
1394
1395 TObject *f;
1396 TF1 *fit = nullptr;
1397 if (fFunctions) {
1398 TIter next(fFunctions);
1399 while ((f = (TObject*) next())) {
1400 if (f->InheritsFrom(TF1::Class())) {
1401 if (f->TestBit(TF1::kNotDraw) == 0) f->Paint("lsame");
1402 fit = (TF1*)f;
1403 } else {
1404 f->Paint();
1405 }
1406 }
1407 }
1408
1409 if (gfit && fit) gfit->PaintStats(fit);
1410}
1411
1412
1413////////////////////////////////////////////////////////////////////////////////
1414/// Divides the active pad and draws all Graphs in the Multigraph separately.
1415
1417{
1418 TIter next(fGraphs);
1419 Int_t neededPads = fGraphs->GetSize();
1420 Int_t existingPads = 0;
1421 TString opt = (TString)option;
1422
1423 TVirtualPad *curPad = gPad;
1424 TObject *obj;
1425 TIter nextPad(curPad->GetListOfPrimitives());
1426
1427 while ((obj = nextPad())) {
1428 if (obj->InheritsFrom(TVirtualPad::Class())) existingPads++;
1429 }
1430 if (existingPads < neededPads) {
1431 curPad->Clear();
1432 Int_t nx = (Int_t)TMath::Sqrt((Double_t)neededPads);
1433 if (nx*nx < neededPads) nx++;
1434 Int_t ny = nx;
1435 if (((nx*ny)-nx) >= neededPads) ny--;
1436 curPad->Divide(nx,ny);
1437 }
1438 Int_t i = 0;
1439 TGraph *g;
1440
1442 obj = 0;
1443
1444 while (lnk) {
1445 g = (TGraph*)lnk->GetObject();
1446 i++;
1447 curPad->cd(i);
1448 TString apopt = lnk->GetOption();
1449 if (strlen(apopt)) {
1450 g->Draw((apopt.Append("A")).Data());
1451 } else {
1452 if (strlen(opt)) g->Draw(opt.Append("A"));
1453 else g->Draw("LA");
1454 }
1455 lnk = (TObjOptLink*)lnk->Next();
1456 }
1457
1458 curPad->cd();
1459}
1460
1461
1462////////////////////////////////////////////////////////////////////////////////
1463/// Paint all the graphs of this multigraph as 3D lines.
1464
1466{
1467 Int_t i, npt=0;
1468 char *l;
1469 Double_t rwxmin=0., rwxmax=0., rwymin=0., rwymax=0.;
1470 TIter next(fGraphs);
1471 TGraph *g;
1472
1473 g = (TGraph*) next();
1474 if (g) {
1475 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1476 npt = g->GetN();
1477 }
1478
1479 if (!fHistogram) {
1480 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1481 }
1482
1483 while ((g = (TGraph*) next())) {
1484 Double_t rx1,ry1,rx2,ry2;
1485 g->ComputeRange(rx1, ry1, rx2, ry2);
1486 if (rx1 < rwxmin) rwxmin = rx1;
1487 if (ry1 < rwymin) rwymin = ry1;
1488 if (rx2 > rwxmax) rwxmax = rx2;
1489 if (ry2 > rwymax) rwymax = ry2;
1490 if (g->GetN() > npt) npt = g->GetN();
1491 }
1492
1493 Int_t ndiv = fGraphs->GetSize();
1494
1495 TH2F* frame = new TH2F("frame","", ndiv, 0., (Double_t)(ndiv), npt, rwxmin, rwxmax);
1496 if (fHistogram) {
1497 frame->SetTitle(fHistogram->GetTitle());
1498 frame->GetYaxis()->SetTitle(fHistogram->GetXaxis()->GetTitle());
1500 frame->GetZaxis()->SetTitle(fHistogram->GetYaxis()->GetTitle());
1501 }
1502
1503 TAxis *Xaxis = frame->GetXaxis();
1504 Xaxis->SetNdivisions(-ndiv);
1505 next.Reset();
1506 for (i=ndiv; i>=1; i--) {
1507 g = (TGraph*) next();
1508 Xaxis->SetBinLabel(i, g->GetTitle());
1509 }
1510
1511 frame->SetStats(kFALSE);
1512 if (fMinimum != -1111) frame->SetMinimum(fMinimum);
1513 else frame->SetMinimum(rwymin);
1514 if (fMaximum != -1111) frame->SetMaximum(fMaximum);
1515 else frame->SetMaximum(rwymax);
1516
1517 l = (char*)strstr(option,"A");
1518 if (l) frame->Paint("lego9,fb,bb");
1519 l = (char*)strstr(option,"BB");
1520 if (!l) frame->Paint("lego9,fb,a,same");
1521
1522 Double_t *x, *y;
1523 Double_t xyz1[3], xyz2[3];
1524
1525 Double_t xl = frame->GetYaxis()->GetBinLowEdge(frame->GetYaxis()->GetFirst());
1526 Double_t xu = frame->GetYaxis()->GetBinUpEdge(frame->GetYaxis()->GetLast());
1527 Double_t yl = frame->GetMinimum();
1528 Double_t yu = frame->GetMaximum();
1529 Double_t xc[2],yc[2];
1530 next.Reset();
1531 Int_t j = ndiv;
1532
1533 while ((g = (TGraph*) next())) {
1534 npt = g->GetN();
1535 x = g->GetX();
1536 y = g->GetY();
1537 gPad->SetLineColor(g->GetLineColor());
1538 gPad->SetLineWidth(g->GetLineWidth());
1539 gPad->SetLineStyle(g->GetLineStyle());
1540 gPad->TAttLine::Modify();
1541 for (i=0; i<npt-1; i++) {
1542 xc[0] = x[i];
1543 xc[1] = x[i+1];
1544 yc[0] = y[i];
1545 yc[1] = y[i+1];
1546 if (gPad->Clip(&xc[0], &yc[0], xl, yl, xu, yu)<2) {
1547 xyz1[0] = j-0.5;
1548 xyz1[1] = xc[0];
1549 xyz1[2] = yc[0];
1550 xyz2[0] = j-0.5;
1551 xyz2[1] = xc[1];
1552 xyz2[2] = yc[1];
1553 gPad->PaintLine3D(xyz1, xyz2);
1554 }
1555 }
1556 j--;
1557 }
1558
1559 l = (char*)strstr(option,"FB");
1560 if (!l) frame->Paint("lego9,bb,a,same");
1561 delete frame;
1562}
1563
1564
1565////////////////////////////////////////////////////////////////////////////////
1566/// Paint all the graphs of this multigraph reverting values along X and/or Y axis.
1567/// New graphs are created.
1568
1570{
1571 auto *h = GetHistogram();
1572 TH1F *hg = nullptr;
1573 TGraph *fg = nullptr;
1574 if (!h)
1575 return;
1576 TString mgopt = option;
1577 mgopt.ToLower();
1578
1579 TIter next(fGraphs);
1580 TGraph *g;
1581 Bool_t first = kTRUE;
1582 TString gopt;
1583 while ((g = (TGraph *)next())) {
1584 gopt = GetGraphDrawOption(g);
1585 gopt.Append(mgopt);
1586 if (first) {
1587 fg = g;
1588 hg = fg->GetHistogram();
1589 fg->SetHistogram(h);
1590 fg->Paint(gopt.Data());
1591 first = kFALSE;
1592 } else {
1593 g->Paint(gopt.ReplaceAll("a", "").Data());
1594 }
1595 }
1596 if (fg)
1597 fg->SetHistogram(hg);
1598}
1599
1600
1601////////////////////////////////////////////////////////////////////////////////
1602/// Print the list of graphs.
1603
1604void TMultiGraph::Print(Option_t *option) const
1605{
1606 TGraph *g;
1607 if (fGraphs) {
1608 TIter next(fGraphs);
1609 while ((g = (TGraph*) next())) {
1610 g->Print(option);
1611 }
1612 }
1613}
1614
1615
1616////////////////////////////////////////////////////////////////////////////////
1617/// Recursively remove this object from a list. Typically implemented
1618/// by classes that can contain multiple references to a same object.
1619
1621{
1622 if (!fGraphs) return;
1623 TObject *objr = fGraphs->Remove(obj);
1624 if (!objr) return;
1625 delete fHistogram; fHistogram = 0;
1626 if (gPad) gPad->Modified();
1627}
1628
1629
1630////////////////////////////////////////////////////////////////////////////////
1631/// Save primitive as a C++ statement(s) on output stream out.
1632
1633void TMultiGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
1634{
1635 char quote = '"';
1636 out<<" "<<std::endl;
1637 if (gROOT->ClassSaved(TMultiGraph::Class())) {
1638 out<<" ";
1639 } else {
1640 out<<" TMultiGraph *";
1641 }
1642 out<<"multigraph = new TMultiGraph();"<<std::endl;
1643 out<<" multigraph->SetName("<<quote<<GetName()<<quote<<");"<<std::endl;
1644 out<<" multigraph->SetTitle("<<quote<<GetTitle()<<quote<<");"<<std::endl;
1645
1646 if (fGraphs) {
1648 TObject *g;
1649
1650 while (lnk) {
1651 g = lnk->GetObject();
1652 g->SavePrimitive(out, Form("multigraph%s",lnk->GetOption()));
1653 lnk = (TObjOptLink*)lnk->Next();
1654 }
1655 }
1656 const char *l = strstr(option,"th2poly");
1657 if (l) {
1658 out<<" "<<l+7<<"->AddBin(multigraph);"<<std::endl;
1659 } else {
1660 out<<" multigraph->Draw(" <<quote<<option<<quote<<");"<<std::endl;
1661 }
1662 TAxis *xaxis = GetXaxis();
1663 TAxis *yaxis = GetYaxis();
1664
1665 if (xaxis) {
1666 out<<" multigraph->GetXaxis()->SetLimits("<<xaxis->GetXmin()<<", "<<xaxis->GetXmax()<<");"<<std::endl;
1667 xaxis->SaveAttributes(out, "multigraph","->GetXaxis()");
1668 }
1669 if (yaxis) yaxis->SaveAttributes(out, "multigraph","->GetYaxis()");
1670 if (fMinimum != -1111) out<<" multigraph->SetMinimum("<<fMinimum<<");"<<std::endl;
1671 if (fMaximum != -1111) out<<" multigraph->SetMaximum("<<fMaximum<<");"<<std::endl;
1672}
1673
1674
1675////////////////////////////////////////////////////////////////////////////////
1676/// Set multigraph maximum.
1677
1679{
1680 fMaximum = maximum;
1681 if (fHistogram) fHistogram->SetMaximum(maximum);
1682}
1683
1684
1685////////////////////////////////////////////////////////////////////////////////
1686/// Set multigraph minimum.
1687
1689{
1690 fMinimum = minimum;
1691 if (fHistogram) fHistogram->SetMinimum(minimum);
1692}
1693
1694
1695////////////////////////////////////////////////////////////////////////////////
1696/// Get iterator over internal graphs list.
1697
1699{
1700 return TIter(fGraphs);
1701}
void Class()
Definition: Class.C:29
#define b(i)
Definition: RSha256.hxx:100
#define f(i)
Definition: RSha256.hxx:104
#define g(i)
Definition: RSha256.hxx:105
#define h(i)
Definition: RSha256.hxx:106
#define e(i)
Definition: RSha256.hxx:103
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
double Axis_t
Definition: RtypesCore.h:85
double Double_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:364
R__EXTERN TEnv * gEnv
Definition: TEnv.h:170
char name[80]
Definition: TGX11.cxx:110
float xmin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
void H1LeastSquareSeqnd(Int_t n, Double_t *a, Int_t idim, Int_t &ifail, Int_t k, Double_t *b)
Extracted from CERN Program library routine DSEQN.
Definition: TH1.cxx:4805
#define gROOT
Definition: TROOT.h:404
char * Form(const char *fmt,...)
void Printf(const char *fmt,...)
#define gPad
Definition: TVirtualPad.h:287
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition: DataRange.h:35
virtual void SetNdivisions(Int_t n=510, Bool_t optim=kTRUE)
Set the number of divisions for this axis.
Definition: TAttAxis.cxx:237
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
Class to manage histogram axis.
Definition: TAxis.h:30
virtual void SetBinLabel(Int_t bin, const char *label)
Set label for bin.
Definition: TAxis.cxx:823
virtual Bool_t GetTimeDisplay() const
Definition: TAxis.h:126
Double_t GetXmax() const
Definition: TAxis.h:134
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:661
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:518
virtual void SetTimeDisplay(Int_t value)
Definition: TAxis.h:161
Int_t GetLast() const
Return last bin on the axis i.e.
Definition: TAxis.cxx:469
virtual void SetLimits(Double_t xmin, Double_t xmax)
Definition: TAxis.h:154
Double_t GetXmin() const
Definition: TAxis.h:133
virtual const char * GetTimeFormat() const
Definition: TAxis.h:127
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:129
virtual void SetTimeFormat(const char *format="")
Change the format used for time plotting.
Definition: TAxis.cxx:1020
virtual void SetRange(Int_t first=0, Int_t last=0)
Set the viewing range for the axis using bin numbers.
Definition: TAxis.cxx:920
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition: TAxis.cxx:528
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition: TAxis.cxx:458
Using a TBrowser one can browse all ROOT objects.
Definition: TBrowser.h:37
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
Definition: TCollection.h:184
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition: TEnv.cxx:491
1-Dim function class
Definition: TF1.h:213
virtual Int_t GetNpar() const
Definition: TF1.h:481
@ kNotDraw
Definition: TF1.h:326
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set limits for parameter ipar.
Definition: TF1.cxx:3518
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:634
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:32
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
@ kResetHisto
fHistogram must be reset in GetHistogram
Definition: TGraph.h:72
virtual void Paint(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:2038
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:2065
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1485
virtual void SetHistogram(TH1F *h)
Definition: TGraph.h:178
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:575
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
virtual void SetDirectory(TDirectory *dir)
By default, when a histogram is created, it is added to the list of histogram objects in the current ...
Definition: TH1.cxx:8780
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition: TH1.cxx:6680
TAxis * GetZaxis()
Definition: TH1.h:322
@ kNoStats
Don't draw stats box.
Definition: TH1.h:164
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:320
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:8388
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:398
TAxis * GetYaxis()
Definition: TH1.h:321
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:399
virtual void Paint(Option_t *option="")
Control routine to paint any kind of histograms.
Definition: TH1.cxx:6157
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:8478
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a line.
Definition: TH1.cxx:2812
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8833
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:251
void Reset()
Definition: TCollection.h:254
Iterator of linked list.
Definition: TList.h:200
Option_t * GetOption() const
Returns the object option stored in the list.
Definition: TList.cxx:1142
A doubly linked list.
Definition: TList.h:44
virtual void Add(TObject *obj)
Definition: TList.h:87
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:822
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:578
virtual TObjLink * FirstLink() const
Definition: TList.h:108
virtual TObject * At(Int_t idx) const
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:357
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
virtual TObject * First() const
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:659
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:36
Double_t fMinimum
Minimum value for plotting along y.
Definition: TMultiGraph.h:43
TH1F * fHistogram
Pointer to histogram used for drawing axis.
Definition: TMultiGraph.h:41
TMultiGraph()
TMultiGraph default constructor.
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.
TList * GetListOfGraphs() const
Definition: TMultiGraph.h:70
TList * fGraphs
Pointer to list of TGraphs.
Definition: TMultiGraph.h:39
virtual void Add(TGraph *graph, Option_t *chopt="")
Add a new graph to the list of graphs.
Double_t fMaximum
Maximum value for plotting along y.
Definition: TMultiGraph.h:42
virtual void SetMinimum(Double_t minimum=-1111)
Set multigraph minimum.
TH1F * GetHistogram()
Returns a pointer to the histogram used to draw the axis.
TF1 * GetFunction(const char *name) const
Return pointer to function with name.
virtual void LeastSquareFit(Int_t m, Double_t *a, Double_t xmin, Double_t xmax)
Least squares lpolynomial fitting without weights.
virtual void InitPolynom(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for a polynom.
void PaintPolyLine3D(Option_t *chopt="")
Paint all the graphs of this multigraph as 3D lines.
virtual void InitExpo(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for an exponential.
virtual void FitPanel()
Display a panel with all histogram fit options.
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to each graph.
virtual void Paint(Option_t *chopt="")
Paint all the graphs of this multigraph.
virtual void Draw(Option_t *chopt="")
Draw this multigraph with its current attributes.
TIter begin() const
Get iterator over internal graphs list.
TMultiGraph & operator=(const TMultiGraph &)
Assignment operator.
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.
void PaintPads(Option_t *chopt="")
Divides the active pad and draws all Graphs in the Multigraph separately.
virtual Option_t * GetGraphDrawOption(const TGraph *gr) const
Return the draw option for the TGraph gr in this TMultiGraph.
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive as a C++ statement(s) on output stream out.
TAxis * GetYaxis()
Get y axis of the graph.
virtual void InitGaus(Double_t xmin, Double_t xmax)
Compute Initial values of parameters for a gaussian.
virtual ~TMultiGraph()
TMultiGraph destructor.
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.
TList * GetListOfFunctions()
Return pointer to list of functions.
virtual void Print(Option_t *chopt="") const
Print the list of graphs.
virtual void SetMaximum(Double_t maximum=-1111)
Set multigraph maximum.
virtual void Browse(TBrowser *b)
Browse multigraph.
void PaintReverse(Option_t *chopt="")
Paint all the graphs of this multigraph reverting values along X and/or Y axis.
TList * fFunctions
Pointer to list of functions (fits and user)
Definition: TMultiGraph.h:40
TAxis * GetXaxis()
Get x axis of the graph.
virtual void RecursiveRemove(TObject *obj)
Recursively remove this object from a list.
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Mother of all ROOT objects.
Definition: TObject.h:37
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:107
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
virtual void Paint(Option_t *option="")
This method must be overridden if a class wants to paint itself.
Definition: TObject.cxx:521
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition: TObject.h:68
@ kMustCleanup
if object destructor must call RecursiveRemove()
Definition: TObject.h:60
Longptr_t ExecPlugin(int nargs, const T &... params)
Int_t LoadPlugin()
Load the plugin library for this handler.
Basic string class.
Definition: TString.h:136
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1150
const char * Data() const
Definition: TString.h:369
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
Bool_t IsNull() const
Definition: TString.h:407
TString & Append(const char *cs)
Definition: TString.h:564
Bool_t IsWhitespace() const
Definition: TString.h:408
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
Abstract Base Class for Fitting.
static TVirtualFitter * GetFitter()
static: return the current Fitter
virtual TObject * GetUserFunc() const
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition: TVirtualPad.h:51
virtual TList * GetListOfPrimitives() const =0
virtual TVirtualPad * cd(Int_t subpadnumber=0)=0
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
void Clear(Option_t *option="") override=0
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
TGraphErrors * gr
Definition: legend1.C:25
TH1F * h1
Definition: legend1.C:5
TF1 * f1
Definition: legend1.C:11
def fit(model, train_loader, val_loader, num_epochs, batch_size, optimizer, criterion, save_best, scheduler)
TFitResultPtr FitObject(TH1 *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
fitting function for a TH1 (called from TH1::Fit)
Definition: HFitImpl.cxx:971
void FitOptionsMake(EFitObjectType type, const char *option, Foption_t &fitOption)
Decode list of options into fitOption.
Definition: HFitImpl.cxx:684
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
static constexpr double mg
Double_t Log(Double_t x)
Definition: TMath.h:760
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
Double_t Log10(Double_t x)
Definition: TMath.h:764
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Definition: first.py:1
Definition: graph.py:1
auto * m
Definition: textangle.C:8
auto * l
Definition: textangle.C:4
auto * a
Definition: textangle.C:12