Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TSPlot.cxx
Go to the documentation of this file.
1// @(#)root/splot:$Id$
2// Author: Muriel Pivk, Anna Kreshuk 10/2005
3
4/**********************************************************************
5 * *
6 * Copyright (c) 2005 ROOT Foundation, CERN/PH-SFT *
7 * *
8 **********************************************************************/
9
10#include "TSPlot.h"
11#include "TVirtualFitter.h"
12#include "TH1.h"
13#include "TTreePlayer.h"
14#include "TTreeFormula.h"
15#include "TTreeFormulaManager.h"
16#include "TSelectorDraw.h"
17#include "TBrowser.h"
18#include "TClass.h"
19#include "TMath.h"
20#include "snprintf.h"
21
22extern void Yields(Int_t &, Double_t *, Double_t &f, Double_t *x, Int_t iflag);
23
25
26/** \class TSPlot
27
28\legacy{TSPlot}
29
30A common method used in High Energy Physics to perform measurements is
31the maximum Likelihood method, exploiting discriminating variables to
32disentangle signal from background. The crucial point for such an
33analysis to be reliable is to use an exhaustive list of sources of
34events combined with an accurate description of all the Probability
35Density Functions (PDF).
36
37To assess the validity of the fit, a convincing quality check
38is to explore further the data sample by examining the distributions of
39control variables. A control variable can be obtained for instance by
40removing one of the discriminating variables before performing again
41the maximum Likelihood fit: this removed variable is a control
42variable. The expected distribution of this control variable, for
43signal, is to be compared to the one extracted, for signal, from the
44data sample. In order to be able to do so, one must be able to unfold
45from the distribution of the whole data sample.
46
47The TSPlot method allows to reconstruct the distributions for
48the control variable, independently for each of the various sources of
49events, without making use of any <em>a priori</em> knowledge on <u>this</u>
50variable. The aim is thus to use the knowledge available for the
51discriminating variables to infer the behaviour of the individual
52sources of events with respect to the control variable.
53
54TSPlot is optimal if the control variable is uncorrelated with the discriminating variables.
55
56
57A detail description of the formalism itself, called \f$\hbox{$_s$}{\cal P}lot\f$, is given
58in [[1](https://arxiv.org/abs/physics/0402083)].
59
60### The method
61
62
63The \f$\hbox{$_s$}{\cal P}lot\f$ technique is developed in the above context of a
64maximum Likelihood method making use of discriminating variables.
65
66One considers a data sample in which are merged several species
67of events. These species represent various signal components and
68background components which all together account for the data sample.
69The different terms of the log-Likelihood are:
70
71 - \f$N\f$ : the total number of events in the data sample,
72 - \f${\rm N}_{\rm s}\f$ : the number of species of events populating the data sample,
73 - \f$N_i\f$ : the number of events expected on the average for the \f$i^{\rm th}\f$ species,
74 - \f${\rm f}_i(y_e)\f$ : the value of the PDFs of the discriminating variables
75 \f$y\f$ for the\f$i^{th}\f$ species and for event\f$e\f$,
76 - \f$x\f$ : the set of control variables which, by definition, do not appear in
77 the expression of the Likelihood function \f${\cal L}\f$.
78
79The extended log-Likelihood reads:
80
81 \f[
82{\cal L}=\sum_{e=1}^{N}\ln \Big\{ \sum_{i=1}^{{\rm N}_{\rm s}}N_i{\rm f}_i(y_e) \Big\} -\sum_{i=1}^{{\rm N}_{\rm s}}N_i \tag{1}
83\f]
84
85From this expression, after maximization of \f${\cal L}\f$ with respect to the \f$N_i\f$ parameters,
86a weight can be computed for every event and each species, in order to obtain later the true distribution
87\f$\hbox{M}_i(x)\f$ of variable \f$x\f$. If \f${\rm n}\f$ is one of the
88 \f${\rm N}_{\rm s}\f$ species present in the data sample, the weight for this species is defined by:
89
90
91\f[
92\fbox{$
93{_s{\cal P}}_{\rm n}(y_e)={\sum_{j=1}^{{\rm N}_{\rm s}} \hbox{V}_{{\rm n}j}{\rm f}_j(y_e)\over\sum_{k=1}^{{\rm N}_{\rm s}}N_k{\rm f}_k(y_e) } $} , \tag{2}
94\f]
95
96
97where \f$\hbox{V}_{{\rm n}j}\f$
98
99is the covariance matrix resulting from the Likelihood maximization.
100This matrix can be used directly from the fit, but this is numerically
101less accurate than the direct computation:
102
103
104\f[
105\hbox{ V}^{-1}_{{\rm n}j}~=~
106{\partial^2(-{\cal L})\over\partial N_{\rm n}\partial N_j}~=~
107\sum_{e=1}^N {{\rm f}_{\rm n}(y_e){\rm f}_j(y_e)\over(\sum_{k=1}^{{\rm N}_{\rm s}}N_k{\rm f}_k(y_e))^2} . \tag{3}
108\f]
109
110
111The distribution of the control variable \f$x\f$ obtained by histogramming the weighted
112events reproduces, on average, the true distribution
113\f${\hbox{ {M}}}_{\rm n}(x)\f$
114
115The class TSPlot allows to reconstruct the true distribution
116\f${\hbox{ {M}}}_{\rm n}(x)\f$
117
118of a control variable \f$x\f$ for each of the \f${\rm N}_{\rm s}\f$ species from
119the sole knowledge of the PDFs of the discriminating variables \f${\rm f}_i(y)\f$.
120The plots obtained thanks to the TSPlot class are called \f$\hbox {$_s$}{\cal P}lots\f$.
121
122
123### Some properties and checks
124
125
126Beside reproducing the true distribution,\f$\hbox {$_s$}{\cal P}lots\f$ bear remarkable properties:
127
128
129 - Each \f$x\f$ - distribution is properly normalized:
130
131\f[
132\sum_{e=1}^{N} {_s{\cal P}}_{\rm n}(y_e)~=~N_{\rm n} ~. \tag{4}
133\f]
134
135
136 - For any event:
137
138\f[
139\sum_{l=1}^{{\rm N}_{\rm s}} {_s{\cal P}}_l(y_e) ~=~1 ~. \tag{5}
140\f]
141
142That is to say that, summing up the \f${\rm N}_{\rm s}\f$ \f$\hbox {$_s$}{\cal P}lots\f$,
143one recovers the data sample distribution in \f$x\f$, and summing up the number of events
144entering in a \f$\hbox{$_s$}{\cal P}lot\f$ for a given species, one recovers the yield of the
145species, as provided by the fit.
146The property 4 is implemented in the TSPlot class as a check.
147
148 - the sum of the statistical uncertainties per bin
149
150\f[
151\sigma[N_{\rm n}\ _s\tilde{\rm M}_{\rm n}(x) {\delta x}]~=~\sqrt{\sum_{e \subset {\delta x}} ({_s{\cal P}}_{\rm n})^2} ~. \tag{6}
152\f]
153
154reproduces the statistical uncertainty on the yield \f$N_{\rm n}\f$, as provided by the fit:
155\f$\sigma[N_{\rm n}]\equiv\sqrt{\hbox{ V}_{{\rm n}{\rm n}}}\f$ .
156Because of that and since the determination of the yields is optimal
157when obtained using a Likelihood fit, one can conclude that the \f$\hbox{$_s$}{\cal P}lot\f$
158technique is itself an optimal method to reconstruct distributions of control variables.
159
160
161### Different steps followed by TSPlot
162
163
164 1. A maximum Likelihood fit is performed to obtain the yields \f$N_i\f$
165 of the various species.The fit relies on discriminating variables \f$y\f$
166 uncorrelated with a control variable \f$x\f$:
167 the later is therefore totally absent from the fit.
168
169 2. The weights \f${_s{\cal P}}\f$ are calculated using Eq.
170 (2)
171 where the covariance matrix is taken from Minuit.
172
173 3. Histograms of \f$x\f$ are filled by weighting the events with \f${_s{\cal P}}\f$ .
174
175 4. Error bars per bin are given by Eq. (6).
176
177
178The \f$\hbox {$_s$}{\cal P}lots\f$ reproduce the true distributions of the species
179in the control variable \f$x\f$, within the above defined statistical uncertainties.
180
181### Illustrations
182
183
184To illustrate the technique, one considers an example derived from the analysis where
185\f$\hbox {$_s$}{\cal P}lots\f$
186have been first used (charmless B decays). One is dealing with a data
187sample in which two species are present: the first is termed signal and
188the second background. A maximum Likelihood fit is performed to obtain
189the two yields \f$N_1\f$ and \f$N_2\f$ . The fit relies on two discriminating
190variables collectively denoted \f$y\f$ which are chosen within three possible
191variables denoted \f${m_{\rm ES}}\f$ , \f$\Delta E\f$ and \f${\cal F}\f$.
192The variable which is not incorporated in \f$y\f$ is used as the control variable
193\f$x\f$ . The six distributions of the three variables are assumed to be the ones
194depicted in Fig. 1.
195
196
197\image html splot_pdfmesNIM.png Figure 1 width=800
198
199Distributions of the three discriminating variables available to perform the Likelihood fit:
200\f${m_{\rm ES}}\f$ , \f$\Delta E\f$ , \f${\cal F}\f$ .
201Among the three variables, two are used to perform the fit while one is
202kept out of the fit to serve the purpose of a control variable. The
203three distributions on the top (resp. bottom) of the figure correspond
204to the signal (resp. background). The unit of the vertical axis is
205chosen such that it indicates the number of entries per bin, if one
206slices the histograms in 25 bins.
207
208A data sample being built through a Monte Carlo simulation based on the
209distributions shown in Fig.
2101,
211one obtains the three distributions of Fig. 2.
212Whereas the distribution of \f$\Delta E\f$ clearly indicates the presence of the signal,
213the distribution of \f${m_{\rm ES}}\f$ and \f${\cal F}\f$ are less obviously populated by signal.
214
215
216\image html splot_genfiTOTNIM.png Figure 2 width=800
217
218Distributions of the three discriminating variables for signal plus
219background. The three distributions are the ones obtained from a data
220sample obtained through a Monte Carlo simulation based on the
221distributions shown in Fig.
2221.
223The data sample consists of 500 signal events and 5000 background events.
224
225
226Choosing \f$\Delta E\f$ and \f${\cal F}\f$ as discriminating variables to determine
227\f$N_1\f$ and \f$N_2\f$ through a maximum Likelihood fit, one builds, for the control
228variable \f${m_{\rm ES}}\f$ which is unknown to the fit, the two \f$\hbox {$_s$}{\cal P}lots\f$
229for signal and background shown in
230Fig. 3.
231One observes that the \f$\hbox{$_s$}{\cal P}lot\f$
232for signal reproduces correctly the PDF even where the latter vanishes,
233although the error bars remain sizeable. This results from the almost
234complete cancellation between positive and negative weights: the sum of
235weights is close to zero while the sum of weights squared is not. The
236occurrence of negative weights occurs through the appearance of the
237covariance matrix, and its negative components, in the definition of
238Eq. (2).
239
240
241A word of caution is in order with respect to the error bars. Whereas
242their sum in quadrature is identical to the statistical uncertainties
243of the yields determined by the fit, and if, in addition, they are
244asymptotically correct, the error bars should be handled with care for
245low statistics and/or for too fine binning. This is because the error
246bars do not incorporate two known properties of the PDFs: PDFs are
247positive definite and can be non-zero in a given x-bin, even if in the
248particular data sample at hand, no event is observed in this bin. The
249latter limitation is not specific to \f$\hbox {$_s$}{\cal P}lots\f$ ,
250rather it is always present when one is willing to infer the PDF at the
251origin of an histogram, when, for some bins, the number of entries does
252not guaranty the applicability of the Gaussian regime. In such
253situations, a satisfactory practice is to attach allowed ranges to the
254histogram to indicate the upper and lower limits of the PDF value which
255are consistent with the actual observation, at a given confidence
256level.
257
258
259\image html splot_mass-bkg-sPlot.png Figure 3 width=600
260
261The \f$\hbox {$_s$}{\cal P}lots\f$ (signal on top, background on bottom)
262obtained for \f${m_{\rm ES}}\f$ are represented as dots with error bars.
263They are obtained from a fit using only information from \f$\Delta E\f$ and
264\f${\cal F}\f$
265
266<p>
267Choosing \f${m_{\rm ES}}\f$ and \f$\Delta E\f$ as discriminating variables to
268determine \f$N_1\f$ and \f$N_2\f$ through a maximum Likelihood fit, one builds,
269for the control variable \f${\cal F}\f$ which is unknown to the fit, the two
270\f$\hbox {$_s$}{\cal P}lots\f$ for signal and background shown in
271Fig. 4.
272In the \f$\hbox{$_s$}{\cal P}lot\f$ for signal one observes that error bars are
273the largest in the \f$x\f$ regions where the background is the largest.
274
275
276\image html splot_fisher-bkg-sPlot.png Figure 4 width=600
277
278The \f$\hbox {$_s$}{\cal P}lots\f$ (signal on top, background on bottom) obtained
279for \f${\cal F}\f$ are represented as dots with error bars. They are obtained
280from a fit using only information from \f${m_{\rm ES}}\f$ and \f$\Delta E\f$
281
282The results above can be obtained by running the tutorial TestSPlot.C
283*/
284
285
286////////////////////////////////////////////////////////////////////////////////
287/// default constructor (used by I/O only)
288
290 fTree(nullptr),
291 fTreename(nullptr),
292 fVarexp(nullptr),
293 fSelection(nullptr)
294{
295 fNx = 0;
296 fNy=0;
297 fNevents = 0;
298 fNSpecies=0;
299 fNumbersOfEvents=nullptr;
300}
301
302////////////////////////////////////////////////////////////////////////////////
303/// Normal TSPlot constructor
304/// - nx : number of control variables
305/// - ny : number of discriminating variables
306/// - ne : total number of events
307/// - ns : number of species
308/// - tree: input data
309
310TSPlot::TSPlot(Int_t nx, Int_t ny, Int_t ne, Int_t ns, TTree *tree) :
311 fTreename(nullptr),
312 fVarexp(nullptr),
313 fSelection(nullptr)
314
315{
316 fNx = nx;
317 fNy=ny;
318 fNevents = ne;
319 fNSpecies=ns;
320
325 fTree = tree;
326 fNumbersOfEvents = nullptr;
327}
328
329////////////////////////////////////////////////////////////////////////////////
330/// Destructor
331
333{
335 delete [] fNumbersOfEvents;
336 if (!fXvarHists.IsEmpty())
338 if (!fYvarHists.IsEmpty())
340 if (!fYpdfHists.IsEmpty())
342}
343
344////////////////////////////////////////////////////////////////////////////////
345/// To browse the histograms
346
348{
349 if (!fSWeightsHists.IsEmpty()) {
350 TIter next(&fSWeightsHists);
351 TH1D* h = nullptr;
352 while ((h = (TH1D*)next()))
353 b->Add(h,h->GetName());
354 }
355
356 if (!fYpdfHists.IsEmpty()) {
357 TIter next(&fYpdfHists);
358 TH1D* h = nullptr;
359 while ((h = (TH1D*)next()))
360 b->Add(h,h->GetName());
361 }
362 if (!fYvarHists.IsEmpty()) {
363 TIter next(&fYvarHists);
364 TH1D* h = nullptr;
365 while ((h = (TH1D*)next()))
366 b->Add(h,h->GetName());
367 }
368 if (!fXvarHists.IsEmpty()) {
369 TIter next(&fXvarHists);
370 TH1D* h = nullptr;
371 while ((h = (TH1D*)next()))
372 b->Add(h,h->GetName());
373 }
374 b->Add(&fSWeights, "sWeights");
375}
376
377////////////////////////////////////////////////////////////////////////////////
378/// Set the initial number of events of each species - used
379/// as initial estimates in minuit
380
382{
383 if (!fNumbersOfEvents)
385 for (Int_t i=0; i<fNSpecies; i++)
386 fNumbersOfEvents[i]=numbers[i];
387}
388
389////////////////////////////////////////////////////////////////////////////////
390/// Calculates the sWeights
391///
392/// The option controls the print level
393/// - "Q" - no print out
394/// - "V" - prints the estimated `#of events` in species - default
395/// - "VV" - as "V" + the minuit printing + sums of weights for control
396
398{
399
400 if (!fNumbersOfEvents){
401 Error("MakeSPlot","Initial numbers of events in species have not been set");
402 return;
403 }
404 Int_t i, j, ispecies;
405
406 TString opt = option;
407 opt.ToUpper();
408 opt.ReplaceAll("VV", "W");
409
410 //make sure that global fitter is minuit
411 char s[]="TFitter";
413 Int_t strdiff=strcmp(TVirtualFitter::GetFitter()->IsA()->GetName(), s);
414 if (strdiff!=0)
416 }
417
418
419 TVirtualFitter *minuit = TVirtualFitter::Fitter(nullptr, 2);
421
422 //now let's do it, excluding different yvars
423 //for iplot = -1 none is excluded
424 for (Int_t iplot=-1; iplot<fNy; iplot++){
425 for (i=0; i<fNevents; i++){
426 for (ispecies=0; ispecies<fNSpecies; ispecies++){
427 fPdfTot(i, ispecies)=1;
428 for (j=0; j<fNy; j++){
429 if (j!=iplot)
430 fPdfTot(i, ispecies)*=fYpdf(i, ispecies*fNy+j);
431 }
432 }
433 }
434 minuit->Clear();
435 minuit->SetFCN(Yields);
436 Double_t arglist[10];
437 //set the print level
438 if (opt.Contains("Q")||opt.Contains("V")){
439 arglist[0]=-1;
440 }
441 if (opt.Contains("W"))
442 arglist[0]=0;
443 minuit->ExecuteCommand("SET PRINT", arglist, 1);
444
445 minuit->SetObjectFit(&fPdfTot); //a tricky way to get fPdfTot matrix to fcn
446 for (ispecies=0; ispecies<fNSpecies; ispecies++)
447 minuit->SetParameter(ispecies, "", fNumbersOfEvents[ispecies], 1, 0, 0);
448
449 minuit->ExecuteCommand("MIGRAD", arglist, 0);
450 for (ispecies=0; ispecies<fNSpecies; ispecies++){
451 fNumbersOfEvents[ispecies]=minuit->GetParameter(ispecies);
452 if (!opt.Contains("Q"))
453 printf("estimated #of events in species %d = %f\n", ispecies, fNumbersOfEvents[ispecies]);
454 }
455 if (!opt.Contains("Q"))
456 printf("\n");
457 Double_t *covmat = minuit->GetCovarianceMatrix();
458 SPlots(covmat, iplot);
459
460 if (opt.Contains("W")){
461 Double_t *sumweight = new Double_t[fNSpecies];
462 for (i=0; i<fNSpecies; i++){
463 sumweight[i]=0;
464 for (j=0; j<fNevents; j++)
465 sumweight[i]+=fSWeights(j, (iplot+1)*fNSpecies + i);
466 printf("checking sum of weights[%d]=%f\n", i, sumweight[i]);
467 }
468 printf("\n");
469 delete [] sumweight;
470 }
471 }
472}
473
474////////////////////////////////////////////////////////////////////////////////
475/// Computes the sWeights from the covariance matrix
476
477void TSPlot::SPlots(Double_t *covmat, Int_t i_excl)
478{
479 Int_t i, ispecies, k;
480 Double_t numerator, denominator;
481 for (i=0; i<fNevents; i++){
482 denominator=0;
483 for (ispecies=0; ispecies<fNSpecies; ispecies++)
484 denominator+=fNumbersOfEvents[ispecies]*fPdfTot(i, ispecies);
485 for (ispecies=0; ispecies<fNSpecies; ispecies++){
486 numerator=0;
487 for (k=0; k<fNSpecies; k++)
488 numerator+=covmat[ispecies*fNSpecies+k]*fPdfTot(i, k);
489 fSWeights(i, (i_excl+1)*fNSpecies + ispecies)=numerator/denominator;
490 }
491 }
492
493}
494
495////////////////////////////////////////////////////////////////////////////////
496/// Returns the matrix of sweights
497
499{
500 if (weights.GetNcols()!=fNSpecies*(fNy+1) || weights.GetNrows()!=fNevents)
501 weights.ResizeTo(fNevents, fNSpecies*(fNy+1));
502 weights = fSWeights;
503}
504
505////////////////////////////////////////////////////////////////////////////////
506/// Returns the matrix of sweights. It is assumed that the array passed in the
507/// argurment is big enough
508
510{
511 for (Int_t i=0; i<fNevents; i++){
512 for (Int_t j=0; j<fNSpecies; j++){
513 weights[i*fNSpecies+j]=fSWeights(i, j);
514 }
515 }
516}
517
518////////////////////////////////////////////////////////////////////////////////
519/// Fills the histograms of x variables (not weighted) with nbins
520
522{
523 Int_t i, j;
524
525 if (!fXvarHists.IsEmpty()){
526 if (((TH1D*)fXvarHists.First())->GetNbinsX()!=nbins)
528 else
529 return;
530 }
531
532 //make the histograms
533 char name[12];
534 for (i=0; i<fNx; i++){
535 snprintf(name,sizeof(name), "x%d", i);
536 TH1D *h = new TH1D(name, name, nbins, fMinmax(0, i), fMinmax(1, i));
537 for (j=0; j<fNevents; j++)
538 h->Fill(fXvar(j, i));
540 }
541
542}
543
544////////////////////////////////////////////////////////////////////////////////
545/// Returns the array of histograms of x variables (not weighted).
546/// If histograms have not already
547/// been filled, they are filled with default binning 100.
548
550{
551 Int_t nbins = 100;
552 if (fXvarHists.IsEmpty())
553 FillXvarHists(nbins);
554 else if (((TH1D*)fXvarHists.First())->GetNbinsX()!=nbins)
555 FillXvarHists(nbins);
556 return &fXvarHists;
557}
558
559////////////////////////////////////////////////////////////////////////////////
560///Returns the histogram of variable ixvar.
561/// If histograms have not already
562/// been filled, they are filled with default binning 100.
563
565{
566 Int_t nbins = 100;
567 if (fXvarHists.IsEmpty())
568 FillXvarHists(nbins);
569 else if (((TH1D*)fXvarHists.First())->GetNbinsX()!=nbins)
570 FillXvarHists(nbins);
571
572 return (TH1D*)fXvarHists.UncheckedAt(ixvar);
573}
574
575////////////////////////////////////////////////////////////////////////////////
576/// Fill the histograms of y variables
577
579{
580 Int_t i, j;
581
582 if (!fYvarHists.IsEmpty()){
583 if (((TH1D*)fYvarHists.First())->GetNbinsX()!=nbins)
585 else
586 return;
587 }
588
589 //make the histograms
590 char name[12];
591 for (i=0; i<fNy; i++){
592 snprintf(name,sizeof(name), "y%d", i);
593 TH1D *h=new TH1D(name, name, nbins, fMinmax(0, fNx+i), fMinmax(1, fNx+i));
594 for (j=0; j<fNevents; j++)
595 h->Fill(fYvar(j, i));
597 }
598}
599
600////////////////////////////////////////////////////////////////////////////////
601/// Returns the array of histograms of y variables. If histograms have not already
602/// been filled, they are filled with default binning 100.
603
605{
606 Int_t nbins = 100;
607 if (fYvarHists.IsEmpty())
608 FillYvarHists(nbins);
609 else if (((TH1D*)fYvarHists.First())->GetNbinsX()!=nbins)
610 FillYvarHists(nbins);
611 return &fYvarHists;
612}
613
614////////////////////////////////////////////////////////////////////////////////
615/// Returns the histogram of variable iyvar.If histograms have not already
616/// been filled, they are filled with default binning 100.
617
619{
620 Int_t nbins = 100;
621 if (fYvarHists.IsEmpty())
622 FillYvarHists(nbins);
623 else if (((TH1D*)fYvarHists.First())->GetNbinsX()!=nbins)
624 FillYvarHists(nbins);
625 return (TH1D*)fYvarHists.UncheckedAt(iyvar);
626}
627
628////////////////////////////////////////////////////////////////////////////////
629/// Fills the histograms of pdf-s of y variables with binning nbins
630
632{
633 Int_t i, j, ispecies;
634
635 if (!fYpdfHists.IsEmpty()){
636 if (((TH1D*)fYpdfHists.First())->GetNbinsX()!=nbins)
638 else
639 return;
640 }
641
642 char name[34];
643 for (ispecies=0; ispecies<fNSpecies; ispecies++){
644 for (i=0; i<fNy; i++){
645 snprintf(name,sizeof(name), "pdf_species%d_y%d", ispecies, i);
646 //TH1D *h = new TH1D(name, name, nbins, ypdfmin[ispecies*fNy+i], ypdfmax[ispecies*fNy+i]);
647 TH1D *h = new TH1D(name, name, nbins, fMinmax(0, fNx+fNy+ispecies*fNy+i), fMinmax(1, fNx+fNy+ispecies*fNy+i));
648 for (j=0; j<fNevents; j++)
649 h->Fill(fYpdf(j, ispecies*fNy+i));
651 }
652 }
653}
654
655////////////////////////////////////////////////////////////////////////////////
656/// Returns the array of histograms of pdf's of y variables with binning nbins.
657/// If histograms have not already
658/// been filled, they are filled with default binning 100.
659
661{
662 Int_t nbins = 100;
663 if (fYpdfHists.IsEmpty())
664 FillYpdfHists(nbins);
665
666 return &fYpdfHists;
667}
668
669////////////////////////////////////////////////////////////////////////////////
670/// Returns the histogram of the pdf of variable iyvar for species `#ispecies`, binning nbins.
671/// If histograms have not already
672/// been filled, they are filled with default binning 100.
673
675{
676 Int_t nbins = 100;
677 if (fYpdfHists.IsEmpty())
678 FillYpdfHists(nbins);
679
680 return (TH1D*)fYpdfHists.UncheckedAt(fNy*ispecies+iyvar);
681}
682
683////////////////////////////////////////////////////////////////////////////////
684/// The order of histograms in the array:
685///
686/// x0_species0, x0_species1,..., x1_species0, x1_species1,..., y0_species0, y0_species1,...
687///
688/// If the histograms have already been filled with a different binning, they are refilled
689/// and all histograms are deleted
690
692{
693 if (fSWeights.GetNoElements()==0){
694 Error("GetSWeightsHists", "SWeights were not computed");
695 return;
696 }
697
698 if (!fSWeightsHists.IsEmpty()){
699 if (((TH1D*)fSWeightsHists.First())->GetNbinsX()!=nbins)
701 else
702 return;
703 }
704
705 char name[30];
706
707 //Fill histograms of x-variables weighted with sWeights
708 for (Int_t ivar=0; ivar<fNx; ivar++){
709 for (Int_t ispecies=0; ispecies<fNSpecies; ispecies++){
710 snprintf(name,30, "x%d_species%d", ivar, ispecies);
711 TH1D *h = new TH1D(name, name, nbins, fMinmax(0, ivar), fMinmax(1, ivar));
712 h->Sumw2();
713 for (Int_t ievent=0; ievent<fNevents; ievent++)
714 h->Fill(fXvar(ievent, ivar), fSWeights(ievent, ispecies));
716 }
717 }
718
719 //Fill histograms of y-variables (excluded from the fit), weighted with sWeights
720 for (Int_t iexcl=0; iexcl<fNy; iexcl++){
721 for(Int_t ispecies=0; ispecies<fNSpecies; ispecies++){
722 snprintf(name,30, "y%d_species%d", iexcl, ispecies);
723 TH1D *h = new TH1D(name, name, nbins, fMinmax(0, fNx+iexcl), fMinmax(1, fNx+iexcl));
724 h->Sumw2();
725 for (Int_t ievent=0; ievent<fNevents; ievent++)
726 h->Fill(fYvar(ievent, iexcl), fSWeights(ievent, fNSpecies*(iexcl+1)+ispecies));
728 }
729 }
730}
731
732////////////////////////////////////////////////////////////////////////////////
733/// Returns an array of all histograms of variables, weighted with sWeights.
734/// If histograms have not been already filled, they are filled with default binning 50
735/// The order of histograms in the array:
736///
737/// x0_species0, x0_species1,..., x1_species0, x1_species1,..., y0_species0, y0_species1,...
738
740{
741 Int_t nbins = 50; //default binning
743 FillSWeightsHists(nbins);
744
745 return &fSWeightsHists;
746}
747
748////////////////////////////////////////////////////////////////////////////////
749/// The Fill...Hist() methods fill the histograms with the real limits on the variables
750/// This method allows to refill the specified histogram with user-set boundaries min and max
751///
752///Parameters:
753///
754/// - type = 1 - histogram of x variable `#%nvar`
755/// - type = 2 - histogram of y variable `#%nvar`
756/// - type = 3 - histogram of y_pdf for y `#%nvar` and species #%nspecies`
757/// - type = 4 - histogram of x variable `#%nvar`, species `#%nspecies`, WITH sWeights
758/// - type = 5 - histogram of y variable `#%nvar`, species `#%nspecies`, WITH sWeights
759
760void TSPlot::RefillHist(Int_t type, Int_t nvar, Int_t nbins, Double_t min, Double_t max, Int_t nspecies)
761{
762 if (type<1 || type>5){
763 Error("RefillHist", "type must lie between 1 and 5");
764 return;
765 }
766 char name[20];
767 Int_t j;
768 TH1D *hremove;
769 if (type==1){
770 hremove = (TH1D*)fXvarHists.RemoveAt(nvar);
771 delete hremove;
772 snprintf(name,20,"x%d",nvar);
773 TH1D *h = new TH1D(name, name, nbins, min, max);
774 for (j=0; j<fNevents;j++)
775 h->Fill(fXvar(j, nvar));
776 fXvarHists.AddAt(h, nvar);
777 }
778 if (type==2){
779 hremove = (TH1D*)fYvarHists.RemoveAt(nvar);
780 delete hremove;
781 snprintf(name,20, "y%d", nvar);
782 TH1D *h = new TH1D(name, name, nbins, min, max);
783 for (j=0; j<fNevents;j++)
784 h->Fill(fYvar(j, nvar));
785 fXvarHists.AddAt(h, nvar);
786 }
787 if (type==3){
788 hremove = (TH1D*)fYpdfHists.RemoveAt(nspecies*fNy+nvar);
789 delete hremove;
790 snprintf(name,20, "pdf_species%d_y%d", nspecies, nvar);
791 TH1D *h=new TH1D(name, name, nbins, min, max);
792 for (j=0; j<fNevents; j++)
793 h->Fill(fYpdf(j, nspecies*fNy+nvar));
794 fYpdfHists.AddAt(h, nspecies*fNy+nvar);
795 }
796 if (type==4){
797 hremove = (TH1D*)fSWeightsHists.RemoveAt(fNSpecies*nvar+nspecies);
798 delete hremove;
799 snprintf(name,20, "x%d_species%d", nvar, nspecies);
800 TH1D *h = new TH1D(name, name, nbins, min, max);
801 h->Sumw2();
802 for (Int_t ievent=0; ievent<fNevents; ievent++)
803 h->Fill(fXvar(ievent, nvar), fSWeights(ievent, nspecies));
804 fSWeightsHists.AddAt(h, fNSpecies*nvar+nspecies);
805 }
806 if (type==5){
807 hremove = (TH1D*)fSWeightsHists.RemoveAt(fNx*fNSpecies + fNSpecies*nvar+nspecies);
808 delete hremove;
809 snprintf(name,20, "y%d_species%d", nvar, nspecies);
810 TH1D *h = new TH1D(name, name, nbins, min, max);
811 h->Sumw2();
812 for (Int_t ievent=0; ievent<fNevents; ievent++)
813 h->Fill(fYvar(ievent, nvar), fSWeights(ievent, nspecies));
814 fSWeightsHists.AddAt(h, fNx*fNSpecies + fNSpecies*nvar+nspecies);
815 }
816}
817////////////////////////////////////////////////////////////////////////////////
818/// Returns the histogram of a variable, weighted with sWeights.
819/// - If histograms have not been already filled, they are filled with default binning 50
820/// - If parameter ixvar!=-1, the histogram of x-variable ixvar is returned for species ispecies
821/// - If parameter ixvar==-1, the histogram of y-variable iyexcl is returned for species ispecies
822/// - If the histogram has already been filled and the binning is different from the parameter nbins
823/// all histograms with old binning will be deleted and refilled.
824
826{
827
828 Int_t nbins = 50; //default binning
830 FillSWeightsHists(nbins);
831
832 if (ixvar==-1)
833 return (TH1D*)fSWeightsHists.UncheckedAt(fNx*fNSpecies + fNSpecies*iyexcl+ispecies);
834 else
835 return (TH1D*)fSWeightsHists.UncheckedAt(fNSpecies*ixvar + ispecies);
836
837}
838
839
840////////////////////////////////////////////////////////////////////////////////
841/// Set the input Tree
842
844{
845 fTree = tree;
846}
847
848////////////////////////////////////////////////////////////////////////////////
849///Specifies the variables from the tree to be used for splot
850///
851///Variables fNx, fNy, fNSpecies and fNEvents should already be set!
852///
853///In the 1st parameter it is assumed that first fNx variables are x(control variables),
854///then fNy y variables (discriminating variables),
855///then fNy*fNSpecies ypdf variables (probability distribution functions of discriminating
856///variables for different species). The order of pdfs should be: species0_y0, species0_y1,...
857///species1_y0, species1_y1,...species[fNSpecies-1]_y0...
858///The 2nd parameter allows to make a cut
859///TTree::Draw method description contains more details on specifying expression and selection
860
861void TSPlot::SetTreeSelection(const char* varexp, const char *selection, Long64_t firstentry)
862{
863 TTreeFormula **var;
864 std::vector<TString> cnames;
865 TList formulaList;
866 TSelectorDraw *selector = (TSelectorDraw*)(((TTreePlayer*)fTree->GetPlayer())->GetSelector());
867
868 Long64_t entry, entryNumber;
869 Int_t i,nch;
870 Int_t ncols;
871 TObjArray *leaves = fTree->GetListOfLeaves();
872
874 if (varexp)
875 fVarexp = new TString(varexp);
876 if (selection)
877 fSelection= new TString(selection);
878
879 nch = varexp ? strlen(varexp) : 0;
880
881
882//*-*- Compile selection expression if there is one
883 TTreeFormula *select = nullptr;
884 if (selection && strlen(selection)) {
885 select = new TTreeFormula("Selection",selection,fTree);
886 if (!select) return;
887 if (!select->GetNdim()) { delete select; return; }
888 formulaList.Add(select);
889 }
890//*-*- if varexp is empty, take first nx + ny + ny*nspecies columns by default
891
892 if (nch == 0) {
893 ncols = fNx + fNy + fNy*fNSpecies;
894 for (i=0;i<ncols;i++) {
895 cnames.push_back( leaves->At(i)->GetName() );
896 }
897//*-*- otherwise select only the specified columns
898 } else {
899 ncols = selector->SplitNames(varexp,cnames);
900 }
901 var = new TTreeFormula* [ncols];
902 Double_t *xvars = new Double_t[ncols];
903
904 fMinmax.ResizeTo(2, ncols);
905 for (i=0; i<ncols; i++){
906 fMinmax(0, i)=1e30;
907 fMinmax(1, i)=-1e30;
908 }
909
910//*-*- Create the TreeFormula objects corresponding to each column
911 for (i=0;i<ncols;i++) {
912 var[i] = new TTreeFormula("Var1",cnames[i].Data(),fTree);
913 formulaList.Add(var[i]);
914 }
915
916//*-*- Create a TreeFormulaManager to coordinate the formulas
917 TTreeFormulaManager *manager = nullptr;
918 if (formulaList.LastIndex()>=0) {
919 manager = new TTreeFormulaManager;
920 for(i=0;i<=formulaList.LastIndex();i++) {
921 manager->Add((TTreeFormula*)formulaList.At(i));
922 }
923 manager->Sync();
924 formulaList.Clear();
925 }
926//*-*- loop on all selected entries
927 // fSelectedRows = 0;
928 Int_t tnumber = -1;
929 Long64_t selectedrows=0;
930 for (entry=firstentry;entry<firstentry+fNevents;entry++) {
931 entryNumber = fTree->GetEntryNumber(entry);
932 if (entryNumber < 0) break;
933 Long64_t localEntry = fTree->LoadTree(entryNumber);
934 if (localEntry < 0) break;
935 if (tnumber != fTree->GetTreeNumber()) {
936 tnumber = fTree->GetTreeNumber();
937 if (manager) manager->UpdateFormulaLeaves();
938 }
939 Int_t ndata = 1;
940 if (manager && manager->GetMultiplicity()) {
941 ndata = manager->GetNdata();
942 }
943
944 for(Int_t inst=0;inst<ndata;inst++) {
945 Bool_t loaded = kFALSE;
946 if (select) {
947 if (select->EvalInstance(inst) == 0) {
948 continue;
949 }
950 }
951
952 if (inst==0) loaded = kTRUE;
953 else if (!loaded) {
954 // EvalInstance(0) always needs to be called so that
955 // the proper branches are loaded.
956 for (i=0;i<ncols;i++) {
957 var[i]->EvalInstance(0);
958 }
959 loaded = kTRUE;
960 }
961
962 for (i=0;i<ncols;i++) {
963 xvars[i] = var[i]->EvalInstance(inst);
964 }
965
966 // curentry = entry-firstentry;
967 //printf("event#%d\n", curentry);
968 //for (i=0; i<ncols; i++)
969 // printf("xvars[%d]=%f\n", i, xvars[i]);
970 //selectedrows++;
971 for (i=0; i<fNx; i++){
972 fXvar(selectedrows, i) = xvars[i];
973 if (fXvar(selectedrows, i) < fMinmax(0, i))
974 fMinmax(0, i)=fXvar(selectedrows, i);
975 if (fXvar(selectedrows, i) > fMinmax(1, i))
976 fMinmax(1, i)=fXvar(selectedrows, i);
977 }
978 for (i=0; i<fNy; i++){
979 fYvar(selectedrows, i) = xvars[i+fNx];
980 //printf("y_in_loop(%d, %d)=%f, xvars[%d]=%f\n", selectedrows, i, fYvar(selectedrows, i), i+fNx, xvars[i+fNx]);
981 if (fYvar(selectedrows, i) < fMinmax(0, i+fNx))
982 fMinmax(0, i+fNx) = fYvar(selectedrows, i);
983 if (fYvar(selectedrows, i) > fMinmax(1, i+fNx))
984 fMinmax(1, i+fNx) = fYvar(selectedrows, i);
985 for (Int_t j=0; j<fNSpecies; j++){
986 fYpdf(selectedrows, j*fNy + i)=xvars[j*fNy + i+fNx+fNy];
987 if (fYpdf(selectedrows, j*fNy+i) < fMinmax(0, j*fNy+i+fNx+fNy))
988 fMinmax(0, j*fNy+i+fNx+fNy) = fYpdf(selectedrows, j*fNy+i);
989 if (fYpdf(selectedrows, j*fNy+i) > fMinmax(1, j*fNy+i+fNx+fNy))
990 fMinmax(1, j*fNy+i+fNx+fNy) = fYpdf(selectedrows, j*fNy+i);
991 }
992 }
993 selectedrows++;
994 }
995 }
996 fNevents=selectedrows;
997 // for (i=0; i<fNevents; i++){
998 // printf("event#%d\n", i);
999 //for (Int_t iy=0; iy<fNy; iy++)
1000 // printf("y[%d]=%f\n", iy, fYvar(i, iy));
1001 //for (Int_t ispecies=0; ispecies<fNSpecies; ispecies++){
1002 // for (Int_t iy=0; iy<fNy; iy++)
1003 // printf("ypdf[sp. %d, y %d]=%f\n", ispecies, iy, fYpdf(i, ispecies*fNy+iy));
1004 // }
1005 //}
1006
1007 delete [] xvars;
1008 delete [] var;
1009}
1010
1011////////////////////////////////////////////////////////////////////////////////
1012/// FCN-function for Minuit
1013
1014void Yields(Int_t &, Double_t *, Double_t &f, Double_t *x, Int_t /*iflag*/)
1015{
1016 Double_t likelihood;
1017 Int_t i, ispecies;
1018
1020 TMatrixD *pdftot = (TMatrixD*)fitter->GetObjectFit();
1021 Int_t nev = pdftot->GetNrows();
1022 Int_t nes = pdftot->GetNcols();
1023 f=0;
1024 for (i=0; i<nev; i++){
1025 likelihood=0;
1026 for (ispecies=0; ispecies<nes; ispecies++)
1027 likelihood+=x[ispecies]*(*pdftot)(i, ispecies);
1028 if (likelihood<0) likelihood=1;
1029 f+=TMath::Log(likelihood);
1030 }
1031 //extended likelihood, equivalent to chi2
1032 Double_t ntot=0;
1033 for (i=0; i<nes; i++)
1034 ntot += x[i];
1035 f = -2*(f-ntot);
1036}
1037
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
constexpr Bool_t kFALSE
Definition RtypesCore.h:101
double Double_t
Definition RtypesCore.h:59
long long Long64_t
Definition RtypesCore.h:80
constexpr Bool_t kTRUE
Definition RtypesCore.h:100
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:377
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 winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char name[80]
Definition TGX11.cxx:110
void Yields(Int_t &, Double_t *, Double_t &f, Double_t *x, Int_t iflag)
FCN-function for Minuit.
Definition TSPlot.cxx:1014
#define snprintf
Definition civetweb.c:1540
virtual Int_t GetNdim() const
Definition TFormula.h:237
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:664
virtual Bool_t Add(TF1 *h1, Double_t c1=1, Option_t *option="")
Performs the operation: this = this + c1*f1 if errors are defined (see TH1::Sumw2),...
Definition TH1.cxx:824
A doubly linked list.
Definition TList.h:38
void Clear(Option_t *option="") override
Remove all objects from the list.
Definition TList.cxx:400
void Add(TObject *obj) override
Definition TList.h:81
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:355
Int_t GetNrows() const
Int_t GetNoElements() const
Int_t GetNcols() const
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
An array of TObjects.
Definition TObjArray.h:31
void AddAt(TObject *obj, Int_t idx) override
Add object at position ids.
void Delete(Option_t *option="") override
Remove all objects from the array AND delete all heap based objects.
TObject * At(Int_t idx) const override
Definition TObjArray.h:164
void AddLast(TObject *obj) override
Add object in the next empty slot in the array.
TObject * UncheckedAt(Int_t i) const
Definition TObjArray.h:84
TObject * First() const override
Return the object in the first slot.
Bool_t IsEmpty() const override
Definition TObjArray.h:65
TObject * RemoveAt(Int_t idx) override
Remove object at index idx.
void Add(TObject *obj) override
Definition TObjArray.h:68
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:439
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:976
<div class="legacybox"><h2>Legacy Code</h2> TSPlot is a legacy interface: there will be no bug fixes ...
Definition TSPlot.h:21
void FillYvarHists(Int_t nbins=100)
Fill the histograms of y variables.
Definition TSPlot.cxx:578
TObjArray fXvarHists
Definition TSPlot.h:30
TObjArray fSWeightsHists
Definition TSPlot.h:33
void RefillHist(Int_t type, Int_t var, Int_t nbins, Double_t min, Double_t max, Int_t nspecies=-1)
The Fill...Hist() methods fill the histograms with the real limits on the variables This method allow...
Definition TSPlot.cxx:760
~TSPlot() override
Destructor.
Definition TSPlot.cxx:332
void Browse(TBrowser *b) override
To browse the histograms.
Definition TSPlot.cxx:347
TMatrixD fSWeights
Definition TSPlot.h:28
Int_t fNevents
Definition TSPlot.h:44
Int_t fNy
Definition TSPlot.h:42
Int_t fNSpecies
Definition TSPlot.h:43
TString * fTreename
Definition TSPlot.h:36
void SPlots(Double_t *covmat, Int_t i_excl)
Computes the sWeights from the covariance matrix.
Definition TSPlot.cxx:477
void GetSWeights(TMatrixD &weights)
Returns the matrix of sweights.
Definition TSPlot.cxx:498
TMatrixD fMinmax
Definition TSPlot.h:27
void FillSWeightsHists(Int_t nbins=50)
The order of histograms in the array:
Definition TSPlot.cxx:691
Int_t fNx
Definition TSPlot.h:41
TString * fSelection
Definition TSPlot.h:38
TH1D * GetSWeightsHist(Int_t ixvar, Int_t ispecies, Int_t iyexcl=-1)
Returns the histogram of a variable, weighted with sWeights.
Definition TSPlot.cxx:825
TMatrixD fPdfTot
Definition TSPlot.h:26
void MakeSPlot(Option_t *option="v")
Calculates the sWeights.
Definition TSPlot.cxx:397
TH1D * GetYpdfHist(Int_t iyvar, Int_t ispecies)
Returns the histogram of the pdf of variable iyvar for species #ispecies, binning nbins.
Definition TSPlot.cxx:674
TObjArray * GetSWeightsHists()
Returns an array of all histograms of variables, weighted with sWeights.
Definition TSPlot.cxx:739
TTree * fTree
Definition TSPlot.h:35
TObjArray * GetYpdfHists()
Returns the array of histograms of pdf's of y variables with binning nbins.
Definition TSPlot.cxx:660
void FillXvarHists(Int_t nbins=100)
Fills the histograms of x variables (not weighted) with nbins.
Definition TSPlot.cxx:521
TObjArray * GetXvarHists()
Returns the array of histograms of x variables (not weighted).
Definition TSPlot.cxx:549
TObjArray fYvarHists
Definition TSPlot.h:31
TString * fVarexp
Definition TSPlot.h:37
void SetTree(TTree *tree)
Set the input Tree.
Definition TSPlot.cxx:843
TObjArray * GetYvarHists()
Returns the array of histograms of y variables.
Definition TSPlot.cxx:604
void FillYpdfHists(Int_t nbins=100)
Fills the histograms of pdf-s of y variables with binning nbins.
Definition TSPlot.cxx:631
void SetTreeSelection(const char *varexp="", const char *selection="", Long64_t firstentry=0)
Specifies the variables from the tree to be used for splot.
Definition TSPlot.cxx:861
TMatrixD fXvar
Definition TSPlot.h:23
TSPlot()
default constructor (used by I/O only)
Definition TSPlot.cxx:289
TClass * IsA() const override
Definition TSPlot.h:90
TH1D * GetYvarHist(Int_t iyvar)
Returns the histogram of variable iyvar.If histograms have not already been filled,...
Definition TSPlot.cxx:618
void SetInitialNumbersOfSpecies(Int_t *numbers)
Set the initial number of events of each species - used as initial estimates in minuit.
Definition TSPlot.cxx:381
TObjArray fYpdfHists
Definition TSPlot.h:32
Double_t * fNumbersOfEvents
Definition TSPlot.h:46
TMatrixD fYpdf
Definition TSPlot.h:25
TH1D * GetXvarHist(Int_t ixvar)
Returns the histogram of variable ixvar.
Definition TSPlot.cxx:564
TMatrixD fYvar
Definition TSPlot.h:24
A specialized TSelector for TTree::Draw.
virtual UInt_t SplitNames(const TString &varexp, std::vector< TString > &names)
Build Index array for names in varexp.
Int_t LastIndex() const
Basic string class.
Definition TString.h:139
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:704
void ToUpper()
Change string to upper case.
Definition TString.cxx:1195
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
Used to coordinate one or more TTreeFormula objects.
virtual void UpdateFormulaLeaves()
This function could be called TTreePlayer::UpdateFormulaLeaves, itself called by TChain::LoadTree whe...
virtual Int_t GetNdata(bool forceLoadDim=false)
Return number of available instances in the formulas.
virtual void Add(TTreeFormula *)
Add a new formula to the list of formulas managed The manager of the formula will be changed and the ...
virtual bool Sync()
Synchronize all the formulae.
virtual Int_t GetMultiplicity() const
Used to pass a selection expression to the Tree drawing routine.
T EvalInstance(Int_t i=0, const char *stringStack[]=nullptr)
Evaluate this treeformula.
Implement some of the functionality of the class TTree requiring access to extra libraries (Histogram...
Definition TTreePlayer.h:37
A TTree represents a columnar dataset.
Definition TTree.h:79
virtual TObjArray * GetListOfLeaves()
Definition TTree.h:489
TVirtualTreePlayer * GetPlayer()
Load the TTreePlayer (if not already done).
Definition TTree.cxx:6305
virtual Long64_t GetEntryNumber(Long64_t entry) const
Return entry number corresponding to entry.
Definition TTree.cxx:5865
virtual Long64_t LoadTree(Long64_t entry)
Set current entry.
Definition TTree.cxx:6473
virtual Int_t GetTreeNumber() const
Definition TTree.h:519
Abstract Base Class for Fitting.
virtual void SetObjectFit(TObject *obj)
virtual TObject * GetObjectFit() const
virtual Double_t * GetCovarianceMatrix() const =0
void Clear(Option_t *option="") override=0
Set name and title to empty strings ("").
virtual void SetFCN(void(*fcn)(Int_t &, Double_t *, Double_t &f, Double_t *, Int_t))
To set the address of the minimization objective function called by the native compiler (see function...
virtual Int_t SetParameter(Int_t ipar, const char *parname, Double_t value, Double_t verr, Double_t vlow, Double_t vhigh)=0
virtual Int_t ExecuteCommand(const char *command, Double_t *args, Int_t nargs)=0
virtual Double_t GetParameter(Int_t ipar) const =0
static TVirtualFitter * GetFitter()
static: return the current Fitter
static TVirtualFitter * Fitter(TObject *obj, Int_t maxpar=25)
Static function returning a pointer to the current fitter.
Double_t x[n]
Definition legend1.C:17
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition TMath.h:756