Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
Reader.cxx
Go to the documentation of this file.
1// @(#)root/tmva $Id$
2// Author: Andreas Hoecker, Joerg Stelzer, Helge Voss, Eckhard von Toerne, Jan Therhaag
3
4/**********************************************************************************
5 * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6 * Package: TMVA *
7 * Class : Reader *
8 * Web : http://tmva.sourceforge.net *
9 * *
10 * Description: *
11 * Reader class to be used in the user application to interpret the trained *
12 * MVAs in an analysis context *
13 * *
14 * Authors (alphabetical order): *
15 * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
16 * Peter Speckmayer <peter.speckmayer@cern.ch> - CERN, Switzerland *
17 * Joerg Stelzer <Joerg.Stelzer@cern.ch> - CERN, Switzerland *
18 * Jan Therhaag <Jan.Therhaag@cern.ch> - U of Bonn, Germany *
19 * Eckhard v. Toerne <evt@uni-bonn.de> - U of Bonn, Germany *
20 * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
21 * Kai Voss <Kai.Voss@cern.ch> - U. of Victoria, Canada *
22 * *
23 * Copyright (c) 2005-2011: *
24 * CERN, Switzerland *
25 * U. of Victoria, Canada *
26 * MPI-K Heidelberg, Germany *
27 * U. of Bonn, Germany *
28 * *
29 * Redistribution and use in source and binary forms, with or without *
30 * modification, are permitted according to the terms listed in LICENSE *
31 * (http://tmva.sourceforge.net/LICENSE) *
32 **********************************************************************************/
33
34/*! \class TMVA::Reader
35\ingroup TMVA
36
37 The Reader class serves to use the MVAs in a specific analysis context.
38 Within an event loop, a vector is filled that corresponds to the variables
39 that were used to train the MVA(s) during the training stage. This vector
40 is transfered to the Reader, who takes care of interpreting the weight
41 file of the MVA of choice, and to return the MVA's output. This is then
42 used by the user for further analysis.
43
44 Usage:
45
46~~~ {.cpp}
47 // ------ before starting the event loop (eg, in the initialisation step)
48
49 //
50 // create TMVA::Reader object
51 //
52 TMVA::Reader *reader = new TMVA::Reader();
53
54 // create a set of variables and declare them to the reader
55 // - the variable names must corresponds in name and type to
56 // those given in the weight file(s) that you use
57 Float_t var1, var2, var3, var4;
58 reader->AddVariable( "var1", &var1 );
59 reader->AddVariable( "var2", &var2 );
60 reader->AddVariable( "var3", &var3 );
61 reader->AddVariable( "var4", &var4 );
62
63 // book the MVA of your choice (prior training of these methods, ie,
64 // existence of the weight files is required)
65 reader->BookMVA( "Fisher method", "weights/Fisher.weights.txt" );
66 reader->BookMVA( "MLP method", "weights/MLP.weights.txt" );
67 // ... etc
68
69 // ------- start your event loop
70
71 for (Long64_t ievt=0; ievt<myTree->GetEntries();ievt++) {
72
73 // fill vector with values of variables computed from those in the tree
74 var1 = myvar1;
75 var2 = myvar2;
76 var3 = myvar3;
77 var4 = myvar4;
78
79 // retrieve the corresponding MVA output
80 double mvaFi = reader->EvaluateMVA( "Fisher method" );
81 double mvaNN = reader->EvaluateMVA( "MLP method" );
82
83 // do something with these ...., e.g., fill them into your ntuple
84
85 } // end of event loop
86
87 delete reader;
88~~~
89*/
90
91#include "TMVA/Reader.h"
92
93#include "TMVA/Config.h"
94#include "TMVA/Configurable.h"
97#include "TMVA/DataSetInfo.h"
98#include "TMVA/DataSetManager.h"
99#include "TMVA/IMethod.h"
100#include "TMVA/MethodBase.h"
101#include "TMVA/MethodCuts.h"
102#include "TMVA/MethodCategory.h"
103#include "TMVA/MsgLogger.h"
104#include "TMVA/Tools.h"
105#include "TMVA/Types.h"
106
107#include "TLeaf.h"
108#include "TString.h"
109#include "TH1D.h"
110#include "TVector.h"
111#include "TXMLEngine.h"
112#include "TMath.h"
113
114#include <cstdlib>
115
116#include <string>
117#include <vector>
118#include <fstream>
119
120////////////////////////////////////////////////////////////////////////////////
121/// constructor
122
123TMVA::Reader::Reader( const TString& theOption, Bool_t verbose )
124: Configurable( theOption ),
125 fDataSetManager( NULL ), // DSMTEST
126 fDataSetInfo(),
127 fVerbose( verbose ),
128 fSilent ( kFALSE ),
129 fColor ( kFALSE ),
130 fCalculateError(kFALSE),
131 fMvaEventError( 0 ),
132 fMvaEventErrorUpper( 0 ),
133 fLogger ( 0 )
134{
137 fLogger = new MsgLogger(this);
140 ParseOptions();
141
142 Init();
143}
144
145////////////////////////////////////////////////////////////////////////////////
146/// constructor
147
148TMVA::Reader::Reader( std::vector<TString>& inputVars, const TString& theOption, Bool_t verbose )
149 : Configurable( theOption ),
150 fDataSetManager( NULL ), // DSMTEST
151 fDataSetInfo(),
152 fVerbose( verbose ),
153 fSilent ( kFALSE ),
154 fColor ( kFALSE ),
155 fCalculateError(kFALSE),
156 fMvaEventError( 0 ),
157 fMvaEventErrorUpper( 0 ), //zjh
158 fLogger ( 0 )
159{
162 fLogger = new MsgLogger(this);
165 ParseOptions();
166
167 // arguments: names of input variables (vector)
168 // verbose flag
169 for (std::vector<TString>::iterator ivar = inputVars.begin(); ivar != inputVars.end(); ++ivar)
170 DataInfo().AddVariable( *ivar );
171
172 Init();
173}
174
175////////////////////////////////////////////////////////////////////////////////
176/// constructor
177
178TMVA::Reader::Reader( std::vector<std::string>& inputVars, const TString& theOption, Bool_t verbose )
179 : Configurable( theOption ),
180 fDataSetManager( NULL ), // DSMTEST
181 fDataSetInfo(),
182 fVerbose( verbose ),
183 fSilent ( kFALSE ),
184 fColor ( kFALSE ),
185 fCalculateError(kFALSE),
186 fMvaEventError( 0 ),
187 fMvaEventErrorUpper( 0 ),
188 fLogger ( 0 )
189{
192 fLogger = new MsgLogger(this);
195 ParseOptions();
196
197 // arguments: names of input variables (vector)
198 // verbose flag
199 for (std::vector<std::string>::iterator ivar = inputVars.begin(); ivar != inputVars.end(); ++ivar)
200 DataInfo().AddVariable( ivar->c_str() );
201
202 Init();
203}
204
205////////////////////////////////////////////////////////////////////////////////
206/// constructor
207
208TMVA::Reader::Reader( const std::string& varNames, const TString& theOption, Bool_t verbose )
209 : Configurable( theOption ),
210 fDataSetManager( NULL ), // DSMTEST
211 fDataSetInfo(),
212 fVerbose( verbose ),
213 fSilent ( kFALSE ),
214 fColor ( kFALSE ),
215 fCalculateError(kFALSE),
216 fMvaEventError( 0 ),
217 fMvaEventErrorUpper( 0 ),
218 fLogger ( 0 )
219{
222 fLogger = new MsgLogger(this);
225 ParseOptions();
226
227 // arguments: names of input variables given in form: "name1:name2:name3"
228 // verbose flag
229 DecodeVarNames(varNames);
230 Init();
231}
232
233////////////////////////////////////////////////////////////////////////////////
234/// constructor
235
236TMVA::Reader::Reader( const TString& varNames, const TString& theOption, Bool_t verbose )
237 : Configurable( theOption ),
238 fDataSetManager( NULL ), // DSMTEST
239 fDataSetInfo(),
240 fVerbose( verbose ),
241 fSilent ( kFALSE ),
242 fColor ( kFALSE ),
243 fCalculateError(kFALSE),
244 fMvaEventError( 0 ),
245 fMvaEventErrorUpper( 0 ),
246 fLogger ( 0 )
247{
250 fLogger = new MsgLogger(this);
253 ParseOptions();
254
255 // arguments: names of input variables given in form: "name1:name2:name3"
256 // verbose flag
257 DecodeVarNames(varNames);
258 Init();
259}
260
261////////////////////////////////////////////////////////////////////////////////
262/// declaration of configuration options
263
265{
266 if (gTools().CheckForSilentOption( GetOptions() )) Log().InhibitOutput(); // make sure is silent if wanted to
267
268 DeclareOptionRef( fVerbose, "V", "Verbose flag" );
269 DeclareOptionRef( fColor, "Color", "Color flag (default True)" );
270 DeclareOptionRef( fSilent, "Silent", "Boolean silent flag (default False)" );
271 DeclareOptionRef( fCalculateError, "Error", "Calculates errors (default False)" );
272}
273
274////////////////////////////////////////////////////////////////////////////////
275/// destructor
276
278{
279 delete fDataSetManager; // DSMTEST
280
281 delete fLogger;
282
283 for (auto it=fMethodMap.begin(); it!=fMethodMap.end(); it++){
284 MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(it->second);
285 delete kl;
286 }
287}
288
289////////////////////////////////////////////////////////////////////////////////
290/// default initialisation (no member variables)
291
293{
294 if (Verbose()) fLogger->SetMinType( kVERBOSE );
295
296 gConfig().SetUseColor( fColor );
297 gConfig().SetSilent ( fSilent );
298}
299
300////////////////////////////////////////////////////////////////////////////////
301/// Add a float variable or expression to the reader
302
303void TMVA::Reader::AddVariable( const TString& expression, Float_t* datalink )
304{
305 DataInfo().AddVariable( expression, "", "", 0, 0, 'F', kFALSE ,(void*)datalink ); // <= should this be F or rather T?
306}
307
308////////////////////////////////////////////////////////////////////////////////
309
310void TMVA::Reader::AddVariable( const TString& expression, Int_t* datalink )
311{
312 Log() << kFATAL << "Reader::AddVariable( const TString& expression, Int_t* datalink ), this function is deprecated, please provide all variables to the reader as floats" << Endl;
313 // Add an integer variable or expression to the reader
314 Log() << kFATAL << "Reader::AddVariable( const TString& expression, Int_t* datalink ), this function is deprecated, please provide all variables to the reader as floats" << Endl;
315 DataInfo().AddVariable(expression, "", "", 0, 0, 'I', kFALSE, (void*)datalink ); // <= should this be F or rather T?
316}
317
318////////////////////////////////////////////////////////////////////////////////
319/// Add a float spectator or expression to the reader
320
321void TMVA::Reader::AddSpectator( const TString& expression, Float_t* datalink )
322{
323 DataInfo().AddSpectator( expression, "", "", 0, 0, 'F', kFALSE ,(void*)datalink );
324}
325
326////////////////////////////////////////////////////////////////////////////////
327/// Add an integer spectator or expression to the reader
328
329void TMVA::Reader::AddSpectator( const TString& expression, Int_t* datalink )
330{
331 DataInfo().AddSpectator(expression, "", "", 0, 0, 'I', kFALSE, (void*)datalink );
332}
333
334////////////////////////////////////////////////////////////////////////////////
335/// read the method type from the file
336
338{
339 std::ifstream fin( filename );
340 if (!fin.good()) { // file not found --> Error
341 Log() << kFATAL << "<BookMVA> fatal error: "
342 << "unable to open input weight file: " << filename << Endl;
343 }
344
345 TString fullMethodName("");
346 if (filename.EndsWith(".xml")) {
347 fin.close();
348 void* doc = gTools().xmlengine().ParseFile(filename,gTools().xmlenginebuffersize());// the default buffer size in TXMLEngine::ParseFile is 100k. Starting with ROOT 5.29 one can set the buffer size, see: http://savannah.cern.ch/bugs/?78864. This might be necessary for large XML files
349 void* rootnode = gTools().xmlengine().DocGetRootElement(doc); // node "MethodSetup"
350 gTools().ReadAttr(rootnode, "Method", fullMethodName);
351 gTools().xmlengine().FreeDoc(doc);
352 }
353 else {
354 char buf[512];
355 fin.getline(buf,512);
356 while (!TString(buf).BeginsWith("Method")) fin.getline(buf,512);
357 fullMethodName = TString(buf);
358 fin.close();
359 }
360 TString methodType = fullMethodName(0,fullMethodName.Index("::"));
361 if (methodType.Contains(" ")) methodType = methodType(methodType.Last(' ')+1,methodType.Length());
362 return methodType;
363}
364
365////////////////////////////////////////////////////////////////////////////////
366/// read method name from weight file
367
368TMVA::IMethod* TMVA::Reader::BookMVA( const TString& methodTag, const TString& weightfile )
369{
370 // assert non-existence
371 if (fMethodMap.find( methodTag ) != fMethodMap.end())
372 Log() << kFATAL << "<BookMVA> method tag \"" << methodTag << "\" already exists!" << Endl;
373
374 TString methodType(GetMethodTypeFromFile(weightfile));
375
376 Log() << kINFO << "Booking \"" << methodTag << "\" of type \"" << methodType << "\" from " << weightfile << "." << Endl;
377
378 MethodBase* method = dynamic_cast<MethodBase*>(this->BookMVA( Types::Instance().GetMethodType(methodType),
379 weightfile ) );
380 if( method && method->GetMethodType() == Types::kCategory ){
381 MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
382 if( !methCat )
383 Log() << kFATAL << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
384 methCat->fDataSetManager = fDataSetManager;
385 }
386
387 return fMethodMap[methodTag] = method;
388}
389
390////////////////////////////////////////////////////////////////////////////////
391/// books MVA method from weightfile
392
394{
395 IMethod *im =
396 ClassifierFactory::Instance().Create(Types::Instance().GetMethodName(methodType).Data(), DataInfo(), weightfile);
397
398 MethodBase *method = (dynamic_cast<MethodBase*>(im));
399
400 if (method==0) return im;
401
402 if( method->GetMethodType() == Types::kCategory ){
403 MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
404 if( !methCat )
405 Log() << kERROR << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
406 methCat->fDataSetManager = fDataSetManager;
407 }
408
409 method->SetupMethod();
410
411 // when reading older weight files, they could include options
412 // that are not supported any longer
414
415 // read weight file
416 method->ReadStateFromFile();
417
418 // check for unused options
419 method->CheckSetup();
420
421 Log() << kINFO << "Booked classifier \"" << method->GetMethodName()
422 << "\" of type: \"" << method->GetMethodTypeName() << "\"" << Endl;
423
424 return method;
425}
426
427////////////////////////////////////////////////////////////////////////////////
428
429TMVA::IMethod* TMVA::Reader::BookMVA( TMVA::Types::EMVA methodType, const char* xmlstr )
430{
431 // books MVA method from weightfile
432 IMethod *im =
433 ClassifierFactory::Instance().Create(Types::Instance().GetMethodName(methodType).Data(), DataInfo(), "");
434
435 MethodBase *method = (dynamic_cast<MethodBase*>(im));
436
437 if(!method) return 0;
438
439 if( method->GetMethodType() == Types::kCategory ){
440 MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
441 if( !methCat )
442 Log() << kFATAL << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
443 methCat->fDataSetManager = fDataSetManager;
444 }
445
446 method->SetupMethod();
447
448 // when reading older weight files, they could include options
449 // that are not supported any longer
451
452 // read weight file
453 method->ReadStateFromXMLString( xmlstr );
454
455 // check for unused options
456 method->CheckSetup();
457
458 Log() << kINFO << "Booked classifier \"" << method->GetMethodName()
459 << "\" of type: \"" << method->GetMethodTypeName() << "\"" << Endl;
460
461 return method;
462}
463
464////////////////////////////////////////////////////////////////////////////////
465/// Evaluate a std::vector<float> of input data for a given method
466/// The parameter aux is obligatory for the cuts method where it represents the efficiency cutoff
467
468Double_t TMVA::Reader::EvaluateMVA( const std::vector<Float_t>& inputVec, const TString& methodTag, Double_t aux )
469{
470 // create a temporary event from the vector.
471 IMethod* imeth = FindMVA( methodTag );
472 MethodBase* meth = dynamic_cast<TMVA::MethodBase*>(imeth);
473 if(meth==0) return 0;
474
475 // Event* tmpEvent=new Event(inputVec, 2); // ToDo resolve magic 2 issue
476 Event* tmpEvent=new Event(inputVec, DataInfo().GetNVariables()); // is this the solution?
477 for (UInt_t i=0; i<inputVec.size(); i++){
478 if (TMath::IsNaN(inputVec[i])) {
479 Log() << kERROR << i << "-th variable of the event is NaN --> return MVA value -999, \n that's all I can do, please fix or remove this event." << Endl;
480 delete tmpEvent;
481 return -999;
482 }
483 }
484
485 if (meth->GetMethodType() == TMVA::Types::kCuts) {
486 TMVA::MethodCuts* mc = dynamic_cast<TMVA::MethodCuts*>(meth);
487 if(mc)
488 mc->SetTestSignalEfficiency( aux );
489 }
490 Double_t val = meth->GetMvaValue( tmpEvent, (fCalculateError?&fMvaEventError:0));
491 delete tmpEvent;
492 return val;
493}
494
495////////////////////////////////////////////////////////////////////////////////
496/// Evaluate a std::vector<double> of input data for a given method
497/// The parameter aux is obligatory for the cuts method where it represents the efficiency cutoff
498
499Double_t TMVA::Reader::EvaluateMVA( const std::vector<Double_t>& inputVec, const TString& methodTag, Double_t aux )
500{
501 // performs a copy to float values which are internally used by all methods
502 if(fTmpEvalVec.size() != inputVec.size())
503 fTmpEvalVec.resize(inputVec.size());
504
505 for (UInt_t idx=0; idx!=inputVec.size(); idx++ )
506 fTmpEvalVec[idx]=inputVec[idx];
507
508 return EvaluateMVA( fTmpEvalVec, methodTag, aux );
509}
510
511////////////////////////////////////////////////////////////////////////////////
512/// evaluates MVA for given set of input variables
513
515{
516 IMethod* method = 0;
517
518 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
519 if (it == fMethodMap.end()) {
520 Log() << kINFO << "<EvaluateMVA> unknown classifier in map; "
521 << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
522 for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
523 Log() << "Check calling string" << kFATAL << Endl;
524 }
525
526 else method = it->second;
527
528 MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
529
530 if(kl==0)
531 Log() << kFATAL << methodTag << " is not a method" << Endl;
532
533 // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
534 // it is not again checked in each of these subsequent calls..
535 const Event* ev = kl->GetEvent();
536 for (UInt_t i=0; i<ev->GetNVariables(); i++){
537 if (TMath::IsNaN(ev->GetValue(i))) {
538 Log() << kERROR << i << "-th variable of the event is NaN --> return MVA value -999, \n that's all I can do, please fix or remove this event." << Endl;
539 return -999;
540 }
541 }
542 return this->EvaluateMVA( kl, aux );
543}
544
545////////////////////////////////////////////////////////////////////////////////
546/// evaluates the MVA
547
549{
550 // the aux value is only needed for MethodCuts: it sets the
551 // required signal efficiency
552 if (method->GetMethodType() == TMVA::Types::kCuts) {
553 TMVA::MethodCuts* mc = dynamic_cast<TMVA::MethodCuts*>(method);
554 if(mc)
555 mc->SetTestSignalEfficiency( aux );
556 }
557
558 return method->GetMvaValue( (fCalculateError?&fMvaEventError:0),
559 (fCalculateError?&fMvaEventErrorUpper:0) );
560}
561
562////////////////////////////////////////////////////////////////////////////////
563/// evaluates MVA for given set of input variables
564
565const std::vector< Float_t >& TMVA::Reader::EvaluateRegression( const TString& methodTag, Double_t aux )
566{
567 IMethod* method = 0;
568
569 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
570 if (it == fMethodMap.end()) {
571 Log() << kINFO << "<EvaluateMVA> unknown method in map; "
572 << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
573 for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
574 Log() << "Check calling string" << kFATAL << Endl;
575 }
576 else method = it->second;
577
578 MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
579
580 if(kl==0)
581 Log() << kFATAL << methodTag << " is not a method" << Endl;
582 // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
583 // it is not again checked in each of these subsequent calls..
584 const Event* ev = kl->GetEvent();
585 for (UInt_t i=0; i<ev->GetNVariables(); i++){
586 if (TMath::IsNaN(ev->GetValue(i))) {
587 Log() << kERROR << i << "-th variable of the event is NaN, \n regression values might evaluate to .. what do I know. \n sorry this warning is all I can do, please fix or remove this event." << Endl;
588 }
589 }
590
591 return this->EvaluateRegression( kl, aux );
592}
593
594////////////////////////////////////////////////////////////////////////////////
595/// evaluates the regression MVA
596/// check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
597/// it is not again checked in each of these subsequent calls.
598
599const std::vector< Float_t >& TMVA::Reader::EvaluateRegression( MethodBase* method, Double_t /*aux*/ )
600{
601 const Event* ev = method->GetEvent();
602 for (UInt_t i=0; i<ev->GetNVariables(); i++){
603 if (TMath::IsNaN(ev->GetValue(i))) {
604 Log() << kERROR << i << "-th variable of the event is NaN, \n regression values might evaluate to .. what do I know. \n sorry this warning is all I can do, please fix or remove this event." << Endl;
605 }
606 }
607 return method->GetRegressionValues();
608}
609
610
611////////////////////////////////////////////////////////////////////////////////
612/// evaluates the regression MVA
613
615{
616 try {
617 return EvaluateRegression(methodTag, aux).at(tgtNumber);
618 }
619 catch (std::out_of_range &) {
620 Log() << kWARNING << "Regression could not be evaluated for target-number " << tgtNumber << Endl;
621 return 0;
622 }
623}
624
625
626
627////////////////////////////////////////////////////////////////////////////////
628/// evaluates MVA for given set of input variables
629
630const std::vector< Float_t >& TMVA::Reader::EvaluateMulticlass( const TString& methodTag, Double_t aux )
631{
632 IMethod* method = 0;
633
634 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
635 if (it == fMethodMap.end()) {
636 Log() << kINFO << "<EvaluateMVA> unknown method in map; "
637 << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
638 for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
639 Log() << "Check calling string" << kFATAL << Endl;
640 }
641 else method = it->second;
642
643 MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
644
645 if(kl==0)
646 Log() << kFATAL << methodTag << " is not a method" << Endl;
647 // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
648 // it is not again checked in each of these subsequent calls..
649
650 const Event* ev = kl->GetEvent();
651 for (UInt_t i=0; i<ev->GetNVariables(); i++){
652 if (TMath::IsNaN(ev->GetValue(i))) {
653 Log() << kERROR << i << "-th variable of the event is NaN, \n regression values might evaluate to .. what do I know. \n sorry this warning is all I can do, please fix or remove this event." << Endl;
654 }
655 }
656
657 return this->EvaluateMulticlass( kl, aux );
658}
659
660////////////////////////////////////////////////////////////////////////////////
661/// evaluates the multiclass MVA
662/// check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
663/// it is not again checked in each of these subsequent calls.
664
665const std::vector< Float_t >& TMVA::Reader::EvaluateMulticlass( MethodBase* method, Double_t /*aux*/ )
666{
667 const Event* ev = method->GetEvent();
668 for (UInt_t i=0; i<ev->GetNVariables(); i++){
669 if (TMath::IsNaN(ev->GetValue(i))) {
670 Log() << kERROR << i << "-th variable of the event is NaN, \n regression values might evaluate to .. what do I know. \n sorry this warning is all I can do, please fix or remove this event." << Endl;
671 }
672 }
673 return method->GetMulticlassValues();
674}
675
676
677////////////////////////////////////////////////////////////////////////////////
678/// evaluates the multiclass MVA
679
681{
682 try {
683 return EvaluateMulticlass(methodTag, aux).at(clsNumber);
684 }
685 catch (std::out_of_range &) {
686 Log() << kWARNING << "Multiclass could not be evaluated for class-number " << clsNumber << Endl;
687 return 0;
688 }
689}
690
691
692////////////////////////////////////////////////////////////////////////////////
693/// return pointer to method with tag "methodTag"
694
696{
697 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
698 if (it != fMethodMap.end()) return it->second;
699 Log() << kERROR << "Method " << methodTag << " not found!" << Endl;
700 return 0;
701}
702
703////////////////////////////////////////////////////////////////////////////////
704/// special function for Cuts to avoid dynamic_casts in ROOT macros,
705/// which are not properly handled by CINT
706
708{
709 return dynamic_cast<MethodCuts*>(FindMVA(methodTag));
710}
711
712////////////////////////////////////////////////////////////////////////////////
713/// evaluates probability of MVA for given set of input variables
714
715Double_t TMVA::Reader::GetProba( const TString& methodTag, Double_t ap_sig, Double_t mvaVal )
716{
717 IMethod* method = 0;
718 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
719 if (it == fMethodMap.end()) {
720 for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "M" << it->first << Endl;
721 Log() << kFATAL << "<EvaluateMVA> unknown classifier in map: " << method << "; "
722 << "you looked for " << methodTag<< " while the available methods are : " << Endl;
723 }
724 else method = it->second;
725
726 MethodBase* kl = dynamic_cast<MethodBase*>(method);
727 if(kl==0) return -1;
728 // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
729 // it is not again checked in each of these subsequent calls..
730 const Event* ev = kl->GetEvent();
731 for (UInt_t i=0; i<ev->GetNVariables(); i++){
732 if (TMath::IsNaN(ev->GetValue(i))) {
733 Log() << kERROR << i << "-th variable of the event is NaN --> return MVA value -999, \n that's all I can do, please fix or remove this event." << Endl;
734 return -999;
735 }
736 }
737
738 if (mvaVal == -9999999) mvaVal = kl->GetMvaValue();
739
740 return kl->GetProba( mvaVal, ap_sig );
741}
742
743////////////////////////////////////////////////////////////////////////////////
744/// evaluates the MVA's rarity
745
747{
748 IMethod* method = 0;
749 std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
750 if (it == fMethodMap.end()) {
751 for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "M" << it->first << Endl;
752 Log() << kFATAL << "<EvaluateMVA> unknown classifier in map: \"" << method << "\"; "
753 << "you looked for \"" << methodTag<< "\" while the available methods are : " << Endl;
754 }
755 else method = it->second;
756
757 MethodBase* kl = dynamic_cast<MethodBase*>(method);
758 if(kl==0) return -1;
759 // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
760 // it is not again checked in each of these subsequent calls..
761 const Event* ev = kl->GetEvent();
762 for (UInt_t i=0; i<ev->GetNVariables(); i++){
763 if (TMath::IsNaN(ev->GetValue(i))) {
764 Log() << kERROR << i << "-th variable of the event is NaN --> return MVA value -999, \n that's all I can do, please fix or remove this event." << Endl;
765 return -999;
766 }
767 }
768
769 if (mvaVal == -9999999) mvaVal = kl->GetMvaValue();
770
771 return kl->GetRarity( mvaVal );
772}
773
774// ---------------------------------------------------------------------------------------
775// ----- methods related to the decoding of the input variable names ---------------------
776// ---------------------------------------------------------------------------------------
777
778////////////////////////////////////////////////////////////////////////////////
779/// decodes "name1:name2:..." form
780
781void TMVA::Reader::DecodeVarNames( const std::string& varNames )
782{
783 size_t ipos = 0, f = 0;
784 while (f != varNames.length()) {
785 f = varNames.find( ':', ipos );
786 if (f > varNames.length()) f = varNames.length();
787 std::string subs = varNames.substr( ipos, f-ipos ); ipos = f+1;
788 DataInfo().AddVariable( subs.c_str() );
789 }
790}
791
792////////////////////////////////////////////////////////////////////////////////
793/// decodes "name1:name2:..." form
794
796{
797 TString format;
798 Int_t n = varNames.Length();
799 TString format_obj;
800
801 for (int i=0; i< n+1 ; i++) {
802 format.Append(varNames(i));
803 if (varNames(i) == ':' || i == n) {
804 format.Chop();
805 format_obj = format;
806 format_obj.ReplaceAll("@","");
807 DataInfo().AddVariable( format_obj );
808 format.Resize(0);
809 }
810 }
811}
#define f(i)
Definition RSha256.hxx:104
const Bool_t kFALSE
Definition RtypesCore.h:101
double Double_t
Definition RtypesCore.h:59
float Float_t
Definition RtypesCore.h:57
IMethod * Create(const std::string &name, const TString &job, const TString &title, DataSetInfo &dsi, const TString &option)
creates the method if needed based on the method name using the creator function the factory has stor...
static ClassifierFactory & Instance()
access to the ClassifierFactory singleton creates the instance if needed
void SetUseColor(Bool_t uc)
Definition Config.h:60
void SetSilent(Bool_t s)
Definition Config.h:63
void SetConfigName(const char *n)
virtual void ParseOptions()
options parser
VariableInfo & AddVariable(const TString &expression, const TString &title="", const TString &unit="", Double_t min=0, Double_t max=0, char varType='F', Bool_t normalized=kTRUE, void *external=0)
add a variable (can be a complex expression) to the set of variables used in the MV analysis
Class that contains all the data information.
DataSetInfo & AddDataSetInfo(DataSetInfo &dsi)
stores a copy of the dataset info object
Float_t GetValue(UInt_t ivar) const
return value of i'th variable
Definition Event.cxx:236
UInt_t GetNVariables() const
accessor to the number of variables
Definition Event.cxx:316
Interface for all concrete MVA method implementations.
Definition IMethod.h:53
Virtual base Class for all MVA method.
Definition MethodBase.h:111
const std::vector< Float_t > & GetRegressionValues(const TMVA::Event *const ev)
Definition MethodBase.h:214
virtual void DeclareCompatibilityOptions()
options that are used ONLY for the READER to ensure backward compatibility they are hence without any...
TString GetMethodTypeName() const
Definition MethodBase.h:332
virtual Double_t GetProba(const Event *ev)
virtual const std::vector< Float_t > & GetMulticlassValues()
Definition MethodBase.h:227
void SetupMethod()
setup of methods
const TString & GetMethodName() const
Definition MethodBase.h:331
const Event * GetEvent() const
Definition MethodBase.h:751
void ReadStateFromXMLString(const char *xmlstr)
for reading from memory
void ReadStateFromFile()
Function to write options and weights to file.
virtual Double_t GetMvaValue(Double_t *errLower=0, Double_t *errUpper=0)=0
virtual Double_t GetRarity(Double_t mvaVal, Types::ESBType reftype=Types::kBackground) const
compute rarity:
Types::EMVA GetMethodType() const
Definition MethodBase.h:333
virtual void CheckSetup()
check may be overridden by derived class (sometimes, eg, fitters are used which can only be implement...
Class for categorizing the phase space.
DataSetManager * fDataSetManager
Multivariate optimisation of signal efficiency for given background efficiency, applying rectangular ...
Definition MethodCuts.h:61
void SetTestSignalEfficiency(Double_t effS)
Definition MethodCuts.h:116
ostringstream derivative to redirect and format output
Definition MsgLogger.h:57
const std::vector< Float_t > & EvaluateRegression(const TString &methodTag, Double_t aux=0)
evaluates MVA for given set of input variables
Definition Reader.cxx:565
virtual const char * GetName() const
Returns name of object.
Definition Reader.h:117
Double_t EvaluateMVA(const std::vector< Float_t > &, const TString &methodTag, Double_t aux=0)
Evaluate a std::vector<float> of input data for a given method The parameter aux is obligatory for th...
Definition Reader.cxx:468
Double_t GetRarity(const TString &methodTag, Double_t mvaVal=-9999999)
evaluates the MVA's rarity
Definition Reader.cxx:746
IMethod * FindMVA(const TString &methodTag)
return pointer to method with tag "methodTag"
Definition Reader.cxx:695
void Init(void)
default initialisation (no member variables)
Definition Reader.cxx:292
Double_t GetProba(const TString &methodTag, Double_t ap_sig=0.5, Double_t mvaVal=-9999999)
evaluates probability of MVA for given set of input variables
Definition Reader.cxx:715
MethodCuts * FindCutsMVA(const TString &methodTag)
special function for Cuts to avoid dynamic_casts in ROOT macros, which are not properly handled by CI...
Definition Reader.cxx:707
TString GetMethodTypeFromFile(const TString &filename)
read the method type from the file
Definition Reader.cxx:337
DataSetManager * fDataSetManager
Definition Reader.h:132
DataSetInfo fDataSetInfo
Definition Reader.h:140
IMethod * BookMVA(const TString &methodTag, const TString &weightfile)
read method name from weight file
Definition Reader.cxx:368
Reader(const TString &theOption="", Bool_t verbose=0)
constructor
Definition Reader.cxx:123
const std::vector< Float_t > & EvaluateMulticlass(const TString &methodTag, Double_t aux=0)
evaluates MVA for given set of input variables
Definition Reader.cxx:630
DataInputHandler fDataInputHandler
Definition Reader.h:142
void DecodeVarNames(const std::string &varNames)
decodes "name1:name2:..." form
Definition Reader.cxx:781
void DeclareOptions()
declaration of configuration options
Definition Reader.cxx:264
void AddSpectator(const TString &expression, Float_t *)
Add a float spectator or expression to the reader.
Definition Reader.cxx:321
void AddVariable(const TString &expression, Float_t *)
Add a float variable or expression to the reader.
Definition Reader.cxx:303
virtual ~Reader(void)
destructor
Definition Reader.cxx:277
MsgLogger * fLogger
Definition Reader.h:165
const DataSetInfo & DataInfo() const
Definition Reader.h:121
TXMLEngine & xmlengine()
Definition Tools.h:262
void ReadAttr(void *node, const char *, T &value)
read attribute from xml
Definition Tools.h:329
static Types & Instance()
the the single instance of "Types" if existing already, or create it (Singleton)
Definition Types.cxx:70
@ kCategory
Definition Types.h:97
Basic string class.
Definition TString.h:136
Ssiz_t Length() const
Definition TString.h:410
Bool_t EndsWith(const char *pat, ECaseCompare cmp=kExact) const
Return true if string ends with the specified string.
Definition TString.cxx:2202
TString & Chop()
Definition TString.h:679
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:692
void Resize(Ssiz_t n)
Resize the string. Truncate or add blanks as necessary.
Definition TString.cxx:1120
Ssiz_t Last(char c) const
Find last occurrence of a character c.
Definition TString.cxx:916
TString & Append(const char *cs)
Definition TString.h:564
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:624
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:639
void FreeDoc(XMLDocPointer_t xmldoc)
frees allocated document data and deletes document itself
XMLNodePointer_t DocGetRootElement(XMLDocPointer_t xmldoc)
returns root node of document
XMLDocPointer_t ParseFile(const char *filename, Int_t maxbuf=100000)
Parses content of file and tries to produce xml structures.
const Int_t n
Definition legend1.C:16
Config & gConfig()
Tools & gTools()
MsgLogger & Endl(MsgLogger &ml)
Definition MsgLogger.h:148
Bool_t IsNaN(Double_t x)
Definition TMath.h:842