Logo ROOT   6.18/05
Reference Guide
Tools.cxx
Go to the documentation of this file.
1// @(#)root/tmva $Id$
2// Author: Andreas Hoecker, Peter Speckmayer, Joerg Stelzer, Helge Voss, Jan Therhaag
3
4/**********************************************************************************
5 * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6 * Package: TMVA *
7 * Class : Tools *
8 * Web : http://tmva.sourceforge.net *
9 * *
10 * Description: *
11 * Implementation (see header for description) *
12 * *
13 * Authors (alphabetical): *
14 * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15 * Peter Speckmayer <Peter.Speckmayer@cern.ch> - CERN, Switzerland *
16 * Jan Therhaag <Jan.Therhaag@cern.ch> - U of Bonn, Germany *
17 * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
18 * Kai Voss <Kai.Voss@cern.ch> - U. of Victoria, Canada *
19 * *
20 * Copyright (c) 2005-2011: *
21 * CERN, Switzerland *
22 * U. of Victoria, Canada *
23 * MPI-K Heidelberg, Germany *
24 * U. of Bonn, Germany *
25 * *
26 * Redistribution and use in source and binary forms, with or without *
27 * modification, are permitted according to the terms listed in LICENSE *
28 * (http://tmva.sourceforge.net/LICENSE) *
29 **********************************************************************************/
30
31/*! \class TMVA::Tools
32\ingroup TMVA
33Global auxiliary applications and data treatment routines.
34*/
35
36#include "TMVA/Tools.h"
37
38#include "TMVA/Config.h"
39#include "TMVA/Event.h"
40#include "TMVA/Version.h"
41#include "TMVA/PDF.h"
42#include "TMVA/MsgLogger.h"
43#include "TMVA/Types.h"
44
45#include "TObjString.h"
46#include "TMath.h"
47#include "TString.h"
48#include "TTree.h"
49#include "TLeaf.h"
50#include "TH1.h"
51#include "TH2.h"
52#include "TList.h"
53#include "TSpline.h"
54#include "TVector.h"
55#include "TMatrixD.h"
56#include "TMatrixDSymEigen.h"
57#include "TVectorD.h"
58#include "TTreeFormula.h"
59#include "TXMLEngine.h"
60#include "TROOT.h"
61#include "TMatrixDSymEigen.h"
62
63#include <algorithm>
64#include <cstdlib>
65
66using namespace std;
67
68#if __cplusplus > 199711L
69std::atomic<TMVA::Tools*> TMVA::Tools::fgTools{0};
70#else
72#endif
73
76#if __cplusplus > 199711L
77 if(!fgTools) {
78 Tools* tmp = new Tools();
79 Tools* expected = 0;
80 if(! fgTools.compare_exchange_strong(expected,tmp)) {
81 //another thread beat us
82 delete tmp;
83 }
84 }
85 return *fgTools;
86#else
87 return fgTools?*(fgTools): *(fgTools = new Tools());
88#endif
89}
91 //NOTE: there is no thread safe way to do this so
92 // one must only call this method ones in an executable
93#if __cplusplus > 199711L
94 if (fgTools != 0) { delete fgTools.load(); fgTools=0; }
95#else
96 if (fgTools != 0) { delete fgTools; fgTools=0; }
97#endif
98}
99
100////////////////////////////////////////////////////////////////////////////////
101/// constructor
102
104 fRegexp("$&|!%^&()'<>?= "),
105 fLogger(new MsgLogger("Tools")),
106 fXMLEngine(new TXMLEngine())
107{
108}
109
110////////////////////////////////////////////////////////////////////////////////
111/// destructor
112
114{
115 delete fLogger;
116 delete fXMLEngine;
117}
118
119////////////////////////////////////////////////////////////////////////////////
120/// normalise to output range: [-1, 1]
121
123{
124 return 2*(x - xmin)/(xmax - xmin) - 1.0;
125}
126
127////////////////////////////////////////////////////////////////////////////////
128/// compute "separation" defined as
129/// \f[
130/// <s2> = \frac{1}{2} \int_{-\infty}^{+\infty} \frac{(S(x) - B(x))^2}{(S(x) + B(x))} dx
131/// \f]
132
134{
135 Double_t separation = 0;
136
137 // sanity checks
138 // signal and background histograms must have same number of bins and
139 // same limits
140 if ((S->GetNbinsX() != B->GetNbinsX()) || (S->GetNbinsX() <= 0)) {
141 Log() << kFATAL << "<GetSeparation> signal and background"
142 << " histograms have different number of bins: "
143 << S->GetNbinsX() << " : " << B->GetNbinsX() << Endl;
144 }
145
146 if (S->GetXaxis()->GetXmin() != B->GetXaxis()->GetXmin() ||
147 S->GetXaxis()->GetXmax() != B->GetXaxis()->GetXmax() ||
148 S->GetXaxis()->GetXmax() <= S->GetXaxis()->GetXmin()) {
149 Log() << kINFO << S->GetXaxis()->GetXmin() << " " << B->GetXaxis()->GetXmin()
150 << " " << S->GetXaxis()->GetXmax() << " " << B->GetXaxis()->GetXmax()
151 << " " << S->GetXaxis()->GetXmax() << " " << S->GetXaxis()->GetXmin() << Endl;
152 Log() << kFATAL << "<GetSeparation> signal and background"
153 << " histograms have different or invalid dimensions:" << Endl;
154 }
155
156 Int_t nstep = S->GetNbinsX();
157 Double_t intBin = (S->GetXaxis()->GetXmax() - S->GetXaxis()->GetXmin())/nstep;
158 Double_t nS = S->GetSumOfWeights()*intBin;
159 Double_t nB = B->GetSumOfWeights()*intBin;
160
161 if (nS > 0 && nB > 0) {
162 for (Int_t bin=0; bin<nstep; bin++) {
163 Double_t s = S->GetBinContent( bin+1 )/Double_t(nS);
164 Double_t b = B->GetBinContent( bin+1 )/Double_t(nB);
165 // separation
166 if (s + b > 0) separation += (s - b)*(s - b)/(s + b);
167 }
168 separation *= (0.5*intBin);
169 }
170 else {
171 Log() << kWARNING << "<GetSeparation> histograms with zero entries: "
172 << nS << " : " << nB << " cannot compute separation"
173 << Endl;
174 separation = 0;
175 }
176
177 return separation;
178}
179
180////////////////////////////////////////////////////////////////////////////////
181/// compute "separation" defined as
182/// \f[
183/// <s2> = \frac{1}{2} \int_{-\infty}^{+\infty} \frac{(S(x) - B(x))^2}{(S(x) + B(x))} dx
184/// \f]
185
186Double_t TMVA::Tools::GetSeparation( const PDF& pdfS, const PDF& pdfB ) const
187{
188 Double_t xmin = pdfS.GetXmin();
189 Double_t xmax = pdfS.GetXmax();
190 // sanity check
191 if (xmin != pdfB.GetXmin() || xmax != pdfB.GetXmax()) {
192 Log() << kFATAL << "<GetSeparation> Mismatch in PDF limits: "
193 << xmin << " " << pdfB.GetXmin() << xmax << " " << pdfB.GetXmax() << Endl;
194 }
195
196 Double_t separation = 0;
197 Int_t nstep = 100;
198 Double_t intBin = (xmax - xmin)/Double_t(nstep);
199 for (Int_t bin=0; bin<nstep; bin++) {
200 Double_t x = (bin + 0.5)*intBin + xmin;
201 Double_t s = pdfS.GetVal( x );
202 Double_t b = pdfB.GetVal( x );
203 // separation
204 if (s + b > 0) separation += (s - b)*(s - b)/(s + b);
205 }
206 separation *= (0.5*intBin);
207
208 return separation;
209}
210
211////////////////////////////////////////////////////////////////////////////////
212/// sanity check
213
214void TMVA::Tools::ComputeStat( const std::vector<TMVA::Event*>& events, std::vector<Float_t>* valVec,
215 Double_t& meanS, Double_t& meanB,
216 Double_t& rmsS, Double_t& rmsB,
218 Int_t signalClass, Bool_t norm )
219{
220 if (0 == valVec)
221 Log() << kFATAL << "<Tools::ComputeStat> value vector is zero pointer" << Endl;
222
223 if ( events.size() != valVec->size() )
224 Log() << kWARNING << "<Tools::ComputeStat> event and value vector have different lengths "
225 << events.size() << "!=" << valVec->size() << Endl;
226
227 Long64_t entries = valVec->size();
228
229 // first fill signal and background in arrays before analysis
230 Double_t* varVecS = new Double_t[entries];
231 Double_t* varVecB = new Double_t[entries];
232 Double_t* wgtVecS = new Double_t[entries];
233 Double_t* wgtVecB = new Double_t[entries];
234 xmin = +DBL_MAX;
235 xmax = -DBL_MAX;
236 Long64_t nEventsS = 0;
237 Long64_t nEventsB = 0;
238 Double_t xmin_ = 0, xmax_ = 0;
239
240 if (norm) {
241 xmin_ = *std::min( valVec->begin(), valVec->end() );
242 xmax_ = *std::max( valVec->begin(), valVec->end() );
243 }
244
245 for (Int_t ievt=0; ievt<entries; ievt++) {
246 Double_t theVar = (*valVec)[ievt];
247 if (norm) theVar = Tools::NormVariable( theVar, xmin_, xmax_ );
248
249 if (Int_t(events[ievt]->GetClass()) == signalClass ){
250 wgtVecS[nEventsS] = events[ievt]->GetWeight(); // this is signal
251 varVecS[nEventsS++] = theVar; // this is signal
252 }
253 else {
254 wgtVecB[nEventsB] = events[ievt]->GetWeight(); // this is signal
255 varVecB[nEventsB++] = theVar; // this is background
256 }
257
258 if (theVar > xmax) xmax = theVar;
259 if (theVar < xmin) xmin = theVar;
260 }
261 // ++nEventsS;
262 // ++nEventsB;
263
264 // basic statistics
265 // !!! TMath::Mean allows for weights, but NOT for negative weights
266 // and TMath::RMS doesn't allow for weights all together...
267 meanS = TMVA::Tools::Mean( nEventsS, varVecS, wgtVecS );
268 meanB = TMVA::Tools::Mean( nEventsB, varVecB, wgtVecB );
269 rmsS = TMVA::Tools::RMS ( nEventsS, varVecS, wgtVecS );
270 rmsB = TMVA::Tools::RMS ( nEventsB, varVecB, wgtVecB );
271
272 delete [] varVecS;
273 delete [] varVecB;
274 delete [] wgtVecS;
275 delete [] wgtVecB;
276}
277
278////////////////////////////////////////////////////////////////////////////////
279/// square-root of symmetric matrix
280/// of course the resulting sqrtMat is also symmetric, but it's easier to
281/// treat it as a general matrix
282
284{
285 Int_t n = symMat->GetNrows();
286
287 // compute eigenvectors
288 TMatrixDSymEigen* eigen = new TMatrixDSymEigen( *symMat );
289
290 // D = ST C S
291 TMatrixD* si = new TMatrixD( eigen->GetEigenVectors() );
292 TMatrixD* s = new TMatrixD( *si ); // copy
293 si->Transpose( *si ); // invert (= transpose)
294
295 // diagonal matrices
296 TMatrixD* d = new TMatrixD( n, n);
297 d->Mult( (*si), (*symMat) ); (*d) *= (*s);
298
299 // sanity check: matrix must be diagonal and positive definit
300 Int_t i, j;
301 Double_t epsilon = 1.0e-8;
302 for (i=0; i<n; i++) {
303 for (j=0; j<n; j++) {
304 if ((i != j && TMath::Abs((*d)(i,j))/TMath::Sqrt((*d)(i,i)*(*d)(j,j)) > epsilon) ||
305 (i == j && (*d)(i,i) < 0)) {
306 //d->Print();
307 Log() << kWARNING << "<GetSQRootMatrix> error in matrix diagonalization; printed S and B" << Endl;
308 }
309 }
310 }
311
312 // make exactly diagonal
313 for (i=0; i<n; i++) for (j=0; j<n; j++) if (j != i) (*d)(i,j) = 0;
314
315 // compute the square-root C' of covariance matrix: C = C'*C'
316 for (i=0; i<n; i++) (*d)(i,i) = TMath::Sqrt((*d)(i,i));
317
318 TMatrixD* sqrtMat = new TMatrixD( n, n );
319 sqrtMat->Mult( (*s), (*d) );
320 (*sqrtMat) *= (*si);
321
322 // invert square-root matrices
323 sqrtMat->Invert();
324
325 delete eigen;
326 delete s;
327 delete si;
328 delete d;
329
330 return sqrtMat;
331}
332
333////////////////////////////////////////////////////////////////////////////////
334/// turns covariance into correlation matrix
335
337{
338
339 if (covMat == 0) return 0;
340 // sanity check
341 Int_t nvar = covMat->GetNrows();
342 if (nvar != covMat->GetNcols())
343 Log() << kFATAL << "<GetCorrelationMatrix> input matrix not quadratic" << Endl;
344
345 Log() << kWARNING;
346 TMatrixD* corrMat = new TMatrixD( nvar, nvar );
347 for (Int_t ivar=0; ivar<nvar; ivar++) {
348 for (Int_t jvar=0; jvar<nvar; jvar++) {
349 if (ivar != jvar) {
350 Double_t d = (*covMat)(ivar, ivar)*(*covMat)(jvar, jvar);
351 if (d > 1E-20) {
352 (*corrMat)(ivar, jvar) = (*covMat)(ivar, jvar)/TMath::Sqrt(d);
353 } else {
354 Log() << "<GetCorrelationMatrix> zero variances for variables "
355 << "(" << ivar << ", " << jvar << ")" << Endl;
356 (*corrMat)(ivar, jvar) = 0;
357 }
358 if (TMath::Abs( (*corrMat)(ivar,jvar)) > 1){
359 Log() << kWARNING
360 << " Element corr("<<ivar<<","<<ivar<<")=" << (*corrMat)(ivar,jvar)
361 << " sigma2="<<d
362 << " cov("<<ivar<<","<<ivar<<")=" <<(*covMat)(ivar, ivar)
363 << " cov("<<jvar<<","<<jvar<<")=" <<(*covMat)(jvar, jvar)
364 << Endl;
365
366 }
367 }
368 else (*corrMat)(ivar, ivar) = 1.0;
369 }
370 }
371 Log() << Endl;
372 return corrMat;
373}
374
375////////////////////////////////////////////////////////////////////////////////
376/// projects variable from tree into normalised histogram
377
378TH1* TMVA::Tools::projNormTH1F( TTree* theTree, const TString& theVarName,
379 const TString& name, Int_t nbins,
380 Double_t xmin, Double_t xmax, const TString& cut )
381{
382 // needed because of ROOT bug (feature) that excludes events that have value == xmax
383 xmax += 0.00001;
384
385 TH1* hist = new TH1F( name, name, nbins, xmin, xmax );
386 hist->Sumw2(); // enable quadratic errors
387 theTree->Project( name, theVarName, cut );
388 NormHist( hist );
389 return hist;
390}
391
392////////////////////////////////////////////////////////////////////////////////
393/// normalises histogram
394
396{
397 if (!theHist) return 0;
398
399 if (theHist->GetSumw2N() == 0) theHist->Sumw2();
400 if (theHist->GetSumOfWeights() != 0) {
401 Double_t w = ( theHist->GetSumOfWeights()
402 *(theHist->GetXaxis()->GetXmax() - theHist->GetXaxis()->GetXmin())/theHist->GetNbinsX() );
403 if (w > 0) theHist->Scale( norm/w );
404 return w;
405 }
406
407 return 1.0;
408}
409
410////////////////////////////////////////////////////////////////////////////////
411/// Parse the string and cut into labels separated by ":"
412
413TList* TMVA::Tools::ParseFormatLine( TString formatString, const char* sep )
414{
415 TList* labelList = new TList();
416 labelList->SetOwner();
417 while (formatString.First(sep)==0) formatString.Remove(0,1); // remove initial separators
418
419 while (formatString.Length()>0) {
420 if (formatString.First(sep) == -1) { // no more separator
421 labelList->Add(new TObjString(formatString.Data()));
422 formatString="";
423 break;
424 }
425
426 Ssiz_t posSep = formatString.First(sep);
427 labelList->Add(new TObjString(TString(formatString(0,posSep)).Data()));
428 formatString.Remove(0,posSep+1);
429
430 while (formatString.First(sep)==0) formatString.Remove(0,1); // remove additional separators
431
432 }
433 return labelList;
434}
435
436////////////////////////////////////////////////////////////////////////////////
437/// parse option string for ANN methods
438/// default settings (should be defined in theOption string)
439///
440/// format and syntax of option string: "3000:N:N+2:N-3:6"
441///
442/// where:
443/// - 3000 - number of training cycles (epochs)
444/// - N - number of nodes in first hidden layer, where N is the number
445/// of discriminating variables used (note that the first ANN
446/// layer necessarily has N nodes, and hence is not given).
447/// - N+2 - number of nodes in 2nd hidden layer (2 nodes more than
448/// number of variables)
449/// - N-3 - number of nodes in 3rd hidden layer (3 nodes less than
450/// number of variables)
451/// - 6 - 6 nodes in last (4th) hidden layer (note that the last ANN
452/// layer in MVA has 2 nodes, each one for signal and background
453/// classes)
454
455vector<Int_t>* TMVA::Tools::ParseANNOptionString( TString theOptions, Int_t nvar,
456 vector<Int_t>* nodes )
457{
458 TList* list = TMVA::Tools::ParseFormatLine( theOptions, ":" );
459
460
461 // sanity check
462 if (list->GetSize() < 1) {
463 Log() << kFATAL << "<ParseANNOptionString> unrecognized option string: " << theOptions << Endl;
464 }
465
466 // add number of cycles
467 nodes->push_back( atoi( ((TObjString*)list->At(0))->GetString() ) );
468
469 Int_t a;
470 if (list->GetSize() > 1) {
471 for (Int_t i=1; i<list->GetSize(); i++) {
472 TString s = ((TObjString*)list->At(i))->GetString();
473 s.ToUpper();
474 if (s(0) == 'N') {
475 if (s.Length() > 1) nodes->push_back( nvar + atoi(&s[1]) );
476 else nodes->push_back( nvar );
477 }
478 else if ((a = atoi( s )) > 0) nodes->push_back( atoi(s ) );
479 else {
480 Log() << kFATAL << "<ParseANNOptionString> unrecognized option string: " << theOptions << Endl;
481 }
482 }
483 }
484
485 return nodes;
486}
487
488////////////////////////////////////////////////////////////////////////////////
489/// check quality of splining by comparing splines and histograms in each bin
490
491Bool_t TMVA::Tools::CheckSplines( const TH1* theHist, const TSpline* theSpline )
492{
493 const Double_t sanityCrit = 0.01; // relative deviation
494
495 Bool_t retval = kTRUE;
496 for (Int_t ibin=1; ibin<=theHist->GetNbinsX(); ibin++) {
497 Double_t x = theHist->GetBinCenter( ibin );
498 Double_t yh = theHist->GetBinContent( ibin ); // the histogram output
499 Double_t ys = theSpline->Eval( x ); // the spline output
500
501 if (ys + yh > 0) {
502 Double_t dev = 0.5*(ys - yh)/(ys + yh);
503 if (TMath::Abs(dev) > sanityCrit) {
504 Log() << kFATAL << "<CheckSplines> Spline failed sanity criterion; "
505 << " relative deviation from histogram: " << dev
506 << " in (bin, value): (" << ibin << ", " << x << ")" << Endl;
507 retval = kFALSE;
508 }
509 }
510 }
511
512 return retval;
513}
514
515////////////////////////////////////////////////////////////////////////////////
516/// computes difference between two vectors
517
518std::vector<Double_t> TMVA::Tools::MVADiff( std::vector<Double_t>& a, std::vector<Double_t>& b )
519{
520 if (a.size() != b.size()) {
521 throw;
522 }
523 vector<Double_t> result(a.size());
524 for (UInt_t i=0; i<a.size();i++) result[i]=a[i]-b[i];
525 return result;
526}
527
528////////////////////////////////////////////////////////////////////////////////
529/// scales double vector
530
531void TMVA::Tools::Scale( std::vector<Double_t>& v, Double_t f )
532{
533 for (UInt_t i=0; i<v.size();i++) v[i]*=f;
534}
535
536////////////////////////////////////////////////////////////////////////////////
537/// scales float vector
538
539void TMVA::Tools::Scale( std::vector<Float_t>& v, Float_t f )
540{
541 for (UInt_t i=0; i<v.size();i++) v[i]*=f;
542}
543
544////////////////////////////////////////////////////////////////////////////////
545/// sort 2D vector (AND in parallel a TString vector) in such a way
546/// that the "first vector is sorted" and the other vectors are reshuffled
547/// in the same way as necessary to have the first vector sorted.
548/// I.e. the correlation between the elements is kept.
549
550void TMVA::Tools::UsefulSortAscending( std::vector<vector<Double_t> >& v, std::vector<TString>* vs ){
551 UInt_t nArrays=v.size();
552 Double_t temp;
553 if (nArrays > 0) {
554 UInt_t sizeofarray=v[0].size();
555 for (UInt_t i=0; i<sizeofarray; i++) {
556 for (UInt_t j=sizeofarray-1; j>i; j--) {
557 if (v[0][j-1] > v[0][j]) {
558 for (UInt_t k=0; k< nArrays; k++) {
559 temp = v[k][j-1]; v[k][j-1] = v[k][j]; v[k][j] = temp;
560 }
561 if (nullptr != vs) {
562 TString temps = (*vs)[j-1]; (*vs)[j-1] = (*vs)[j]; (*vs)[j] = temps;
563 }
564 }
565 }
566 }
567 }
568}
569
570////////////////////////////////////////////////////////////////////////////////
571/// sort 2D vector (AND in parallel a TString vector) in such a way
572/// that the "first vector is sorted" and the other vectors are reshuffled
573/// in the same way as necessary to have the first vector sorted.
574/// I.e. the correlation between the elements is kept.
575
576void TMVA::Tools::UsefulSortDescending( std::vector<std::vector<Double_t> >& v, std::vector<TString>* vs )
577{
578 UInt_t nArrays=v.size();
579 Double_t temp;
580 if (nArrays > 0) {
581 UInt_t sizeofarray=v[0].size();
582 for (UInt_t i=0; i<sizeofarray; i++) {
583 for (UInt_t j=sizeofarray-1; j>i; j--) {
584 if (v[0][j-1] < v[0][j]) {
585 for (UInt_t k=0; k< nArrays; k++) {
586 temp = v[k][j-1]; v[k][j-1] = v[k][j]; v[k][j] = temp;
587 }
588 if (nullptr != vs) {
589 TString temps = (*vs)[j-1]; (*vs)[j-1] = (*vs)[j]; (*vs)[j] = temps;
590 }
591 }
592 }
593 }
594 }
595}
596
597////////////////////////////////////////////////////////////////////////////////
598/// Mutual Information method for non-linear correlations estimates in 2D histogram
599/// Author: Moritz Backes, Geneva (2009)
600
602{
603 Double_t hi = h_.Integral();
604 if (hi == 0) return -1;
605
606 // copy histogram and rebin to speed up procedure
607 TH2F h( h_ );
608 h.RebinX(2);
609 h.RebinY(2);
610
611 Double_t mutualInfo = 0.;
612 Int_t maxBinX = h.GetNbinsX();
613 Int_t maxBinY = h.GetNbinsY();
614 for (Int_t x = 1; x <= maxBinX; x++) {
615 for (Int_t y = 1; y <= maxBinY; y++) {
616 Double_t p_xy = h.GetBinContent(x,y)/hi;
617 Double_t p_x = h.Integral(x,x,1,maxBinY)/hi;
618 Double_t p_y = h.Integral(1,maxBinX,y,y)/hi;
619 if (p_x > 0. && p_y > 0. && p_xy > 0.){
620 mutualInfo += p_xy*TMath::Log(p_xy / (p_x * p_y));
621 }
622 }
623 }
624
625 return mutualInfo;
626}
627
628////////////////////////////////////////////////////////////////////////////////
629/// Compute Correlation Ratio of 2D histogram to estimate functional dependency between two variables
630/// Author: Moritz Backes, Geneva (2009)
631
633{
634 Double_t hi = h_.Integral();
635 if (hi == 0.) return -1;
636
637 // copy histogram and rebin to speed up procedure
638 TH2F h( h_ );
639 h.RebinX(2);
640 h.RebinY(2);
641
642 Double_t corrRatio = 0.;
643 Double_t y_mean = h.ProjectionY()->GetMean();
644 for (Int_t ix=1; ix<=h.GetNbinsX(); ix++) {
645 corrRatio += (h.Integral(ix,ix,1,h.GetNbinsY())/hi)*pow((GetYMean_binX(h,ix)-y_mean),2);
646 }
647 corrRatio /= pow(h.ProjectionY()->GetRMS(),2);
648 return corrRatio;
649}
650
651////////////////////////////////////////////////////////////////////////////////
652/// Compute the mean in Y for a given bin X of a 2D histogram
653
655{
656 if (h.Integral(bin_x,bin_x,1,h.GetNbinsY()) == 0.) {return 0;}
657 Double_t y_bin_mean = 0.;
658 TH1* py = h.ProjectionY();
659 for (Int_t y = 1; y <= h.GetNbinsY(); y++){
660 y_bin_mean += h.GetBinContent(bin_x,y)*py->GetBinCenter(y);
661 }
662 y_bin_mean /= h.Integral(bin_x,bin_x,1,h.GetNbinsY());
663 return y_bin_mean;
664}
665
666////////////////////////////////////////////////////////////////////////////////
667/// Transpose quadratic histogram
668
670{
671 // sanity check
672 if (h.GetNbinsX() != h.GetNbinsY()) {
673 Log() << kFATAL << "<TransposeHist> cannot transpose non-quadratic histogram" << Endl;
674 }
675
676 TH2F *transposedHisto = new TH2F( h );
677 for (Int_t ix=1; ix <= h.GetNbinsX(); ix++){
678 for (Int_t iy=1; iy <= h.GetNbinsY(); iy++){
679 transposedHisto->SetBinContent(iy,ix,h.GetBinContent(ix,iy));
680 }
681 }
682
683 // copy stats (thanks to Swagato Banerjee for pointing out the missing stats information)
684 Double_t stats_old[7];
685 Double_t stats_new[7];
686
687 h.GetStats(stats_old);
688 stats_new[0] = stats_old[0];
689 stats_new[1] = stats_old[1];
690 stats_new[2] = stats_old[4];
691 stats_new[3] = stats_old[5];
692 stats_new[4] = stats_old[2];
693 stats_new[5] = stats_old[3];
694 stats_new[6] = stats_old[6];
695 transposedHisto->PutStats(stats_new);
696
697 return transposedHisto; // ownership returned
698}
699
700////////////////////////////////////////////////////////////////////////////////
701/// check for "silence" option in configuration option string
702
704{
705 Bool_t isSilent = kFALSE;
706
707 TString s( cs );
708 s.ToLower();
709 s.ReplaceAll(" ","");
710 if (s.Contains("silent") && !s.Contains("silent=f")) {
711 if (!s.Contains("!silent") || s.Contains("silent=t")) isSilent = kTRUE;
712 }
713
714 return isSilent;
715}
716
717////////////////////////////////////////////////////////////////////////////////
718/// check if verbosity "V" set in option
719
721{
722 Bool_t isVerbose = kFALSE;
723
724 TString s( cs );
725 s.ToLower();
726 s.ReplaceAll(" ","");
727 std::vector<TString> v = SplitString( s, ':' );
728 for (std::vector<TString>::iterator it = v.begin(); it != v.end(); ++it) {
729 if ((*it == "v" || *it == "verbose") && !it->Contains("!")) isVerbose = kTRUE;
730 }
731
732 return isVerbose;
733}
734
735////////////////////////////////////////////////////////////////////////////////
736/// sort vector
737
738void TMVA::Tools::UsefulSortDescending( std::vector<Double_t>& v )
739{
740 vector< vector<Double_t> > vtemp;
741 vtemp.push_back(v);
742 UsefulSortDescending(vtemp);
743 v = vtemp[0];
744}
745
746////////////////////////////////////////////////////////////////////////////////
747/// sort vector
748
749void TMVA::Tools::UsefulSortAscending( std::vector<Double_t>& v )
750{
751 vector<vector<Double_t> > vtemp;
752 vtemp.push_back(v);
753 UsefulSortAscending(vtemp);
754 v = vtemp[0];
755}
756
757////////////////////////////////////////////////////////////////////////////////
758/// find index of maximum entry in vector
759
760Int_t TMVA::Tools::GetIndexMaxElement( std::vector<Double_t>& v )
761{
762 if (v.empty()) return -1;
763
764 Int_t pos=0; Double_t mx=v[0];
765 for (UInt_t i=0; i<v.size(); i++){
766 if (v[i] > mx){
767 mx=v[i];
768 pos=i;
769 }
770 }
771 return pos;
772}
773
774////////////////////////////////////////////////////////////////////////////////
775/// find index of minimum entry in vector
776
777Int_t TMVA::Tools::GetIndexMinElement( std::vector<Double_t>& v )
778{
779 if (v.empty()) return -1;
780
781 Int_t pos=0; Double_t mn=v[0];
782 for (UInt_t i=0; i<v.size(); i++){
783 if (v[i] < mn){
784 mn=v[i];
785 pos=i;
786 }
787 }
788 return pos;
789}
790
791
792////////////////////////////////////////////////////////////////////////////////
793/// check if regular expression
794/// helper function to search for "$!%^&()'<>?= " in a string
795
797{
798 Bool_t regular = kFALSE;
799 for (Int_t i = 0; i < Tools::fRegexp.Length(); i++)
800 if (s.Contains( Tools::fRegexp[i] )) { regular = kTRUE; break; }
801
802 return regular;
803}
804
805////////////////////////////////////////////////////////////////////////////////
806/// replace regular expressions
807/// helper function to remove all occurrences "$!%^&()'<>?= " from a string
808/// and replace all ::,$,*,/,+,- with _M_,_S_,_T_,_D_,_P_,_M_ respectively
809
811{
812 TString snew = s;
813 for (Int_t i = 0; i < Tools::fRegexp.Length(); i++)
814 snew.ReplaceAll( Tools::fRegexp[i], r );
815
816 snew.ReplaceAll( "::", r );
817 snew.ReplaceAll( "$", "_S_" );
818 snew.ReplaceAll( "&", "_A_" );
819 snew.ReplaceAll( "%", "_MOD_" );
820 snew.ReplaceAll( "|", "_O_" );
821 snew.ReplaceAll( "*", "_T_" );
822 snew.ReplaceAll( "/", "_D_" );
823 snew.ReplaceAll( "+", "_P_" );
824 snew.ReplaceAll( "-", "_M_" );
825 snew.ReplaceAll( " ", "_" );
826 snew.ReplaceAll( "[", "_" );
827 snew.ReplaceAll( "]", "_" );
828 snew.ReplaceAll( "=", "_E_" );
829 snew.ReplaceAll( ">", "_GT_" );
830 snew.ReplaceAll( "<", "_LT_" );
831 snew.ReplaceAll( "(", "_" );
832 snew.ReplaceAll( ")", "_" );
833
834 return snew;
835}
836
837////////////////////////////////////////////////////////////////////////////////
838/// human readable color strings
839
841{
842 static const TString gClr_none = "" ;
843 static const TString gClr_white = "\033[1;37m"; // white
844 static const TString gClr_black = "\033[30m"; // black
845 static const TString gClr_blue = "\033[34m"; // blue
846 static const TString gClr_red = "\033[1;31m" ; // red
847 static const TString gClr_yellow = "\033[1;33m"; // yellow
848 static const TString gClr_darkred = "\033[31m"; // dark red
849 static const TString gClr_darkgreen = "\033[32m"; // dark green
850 static const TString gClr_darkyellow = "\033[33m"; // dark yellow
851
852 static const TString gClr_bold = "\033[1m" ; // bold
853 static const TString gClr_black_b = "\033[30m" ; // bold black
854 static const TString gClr_lblue_b = "\033[1;34m" ; // bold light blue
855 static const TString gClr_cyan_b = "\033[0;36m" ; // bold cyan
856 static const TString gClr_lgreen_b = "\033[1;32m"; // bold light green
857
858 static const TString gClr_blue_bg = "\033[44m"; // blue background
859 static const TString gClr_red_bg = "\033[1;41m"; // white on red background
860 static const TString gClr_whiteonblue = "\033[1;44m"; // white on blue background
861 static const TString gClr_whiteongreen = "\033[1;42m"; // white on green background
862 static const TString gClr_grey_bg = "\033[47m"; // grey background
863
864 static const TString gClr_reset = "\033[0m"; // reset
865
866 if (!gConfig().UseColor()) return gClr_none;
867
868 if (c == "white" ) return gClr_white;
869 if (c == "blue" ) return gClr_blue;
870 if (c == "black" ) return gClr_black;
871 if (c == "lightblue") return gClr_cyan_b;
872 if (c == "yellow") return gClr_yellow;
873 if (c == "red" ) return gClr_red;
874 if (c == "dred" ) return gClr_darkred;
875 if (c == "dgreen") return gClr_darkgreen;
876 if (c == "lgreenb") return gClr_lgreen_b;
877 if (c == "dyellow") return gClr_darkyellow;
878
879 if (c == "bold") return gClr_bold;
880 if (c == "bblack") return gClr_black_b;
881
882 if (c == "blue_bgd") return gClr_blue_bg;
883 if (c == "red_bgd" ) return gClr_red_bg;
884
885 if (c == "white_on_blue" ) return gClr_whiteonblue;
886 if (c == "white_on_green") return gClr_whiteongreen;
887
888 if (c == "reset") return gClr_reset;
889
890 std::cout << "Unknown color " << c << std::endl;
891 exit(1);
892
893 return gClr_none;
894}
895
896////////////////////////////////////////////////////////////////////////////////
897/// formatted output of simple table
898
899void TMVA::Tools::FormattedOutput( const std::vector<Double_t>& values, const std::vector<TString>& V,
900 const TString titleVars, const TString titleValues, MsgLogger& logger,
901 TString format )
902{
903 // sanity check
904 UInt_t nvar = V.size();
905 if ((UInt_t)values.size() != nvar) {
906 logger << kFATAL << "<FormattedOutput> fatal error with dimensions: "
907 << values.size() << " OR " << " != " << nvar << Endl;
908 }
909
910 // find maximum length in V (and column title)
911 UInt_t maxL = 7;
912 std::vector<UInt_t> vLengths;
913 for (UInt_t ivar=0; ivar<nvar; ivar++) maxL = TMath::Max( (UInt_t)V[ivar].Length(), maxL );
914 maxL = TMath::Max( (UInt_t)titleVars.Length(), maxL );
915
916 // column length
917 UInt_t maxV = 7;
918 maxV = TMath::Max( (UInt_t)titleValues.Length() + 1, maxL );
919
920 // full column length
921 UInt_t clen = maxL + maxV + 3;
922
923 // bar line
924 for (UInt_t i=0; i<clen; i++) logger << "-";
925 logger << Endl;
926
927 // title bar
928 logger << setw(maxL) << titleVars << ":";
929 logger << setw(maxV+1) << titleValues << ":";
930 logger << Endl;
931 for (UInt_t i=0; i<clen; i++) logger << "-";
932 logger << Endl;
933
934 // the numbers
935 for (UInt_t irow=0; irow<nvar; irow++) {
936 logger << setw(maxL) << V[irow] << ":";
937 logger << setw(maxV+1) << Form( format.Data(), values[irow] );
938 logger << Endl;
939 }
940
941 // bar line
942 for (UInt_t i=0; i<clen; i++) logger << "-";
943 logger << Endl;
944}
945
946////////////////////////////////////////////////////////////////////////////////
947/// formatted output of matrix (with labels)
948
949void TMVA::Tools::FormattedOutput( const TMatrixD& M, const std::vector<TString>& V, MsgLogger& logger )
950{
951 // sanity check: matrix must be quadratic
952 UInt_t nvar = V.size();
953 if ((UInt_t)M.GetNcols() != nvar || (UInt_t)M.GetNrows() != nvar) {
954 logger << kFATAL << "<FormattedOutput> fatal error with dimensions: "
955 << M.GetNcols() << " OR " << M.GetNrows() << " != " << nvar << " ==> abort" << Endl;
956 }
957
958 // get length of each variable, and maximum length
959 UInt_t minL = 7;
960 UInt_t maxL = minL;
961 std::vector<UInt_t> vLengths;
962 for (UInt_t ivar=0; ivar<nvar; ivar++) {
963 vLengths.push_back(TMath::Max( (UInt_t)V[ivar].Length(), minL ));
964 maxL = TMath::Max( vLengths.back(), maxL );
965 }
966
967 // count column length
968 UInt_t clen = maxL+1;
969 for (UInt_t icol=0; icol<nvar; icol++) clen += vLengths[icol]+1;
970
971 // bar line
972 for (UInt_t i=0; i<clen; i++) logger << "-";
973 logger << Endl;
974
975 // title bar
976 logger << setw(maxL+1) << " ";
977 for (UInt_t icol=0; icol<nvar; icol++) logger << setw(vLengths[icol]+1) << V[icol];
978 logger << Endl;
979
980 // the numbers
981 for (UInt_t irow=0; irow<nvar; irow++) {
982 logger << setw(maxL) << V[irow] << ":";
983 for (UInt_t icol=0; icol<nvar; icol++) {
984 logger << setw(vLengths[icol]+1) << Form( "%+1.3f", M(irow,icol) );
985 }
986 logger << Endl;
987 }
988
989 // bar line
990 for (UInt_t i=0; i<clen; i++) logger << "-";
991 logger << Endl;
992}
993
994////////////////////////////////////////////////////////////////////////////////
995/// formatted output of matrix (with labels)
996
998 const std::vector<TString>& vert, const std::vector<TString>& horiz,
999 MsgLogger& logger )
1000{
1001 // sanity check: matrix must be quadratic
1002 UInt_t nvvar = vert.size();
1003 UInt_t nhvar = horiz.size();
1004
1005 // get length of each variable, and maximum length
1006 UInt_t minL = 7;
1007 UInt_t maxL = minL;
1008 std::vector<UInt_t> vLengths;
1009 for (UInt_t ivar=0; ivar<nvvar; ivar++) {
1010 vLengths.push_back(TMath::Max( (UInt_t)vert[ivar].Length(), minL ));
1011 maxL = TMath::Max( vLengths.back(), maxL );
1012 }
1013
1014 // count column length
1015 UInt_t minLh = 7;
1016 UInt_t maxLh = minLh;
1017 std::vector<UInt_t> hLengths;
1018 for (UInt_t ivar=0; ivar<nhvar; ivar++) {
1019 hLengths.push_back(TMath::Max( (UInt_t)horiz[ivar].Length(), minL ));
1020 maxLh = TMath::Max( hLengths.back(), maxLh );
1021 }
1022
1023 UInt_t clen = maxLh+1;
1024 for (UInt_t icol=0; icol<nhvar; icol++) clen += hLengths[icol]+1;
1025
1026 // bar line
1027 for (UInt_t i=0; i<clen; i++) logger << "-";
1028 logger << Endl;
1029
1030 // title bar
1031 logger << setw(maxL+1) << " ";
1032 for (UInt_t icol=0; icol<nhvar; icol++) logger << setw(hLengths[icol]+1) << horiz[icol];
1033 logger << Endl;
1034
1035 // the numbers
1036 for (UInt_t irow=0; irow<nvvar; irow++) {
1037 logger << setw(maxL) << vert[irow] << ":";
1038 for (UInt_t icol=0; icol<nhvar; icol++) {
1039 logger << setw(hLengths[icol]+1) << Form( "%+1.3f", M(irow,icol) );
1040 }
1041 logger << Endl;
1042 }
1043
1044 // bar line
1045 for (UInt_t i=0; i<clen; i++) logger << "-";
1046 logger << Endl;
1047}
1048
1049////////////////////////////////////////////////////////////////////////////////
1050/// histogramming utility
1051
1053{
1054 return ( unit == "" ? title : ( title + " [" + unit + "]" ) );
1055}
1056
1057////////////////////////////////////////////////////////////////////////////////
1058/// histogramming utility
1059
1060TString TMVA::Tools::GetYTitleWithUnit( const TH1& h, const TString& unit, Bool_t normalised )
1061{
1062 TString retval = ( normalised ? "(1/N) " : "" );
1063 retval += Form( "dN_{ }/^{ }%.3g %s", h.GetXaxis()->GetBinWidth(1), unit.Data() );
1064 return retval;
1065}
1066
1067////////////////////////////////////////////////////////////////////////////////
1068/// writes a float value with the available precision to a stream
1069
1071{
1072 os << val << " :: ";
1073 void * c = &val;
1074 for (int i=0; i<4; i++) {
1075 Int_t ic = *((char*)c+i)-'\0';
1076 if (ic<0) ic+=256;
1077 os << ic << " ";
1078 }
1079 os << ":: ";
1080}
1081
1082////////////////////////////////////////////////////////////////////////////////
1083/// reads a float value with the available precision from a stream
1084
1086{
1087 Float_t a = 0;
1088 is >> a;
1089 TString dn;
1090 is >> dn;
1091 Int_t c[4];
1092 void * ap = &a;
1093 for (int i=0; i<4; i++) {
1094 is >> c[i];
1095 *((char*)ap+i) = '\0'+c[i];
1096 }
1097 is >> dn;
1098 val = a;
1099}
1100
1101// XML file reading/writing helper functions
1102
1103////////////////////////////////////////////////////////////////////////////////
1104/// add attribute from xml
1105
1106Bool_t TMVA::Tools::HasAttr( void* node, const char* attrname )
1107{
1108 return xmlengine().HasAttr(node, attrname);
1109}
1110
1111////////////////////////////////////////////////////////////////////////////////
1112/// add attribute from xml
1113
1114void TMVA::Tools::ReadAttr( void* node, const char* attrname, TString& value )
1115{
1116 if (!HasAttr(node, attrname)) {
1117 const char * nodename = xmlengine().GetNodeName(node);
1118 Log() << kFATAL << "Trying to read non-existing attribute '" << attrname << "' from xml node '" << nodename << "'" << Endl;
1119 }
1120 const char* val = xmlengine().GetAttr(node, attrname);
1121 value = TString(val);
1122}
1123
1124////////////////////////////////////////////////////////////////////////////////
1125/// add attribute to node
1126
1127void TMVA::Tools::AddAttr( void* node, const char* attrname, const char* value )
1128{
1129 if( node == 0 ) return;
1130 gTools().xmlengine().NewAttr(node, 0, attrname, value );
1131}
1132
1133////////////////////////////////////////////////////////////////////////////////
1134/// add child node
1135
1136void* TMVA::Tools::AddChild( void* parent, const char* childname, const char* content, bool isRootNode )
1137{
1138 if( !isRootNode && parent == 0 ) return 0;
1139 return gTools().xmlengine().NewChild(parent, 0, childname, content);
1140}
1141
1142////////////////////////////////////////////////////////////////////////////////
1143
1144Bool_t TMVA::Tools::AddComment( void* node, const char* comment ) {
1145 if( node == 0 ) return kFALSE;
1146 return gTools().xmlengine().AddComment(node, comment);
1147}
1148
1149////////////////////////////////////////////////////////////////////////////////
1150/// get parent node
1151
1152void* TMVA::Tools::GetParent( void* child)
1153{
1154 void* par = xmlengine().GetParent(child);
1155
1156 return par;
1157}
1158
1159////////////////////////////////////////////////////////////////////////////////
1160/// get child node
1161
1162void* TMVA::Tools::GetChild( void* parent, const char* childname )
1163{
1164 void* ch = xmlengine().GetChild(parent);
1165 if (childname != 0) {
1166 while (ch!=0 && strcmp(xmlengine().GetNodeName(ch),childname) != 0) ch = xmlengine().GetNext(ch);
1167 }
1168 return ch;
1169}
1170
1171////////////////////////////////////////////////////////////////////////////////
1172/// XML helpers
1173
1174void* TMVA::Tools::GetNextChild( void* prevchild, const char* childname )
1175{
1176 void* ch = xmlengine().GetNext(prevchild);
1177 if (childname != 0) {
1178 while (ch!=0 && strcmp(xmlengine().GetNodeName(ch),childname)!=0) ch = xmlengine().GetNext(ch);
1179 }
1180 return ch;
1181}
1182
1183////////////////////////////////////////////////////////////////////////////////
1184/// XML helpers
1185
1186const char* TMVA::Tools::GetContent( void* node )
1187{
1188 return xmlengine().GetNodeContent(node);
1189}
1190
1191////////////////////////////////////////////////////////////////////////////////
1192/// XML helpers
1193
1194const char* TMVA::Tools::GetName( void* node )
1195{
1196 return xmlengine().GetNodeName(node);
1197}
1198
1199////////////////////////////////////////////////////////////////////////////////
1200/// XML helpers
1201
1202Bool_t TMVA::Tools::AddRawLine( void* node, const char * raw )
1203{
1204 return xmlengine().AddRawLine( node, raw );
1205}
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// splits the option string at 'separator' and fills the list
1209/// 'splitV' with the primitive strings
1210
1211std::vector<TString> TMVA::Tools::SplitString(const TString& theOpt, const char separator ) const
1212{
1213 std::vector<TString> splitV;
1214 TString splitOpt(theOpt);
1215 splitOpt.ReplaceAll("\n"," ");
1216 splitOpt = splitOpt.Strip(TString::kBoth,separator);
1217 while (splitOpt.Length()>0) {
1218 if ( !splitOpt.Contains(separator) ) {
1219 splitV.push_back(splitOpt);
1220 break;
1221 }
1222 else {
1223 TString toSave = splitOpt(0,splitOpt.First(separator));
1224 splitV.push_back(toSave);
1225 splitOpt = splitOpt(splitOpt.First(separator),splitOpt.Length());
1226 }
1227 splitOpt = splitOpt.Strip(TString::kLeading,separator);
1228 }
1229 return splitV;
1230}
1231
1232////////////////////////////////////////////////////////////////////////////////
1233/// string tools
1234
1236{
1237 std::stringstream s;
1238 s << i;
1239 return TString(s.str().c_str());
1240}
1241
1242////////////////////////////////////////////////////////////////////////////////
1243/// string tools
1244
1246{
1247 std::stringstream s;
1248 s << Form( "%5.8e", d );
1249 return TString(s.str().c_str());
1250}
1251
1252////////////////////////////////////////////////////////////////////////////////
1253/// XML helpers
1254
1255void TMVA::Tools::WriteTMatrixDToXML( void* node, const char* name, TMatrixD* mat )
1256{
1257 void* matnode = xmlengine().NewChild(node, 0, name);
1258 xmlengine().NewAttr(matnode,0,"Rows", StringFromInt(mat->GetNrows()) );
1259 xmlengine().NewAttr(matnode,0,"Columns", StringFromInt(mat->GetNcols()) );
1260 std::stringstream s;
1261 for (Int_t row = 0; row<mat->GetNrows(); row++) {
1262 for (Int_t col = 0; col<mat->GetNcols(); col++) {
1263 s << Form( "%5.15e ", (*mat)[row][col] );
1264 }
1265 }
1266 xmlengine().AddRawLine( matnode, s.str().c_str() );
1267}
1268
1269////////////////////////////////////////////////////////////////////////////////
1270
1271void TMVA::Tools::WriteTVectorDToXML( void* node, const char* name, TVectorD* vec )
1272{
1273 TMatrixD mat(1,vec->GetNoElements(),&((*vec)[0]));
1274 WriteTMatrixDToXML( node, name, &mat );
1275}
1276
1277////////////////////////////////////////////////////////////////////////////////
1278
1279void TMVA::Tools::ReadTVectorDFromXML( void* node, const char* name, TVectorD* vec )
1280{
1281 TMatrixD mat(1,vec->GetNoElements(),&((*vec)[0]));
1282 ReadTMatrixDFromXML( node, name, &mat );
1283 for (int i=0;i<vec->GetNoElements();++i) (*vec)[i] = mat[0][i];
1284}
1285
1286////////////////////////////////////////////////////////////////////////////////
1287
1288void TMVA::Tools::ReadTMatrixDFromXML( void* node, const char* name, TMatrixD* mat )
1289{
1290 if (strcmp(xmlengine().GetNodeName(node),name)!=0){
1291 Log() << kWARNING << "Possible Error: Name of matrix in weight file"
1292 << " does not match name of matrix passed as argument!" << Endl;
1293 }
1294 Int_t nrows, ncols;
1295 ReadAttr( node, "Rows", nrows );
1296 ReadAttr( node, "Columns", ncols );
1297 if (mat->GetNrows() != nrows || mat->GetNcols() != ncols){
1298 Log() << kWARNING << "Possible Error: Dimension of matrix in weight file"
1299 << " does not match dimension of matrix passed as argument!" << Endl;
1300 mat->ResizeTo(nrows,ncols);
1301 }
1302 const char* content = xmlengine().GetNodeContent(node);
1303 std::stringstream s(content);
1304 for (Int_t row = 0; row<nrows; row++) {
1305 for (Int_t col = 0; col<ncols; col++) {
1306 s >> (*mat)[row][col];
1307 }
1308 }
1309}
1310
1311////////////////////////////////////////////////////////////////////////////////
1312/// direct output, eg, when starting ROOT session -> no use of Logger here
1313
1315{
1316 std::cout << std::endl;
1317 std::cout << Color("bold") << "TMVA -- Toolkit for Multivariate Data Analysis" << Color("reset") << std::endl;
1318 std::cout << " " << "Version " << TMVA_RELEASE << ", " << TMVA_RELEASE_DATE << std::endl;
1319 std::cout << " " << "Copyright (C) 2005-2010 CERN, MPI-K Heidelberg, Us of Bonn and Victoria" << std::endl;
1320 std::cout << " " << "Home page: http://tmva.sf.net" << std::endl;
1321 std::cout << " " << "Citation info: http://tmva.sf.net/citeTMVA.html" << std::endl;
1322 std::cout << " " << "License: http://tmva.sf.net/LICENSE" << std::endl << std::endl;
1323}
1324
1325////////////////////////////////////////////////////////////////////////////////
1326/// prints the TMVA release number and date
1327
1329{
1330 logger << "___________TMVA Version " << TMVA_RELEASE << ", " << TMVA_RELEASE_DATE
1331 << "" << Endl;
1332}
1333
1334////////////////////////////////////////////////////////////////////////////////
1335/// prints the ROOT release number and date
1336
1338{
1339 static const char * const months[] = { "Jan","Feb","Mar","Apr","May",
1340 "Jun","Jul","Aug","Sep","Oct",
1341 "Nov","Dec" };
1342 Int_t idatqq = gROOT->GetVersionDate();
1343 Int_t iday = idatqq%100;
1344 Int_t imonth = (idatqq/100)%100;
1345 Int_t iyear = (idatqq/10000);
1346 TString versionDate = Form("%s %d, %4d",months[imonth-1],iday,iyear);
1347
1348 logger << kHEADER ;
1349 logger << "You are running ROOT Version: " << gROOT->GetVersion() << ", " << versionDate << Endl;
1350}
1351
1352////////////////////////////////////////////////////////////////////////////////
1353/// various kinds of welcome messages
1354/// ASCII text generated by this site: http://www.network-science.de/ascii/
1355
1357{
1358 switch (msgType) {
1359
1360 case kStandardWelcomeMsg:
1361 logger << Color("white") << "TMVA -- Toolkit for Multivariate Analysis" << Color("reset") << Endl;
1362 logger << "Copyright (C) 2005-2006 CERN, LAPP & MPI-K Heidelberg and Victoria U." << Endl;
1363 logger << "Home page http://tmva.sourceforge.net" << Endl;
1364 logger << "All rights reserved, please read http://tmva.sf.net/license.txt" << Endl << Endl;
1365 break;
1366
1367 case kIsometricWelcomeMsg:
1368 logger << " ___ ___ ___ ___ " << Endl;
1369 logger << " /\\ \\ /\\__\\ /\\__\\ /\\ \\ " << Endl;
1370 logger << " \\:\\ \\ /::| | /:/ / /::\\ \\ " << Endl;
1371 logger << " \\:\\ \\ /:|:| | /:/ / /:/\\:\\ \\ " << Endl;
1372 logger << " /::\\ \\ /:/|:|__|__ /:/__/ ___ /::\\~\\:\\ \\ " << Endl;
1373 logger << " /:/\\:\\__\\ /:/ |::::\\__\\ |:| | /\\__\\ /:/\\:\\ \\:\\__\\ " << Endl;
1374 logger << " /:/ \\/__/ \\/__/~~/:/ / |:| |/:/ / \\/__\\:\\/:/ / " << Endl;
1375 logger << "/:/ / /:/ / |:|__/:/ / \\::/ / " << Endl;
1376 logger << "\\/__/ /:/ / \\::::/__/ /:/ / " << Endl;
1377 logger << " /:/ / ~~~~ /:/ / " << Endl;
1378 logger << " \\/__/ \\/__/ " << Endl << Endl;
1379 break;
1380
1381 case kBlockWelcomeMsg:
1382 logger << Endl;
1383 logger << "_|_|_|_|_| _| _| _| _| _|_| " << Endl;
1384 logger << " _| _|_| _|_| _| _| _| _| " << Endl;
1385 logger << " _| _| _| _| _| _| _|_|_|_| " << Endl;
1386 logger << " _| _| _| _| _| _| _| " << Endl;
1387 logger << " _| _| _| _| _| _| " << Endl << Endl;
1388 break;
1389
1390 case kLeanWelcomeMsg:
1391 logger << Endl;
1392 logger << "_/_/_/_/_/ _/ _/ _/ _/ _/_/ " << Endl;
1393 logger << " _/ _/_/ _/_/ _/ _/ _/ _/ " << Endl;
1394 logger << " _/ _/ _/ _/ _/ _/ _/_/_/_/ " << Endl;
1395 logger << " _/ _/ _/ _/ _/ _/ _/ " << Endl;
1396 logger << "_/ _/ _/ _/ _/ _/ " << Endl << Endl;
1397 break;
1398
1399 case kLogoWelcomeMsg:
1400 logger << Endl;
1401 logger << "_/_/_/_/_/ _| _| _| _| _|_| " << Endl;
1402 logger << " _/ _|_| _|_| _| _| _| _| " << Endl;
1403 logger << " _/ _| _| _| _| _| _|_|_|_| " << Endl;
1404 logger << " _/ _| _| _| _| _| _| " << Endl;
1405 logger << "_/ _| _| _| _| _| " << Endl << Endl;
1406 break;
1407
1408 case kSmall1WelcomeMsg:
1409 logger << " _____ __ ____ ___ " << Endl;
1410 logger << "|_ _| \\/ \\ \\ / /_\\ " << Endl;
1411 logger << " | | | |\\/| |\\ V / _ \\ " << Endl;
1412 logger << " |_| |_| |_| \\_/_/ \\_\\" << Endl << Endl;
1413 break;
1414
1415 case kSmall2WelcomeMsg:
1416 logger << " _____ __ ____ ___ " << Endl;
1417 logger << "|_ _| \\/ \\ \\ / / \\ " << Endl;
1418 logger << " | | | |\\/| |\\ \\ / / _ \\ " << Endl;
1419 logger << " | | | | | | \\ V / ___ \\ " << Endl;
1420 logger << " |_| |_| |_| \\_/_/ \\_\\ " << Endl << Endl;
1421 break;
1422
1423 case kOriginalWelcomeMsgColor:
1424 logger << kINFO << "" << Color("red")
1425 << "_______________________________________" << Color("reset") << Endl;
1426 logger << kINFO << "" << Color("blue")
1427 << Color("red_bgd") << Color("bwhite") << " // " << Color("reset")
1428 << Color("white") << Color("blue_bgd")
1429 << "|\\ /|| \\ // /\\\\\\\\\\\\\\\\\\\\\\\\ \\ \\ \\ " << Color("reset") << Endl;
1430 logger << kINFO << ""<< Color("blue")
1431 << Color("red_bgd") << Color("white") << "// " << Color("reset")
1432 << Color("white") << Color("blue_bgd")
1433 << "| \\/ || \\// /--\\\\\\\\\\\\\\\\\\\\\\\\ \\ \\ \\" << Color("reset") << Endl;
1434 break;
1435
1436 case kOriginalWelcomeMsgBW:
1437 logger << kINFO << ""
1438 << "_______________________________________" << Endl;
1439 logger << kINFO << " // "
1440 << "|\\ /|| \\ // /\\\\\\\\\\\\\\\\\\\\\\\\ \\ \\ \\ " << Endl;
1441 logger << kINFO << "// "
1442 << "| \\/ || \\// /--\\\\\\\\\\\\\\\\\\\\\\\\ \\ \\ \\" << Endl;
1443 break;
1444
1445 default:
1446 logger << kFATAL << "unknown message type: " << msgType << Endl;
1447 }
1448}
1449
1450////////////////////////////////////////////////////////////////////////////////
1451/// kinds of TMVA citation
1452
1454{
1455 switch (citType) {
1456
1457 case kPlainText:
1458 logger << "A. Hoecker, P. Speckmayer, J. Stelzer, J. Therhaag, E. von Toerne, H. Voss" << Endl;
1459 logger << "\"TMVA - Toolkit for Multivariate Data Analysis\" PoS ACAT:040,2007. e-Print: physics/0703039" << Endl;
1460 break;
1461
1462 case kBibTeX:
1463 logger << "@Article{TMVA2007," << Endl;
1464 logger << " author = \"Hoecker, Andreas and Speckmayer, Peter and Stelzer, Joerg " << Endl;
1465 logger << " and Therhaag, Jan and von Toerne, Eckhard and Voss, Helge\"," << Endl;
1466 logger << " title = \"{TMVA: Toolkit for multivariate data analysis}\"," << Endl;
1467 logger << " journal = \"PoS\"," << Endl;
1468 logger << " volume = \"ACAT\"," << Endl;
1469 logger << " year = \"2007\"," << Endl;
1470 logger << " pages = \"040\"," << Endl;
1471 logger << " eprint = \"physics/0703039\"," << Endl;
1472 logger << " archivePrefix = \"arXiv\"," << Endl;
1473 logger << " SLACcitation = \"%%CITATION = PHYSICS/0703039;%%\"" << Endl;
1474 logger << "}" << Endl;
1475 break;
1476
1477 case kLaTeX:
1478 logger << "%\\cite{TMVA2007}" << Endl;
1479 logger << "\\bibitem{TMVA2007}" << Endl;
1480 logger << " A.~Hoecker, P.~Speckmayer, J.~Stelzer, J.~Therhaag, E.~von Toerne, H.~Voss" << Endl;
1481 logger << " %``TMVA: Toolkit for multivariate data analysis,''" << Endl;
1482 logger << " PoS A {\\bf CAT} (2007) 040" << Endl;
1483 logger << " [arXiv:physics/0703039]." << Endl;
1484 logger << " %%CITATION = POSCI,ACAT,040;%%" << Endl;
1485 break;
1486
1487 case kHtmlLink:
1488 // logger << kINFO << " " << Endl;
1489 logger << kHEADER << gTools().Color("bold")
1490 << "Thank you for using TMVA!" << gTools().Color("reset") << Endl;
1491 logger << kINFO << gTools().Color("bold")
1492 << "For citation information, please visit: http://tmva.sf.net/citeTMVA.html"
1493 << gTools().Color("reset") << Endl;
1494 }
1495}
1496
1497////////////////////////////////////////////////////////////////////////////////
1498
1500{
1501 return !(h.GetXaxis()->GetXbins()->fN);
1502}
1503
1504////////////////////////////////////////////////////////////////////////////////
1505
1506std::vector<TMatrixDSym*>*
1507TMVA::Tools::CalcCovarianceMatrices( const std::vector<const Event*>& events, Int_t maxCls, VariableTransformBase* transformBase )
1508{
1509 std::vector<Event*> eventVector;
1510 for (std::vector<const Event*>::const_iterator it = events.begin(), itEnd = events.end(); it != itEnd; ++it)
1511 {
1512 eventVector.push_back (new Event(*(*it)));
1513 }
1514 std::vector<TMatrixDSym*>* returnValue = CalcCovarianceMatrices (eventVector, maxCls, transformBase);
1515 for (std::vector<Event*>::const_iterator it = eventVector.begin(), itEnd = eventVector.end(); it != itEnd; ++it)
1516 {
1517 delete (*it);
1518 }
1519 return returnValue;
1520}
1521
1522////////////////////////////////////////////////////////////////////////////////
1523/// compute covariance matrices
1524
1525std::vector<TMatrixDSym*>*
1526TMVA::Tools::CalcCovarianceMatrices( const std::vector<Event*>& events, Int_t maxCls, VariableTransformBase* transformBase )
1527{
1528 if (events.empty()) {
1529 Log() << kWARNING << " Asked to calculate a covariance matrix for an empty event vectors.. sorry cannot do that -> return NULL"<<Endl;
1530 return 0;
1531 }
1532
1533 UInt_t nvars=0, ntgts=0, nspcts=0;
1534 if (transformBase)
1535 transformBase->CountVariableTypes( nvars, ntgts, nspcts );
1536 else {
1537 nvars =events.at(0)->GetNVariables ();
1538 ntgts =events.at(0)->GetNTargets ();
1539 nspcts=events.at(0)->GetNSpectators();
1540 }
1541
1542
1543 // init matrices
1544 Int_t matNum = maxCls;
1545 if (maxCls > 1 ) matNum++; // if more than one classes, then produce one matrix for all events as well (beside the matrices for each class)
1546
1547 std::vector<TVectorD*>* vec = new std::vector<TVectorD*>(matNum);
1548 std::vector<TMatrixD*>* mat2 = new std::vector<TMatrixD*>(matNum);
1549 std::vector<Double_t> count(matNum);
1550 count.assign(matNum,0);
1551
1552 Int_t cls = 0;
1553 TVectorD* v;
1554 TMatrixD* m;
1555 UInt_t ivar=0, jvar=0;
1556 for (cls = 0; cls < matNum ; cls++) {
1557 vec->at(cls) = new TVectorD(nvars);
1558 mat2->at(cls) = new TMatrixD(nvars,nvars);
1559 v = vec->at(cls);
1560 m = mat2->at(cls);
1561
1562 for (ivar=0; ivar<nvars; ivar++) {
1563 (*v)(ivar) = 0;
1564 for (jvar=0; jvar<nvars; jvar++) {
1565 (*m)(ivar, jvar) = 0;
1566 }
1567 }
1568 }
1569
1570 // perform event loop
1571 for (UInt_t i=0; i<events.size(); i++) {
1572
1573 // fill the event
1574 const Event * ev = events[i];
1575 cls = ev->GetClass();
1576 Double_t weight = ev->GetWeight();
1577
1578 std::vector<Float_t> input;
1579 std::vector<Char_t> mask; // entries with kTRUE must not be transformed
1580 // Bool_t hasMaskedEntries = kFALSE;
1581 if (transformBase) {
1582 /* hasMaskedEntries = */ transformBase->GetInput (ev, input, mask);
1583 } else {
1584 for (ivar=0; ivar<nvars; ++ivar) {
1585 input.push_back (ev->GetValue(ivar));
1586 }
1587 }
1588
1589 if (maxCls > 1) {
1590 v = vec->at(matNum-1);
1591 m = mat2->at(matNum-1);
1592
1593 count.at(matNum-1)+=weight; // count used events
1594 for (ivar=0; ivar<nvars; ivar++) {
1595
1596 Double_t xi = input.at (ivar);
1597 (*v)(ivar) += xi*weight;
1598 (*m)(ivar, ivar) += (xi*xi*weight);
1599
1600 for (jvar=ivar+1; jvar<nvars; jvar++) {
1601 Double_t xj = input.at (jvar);
1602 (*m)(ivar, jvar) += (xi*xj*weight);
1603 (*m)(jvar, ivar) = (*m)(ivar, jvar); // symmetric matrix
1604 }
1605 }
1606 }
1607
1608 count.at(cls)+=weight; // count used events
1609 v = vec->at(cls);
1610 m = mat2->at(cls);
1611 for (ivar=0; ivar<nvars; ivar++) {
1612 Double_t xi = input.at (ivar);
1613 (*v)(ivar) += xi*weight;
1614 (*m)(ivar, ivar) += (xi*xi*weight);
1615
1616 for (jvar=ivar+1; jvar<nvars; jvar++) {
1617 Double_t xj = input.at (jvar);
1618 (*m)(ivar, jvar) += (xi*xj*weight);
1619 (*m)(jvar, ivar) = (*m)(ivar, jvar); // symmetric matrix
1620 }
1621 }
1622 }
1623
1624 // variance-covariance
1625 std::vector<TMatrixDSym*>* mat = new std::vector<TMatrixDSym*>(matNum);
1626 for (cls = 0; cls < matNum; cls++) {
1627 v = vec->at(cls);
1628 m = mat2->at(cls);
1629
1630 mat->at(cls) = new TMatrixDSym(nvars);
1631
1632 Double_t n = count.at(cls);
1633 for (ivar=0; ivar<nvars; ivar++) {
1634 for (jvar=0; jvar<nvars; jvar++) {
1635 (*(mat->at(cls)))(ivar, jvar) = (*m)(ivar, jvar)/n - (*v)(ivar)*(*v)(jvar)/(n*n);
1636 }
1637 }
1638 delete v;
1639 delete m;
1640 }
1641
1642 delete mat2;
1643 delete vec;
1644
1645 return mat;
1646}
1647
1648////////////////////////////////////////////////////////////////////////////////
1649/// Return the weighted mean of an array defined by the first and
1650/// last iterators. The w iterator should point to the first element
1651/// of a vector of weights of the same size as the main array.
1652
1653template <typename Iterator, typename WeightIterator>
1654Double_t TMVA::Tools::Mean ( Iterator first, Iterator last, WeightIterator w)
1655{
1656 Double_t sum = 0;
1657 Double_t sumw = 0;
1658 int i = 0;
1659 if (w==NULL)
1660 {
1661 while ( first != last )
1662 {
1663 // if ( *w < 0) {
1664 // ::Error("TMVA::Tools::Mean","w[%d] = %.4e < 0 ?!",i,*w);
1665 // return 0;
1666 // } // SURE, why wouldn't you allow for negative event weights here ?? :)
1667 sum += (*first);
1668 sumw += 1.0 ;
1669 ++first;
1670 ++i;
1671 }
1672 if (sumw <= 0) {
1673 ::Error("TMVA::Tools::Mean","sum of weights <= 0 ?! that's a bit too much of negative event weights :) ");
1674 return 0;
1675 }
1676 }
1677 else
1678 {
1679 while ( first != last )
1680 {
1681 // if ( *w < 0) {
1682 // ::Error("TMVA::Tools::Mean","w[%d] = %.4e < 0 ?!",i,*w);
1683 // return 0;
1684 // } // SURE, why wouldn't you allow for negative event weights here ?? :)
1685 sum += (*w) * (*first);
1686 sumw += (*w) ;
1687 ++w;
1688 ++first;
1689 ++i;
1690 }
1691 if (sumw <= 0) {
1692 ::Error("TMVA::Tools::Mean","sum of weights <= 0 ?! that's a bit too much of negative event weights :) ");
1693 return 0;
1694 }
1695 }
1696 return sum/sumw;
1697}
1698
1699////////////////////////////////////////////////////////////////////////////////
1700/// Return the weighted mean of an array a with length n.
1701
1702template <typename T>
1704{
1705 if (w) {
1706 return TMVA::Tools::Mean(a, a+n, w);
1707 } else {
1708 return TMath::Mean(a, a+n);
1709 }
1710}
1711
1712////////////////////////////////////////////////////////////////////////////////
1713/// Return the Standard Deviation of an array defined by the iterators.
1714/// Note that this function returns the sigma(standard deviation) and
1715/// not the root mean square of the array.
1716
1717template <typename Iterator, typename WeightIterator>
1718Double_t TMVA::Tools::RMS(Iterator first, Iterator last, WeightIterator w)
1719{
1720
1721 Double_t sum = 0;
1722 Double_t sum2 = 0;
1723 Double_t sumw = 0;
1724
1725 Double_t adouble;
1726 if (w==NULL)
1727 {
1728 while ( first != last ) {
1729 adouble=Double_t(*first);
1730 sum += adouble;
1731 sum2 += adouble*adouble;
1732 sumw += 1.0;
1733 ++first;
1734 }
1735 }
1736 else
1737 {
1738 while ( first != last ) {
1739 adouble=Double_t(*first);
1740 sum += adouble * (*w);
1741 sum2 += adouble*adouble * (*w);
1742 sumw += (*w);
1743 ++first;
1744 ++w;
1745 }
1746 }
1747 Double_t norm = 1./sumw;
1748 Double_t mean = sum*norm;
1749 Double_t rms = TMath::Sqrt(TMath::Abs(sum2*norm -mean*mean));
1750 return rms;
1751}
1752
1753////////////////////////////////////////////////////////////////////////////////
1754/// Return the Standard Deviation of an array a with length n.
1755/// Note that this function returns the sigma(standard deviation) and
1756/// not the root mean square of the array.
1757
1758template <typename T>
1760{
1761
1762 if (w) {
1763 return TMVA::Tools::RMS(a, a+n, w);
1764 } else {
1765 return TMath::RMS(a, a+n);
1766 }
1767}
1768
1769////////////////////////////////////////////////////////////////////////////////
1770/// get the cumulative distribution of a histogram
1771
1773{
1774 TH1* cumulativeDist= (TH1*) h->Clone(Form("%sCumul",h->GetTitle()));
1775 //cumulativeDist->Smooth(5); // with this, I get less beautiful ROC curves, hence out!
1776
1777 Float_t partialSum = 0;
1778 Float_t inverseSum = 0.;
1779
1780 Float_t val;
1781 for (Int_t ibinEnd=1, ibin=cumulativeDist->GetNbinsX(); ibin >=ibinEnd ; ibin--){
1782 val = cumulativeDist->GetBinContent(ibin);
1783 if (val>0) inverseSum += val;
1784 }
1785 inverseSum = 1/inverseSum; // as I learned multiplications are much faster than division, and later I need one per bin. Well, not that it would really matter here I guess :)
1786
1787 for (Int_t ibinEnd=1, ibin=cumulativeDist->GetNbinsX(); ibin >=ibinEnd ; ibin--){
1788 val = cumulativeDist->GetBinContent(ibin);
1789 if (val>0) partialSum += val;
1790 cumulativeDist->SetBinContent(ibin,partialSum*inverseSum);
1791 }
1792 return cumulativeDist;
1793}
1794
1795void TMVA::Tools::ReadAttr(void *node, const char *attrname, float &value)
1796{
1797 // read attribute from xml
1798 const char *val = xmlengine().GetAttr(node, attrname);
1799 if (val == nullptr) {
1800 const char *nodename = xmlengine().GetNodeName(node);
1801 Log() << kFATAL << "Trying to read non-existing attribute '" << attrname << "' from xml node '" << nodename << "'"
1802 << Endl;
1803 } else
1804 value = atof(val);
1805}
1806
1807void TMVA::Tools::ReadAttr(void *node, const char *attrname, int &value)
1808{
1809 // read attribute from xml
1810 const char *val = xmlengine().GetAttr(node, attrname);
1811 if (val == nullptr) {
1812 const char *nodename = xmlengine().GetNodeName(node);
1813 Log() << kFATAL << "Trying to read non-existing attribute '" << attrname << "' from xml node '" << nodename << "'"
1814 << Endl;
1815 } else
1816 value = atoi(val);
1817}
1818
1819void TMVA::Tools::ReadAttr(void *node, const char *attrname, short &value)
1820{
1821 // read attribute from xml
1822 const char *val = xmlengine().GetAttr(node, attrname);
1823 if (val == nullptr) {
1824 const char *nodename = xmlengine().GetNodeName(node);
1825 Log() << kFATAL << "Trying to read non-existing attribute '" << attrname << "' from xml node '" << nodename << "'"
1826 << Endl;
1827 } else
1828 value = atoi(val);
1829}
SVector< double, 2 > v
Definition: Dict.h:5
ROOT::R::TRInterface & r
Definition: Object.C:4
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define f(i)
Definition: RSha256.hxx:104
#define c(i)
Definition: RSha256.hxx:101
#define h(i)
Definition: RSha256.hxx:106
int Int_t
Definition: RtypesCore.h:41
int Ssiz_t
Definition: RtypesCore.h:63
unsigned int UInt_t
Definition: RtypesCore.h:42
const Bool_t kFALSE
Definition: RtypesCore.h:88
long Long_t
Definition: RtypesCore.h:50
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
long long Long64_t
Definition: RtypesCore.h:69
float Float_t
Definition: RtypesCore.h:53
const Bool_t kTRUE
Definition: RtypesCore.h:87
void Error(const char *location, const char *msgfmt,...)
char name[80]
Definition: TGX11.cxx:109
float xmin
Definition: THbookFile.cxx:93
float xmax
Definition: THbookFile.cxx:93
float type_of_call hi(const int &, const int &)
double pow(double, double)
TMatrixTSym< Double_t > TMatrixDSym
TMatrixT< Double_t > TMatrixD
Definition: TMatrixDfwd.h:22
#define gROOT
Definition: TROOT.h:414
char * Form(const char *fmt,...)
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:22
#define TMVA_RELEASE
Definition: Version.h:44
#define TMVA_RELEASE_DATE
Definition: Version.h:45
Double_t GetXmax() const
Definition: TAxis.h:134
Double_t GetXmin() const
Definition: TAxis.h:133
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
Definition: TCollection.h:182
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:571
The TH1 histogram class.
Definition: TH1.h:56
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition: TH1.cxx:8554
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition: TH1.h:316
virtual Int_t GetNbinsX() const
Definition: TH1.h:292
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8635
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4882
virtual void Scale(Double_t c1=1, Option_t *option="")
Multiply this histogram by a constant c1.
Definition: TH1.cxx:6218
virtual Int_t GetSumw2N() const
Definition: TH1.h:310
virtual Double_t GetSumOfWeights() const
Return the sum of weights excluding under/overflows.
Definition: TH1.cxx:7389
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition: TH1.cxx:8433
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:248
Service class for 2-Dim histogram classes.
Definition: TH2.h:30
virtual void PutStats(Double_t *stats)
Replace current statistics with the values in array stats.
Definition: TH2.cxx:2306
virtual Double_t Integral(Option_t *option="") const
Return integral of bin contents.
Definition: TH2.cxx:1142
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content.
Definition: TH2.cxx:2440
A doubly linked list.
Definition: TList.h:44
virtual void Add(TObject *obj)
Definition: TList.h:87
virtual TObject * At(Int_t idx) const
Returns the object at position idx. Returns 0 if idx is out of range.
Definition: TList.cxx:354
Float_t GetValue(UInt_t ivar) const
return value of i'th variable
Definition: Event.cxx:237
Double_t GetWeight() const
return the event weight - depending on whether the flag IgnoreNegWeightsInTraining is or not.
Definition: Event.cxx:382
UInt_t GetClass() const
Definition: Event.h:87
ostringstream derivative to redirect and format output
Definition: MsgLogger.h:59
PDF wrapper for histograms; uses user-defined spline interpolation.
Definition: PDF.h:63
Double_t GetXmin() const
Definition: PDF.h:104
Double_t GetXmax() const
Definition: PDF.h:105
Double_t GetVal(Double_t x) const
returns value PDF(x)
Definition: PDF.cxx:704
Global auxiliary applications and data treatment routines.
Definition: Tools.h:80
void ComputeStat(const std::vector< TMVA::Event * > &, std::vector< Float_t > *, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Int_t signalClass, Bool_t norm=kFALSE)
sanity check
Definition: Tools.cxx:214
void * GetParent(void *child)
get parent node
Definition: Tools.cxx:1152
void FormattedOutput(const std::vector< Double_t > &, const std::vector< TString > &, const TString titleVars, const TString titleValues, MsgLogger &logger, TString format="%+1.3f")
formatted output of simple table
Definition: Tools.cxx:899
void UsefulSortDescending(std::vector< std::vector< Double_t > > &, std::vector< TString > *vs=0)
sort 2D vector (AND in parallel a TString vector) in such a way that the "first vector is sorted" and...
Definition: Tools.cxx:576
Bool_t HistoHasEquidistantBins(const TH1 &h)
Definition: Tools.cxx:1499
TList * ParseFormatLine(TString theString, const char *sep=":")
Parse the string and cut into labels separated by ":".
Definition: Tools.cxx:413
Double_t NormVariable(Double_t x, Double_t xmin, Double_t xmax)
normalise to output range: [-1, 1]
Definition: Tools.cxx:122
void WriteFloatArbitraryPrecision(Float_t val, std::ostream &os)
writes a float value with the available precision to a stream
Definition: Tools.cxx:1070
void * GetNextChild(void *prevchild, const char *childname=0)
XML helpers.
Definition: Tools.cxx:1174
TString GetYTitleWithUnit(const TH1 &h, const TString &unit, Bool_t normalised)
histogramming utility
Definition: Tools.cxx:1060
Double_t GetSeparation(TH1 *S, TH1 *B) const
compute "separation" defined as
Definition: Tools.cxx:133
static Tools & Instance()
Definition: Tools.cxx:75
TMatrixD * GetSQRootMatrix(TMatrixDSym *symMat)
square-root of symmetric matrix of course the resulting sqrtMat is also symmetric,...
Definition: Tools.cxx:283
std::vector< Double_t > MVADiff(std::vector< Double_t > &, std::vector< Double_t > &)
computes difference between two vectors
Definition: Tools.cxx:518
Int_t GetIndexMinElement(std::vector< Double_t > &)
find index of minimum entry in vector
Definition: Tools.cxx:777
TH1 * projNormTH1F(TTree *theTree, const TString &theVarName, const TString &name, Int_t nbins, Double_t xmin, Double_t xmax, const TString &cut)
projects variable from tree into normalised histogram
Definition: Tools.cxx:378
void ROOTVersionMessage(MsgLogger &logger)
prints the ROOT release number and date
Definition: Tools.cxx:1337
TString ReplaceRegularExpressions(const TString &s, const TString &replace="+")
replace regular expressions helper function to remove all occurrences "$!%^&()'<>?...
Definition: Tools.cxx:810
void * AddChild(void *parent, const char *childname, const char *content=0, bool isRootNode=false)
add child node
Definition: Tools.cxx:1136
void ReadTVectorDFromXML(void *node, const char *name, TVectorD *vec)
Definition: Tools.cxx:1279
TH1 * GetCumulativeDist(TH1 *h)
get the cumulative distribution of a histogram
Definition: Tools.cxx:1772
void ReadFloatArbitraryPrecision(Float_t &val, std::istream &is)
reads a float value with the available precision from a stream
Definition: Tools.cxx:1085
Bool_t AddRawLine(void *node, const char *raw)
XML helpers.
Definition: Tools.cxx:1202
~Tools()
destructor
Definition: Tools.cxx:113
Bool_t ContainsRegularExpression(const TString &s)
check if regular expression helper function to search for "$!%^&()'<>?= " in a string
Definition: Tools.cxx:796
static void DestroyInstance()
Definition: Tools.cxx:90
Double_t Mean(Long64_t n, const T *a, const Double_t *w=0)
Return the weighted mean of an array a with length n.
Definition: Tools.cxx:1703
Double_t GetMutualInformation(const TH2F &)
Mutual Information method for non-linear correlations estimates in 2D histogram Author: Moritz Backes...
Definition: Tools.cxx:601
std::vector< TString > SplitString(const TString &theOpt, const char separator) const
splits the option string at 'separator' and fills the list 'splitV' with the primitive strings
Definition: Tools.cxx:1211
const TString & Color(const TString &)
human readable color strings
Definition: Tools.cxx:840
TString StringFromDouble(Double_t d)
string tools
Definition: Tools.cxx:1245
Bool_t CheckForVerboseOption(const TString &) const
check if verbosity "V" set in option
Definition: Tools.cxx:720
Double_t RMS(Long64_t n, const T *a, const Double_t *w=0)
Return the Standard Deviation of an array a with length n.
Definition: Tools.cxx:1759
const char * GetContent(void *node)
XML helpers.
Definition: Tools.cxx:1186
void * GetChild(void *parent, const char *childname=0)
get child node
Definition: Tools.cxx:1162
void WriteTVectorDToXML(void *node, const char *name, TVectorD *vec)
Definition: Tools.cxx:1271
TXMLEngine & xmlengine()
Definition: Tools.h:270
void Scale(std::vector< Double_t > &, Double_t)
scales double vector
Definition: Tools.cxx:531
Bool_t CheckSplines(const TH1 *, const TSpline *)
check quality of splining by comparing splines and histograms in each bin
Definition: Tools.cxx:491
const char * GetName(void *node)
XML helpers.
Definition: Tools.cxx:1194
const TMatrixD * GetCorrelationMatrix(const TMatrixD *covMat)
turns covariance into correlation matrix
Definition: Tools.cxx:336
void WriteTMatrixDToXML(void *node, const char *name, TMatrixD *mat)
XML helpers.
Definition: Tools.cxx:1255
Double_t GetCorrelationRatio(const TH2F &)
Compute Correlation Ratio of 2D histogram to estimate functional dependency between two variables Aut...
Definition: Tools.cxx:632
void ReadAttr(void *node, const char *, T &value)
read attribute from xml
Definition: Tools.h:337
Bool_t AddComment(void *node, const char *comment)
Definition: Tools.cxx:1144
void ReadTMatrixDFromXML(void *node, const char *name, TMatrixD *mat)
Definition: Tools.cxx:1288
TString GetXTitleWithUnit(const TString &title, const TString &unit)
histogramming utility
Definition: Tools.cxx:1052
void AddAttr(void *node, const char *, const T &value, Int_t precision=16)
add attribute to xml
Definition: Tools.h:355
TString StringFromInt(Long_t i)
string tools
Definition: Tools.cxx:1235
Tools()
constructor
Definition: Tools.cxx:103
Double_t NormHist(TH1 *theHist, Double_t norm=1.0)
normalises histogram
Definition: Tools.cxx:395
std::vector< TMatrixDSym * > * CalcCovarianceMatrices(const std::vector< Event * > &events, Int_t maxCls, VariableTransformBase *transformBase=0)
compute covariance matrices
Definition: Tools.cxx:1526
Double_t GetYMean_binX(const TH2 &, Int_t bin_x)
Compute the mean in Y for a given bin X of a 2D histogram.
Definition: Tools.cxx:654
void TMVACitation(MsgLogger &logger, ECitation citType=kPlainText)
kinds of TMVA citation
Definition: Tools.cxx:1453
void TMVAVersionMessage(MsgLogger &logger)
prints the TMVA release number and date
Definition: Tools.cxx:1328
void TMVAWelcomeMessage()
direct output, eg, when starting ROOT session -> no use of Logger here
Definition: Tools.cxx:1314
std::vector< Int_t > * ParseANNOptionString(TString theOptions, Int_t nvar, std::vector< Int_t > *nodes)
parse option string for ANN methods default settings (should be defined in theOption string)
Definition: Tools.cxx:455
Int_t GetIndexMaxElement(std::vector< Double_t > &)
find index of maximum entry in vector
Definition: Tools.cxx:760
TH2F * TransposeHist(const TH2F &)
Transpose quadratic histogram.
Definition: Tools.cxx:669
static Tools * fgTools
Definition: Tools.h:236
EWelcomeMessage
Definition: Tools.h:202
Bool_t HasAttr(void *node, const char *attrname)
add attribute from xml
Definition: Tools.cxx:1106
void UsefulSortAscending(std::vector< std::vector< Double_t > > &, std::vector< TString > *vs=0)
sort 2D vector (AND in parallel a TString vector) in such a way that the "first vector is sorted" and...
Definition: Tools.cxx:550
const TString fRegexp
Definition: Tools.h:230
Bool_t CheckForSilentOption(const TString &) const
check for "silence" option in configuration option string
Definition: Tools.cxx:703
Linear interpolation class.
virtual Bool_t GetInput(const Event *event, std::vector< Float_t > &input, std::vector< Char_t > &mask, Bool_t backTransform=kFALSE) const
select the values from the event
virtual void CountVariableTypes(UInt_t &nvars, UInt_t &ntgts, UInt_t &nspcts) const
count variables, targets and spectators
TMatrixDSymEigen.
const TMatrixD & GetEigenVectors() const
Int_t GetNrows() const
Definition: TMatrixTBase.h:124
Int_t GetNcols() const
Definition: TMatrixTBase.h:127
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:1210
TMatrixT< Element > & Invert(Double_t *det=0)
Invert the matrix and calculate its determinant.
Definition: TMatrixT.cxx:1396
TMatrixT< Element > & Transpose(const TMatrixT< Element > &source)
Transpose matrix source.
Definition: TMatrixT.cxx:1469
void Mult(const TMatrixT< Element > &a, const TMatrixT< Element > &b)
General matrix multiplication. Create a matrix C such that C = A * B.
Definition: TMatrixT.cxx:648
Collectable string class.
Definition: TObjString.h:28
Base class for spline implementation containing the Draw/Paint methods.
Definition: TSpline.h:22
virtual Double_t Eval(Double_t x) const =0
Basic string class.
Definition: TString.h:131
Ssiz_t Length() const
Definition: TString.h:405
TSubString Strip(EStripType s=kTrailing, char c=' ') const
Return a substring of self stripped at beginning and/or end.
Definition: TString.cxx:1106
Ssiz_t First(char c) const
Find first occurrence of a character c.
Definition: TString.cxx:499
const char * Data() const
Definition: TString.h:364
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:687
@ kLeading
Definition: TString.h:262
@ kBoth
Definition: TString.h:262
TString & Remove(Ssiz_t pos)
Definition: TString.h:668
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:619
A TTree represents a columnar dataset.
Definition: TTree.h:71
virtual Long64_t Project(const char *hname, const char *varexp, const char *selection="", Option_t *option="", Long64_t nentries=kMaxEntries, Long64_t firstentry=0)
Make a projection of a tree using selections.
Definition: TTree.cxx:7235
Int_t GetNoElements() const
Definition: TVectorT.h:76
XMLAttrPointer_t NewAttr(XMLNodePointer_t xmlnode, XMLNsPointer_t, const char *name, const char *value)
creates new attribute for xmlnode, namespaces are not supported for attributes
Definition: TXMLEngine.cxx:580
Bool_t AddComment(XMLNodePointer_t parent, const char *comment)
Adds comment line to the node.
Definition: TXMLEngine.cxx:872
XMLNodePointer_t NewChild(XMLNodePointer_t parent, XMLNsPointer_t ns, const char *name, const char *content=0)
create new child element for parent node
Definition: TXMLEngine.cxx:709
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
static double B[]
double T(double x)
Definition: ChebyshevPol.h:34
TClass * GetClass(T *)
Definition: TClass.h:608
RooArgSet S(const RooAbsArg &v1)
RooCmdArg Color(Color_t color)
static constexpr double s
Config & gConfig()
Tools & gTools()
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:158
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:97
Double_t Log(Double_t x)
Definition: TMath.h:748
Double_t Mean(Long64_t n, const T *a, const Double_t *w=0)
Return the weighted mean of an array a with length n.
Definition: TMath.h:1061
Double_t Sqrt(Double_t x)
Definition: TMath.h:679
Double_t RMS(Long64_t n, const T *a, const Double_t *w=0)
Return the Standard Deviation of an array a with length n.
Definition: TMath.h:1155
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Definition: first.py:1
auto * m
Definition: textangle.C:8
auto * a
Definition: textangle.C:12
static long int sum(long int i)
Definition: Factory.cxx:2258
REAL epsilon
Definition: triangle.c:617