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