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