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}
372
373
374////////////////////////////////////////////////////////////////////////////////
375/// Constructor with name and title.
376
377TMultiGraph::TMultiGraph(const char *name, const char *title)
378 : TNamed(name,title)
379{
380}
381
382////////////////////////////////////////////////////////////////////////////////
383/// TMultiGraph destructor.
384
386{
387 if (!fGraphs) return;
388 TObject *g;
389 TIter next(fGraphs);
390 while ((g = next())) {
391 g->ResetBit(kMustCleanup);
392 }
393 fGraphs->Delete();
394 delete fGraphs;
395 fGraphs = nullptr;
396 delete fHistogram;
397 fHistogram = nullptr;
398 if (fFunctions) {
400 //special logic to support the case where the same object is
401 //added multiple times in fFunctions.
402 //This case happens when the same object is added with different
403 //drawing modes
404 TObject *obj;
405 while ((obj = fFunctions->First())) {
406 while (fFunctions->Remove(obj)) { }
407 delete obj;
408 }
409 delete fFunctions;
410 fFunctions = nullptr;
411 }
412}
413
414
415////////////////////////////////////////////////////////////////////////////////
416/// Add a new graph to the list of graphs.
417/// Note that the graph is now owned by the TMultigraph.
418/// Deleting the TMultiGraph object will automatically delete the graphs.
419/// You should not delete the graphs when the TMultigraph is still active.
420
422{
423 if (!fGraphs) fGraphs = new TList();
424 graph->SetBit(kMustCleanup);
425 fGraphs->Add(graph,chopt);
426}
427
428
429////////////////////////////////////////////////////////////////////////////////
430/// Add all the graphs in "multigraph" to the list of graphs.
431///
432/// - If "chopt" is defined all the graphs in "multigraph" will be added with
433/// the "chopt" option.
434/// - If "chopt" is undefined each graph will be added with the option it had
435/// in "multigraph".
436
437void TMultiGraph::Add(TMultiGraph *multigraph, Option_t *chopt)
438{
439 TList *graphlist = multigraph->GetListOfGraphs();
440 if (!graphlist) return;
441
442 if (!fGraphs) fGraphs = new TList();
443
444 TObjOptLink *lnk = (TObjOptLink*)graphlist->FirstLink();
445 TObject *obj = nullptr;
446
447 while (lnk) {
448 obj = lnk->GetObject();
449 if (!strlen(chopt)) fGraphs->Add(obj,lnk->GetOption());
450 else fGraphs->Add(obj,chopt);
451 lnk = (TObjOptLink*)lnk->Next();
452 }
453}
454
455
456////////////////////////////////////////////////////////////////////////////////
457/// Browse multigraph.
458
460{
461 TString opt = gEnv->GetValue("TGraph.BrowseOption", "");
462 if (opt.IsNull()) {
463 opt = b ? b->GetDrawOption() : "alp";
464 opt = (opt == "") ? "alp" : opt.Data();
465 }
466 Draw(opt.Data());
467 gPad->Update();
468}
469
470
471////////////////////////////////////////////////////////////////////////////////
472/// Compute distance from point px,py to each graph.
473
475{
476 // Are we on the axis?
477 const Int_t kMaxDiff = 10;
478 Int_t distance = 9999;
479 if (fHistogram) {
480 distance = fHistogram->DistancetoPrimitive(px,py);
481 if (distance <= 0) return distance;
482 }
483
484 // Loop on the list of graphs
485 if (!fGraphs) return distance;
486 TGraph *g;
487 TIter next(fGraphs);
488 while ((g = (TGraph*) next())) {
489 Int_t dist = g->DistancetoPrimitive(px,py);
490 if (dist <= 0) return 0;
491 if (dist < kMaxDiff) {gPad->SetSelected(g); return dist;}
492 }
493 return distance;
494}
495
496
497////////////////////////////////////////////////////////////////////////////////
498/// Draw this multigraph with its current attributes.
499///
500/// Options to draw a graph are described in TGraphPainter.
501///
502/// The drawing option for each TGraph may be specified as an optional
503/// second argument of the Add function. You can use GetGraphDrawOption
504/// to return this option.
505///
506/// If a draw option is specified, it will be used to draw the graph,
507/// otherwise the graph will be drawn with the option specified in
508/// TMultiGraph::Draw. Use GetDrawOption to return the option specified
509/// when drawing the TMultiGraph.
510
512{
513 TString opt = option;
514 opt.ToLower();
515
516 if (gPad) {
517 if (!gPad->IsEditable()) gROOT->MakeDefCanvas();
518 if (opt.Contains("a")) gPad->Clear();
519 }
521}
522
523
524////////////////////////////////////////////////////////////////////////////////
525/// Fit this graph with function with name fname.
526///
527/// interface to TF1::Fit(TF1 *f1...
528
530{
531 char *linear = (char*)strstr(fname, "++");
532 if (linear) {
533 TF1 f1(fname, fname, xmin, xmax);
534 return Fit(&f1,option,"",xmin,xmax);
535 }
536 TF1 * f1 = (TF1*)gROOT->GetFunction(fname);
537 if (!f1) { Printf("Unknown function: %s",fname); return -1; }
538
539 return Fit(f1,option,"",xmin,xmax);
540}
541
542
543////////////////////////////////////////////////////////////////////////////////
544/// Fit this multigraph with function f1.
545///
546/// In this function all graphs of the multigraph are fitted simultaneously
547///
548/// f1 is an already predefined function created by TF1.
549/// Predefined functions such as gaus, expo and poln are automatically
550/// created by ROOT.
551///
552/// The list of fit options is given in parameter `option`which may takes the
553/// following values:
554///
555/// - "W" Ignore all the point errors
556/// - "U" Use a User specified fitting algorithm (via SetFCN)
557/// - "Q" Quiet mode (minimum printing)
558/// - "V" Verbose mode (default is between Q and V)
559/// - "B" Use this option when you want to fix one or more parameters
560/// and the fitting function is like "gaus","expo","poln","landau".
561/// - "R" Use the Range specified in the function range
562/// - "N" Do not store the graphics function, do not draw
563/// - "0" Do not plot the result of the fit. By default the fitted function
564/// is drawn unless the option"N" above is specified.
565/// - "+" Add this new fitted function to the list of fitted functions
566/// (by default, any previous function is deleted)
567/// - "C" In case of linear fitting, not calculate the chisquare (saves time)
568/// - "F" If fitting a polN, switch to minuit fitter
569/// - "ROB" In case of linear fitting, compute the LTS regression
570/// coefficients (robust(resistant) regression), using
571/// the default fraction of good points
572/// - "ROB=0.x" - compute the LTS regression coefficients, using
573/// 0.x as a fraction of good points
574///
575/// When the fit is drawn (by default), the parameter goption may be used
576/// to specify a list of graphics options. See TGraph::Paint for a complete
577/// list of these options.
578///
579/// In order to use the Range option, one must first create a function
580/// with the expression to be fitted. For example, if your graph
581/// has a defined range between -4 and 4 and you want to fit a gaussian
582/// only in the interval 1 to 3, you can do:
583/// ~~~ {.cpp}
584/// TF1 *f1 = new TF1("f1","gaus",1,3);
585/// graph->Fit("f1","R");
586/// ~~~
587///
588/// ### Who is calling this function ?
589///
590/// Note that this function is called when calling TGraphErrors::Fit
591/// or TGraphAsymmErrors::Fit ot TGraphBentErrors::Fit
592/// see the discussion below on the errors calculation.
593///
594/// ### Setting initial conditions
595///
596/// Parameters must be initialized before invoking the Fit function.
597/// The setting of the parameter initial values is automatic for the
598/// predefined functions : poln, expo, gaus, landau. One can however disable
599/// this automatic computation by specifying the option "B".
600/// You can specify boundary limits for some or all parameters via
601/// ~~~ {.cpp}
602/// f1->SetParLimits(p_number, parmin, parmax);
603/// ~~~
604/// if `parmin>=parmax`, the parameter is fixed
605/// Note that you are not forced to fix the limits for all parameters.
606/// For example, if you fit a function with 6 parameters, you can do:
607/// ~~~ {.cpp}
608/// func->SetParameters(0,3.1,1.e-6,0.1,-8,100);
609/// func->SetParLimits(4,-10,-4);
610/// func->SetParLimits(5, 1,1);
611/// ~~~
612/// With this setup, parameters 0->3 can vary freely
613/// Parameter 4 has boundaries [-10,-4] with initial value -8
614/// Parameter 5 is fixed to 100.
615///
616/// ### Fit range
617///
618/// The fit range can be specified in two ways:
619///
620/// - specify rxmax > rxmin (default is rxmin=rxmax=0)
621/// - specify the option "R". In this case, the function will be taken
622/// instead of the full graph range.
623///
624/// ### Changing the fitting function
625///
626/// By default a chi2 fitting function is used for fitting the TGraphs's.
627/// The function is implemented in `FitUtil::EvaluateChi2`.
628/// In case of TGraphErrors an effective chi2 is used
629/// (see TGraphErrors fit in TGraph::Fit) and is implemented in
630/// `FitUtil::EvaluateChi2Effective`
631/// To specify a User defined fitting function, specify option "U" and
632/// call the following function:
633/// ~~~ {.cpp}
634/// TVirtualFitter::Fitter(mygraph)->SetFCN(MyFittingFunction)
635/// ~~~
636/// where MyFittingFunction is of type:
637/// ~~~ {.cpp}
638/// extern void MyFittingFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *u, Int_t flag);
639/// ~~~
640///
641/// ### Access to the fit result
642///
643/// The function returns a TFitResultPtr which can hold a pointer to a TFitResult object.
644/// By default the TFitResultPtr contains only the status of the fit and it converts
645/// automatically to an integer. If the option "S" is instead used, TFitResultPtr contains
646/// the TFitResult and behaves as a smart pointer to it. For example one can do:
647/// ~~~ {.cpp}
648/// TFitResultPtr r = graph->Fit("myFunc","S");
649/// TMatrixDSym cov = r->GetCovarianceMatrix(); // to access the covariance matrix
650/// Double_t par0 = r->Parameter(0); // retrieve the value for the parameter 0
651/// Double_t err0 = r->ParError(0); // retrieve the error for the parameter 0
652/// r->Print("V"); // print full information of fit including covariance matrix
653/// r->Write(); // store the result in a file
654/// ~~~
655///
656/// The fit parameters, error and chi2 (but not covariance matrix) can be retrieved also
657/// from the fitted function.
658///
659/// ### Associated functions
660///
661/// One or more object (typically a TF1*) can be added to the list
662/// of functions (fFunctions) associated to each graph.
663/// When TGraph::Fit is invoked, the fitted function is added to this list.
664/// Given a graph gr, one can retrieve an associated function
665/// with:
666/// ~~~ {.cpp}
667/// TF1 *myfunc = gr->GetFunction("myfunc");
668/// ~~~
669///
670/// If the graph is made persistent, the list of
671/// associated functions is also persistent. Given a pointer (see above)
672/// to an associated function myfunc, one can retrieve the function/fit
673/// parameters with calls such as:
674/// ~~~ {.cpp}
675/// Double_t chi2 = myfunc->GetChisquare();
676/// Double_t par0 = myfunc->GetParameter(0); //value of 1st parameter
677/// Double_t err0 = myfunc->GetParError(0); //error on first parameter
678/// ~~~
679///
680/// ### Fit Statistics
681///
682/// You can change the statistics box to display the fit parameters with
683/// the TStyle::SetOptFit(mode) method. This mode has four digits.
684/// mode = pcev (default = 0111)
685///
686/// - v = 1; print name/values of parameters
687/// - e = 1; print errors (if e=1, v must be 1)
688/// - c = 1; print Chisquare/Number of degrees of freedom
689/// - p = 1; print Probability
690///
691/// For example: `gStyle->SetOptFit(1011);`
692/// prints the fit probability, parameter names/values, and errors.
693/// You can change the position of the statistics box with these lines
694/// (where g is a pointer to the TGraph):
695///
696/// ~~~ {.cpp}
697/// Root > TPaveStats *st = (TPaveStats*)g->GetListOfFunctions()->FindObject("stats")
698/// Root > st->SetX1NDC(newx1); //new x start position
699/// Root > st->SetX2NDC(newx2); //new x end position
700/// ~~~
701
703{
704 // internal multigraph fitting methods
705 Foption_t fitOption;
707
708 // create range and minimizer options with default values
709 ROOT::Fit::DataRange range(rxmin,rxmax);
711 return ROOT::Fit::FitObject(this, f1 , fitOption , minOption, goption, range);
712
713}
714
715////////////////////////////////////////////////////////////////////////////////
716/// Display a panel with all histogram fit options.
717/// See class TFitPanel for example
718
720{
721 if (!gPad)
722 gROOT->MakeDefCanvas();
723
724 if (!gPad) {
725 Error("FitPanel", "Unable to create a default canvas");
726 return;
727 }
728
729 // use plugin manager to create instance of TFitEditor
730 TPluginHandler *handler = gROOT->GetPluginManager()->FindHandler("TFitEditor");
731 if (handler && handler->LoadPlugin() != -1) {
732 if (handler->ExecPlugin(2, gPad, this) == 0)
733 Error("FitPanel", "Unable to crate the FitPanel");
734 }
735 else
736 Error("FitPanel", "Unable to find the FitPanel plug-in");
737}
738
739////////////////////////////////////////////////////////////////////////////////
740/// Return the draw option for the TGraph `gr` in this TMultiGraph.
741/// The return option is the one specified when calling TMultiGraph::Add(gr,option).
742
744{
745 if (!fGraphs || !gr) return "";
746 TListIter next(fGraphs);
747 TObject *obj;
748 while ((obj = next())) {
749 if (obj == (TObject*)gr) return next.GetOption();
750 }
751 return "";
752}
753
754
755////////////////////////////////////////////////////////////////////////////////
756/// Compute Initial values of parameters for a gaussian.
757
759{
760 Double_t allcha, sumx, sumx2, x, val, rms, mean;
761 Int_t bin;
762 const Double_t sqrtpi = 2.506628;
763
764 // Compute mean value and RMS of the graph in the given range
765 Int_t np = 0;
766 allcha = sumx = sumx2 = 0;
767 TGraph *g;
768 TIter next(fGraphs);
769 Double_t *px, *py;
770 Int_t npp; //number of points in each graph
771 while ((g = (TGraph*) next())) {
772 px=g->GetX();
773 py=g->GetY();
774 npp=g->GetN();
775 for (bin=0; bin<npp; bin++) {
776 x=px[bin];
777 if (x<xmin || x>xmax) continue;
778 np++;
779 val=py[bin];
780 sumx+=val*x;
781 sumx2+=val*x*x;
782 allcha+=val;
783 }
784 }
785 if (np == 0 || allcha == 0) return;
786 mean = sumx/allcha;
787 rms = TMath::Sqrt(sumx2/allcha - mean*mean);
788
789 Double_t binwidx = TMath::Abs((xmax-xmin)/np);
790 if (rms == 0) rms = 1;
792 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
793 f1->SetParameter(0,binwidx*allcha/(sqrtpi*rms));
794 f1->SetParameter(1,mean);
795 f1->SetParameter(2,rms);
796 f1->SetParLimits(2,0,10*rms);
797}
798
799
800////////////////////////////////////////////////////////////////////////////////
801/// Compute Initial values of parameters for an exponential.
802
804{
805 Double_t constant, slope;
806 Int_t ifail;
807
808 LeastSquareLinearFit(-1, constant, slope, ifail, xmin, xmax);
809
811 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
812 f1->SetParameter(0,constant);
813 f1->SetParameter(1,slope);
814}
815
816
817////////////////////////////////////////////////////////////////////////////////
818/// Compute Initial values of parameters for a polynom.
819
821{
822 Double_t fitpar[25];
823
825 TF1 *f1 = (TF1*)grFitter->GetUserFunc();
826 Int_t npar = f1->GetNpar();
827
828 LeastSquareFit(npar, fitpar, xmin, xmax);
829
830 for (Int_t i=0;i<npar;i++) f1->SetParameter(i, fitpar[i]);
831}
832
833
834////////////////////////////////////////////////////////////////////////////////
835/// Least squares lpolynomial fitting without weights.
836///
837/// - m number of parameters
838/// - a array of parameters
839/// - first 1st point number to fit (default =0)
840/// - last last point number to fit (default=fNpoints-1)
841///
842/// based on CERNLIB routine LSQ: Translated to C++ by Rene Brun
843
845{
846 const Double_t zero = 0.;
847 const Double_t one = 1.;
848 const Int_t idim = 20;
849
850 Double_t b[400] /* was [20][20] */;
851 Int_t i, k, l, ifail, bin;
852 Double_t power;
853 Double_t da[20], xk, yk;
854
855
856 //count the total number of points to fit
857 TGraph *g;
858 TIter next(fGraphs);
859 Double_t *px, *py;
860 Int_t n=0;
861 Int_t npp;
862 while ((g = (TGraph*) next())) {
863 px=g->GetX();
864 npp=g->GetN();
865 for (bin=0; bin<npp; bin++) {
866 xk=px[bin];
867 if (xk < xmin || xk > xmax) continue;
868 n++;
869 }
870 }
871 if (m <= 2) {
872 LeastSquareLinearFit(n, a[0], a[1], ifail, xmin, xmax);
873 return;
874 }
875 if (m > idim || m > n) return;
876 da[0] = zero;
877 for (l = 2; l <= m; ++l) {
878 b[l-1] = zero;
879 b[m + l*20 - 21] = zero;
880 da[l-1] = zero;
881 }
882 Int_t np = 0;
883
884 next.Reset();
885 while ((g = (TGraph*) next())) {
886 px=g->GetX();
887 py=g->GetY();
888 npp=g->GetN();
889
890 for (k = 0; k <= npp; ++k) {
891 xk = px[k];
892 if (xk < xmin || xk > xmax) continue;
893 np++;
894 yk = py[k];
895 power = one;
896 da[0] += yk;
897 for (l = 2; l <= m; ++l) {
898 power *= xk;
899 b[l-1] += power;
900 da[l-1] += power*yk;
901 }
902 for (l = 2; l <= m; ++l) {
903 power *= xk;
904 b[m + l*20 - 21] += power;
905 }
906 }
907 }
908 b[0] = Double_t(np);
909 for (i = 3; i <= m; ++i) {
910 for (k = i; k <= m; ++k) {
911 b[k - 1 + (i-1)*20 - 21] = b[k + (i-2)*20 - 21];
912 }
913 }
914 H1LeastSquareSeqnd(m, b, idim, ifail, 1, da);
915
916 if (ifail < 0) {
917 //a[0] = fY[0];
918 py=((TGraph *)fGraphs->First())->GetY();
919 a[0]=py[0];
920 for (i=1; i<m; ++i) a[i] = 0;
921 return;
922 }
923 for (i=0; i<m; ++i) a[i] = da[i];
924}
925
926
927////////////////////////////////////////////////////////////////////////////////
928/// Least square linear fit without weights.
929///
930/// Fit a straight line (a0 + a1*x) to the data in this graph.
931///
932/// - ndata: number of points to fit
933/// - first: first point number to fit
934/// - last: last point to fit O(ndata should be last-first
935/// - ifail: return parameter indicating the status of the fit (ifail=0, fit is OK)
936///
937/// extracted from CERNLIB LLSQ: Translated to C++ by Rene Brun
938
941{
942 Double_t xbar, ybar, x2bar;
943 Int_t i;
944 Double_t xybar;
945 Double_t fn, xk, yk;
946 Double_t det;
947
948 ifail = -2;
949 xbar = ybar = x2bar = xybar = 0;
950 Int_t np = 0;
951 TGraph *g;
952 TIter next(fGraphs);
953 Double_t *px, *py;
954 Int_t npp;
955 while ((g = (TGraph*) next())) {
956 px=g->GetX();
957 py=g->GetY();
958 npp=g->GetN();
959 for (i = 0; i < npp; ++i) {
960 xk = px[i];
961 if (xk < xmin || xk > xmax) continue;
962 np++;
963 yk = py[i];
964 if (ndata < 0) {
965 if (yk <= 0) yk = 1e-9;
966 yk = TMath::Log(yk);
967 }
968 xbar += xk;
969 ybar += yk;
970 x2bar += xk*xk;
971 xybar += xk*yk;
972 }
973 }
974 fn = Double_t(np);
975 det = fn*x2bar - xbar*xbar;
976 ifail = -1;
977 if (det <= 0) {
978 if (fn > 0) a0 = ybar/fn;
979 else a0 = 0;
980 a1 = 0;
981 return;
982 }
983 ifail = 0;
984 a0 = (x2bar*ybar - xbar*xybar) / det;
985 a1 = (fn*xybar - xbar*ybar) / det;
986}
987
988
989////////////////////////////////////////////////////////////////////////////////
990/// Return 1 if the point (x,y) is inside one of the graphs 0 otherwise.
991
993{
994 Int_t in = 0;
995 if (!fGraphs) return in;
996 TGraph *g;
997 TIter next(fGraphs);
998 while ((g = (TGraph*) next())) {
999 in = g->IsInside(x, y);
1000 if (in) return in;
1001 }
1002 return in;
1003}
1004
1005
1006////////////////////////////////////////////////////////////////////////////////
1007/// Returns a pointer to the histogram used to draw the axis.
1008/// Takes into account following cases.
1009///
1010/// 1. if `fHistogram` exists it is returned
1011/// 2. if `fHistogram` doesn't exists and `gPad` exists `gPad` is updated. That
1012/// may trigger the creation of `fHistogram`. If `fHistogram` still does not
1013/// exit but `hframe` does (if user called `TPad::DrawFrame`) the pointer to
1014/// `hframe` histogram is returned
1015/// 3. after the two previous steps, if `fHistogram` still doesn't exist, then
1016/// it is created.
1017
1019{
1020 if (fHistogram) return fHistogram;
1021
1022 if (gPad) {
1023 gPad->Modified();
1024 gPad->Update();
1025 if (fHistogram) return fHistogram;
1026 TH1F *h1 = (TH1F*)gPad->FindObject("hframe");
1027 if (h1) return h1;
1028 }
1029
1030 Bool_t initialrangeset = kFALSE;
1031 Double_t rwxmin = 0.,rwxmax = 0.,rwymin = 0.,rwymax = 0.;
1032 TGraph *g;
1033 Int_t npt = 100 ;
1034 TIter next(fGraphs);
1035 while ((g = (TGraph*) next())) {
1036 if (g->GetN() <= 0) continue;
1037 if (initialrangeset) {
1038 Double_t rx1,ry1,rx2,ry2;
1039 g->ComputeRange(rx1, ry1, rx2, ry2);
1040 if (rx1 < rwxmin) rwxmin = rx1;
1041 if (ry1 < rwymin) rwymin = ry1;
1042 if (rx2 > rwxmax) rwxmax = rx2;
1043 if (ry2 > rwymax) rwymax = ry2;
1044 } else {
1045 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1046 initialrangeset = kTRUE;
1047 }
1048 if (g->GetN() > npt) npt = g->GetN();
1049 }
1050 if (rwxmin == rwxmax) rwxmax += 1.;
1051 if (rwymin == rwymax) rwymax += 1.;
1052 double dx = 0.05*(rwxmax-rwxmin);
1053 double dy = 0.05*(rwymax-rwymin);
1054 rwxmin = rwxmin - dx;
1055 rwxmax = rwxmax + dx;
1056 if (gPad && gPad->GetLogy()) {
1057 if (rwymin <= 0) rwymin = 0.001*rwymax;
1058 double r = rwymax/rwymin;
1059 rwymin = rwymin/(1+0.5*TMath::Log10(r));
1060 rwymax = rwymax*(1+0.2*TMath::Log10(r));
1061 } else {
1062 rwymin = rwymin - dy;
1063 rwymax = rwymax + dy;
1064 }
1065 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1066 if (!fHistogram) return 0;
1067 fHistogram->SetMinimum(rwymin);
1069 fHistogram->SetMaximum(rwymax);
1070 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1072 return fHistogram;
1073}
1074
1075
1076////////////////////////////////////////////////////////////////////////////////
1077/// Return pointer to function with name.
1078///
1079/// Functions such as TGraph::Fit store the fitted function in the list of
1080/// functions of this graph.
1081
1083{
1084 if (!fFunctions) return nullptr;
1085 return (TF1*)fFunctions->FindObject(name);
1086}
1087
1088////////////////////////////////////////////////////////////////////////////////
1089/// Return pointer to list of functions.
1090/// If pointer is null create the list
1091
1093{
1094 if (!fFunctions) fFunctions = new TList();
1095 return fFunctions;
1096}
1097
1098
1099////////////////////////////////////////////////////////////////////////////////
1100/// Get x axis of the graph.
1101/// This method returns a valid axis only after the TMultigraph has been drawn.
1102
1104{
1105 TH1 *h = GetHistogram();
1106 if (!h) return nullptr;
1107 return h->GetXaxis();
1108}
1109
1110
1111////////////////////////////////////////////////////////////////////////////////
1112/// Get y axis of the graph.
1113/// This method returns a valid axis only after the TMultigraph has been drawn.
1114
1116{
1117 TH1 *h = GetHistogram();
1118 if (!h) return nullptr;
1119 return h->GetYaxis();
1120}
1121
1122
1123////////////////////////////////////////////////////////////////////////////////
1124/// Paint all the graphs of this multigraph.
1125
1127{
1128 const TPickerStackGuard pushGuard(this);
1129
1130 if (!fGraphs) return;
1131 if (fGraphs->GetSize() == 0) return;
1132
1133 char option[128];
1134 strlcpy(option,choptin,128);
1135 Int_t nch = choptin ? strlen(choptin) : 0;
1136 for (Int_t i=0;i<nch;i++) option[i] = toupper(option[i]);
1137
1138 // Automatic color
1139 char *l1 = strstr(option,"PFC"); // Automatic Fill Color
1140 char *l2 = strstr(option,"PLC"); // Automatic Line Color
1141 char *l3 = strstr(option,"PMC"); // Automatic Marker Color
1142 if (l1 || l2 || l3) {
1143 TString opt1 = option; opt1.ToLower();
1144 if (l1) memcpy(l1," ",3);
1145 if (l2) memcpy(l2," ",3);
1146 if (l3) memcpy(l3," ",3);
1148 TGraph* gAti;
1149 Int_t ngraphs = fGraphs->GetSize();
1150 Int_t ic;
1151 gPad->IncrementPaletteColor(ngraphs, opt1);
1152 for (Int_t i=0;i<ngraphs;i++) {
1153 ic = gPad->NextPaletteColor();
1154 gAti = (TGraph*)(fGraphs->At(i));
1155 if (l1) gAti->SetFillColor(ic);
1156 if (l2) gAti->SetLineColor(ic);
1157 if (l3) gAti->SetMarkerColor(ic);
1158 lnk = (TObjOptLink*)lnk->Next();
1159 }
1160 }
1161
1162 TString chopt = option;
1163
1164 char *l = (char*)strstr(chopt.Data(),"3D");
1165 if (l) {
1166 l = (char*)strstr(chopt.Data(),"L");
1167 if (l) PaintPolyLine3D(chopt.Data());
1168 return;
1169 }
1170
1171 l = (char*)strstr(chopt.Data(),"PADS");
1172 if (l) {
1173 chopt.ReplaceAll("PADS","");
1174 PaintPads(chopt.Data());
1175 return;
1176 }
1177
1178 char *lrx = (char *)strstr(chopt.Data(), "RX"); // Reverse graphs along X axis
1179 char *lry = (char *)strstr(chopt.Data(), "RY"); // Reverse graphs along Y axis
1180 if (lrx || lry) {
1181 PaintReverse(chopt.Data());
1182 return;
1183 }
1184
1185 TGraph *g;
1186
1187 l = (char*)strstr(chopt.Data(),"A");
1188 if (l) {
1189 *l = ' ';
1190 TIter next(fGraphs);
1191 Int_t npt = 100;
1192 Double_t maximum, minimum, rwxmin, rwxmax, rwymin, rwymax, uxmin, uxmax, dx, dy;
1193 rwxmin = gPad->GetUxmin();
1194 rwxmax = gPad->GetUxmax();
1195 rwymin = gPad->GetUymin();
1196 rwymax = gPad->GetUymax();
1197 std::string xtitle, ytitle, timeformat;
1198 Int_t firstx = 0;
1199 Int_t lastx = 0;
1200 Bool_t timedisplay = kFALSE;
1201
1202 if (fHistogram) {
1203 //cleanup in case of a previous unzoom and in case one of the TGraph has changed
1205 Int_t ngraphs = fGraphs->GetSize();
1206 Bool_t reset_hist = kFALSE;
1207 for (Int_t i=0;i<ngraphs;i++) {
1208 TGraph* gAti = (TGraph*)(fGraphs->At(i));
1209 if(gAti->TestBit(TGraph::kResetHisto)) {reset_hist = kTRUE; break;}
1210 lnk = (TObjOptLink*)lnk->Next();
1211 }
1212 if (fHistogram->GetMinimum() >= fHistogram->GetMaximum() || reset_hist) {
1213 firstx = fHistogram->GetXaxis()->GetFirst();
1214 lastx = fHistogram->GetXaxis()->GetLast();
1215 timedisplay = fHistogram->GetXaxis()->GetTimeDisplay();
1216 if (strlen(fHistogram->GetXaxis()->GetTitle()) > 0)
1217 xtitle = fHistogram->GetXaxis()->GetTitle();
1218 if (strlen(fHistogram->GetYaxis()->GetTitle()) > 0)
1219 ytitle = fHistogram->GetYaxis()->GetTitle();
1220 if (strlen(fHistogram->GetXaxis()->GetTimeFormat()) > 0)
1221 timeformat = fHistogram->GetXaxis()->GetTimeFormat();
1222 delete fHistogram;
1223 fHistogram = nullptr;
1224 }
1225 }
1226 if (fHistogram) {
1227 minimum = fHistogram->GetYaxis()->GetXmin();
1228 maximum = fHistogram->GetYaxis()->GetXmax();
1229 uxmin = gPad->PadtoX(rwxmin);
1230 uxmax = gPad->PadtoX(rwxmax);
1231 } else {
1232 Bool_t initialrangeset = kFALSE;
1233 while ((g = (TGraph*) next())) {
1234 if (g->GetN() <= 0) continue;
1235 if (initialrangeset) {
1236 Double_t rx1,ry1,rx2,ry2;
1237 g->ComputeRange(rx1, ry1, rx2, ry2);
1238 if (rx1 < rwxmin) rwxmin = rx1;
1239 if (ry1 < rwymin) rwymin = ry1;
1240 if (rx2 > rwxmax) rwxmax = rx2;
1241 if (ry2 > rwymax) rwymax = ry2;
1242 } else {
1243 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1244 initialrangeset = kTRUE;
1245 }
1246 if (g->GetN() > npt) npt = g->GetN();
1247 }
1248 if (rwxmin == rwxmax) rwxmax += 1.;
1249 if (rwymin == rwymax) rwymax += 1.;
1250 dx = 0.05*(rwxmax-rwxmin);
1251 dy = 0.05*(rwymax-rwymin);
1252 uxmin = rwxmin - dx;
1253 uxmax = rwxmax + dx;
1254 if (gPad->GetLogy()) {
1255 if (rwymin <= 0) rwymin = 0.001*rwymax;
1256 minimum = rwymin/(1+0.5*TMath::Log10(rwymax/rwymin));
1257 maximum = rwymax*(1+0.2*TMath::Log10(rwymax/rwymin));
1258 } else {
1259 minimum = rwymin - dy;
1260 maximum = rwymax + dy;
1261 }
1262 if (minimum < 0 && rwymin >= 0) minimum = 0;
1263 if (maximum > 0 && rwymax <= 0) maximum = 0;
1264 }
1265
1266 if (fMinimum != -1111) rwymin = minimum = fMinimum;
1267 if (fMaximum != -1111) rwymax = maximum = fMaximum;
1268 if (uxmin < 0 && rwxmin >= 0) {
1269 if (gPad->GetLogx()) uxmin = 0.9*rwxmin;
1270 //else uxmin = 0;
1271 }
1272 if (uxmax > 0 && rwxmax <= 0) {
1273 if (gPad->GetLogx()) uxmax = 1.1*rwxmax;
1274 //else uxmax = 0;
1275 }
1276 if (minimum < 0 && rwymin >= 0) {
1277 if (gPad->GetLogy()) minimum = 0.9*rwymin;
1278 //else minimum = 0;
1279 }
1280 if (maximum > 0 && rwymax <= 0) {
1281 if (gPad->GetLogy()) maximum = 1.1*rwymax;
1282 //else maximum = 0;
1283 }
1284 if (minimum <= 0 && gPad->GetLogy()) minimum = 0.001*maximum;
1285 if (uxmin <= 0 && gPad->GetLogx()) {
1286 if (uxmax > 1000) uxmin = 1;
1287 else uxmin = 0.001*uxmax;
1288 }
1289 rwymin = minimum;
1290 rwymax = maximum;
1291 if (fHistogram) {
1292 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1293 }
1294
1295 // Create a temporary histogram to draw the axis
1296 if (!fHistogram) {
1297 // the graph is created with at least as many channels as there are points
1298 // to permit zooming on the full range
1299 rwxmin = uxmin;
1300 rwxmax = uxmax;
1301 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1302 if (!fHistogram) return;
1303 fHistogram->SetMinimum(rwymin);
1305 fHistogram->SetMaximum(rwymax);
1306 fHistogram->GetYaxis()->SetLimits(rwymin,rwymax);
1308 if (!xtitle.empty()) fHistogram->GetXaxis()->SetTitle(xtitle.c_str());
1309 if (!ytitle.empty()) fHistogram->GetYaxis()->SetTitle(ytitle.c_str());
1310 if (firstx != lastx) fHistogram->GetXaxis()->SetRange(firstx,lastx);
1311 if (timedisplay) {fHistogram->GetXaxis()->SetTimeDisplay(timedisplay);}
1312 if (!timeformat.empty()) fHistogram->GetXaxis()->SetTimeFormat(timeformat.c_str());
1313 }
1314 TString chopth = "0";
1315 if ((char*)strstr(chopt.Data(),"X+")) chopth.Append("X+");
1316 if ((char*)strstr(chopt.Data(),"Y+")) chopth.Append("Y+");
1317 if ((char*)strstr(chopt.Data(),"I")) chopth.Append("A");
1318 fHistogram->Paint(chopth.Data());
1319 }
1320
1321 TGraph *gfit = nullptr;
1322 if (fGraphs) {
1324 TObject *obj = 0;
1325
1326 chopt.ReplaceAll("A","");
1327
1328 while (lnk) {
1329
1330 obj = lnk->GetObject();
1331
1332 gPad->PushSelectableObject(obj);
1333
1334 if (!gPad->PadInHighlightMode() || (gPad->PadInHighlightMode() && obj == gPad->GetSelected())) {
1335 TString opt = lnk->GetOption();
1336 if (!opt.IsWhitespace())
1337 obj->Paint(opt.ReplaceAll("A","").Data());
1338 else {
1339 if (!chopt.IsWhitespace()) obj->Paint(chopt.Data());
1340 else obj->Paint("L");
1341 }
1342 }
1343
1344 lnk = (TObjOptLink*)lnk->Next();
1345 }
1346
1347 gfit = (TGraph*)obj; // pick one TGraph in the list to paint the fit parameters.
1348 }
1349
1350 TObject *f;
1351 TF1 *fit = nullptr;
1352 if (fFunctions) {
1353 TIter next(fFunctions);
1354 while ((f = (TObject*) next())) {
1355 if (f->InheritsFrom(TF1::Class())) {
1356 if (f->TestBit(TF1::kNotDraw) == 0) f->Paint("lsame");
1357 fit = (TF1*)f;
1358 } else {
1359 f->Paint();
1360 }
1361 }
1362 }
1363
1364 if (gfit && fit) gfit->PaintStats(fit);
1365}
1366
1367
1368////////////////////////////////////////////////////////////////////////////////
1369/// Divides the active pad and draws all Graphs in the Multigraph separately.
1370
1372{
1373 if (!gPad) return;
1374
1375 Int_t neededPads = fGraphs->GetSize();
1376 Int_t existingPads = 0;
1377
1378 TVirtualPad *curPad = gPad;
1379 TIter nextPad(curPad->GetListOfPrimitives());
1380
1381 while (auto obj = nextPad()) {
1382 if (obj->InheritsFrom(TVirtualPad::Class()))
1383 existingPads++;
1384 }
1385 if (existingPads < neededPads) {
1386 curPad->Clear();
1387 Int_t nx = (Int_t)TMath::Sqrt((Double_t)neededPads);
1388 if (nx*nx < neededPads) nx++;
1389 Int_t ny = nx;
1390 if (((nx*ny)-nx) >= neededPads) ny--;
1391 curPad->Divide(nx,ny);
1392 }
1393 Int_t i = 0;
1394
1395 TIter nextGraph(fGraphs);
1396 while (auto g = (TGraph *) nextGraph()) {
1397 curPad->cd(++i);
1398 TString apopt = nextGraph.GetOption();
1399 if ((apopt.Length() == 0) && option) apopt = option;
1400 if (apopt.Length() == 0) apopt = "L";
1401 g->Draw(apopt.Append("A").Data());
1402 }
1403
1404 curPad->cd();
1405}
1406
1407
1408////////////////////////////////////////////////////////////////////////////////
1409/// Paint all the graphs of this multigraph as 3D lines.
1410
1412{
1413 Int_t i, npt = 0;
1414 Double_t rwxmin=0., rwxmax=0., rwymin=0., rwymax=0.;
1415 TIter next(fGraphs);
1416
1417 TGraph *g = (TGraph*) next();
1418 if (g) {
1419 g->ComputeRange(rwxmin, rwymin, rwxmax, rwymax);
1420 npt = g->GetN();
1421 }
1422
1423 if (!fHistogram)
1424 fHistogram = new TH1F(GetName(),GetTitle(),npt,rwxmin,rwxmax);
1425
1426 while ((g = (TGraph*) next())) {
1427 Double_t rx1,ry1,rx2,ry2;
1428 g->ComputeRange(rx1, ry1, rx2, ry2);
1429 if (rx1 < rwxmin) rwxmin = rx1;
1430 if (ry1 < rwymin) rwymin = ry1;
1431 if (rx2 > rwxmax) rwxmax = rx2;
1432 if (ry2 > rwymax) rwymax = ry2;
1433 if (g->GetN() > npt) npt = g->GetN();
1434 }
1435
1436 Int_t ndiv = fGraphs->GetSize();
1437
1438 TH2F* frame = new TH2F("frame","", ndiv, 0., (Double_t)(ndiv), npt, rwxmin, rwxmax);
1439 if (fHistogram) {
1440 frame->SetTitle(fHistogram->GetTitle());
1441 frame->GetYaxis()->SetTitle(fHistogram->GetXaxis()->GetTitle());
1443 frame->GetZaxis()->SetTitle(fHistogram->GetYaxis()->GetTitle());
1444 }
1445
1446 TAxis *Xaxis = frame->GetXaxis();
1447 Xaxis->SetNdivisions(-ndiv);
1448 next.Reset();
1449 for (i=ndiv; i>=1; i--) {
1450 g = (TGraph*) next();
1451 Xaxis->SetBinLabel(i, g->GetTitle());
1452 }
1453
1454 frame->SetStats(kFALSE);
1455 if (fMinimum != -1111) frame->SetMinimum(fMinimum);
1456 else frame->SetMinimum(rwymin);
1457 if (fMaximum != -1111) frame->SetMaximum(fMaximum);
1458 else frame->SetMaximum(rwymax);
1459
1460 if (strstr(option,"A"))
1461 frame->Paint("lego9,fb,bb");
1462
1463 if (!strstr(option,"BB"))
1464 frame->Paint("lego9,fb,a,same");
1465
1466 Double_t xyz1[3], xyz2[3];
1467
1468 Double_t xl = frame->GetYaxis()->GetBinLowEdge(frame->GetYaxis()->GetFirst());
1469 Double_t xu = frame->GetYaxis()->GetBinUpEdge(frame->GetYaxis()->GetLast());
1470 Double_t yl = frame->GetMinimum();
1471 Double_t yu = frame->GetMaximum();
1472 Double_t xc[2],yc[2];
1473 next.Reset();
1474 Int_t j = ndiv;
1475
1476 while ((g = (TGraph*) next())) {
1477 npt = g->GetN();
1478 auto x = g->GetX();
1479 auto y = g->GetY();
1480 gPad->SetLineColor(g->GetLineColor());
1481 gPad->SetLineWidth(g->GetLineWidth());
1482 gPad->SetLineStyle(g->GetLineStyle());
1483 gPad->TAttLine::Modify();
1484 for (i=0; i<npt-1; i++) {
1485 xc[0] = x[i];
1486 xc[1] = x[i+1];
1487 yc[0] = y[i];
1488 yc[1] = y[i+1];
1489 if (gPad->Clip(&xc[0], &yc[0], xl, yl, xu, yu)<2) {
1490 xyz1[0] = j-0.5;
1491 xyz1[1] = xc[0];
1492 xyz1[2] = yc[0];
1493 xyz2[0] = j-0.5;
1494 xyz2[1] = xc[1];
1495 xyz2[2] = yc[1];
1496 gPad->PaintLine3D(xyz1, xyz2);
1497 }
1498 }
1499 j--;
1500 }
1501
1502 if (!strstr(option,"FB"))
1503 frame->Paint("lego9,bb,a,same");
1504 delete frame;
1505}
1506
1507
1508////////////////////////////////////////////////////////////////////////////////
1509/// Paint all the graphs of this multigraph reverting values along X and/or Y axis.
1510/// New graphs are created.
1511
1513{
1514 auto *h = GetHistogram();
1515 TH1F *hg = nullptr;
1516 TGraph *fg = nullptr;
1517 if (!h)
1518 return;
1519 TString mgopt = option;
1520 mgopt.ToLower();
1521
1522 TIter next(fGraphs);
1523 TGraph *g;
1524 Bool_t first = kTRUE;
1525 TString gopt;
1526 while ((g = (TGraph *)next())) {
1527 gopt = GetGraphDrawOption(g);
1528 gopt.Append(mgopt);
1529 if (first) {
1530 fg = g;
1531 hg = fg->GetHistogram();
1532 fg->SetHistogram(h);
1533 fg->Paint(gopt.Data());
1534 first = kFALSE;
1535 } else {
1536 g->Paint(gopt.ReplaceAll("a", "").Data());
1537 }
1538 }
1539 if (fg)
1540 fg->SetHistogram(hg);
1541}
1542
1543
1544////////////////////////////////////////////////////////////////////////////////
1545/// Print the list of graphs.
1546
1548{
1549 TGraph *g;
1550 if (fGraphs) {
1551 TIter next(fGraphs);
1552 while ((g = (TGraph*) next())) {
1553 g->Print(option);
1554 }
1555 }
1556}
1557
1558
1559////////////////////////////////////////////////////////////////////////////////
1560/// Recursively remove this object from a list. Typically implemented
1561/// by classes that can contain multiple references to a same object.
1562
1564{
1565 if (obj == fHistogram) {
1566 fHistogram = nullptr;
1567 return;
1568 }
1569
1570 if (fFunctions) {
1571 auto f = fFunctions->Remove(obj);
1572 if (f) return;
1573 }
1574
1575 if (!fGraphs) return;
1576 auto objr = fGraphs->Remove(obj);
1577 if (!objr) return;
1578
1579 delete fHistogram; fHistogram = nullptr;
1580 if (gPad) gPad->Modified();
1581}
1582
1583
1584////////////////////////////////////////////////////////////////////////////////
1585/// Save primitive as a C++ statement(s) on output stream out.
1586
1587void TMultiGraph::SavePrimitive(std::ostream &out, Option_t *option /*= ""*/)
1588{
1589 char quote = '"';
1590 out<<" "<<std::endl;
1591 if (gROOT->ClassSaved(TMultiGraph::Class())) {
1592 out<<" ";
1593 } else {
1594 out<<" TMultiGraph *";
1595 }
1596 out<<"multigraph = new TMultiGraph();"<<std::endl;
1597 out<<" multigraph->SetName("<<quote<<GetName()<<quote<<");"<<std::endl;
1598 out<<" multigraph->SetTitle("<<quote<<GetTitle()<<quote<<");"<<std::endl;
1599
1600 if (fGraphs) {
1602 TObject *g;
1603
1604 while (lnk) {
1605 g = lnk->GetObject();
1606 g->SavePrimitive(out, TString::Format("multigraph%s",lnk->GetOption()).Data());
1607 lnk = (TObjOptLink*)lnk->Next();
1608 }
1609 }
1610 const char *l = strstr(option,"th2poly");
1611 if (l) {
1612 out<<" "<<l+7<<"->AddBin(multigraph);"<<std::endl;
1613 } else {
1614 out<<" multigraph->Draw(" <<quote<<option<<quote<<");"<<std::endl;
1615 }
1616 TAxis *xaxis = GetXaxis();
1617 TAxis *yaxis = GetYaxis();
1618
1619 if (xaxis) {
1620 out<<" multigraph->GetXaxis()->SetLimits("<<xaxis->GetXmin()<<", "<<xaxis->GetXmax()<<");"<<std::endl;
1621 xaxis->SaveAttributes(out, "multigraph","->GetXaxis()");
1622 }
1623 if (yaxis) yaxis->SaveAttributes(out, "multigraph","->GetYaxis()");
1624 if (fMinimum != -1111) out<<" multigraph->SetMinimum("<<fMinimum<<");"<<std::endl;
1625 if (fMaximum != -1111) out<<" multigraph->SetMaximum("<<fMaximum<<");"<<std::endl;
1626}
1627
1628
1629////////////////////////////////////////////////////////////////////////////////
1630/// Set multigraph maximum.
1631
1633{
1634 fMaximum = maximum;
1635 if (fHistogram) fHistogram->SetMaximum(maximum);
1636}
1637
1638
1639////////////////////////////////////////////////////////////////////////////////
1640/// Set multigraph minimum.
1641
1643{
1644 fMinimum = minimum;
1645 if (fHistogram) fHistogram->SetMinimum(minimum);
1646}
1647
1648
1649////////////////////////////////////////////////////////////////////////////////
1650/// Get iterator over internal graphs list.
1651
1653{
1654 return TIter(fGraphs);
1655}
#define f(i)
Definition: RSha256.hxx:104
#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
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:375
R__EXTERN TEnv * gEnv
Definition: TEnv.h:170
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t b
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t np
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t g
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:4836
#define gROOT
Definition: TROOT.h:404
void Printf(const char *fmt,...)
Formats a string in a circular formatting buffer and prints the string.
Definition: TString.cxx:2470
#define gPad
Definition: TVirtualPad.h:288
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:236
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:851
virtual Bool_t GetTimeDisplay() const
Definition: TAxis.h:127
const char * GetTitle() const override
Returns title of object.
Definition: TAxis.h:130
Double_t GetXmax() const
Definition: TAxis.h:135
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:162
void SaveAttributes(std::ostream &out, const char *name, const char *subname) override
Save axis attributes as C++ statement(s) on output stream out.
Definition: TAxis.cxx:689
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:155
Double_t GetXmin() const
Definition: TAxis.h:134
virtual const char * GetTimeFormat() const
Definition: TAxis.h:128
virtual void SetTimeFormat(const char *format="")
Change the format used for time plotting.
Definition: TAxis.cxx:1052
virtual void SetRange(Int_t first=0, Int_t last=0)
Set the viewing range for the axis using bin numbers.
Definition: TAxis.cxx:952
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
static TClass * Class()
virtual Int_t GetNpar() const
Definition: TF1.h:486
@ kNotDraw
Definition: TF1.h:326
virtual void SetParLimits(Int_t ipar, Double_t parmin, Double_t parmax)
Set lower and upper limits for parameter ipar.
Definition: TF1.cxx:3536
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:639
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
void Paint(Option_t *chopt="") override
Draw this graph with its current attributes.
Definition: TGraph.cxx:1957
@ kResetHisto
fHistogram must be reset in GetHistogram
Definition: TGraph.h:72
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:1984
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition: TGraph.cxx:1404
virtual void SetHistogram(TH1F *h)
Definition: TGraph.h:180
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:577
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:8804
TAxis * GetZaxis()
Definition: TH1.h:324
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to a line.
Definition: TH1.cxx:2798
void SetTitle(const char *title) override
Change (i.e.
Definition: TH1.cxx:6700
@ kNoStats
Don't draw stats box.
Definition: TH1.h:163
TAxis * GetXaxis()
Definition: TH1.h:322
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:8412
virtual void SetMaximum(Double_t maximum=-1111)
Definition: TH1.h:400
TAxis * GetYaxis()
Definition: TH1.h:323
virtual void SetMinimum(Double_t minimum=-1111)
Definition: TH1.h:401
void Paint(Option_t *option="") override
Control routine to paint any kind of histograms.
Definition: TH1.cxx:6188
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:8502
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8857
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:257
Option_t * GetOption() const
Definition: TCollection.h:253
void Reset()
Definition: TCollection.h:254
Iterator of linked list.
Definition: TList.h:191
Option_t * GetOption() const override
Returns the object option stored in the list.
Definition: TList.cxx:1142
A doubly linked list.
Definition: TList.h:38
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition: TList.cxx:578
void Add(TObject *obj) override
Definition: TList.h:81
TObject * Remove(TObject *obj) override
Remove object from the list.
Definition: TList.cxx:822
TObject * First() const override
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:659
virtual TObjLink * FirstLink() const
Definition: TList.h:102
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:357
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:34
Double_t fMinimum
Minimum value for plotting along y.
Definition: TMultiGraph.h:41
TH1F * fHistogram
Pointer to histogram used for drawing axis.
Definition: TMultiGraph.h:39
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.
void Paint(Option_t *chopt="") override
Paint all the graphs of this multigraph.
TList * GetListOfGraphs() const
Definition: TMultiGraph.h:68
TList * fGraphs
Pointer to list of TGraphs.
Definition: TMultiGraph.h:37
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:40
virtual void SetMinimum(Double_t minimum=-1111)
Set multigraph minimum.
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Compute distance from point px,py to each graph.
void Browse(TBrowser *b) override
Browse multigraph.
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.
void Draw(Option_t *chopt="") override
Draw this multigraph with its current attributes.
static TClass * Class()
virtual void FitPanel()
Display a panel with all histogram fit options.
void RecursiveRemove(TObject *obj) override
Recursively remove this object from a list.
void Print(Option_t *chopt="") const override
Print the list of graphs.
TIter begin() const
Get iterator over internal graphs list.
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.
~TMultiGraph() override
TMultiGraph destructor.
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.
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 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 SetMaximum(Double_t maximum=-1111)
Set multigraph maximum.
void PaintReverse(Option_t *chopt="")
Paint all the graphs of this multigraph reverting values along X and/or Y axis.
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save primitive as a C++ statement(s) on output stream out.
TList * fFunctions
Pointer to list of functions (fits and user)
Definition: TMultiGraph.h:38
TAxis * GetXaxis()
Get x axis of the graph.
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
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
Mother of all ROOT objects.
Definition: TObject.h:41
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:201
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:184
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:775
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:970
virtual void Paint(Option_t *option="")
This method must be overridden if a class wants to paint itself.
Definition: TObject.cxx:602
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition: TObject.h:72
@ kMustCleanup
if object destructor must call RecursiveRemove()
Definition: TObject.h:64
Longptr_t ExecPlugin(int nargs, const T &... params)
Int_t LoadPlugin()
Load the plugin library for this handler.
SCoord_t GetY() const
Definition: TPoint.h:47
SCoord_t GetX() const
Definition: TPoint.h:46
Basic string class.
Definition: TString.h:136
Ssiz_t Length() const
Definition: TString.h:410
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1159
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
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition: TString.cxx:2345
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
static TClass * Class()
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:965
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
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition: TMath.h:754
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition: TMath.h:660
Double_t Log10(Double_t x)
Returns the common (base-10) logarithm of x.
Definition: TMath.h:760
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123
Definition: first.py:1
Definition: graph.py:1
TMarker m
Definition: textangle.C:8
TLine l
Definition: textangle.C:4
TArc a
Definition: textangle.C:12