Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
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"
95 #include "TMVA/ClassifierFactory.h"
96 #include "TMVA/DataInputHandler.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 "TTree.h"
108 #include "TLeaf.h"
109 #include "TString.h"
110 #include "TClass.h"
111 #include "TH1D.h"
112 #include "TKey.h"
113 #include "TVector.h"
114 #include "TXMLEngine.h"
115 #include "TMath.h"
116 
117 #include <cstdlib>
118 
119 #include <string>
120 #include <vector>
121 #include <fstream>
122 
123 #include <iostream>
124 
125 ////////////////////////////////////////////////////////////////////////////////
126 /// constructor
127 
128 TMVA::Reader::Reader( const TString& theOption, Bool_t verbose )
129 : Configurable( theOption ),
130  fDataSetManager( NULL ), // DSMTEST
131  fDataSetInfo(),
132  fVerbose( verbose ),
133  fSilent ( kFALSE ),
134  fColor ( kFALSE ),
135  fCalculateError(kFALSE),
136  fMvaEventError( 0 ),
137  fMvaEventErrorUpper( 0 ),
138  fLogger ( 0 )
139 {
140  fDataSetManager = new DataSetManager( fDataInputHandler );
141  fDataSetManager->AddDataSetInfo(fDataSetInfo);
142  fLogger = new MsgLogger(this);
143  SetConfigName( GetName() );
144  DeclareOptions();
145  ParseOptions();
146 
147  Init();
148 }
149 
150 ////////////////////////////////////////////////////////////////////////////////
151 /// constructor
152 
153 TMVA::Reader::Reader( std::vector<TString>& inputVars, const TString& theOption, Bool_t verbose )
154  : Configurable( theOption ),
155  fDataSetManager( NULL ), // DSMTEST
156  fDataSetInfo(),
157  fVerbose( verbose ),
158  fSilent ( kFALSE ),
159  fColor ( kFALSE ),
160  fCalculateError(kFALSE),
161  fMvaEventError( 0 ),
162  fMvaEventErrorUpper( 0 ), //zjh
163  fLogger ( 0 )
164 {
165  fDataSetManager = new DataSetManager( fDataInputHandler );
166  fDataSetManager->AddDataSetInfo(fDataSetInfo);
167  fLogger = new MsgLogger(this);
168  SetConfigName( GetName() );
169  DeclareOptions();
170  ParseOptions();
171 
172  // arguments: names of input variables (vector)
173  // verbose flag
174  for (std::vector<TString>::iterator ivar = inputVars.begin(); ivar != inputVars.end(); ++ivar)
175  DataInfo().AddVariable( *ivar );
176 
177  Init();
178 }
179 
180 ////////////////////////////////////////////////////////////////////////////////
181 /// constructor
182 
183 TMVA::Reader::Reader( std::vector<std::string>& inputVars, const TString& theOption, Bool_t verbose )
184  : Configurable( theOption ),
185  fDataSetManager( NULL ), // DSMTEST
186  fDataSetInfo(),
187  fVerbose( verbose ),
188  fSilent ( kFALSE ),
189  fColor ( kFALSE ),
190  fCalculateError(kFALSE),
191  fMvaEventError( 0 ),
192  fMvaEventErrorUpper( 0 ),
193  fLogger ( 0 )
194 {
195  fDataSetManager = new DataSetManager( fDataInputHandler );
196  fDataSetManager->AddDataSetInfo(fDataSetInfo);
197  fLogger = new MsgLogger(this);
198  SetConfigName( GetName() );
199  DeclareOptions();
200  ParseOptions();
201 
202  // arguments: names of input variables (vector)
203  // verbose flag
204  for (std::vector<std::string>::iterator ivar = inputVars.begin(); ivar != inputVars.end(); ++ivar)
205  DataInfo().AddVariable( ivar->c_str() );
206 
207  Init();
208 }
209 
210 ////////////////////////////////////////////////////////////////////////////////
211 /// constructor
212 
213 TMVA::Reader::Reader( const std::string& varNames, const TString& theOption, Bool_t verbose )
214  : Configurable( theOption ),
215  fDataSetManager( NULL ), // DSMTEST
216  fDataSetInfo(),
217  fVerbose( verbose ),
218  fSilent ( kFALSE ),
219  fColor ( kFALSE ),
220  fCalculateError(kFALSE),
221  fMvaEventError( 0 ),
222  fMvaEventErrorUpper( 0 ),
223  fLogger ( 0 )
224 {
225  fDataSetManager = new DataSetManager( fDataInputHandler );
226  fDataSetManager->AddDataSetInfo(fDataSetInfo);
227  fLogger = new MsgLogger(this);
228  SetConfigName( GetName() );
229  DeclareOptions();
230  ParseOptions();
231 
232  // arguments: names of input variables given in form: "name1:name2:name3"
233  // verbose flag
234  DecodeVarNames(varNames);
235  Init();
236 }
237 
238 ////////////////////////////////////////////////////////////////////////////////
239 /// constructor
240 
241 TMVA::Reader::Reader( const TString& varNames, const TString& theOption, Bool_t verbose )
242  : Configurable( theOption ),
243  fDataSetManager( NULL ), // DSMTEST
244  fDataSetInfo(),
245  fVerbose( verbose ),
246  fSilent ( kFALSE ),
247  fColor ( kFALSE ),
248  fCalculateError(kFALSE),
249  fMvaEventError( 0 ),
250  fMvaEventErrorUpper( 0 ),
251  fLogger ( 0 )
252 {
253  fDataSetManager = new DataSetManager( fDataInputHandler );
254  fDataSetManager->AddDataSetInfo(fDataSetInfo);
255  fLogger = new MsgLogger(this);
256  SetConfigName( GetName() );
257  DeclareOptions();
258  ParseOptions();
259 
260  // arguments: names of input variables given in form: "name1:name2:name3"
261  // verbose flag
262  DecodeVarNames(varNames);
263  Init();
264 }
265 
266 ////////////////////////////////////////////////////////////////////////////////
267 /// declaration of configuration options
268 
269 void TMVA::Reader::DeclareOptions()
270 {
271  if (gTools().CheckForSilentOption( GetOptions() )) Log().InhibitOutput(); // make sure is silent if wanted to
272 
273  DeclareOptionRef( fVerbose, "V", "Verbose flag" );
274  DeclareOptionRef( fColor, "Color", "Color flag (default True)" );
275  DeclareOptionRef( fSilent, "Silent", "Boolean silent flag (default False)" );
276  DeclareOptionRef( fCalculateError, "Error", "Calculates errors (default False)" );
277 }
278 
279 ////////////////////////////////////////////////////////////////////////////////
280 /// destructor
281 
282 TMVA::Reader::~Reader( void )
283 {
284  delete fDataSetManager; // DSMTEST
285 
286  delete fLogger;
287 
288  for (auto it=fMethodMap.begin(); it!=fMethodMap.end(); it++){
289  MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(it->second);
290  delete kl;
291  }
292 }
293 
294 ////////////////////////////////////////////////////////////////////////////////
295 /// default initialisation (no member variables)
296 
297 void TMVA::Reader::Init( void )
298 {
299  if (Verbose()) fLogger->SetMinType( kVERBOSE );
300 
301  gConfig().SetUseColor( fColor );
302  gConfig().SetSilent ( fSilent );
303 }
304 
305 ////////////////////////////////////////////////////////////////////////////////
306 /// Add a float variable or expression to the reader
307 
308 void TMVA::Reader::AddVariable( const TString& expression, Float_t* datalink )
309 {
310  DataInfo().AddVariable( expression, "", "", 0, 0, 'F', kFALSE ,(void*)datalink ); // <= should this be F or rather T?
311 }
312 
313 ////////////////////////////////////////////////////////////////////////////////
314 
315 void TMVA::Reader::AddVariable( const TString& expression, Int_t* datalink )
316 {
317  Log() << kFATAL << "Reader::AddVariable( const TString& expression, Int_t* datalink ), this function is deprecated, please provide all variables to the reader as floats" << Endl;
318  // Add an integer variable or expression to the reader
319  Log() << kFATAL << "Reader::AddVariable( const TString& expression, Int_t* datalink ), this function is deprecated, please provide all variables to the reader as floats" << Endl;
320  DataInfo().AddVariable(expression, "", "", 0, 0, 'I', kFALSE, (void*)datalink ); // <= should this be F or rather T?
321 }
322 
323 ////////////////////////////////////////////////////////////////////////////////
324 /// Add a float spectator or expression to the reader
325 
326 void TMVA::Reader::AddSpectator( const TString& expression, Float_t* datalink )
327 {
328  DataInfo().AddSpectator( expression, "", "", 0, 0, 'F', kFALSE ,(void*)datalink );
329 }
330 
331 ////////////////////////////////////////////////////////////////////////////////
332 /// Add an integer spectator or expression to the reader
333 
334 void TMVA::Reader::AddSpectator( const TString& expression, Int_t* datalink )
335 {
336  DataInfo().AddSpectator(expression, "", "", 0, 0, 'I', kFALSE, (void*)datalink );
337 }
338 
339 ////////////////////////////////////////////////////////////////////////////////
340 /// read the method type from the file
341 
342 TString TMVA::Reader::GetMethodTypeFromFile( const TString& filename )
343 {
344  std::ifstream fin( filename );
345  if (!fin.good()) { // file not found --> Error
346  Log() << kFATAL << "<BookMVA> fatal error: "
347  << "unable to open input weight file: " << filename << Endl;
348  }
349 
350  TString fullMethodName("");
351  if (filename.EndsWith(".xml")) {
352  fin.close();
353  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
354  void* rootnode = gTools().xmlengine().DocGetRootElement(doc); // node "MethodSetup"
355  gTools().ReadAttr(rootnode, "Method", fullMethodName);
356  gTools().xmlengine().FreeDoc(doc);
357  }
358  else {
359  char buf[512];
360  fin.getline(buf,512);
361  while (!TString(buf).BeginsWith("Method")) fin.getline(buf,512);
362  fullMethodName = TString(buf);
363  fin.close();
364  }
365  TString methodType = fullMethodName(0,fullMethodName.Index("::"));
366  if (methodType.Contains(" ")) methodType = methodType(methodType.Last(' ')+1,methodType.Length());
367  return methodType;
368 }
369 
370 ////////////////////////////////////////////////////////////////////////////////
371 /// read method name from weight file
372 
373 TMVA::IMethod* TMVA::Reader::BookMVA( const TString& methodTag, const TString& weightfile )
374 {
375  // assert non-existence
376  if (fMethodMap.find( methodTag ) != fMethodMap.end())
377  Log() << kFATAL << "<BookMVA> method tag \"" << methodTag << "\" already exists!" << Endl;
378 
379  TString methodType(GetMethodTypeFromFile(weightfile));
380 
381  Log() << kINFO << "Booking \"" << methodTag << "\" of type \"" << methodType << "\" from " << weightfile << "." << Endl;
382 
383  MethodBase* method = dynamic_cast<MethodBase*>(this->BookMVA( Types::Instance().GetMethodType(methodType),
384  weightfile ) );
385  if( method && method->GetMethodType() == Types::kCategory ){
386  MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
387  if( !methCat )
388  Log() << kFATAL << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
389  methCat->fDataSetManager = fDataSetManager;
390  }
391 
392  return fMethodMap[methodTag] = method;
393 }
394 
395 ////////////////////////////////////////////////////////////////////////////////
396 /// books MVA method from weightfile
397 
398 TMVA::IMethod* TMVA::Reader::BookMVA( TMVA::Types::EMVA methodType, const TString& weightfile )
399 {
400  IMethod *im =
401  ClassifierFactory::Instance().Create(Types::Instance().GetMethodName(methodType).Data(), DataInfo(), weightfile);
402 
403  MethodBase *method = (dynamic_cast<MethodBase*>(im));
404 
405  if (method==0) return im;
406 
407  if( method->GetMethodType() == Types::kCategory ){
408  MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
409  if( !methCat )
410  Log() << kERROR << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
411  methCat->fDataSetManager = fDataSetManager;
412  }
413 
414  method->SetupMethod();
415 
416  // when reading older weight files, they could include options
417  // that are not supported any longer
418  method->DeclareCompatibilityOptions();
419 
420  // read weight file
421  method->ReadStateFromFile();
422 
423  // check for unused options
424  method->CheckSetup();
425 
426  Log() << kINFO << "Booked classifier \"" << method->GetMethodName()
427  << "\" of type: \"" << method->GetMethodTypeName() << "\"" << Endl;
428 
429  return method;
430 }
431 
432 ////////////////////////////////////////////////////////////////////////////////
433 
434 TMVA::IMethod* TMVA::Reader::BookMVA( TMVA::Types::EMVA methodType, const char* xmlstr )
435 {
436  // books MVA method from weightfile
437  IMethod *im =
438  ClassifierFactory::Instance().Create(Types::Instance().GetMethodName(methodType).Data(), DataInfo(), "");
439 
440  MethodBase *method = (dynamic_cast<MethodBase*>(im));
441 
442  if(!method) return 0;
443 
444  if( method->GetMethodType() == Types::kCategory ){
445  MethodCategory *methCat = (dynamic_cast<MethodCategory*>(method));
446  if( !methCat )
447  Log() << kFATAL << "Method with type kCategory cannot be casted to MethodCategory. /Reader" << Endl;
448  methCat->fDataSetManager = fDataSetManager;
449  }
450 
451  method->SetupMethod();
452 
453  // when reading older weight files, they could include options
454  // that are not supported any longer
455  method->DeclareCompatibilityOptions();
456 
457  // read weight file
458  method->ReadStateFromXMLString( xmlstr );
459 
460  // check for unused options
461  method->CheckSetup();
462 
463  Log() << kINFO << "Booked classifier \"" << method->GetMethodName()
464  << "\" of type: \"" << method->GetMethodTypeName() << "\"" << Endl;
465 
466  return method;
467 }
468 
469 ////////////////////////////////////////////////////////////////////////////////
470 /// Evaluate a std::vector<float> of input data for a given method
471 /// The parameter aux is obligatory for the cuts method where it represents the efficiency cutoff
472 
473 Double_t TMVA::Reader::EvaluateMVA( const std::vector<Float_t>& inputVec, const TString& methodTag, Double_t aux )
474 {
475  // create a temporary event from the vector.
476  IMethod* imeth = FindMVA( methodTag );
477  MethodBase* meth = dynamic_cast<TMVA::MethodBase*>(imeth);
478  if(meth==0) return 0;
479 
480  // Event* tmpEvent=new Event(inputVec, 2); // ToDo resolve magic 2 issue
481  Event* tmpEvent=new Event(inputVec, DataInfo().GetNVariables()); // is this the solution?
482  for (UInt_t i=0; i<inputVec.size(); i++){
483  if (TMath::IsNaN(inputVec[i])) {
484  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;
485  delete tmpEvent;
486  return -999;
487  }
488  }
489 
490  if (meth->GetMethodType() == TMVA::Types::kCuts) {
491  TMVA::MethodCuts* mc = dynamic_cast<TMVA::MethodCuts*>(meth);
492  if(mc)
493  mc->SetTestSignalEfficiency( aux );
494  }
495  Double_t val = meth->GetMvaValue( tmpEvent, (fCalculateError?&fMvaEventError:0));
496  delete tmpEvent;
497  return val;
498 }
499 
500 ////////////////////////////////////////////////////////////////////////////////
501 /// Evaluate a std::vector<double> of input data for a given method
502 /// The parameter aux is obligatory for the cuts method where it represents the efficiency cutoff
503 
504 Double_t TMVA::Reader::EvaluateMVA( const std::vector<Double_t>& inputVec, const TString& methodTag, Double_t aux )
505 {
506  // performs a copy to float values which are internally used by all methods
507  if(fTmpEvalVec.size() != inputVec.size())
508  fTmpEvalVec.resize(inputVec.size());
509 
510  for (UInt_t idx=0; idx!=inputVec.size(); idx++ )
511  fTmpEvalVec[idx]=inputVec[idx];
512 
513  return EvaluateMVA( fTmpEvalVec, methodTag, aux );
514 }
515 
516 ////////////////////////////////////////////////////////////////////////////////
517 /// evaluates MVA for given set of input variables
518 
519 Double_t TMVA::Reader::EvaluateMVA( const TString& methodTag, Double_t aux )
520 {
521  IMethod* method = 0;
522 
523  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
524  if (it == fMethodMap.end()) {
525  Log() << kINFO << "<EvaluateMVA> unknown classifier in map; "
526  << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
527  for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
528  Log() << "Check calling string" << kFATAL << Endl;
529  }
530 
531  else method = it->second;
532 
533  MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
534 
535  if(kl==0)
536  Log() << kFATAL << methodTag << " is not a method" << Endl;
537 
538  // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
539  // it is not again checked in each of these subsequent calls..
540  const Event* ev = kl->GetEvent();
541  for (UInt_t i=0; i<ev->GetNVariables(); i++){
542  if (TMath::IsNaN(ev->GetValue(i))) {
543  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;
544  return -999;
545  }
546  }
547  return this->EvaluateMVA( kl, aux );
548 }
549 
550 ////////////////////////////////////////////////////////////////////////////////
551 /// evaluates the MVA
552 
553 Double_t TMVA::Reader::EvaluateMVA( MethodBase* method, Double_t aux )
554 {
555  // the aux value is only needed for MethodCuts: it sets the
556  // required signal efficiency
557  if (method->GetMethodType() == TMVA::Types::kCuts) {
558  TMVA::MethodCuts* mc = dynamic_cast<TMVA::MethodCuts*>(method);
559  if(mc)
560  mc->SetTestSignalEfficiency( aux );
561  }
562 
563  return method->GetMvaValue( (fCalculateError?&fMvaEventError:0),
564  (fCalculateError?&fMvaEventErrorUpper:0) );
565 }
566 
567 ////////////////////////////////////////////////////////////////////////////////
568 /// evaluates MVA for given set of input variables
569 
570 const std::vector< Float_t >& TMVA::Reader::EvaluateRegression( const TString& methodTag, Double_t aux )
571 {
572  IMethod* method = 0;
573 
574  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
575  if (it == fMethodMap.end()) {
576  Log() << kINFO << "<EvaluateMVA> unknown method in map; "
577  << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
578  for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
579  Log() << "Check calling string" << kFATAL << Endl;
580  }
581  else method = it->second;
582 
583  MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
584 
585  if(kl==0)
586  Log() << kFATAL << methodTag << " is not a method" << Endl;
587  // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
588  // it is not again checked in each of these subsequent calls..
589  const Event* ev = kl->GetEvent();
590  for (UInt_t i=0; i<ev->GetNVariables(); i++){
591  if (TMath::IsNaN(ev->GetValue(i))) {
592  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;
593  }
594  }
595 
596  return this->EvaluateRegression( kl, aux );
597 }
598 
599 ////////////////////////////////////////////////////////////////////////////////
600 /// evaluates the regression MVA
601 /// check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
602 /// it is not again checked in each of these subsequent calls..
603 
604 const std::vector< Float_t >& TMVA::Reader::EvaluateRegression( MethodBase* method, Double_t /*aux*/ )
605 {
606  const Event* ev = method->GetEvent();
607  for (UInt_t i=0; i<ev->GetNVariables(); i++){
608  if (TMath::IsNaN(ev->GetValue(i))) {
609  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;
610  }
611  }
612  return method->GetRegressionValues();
613 }
614 
615 
616 ////////////////////////////////////////////////////////////////////////////////
617 /// evaluates the regression MVA
618 
619 Float_t TMVA::Reader::EvaluateRegression( UInt_t tgtNumber, const TString& methodTag, Double_t aux )
620 {
621  try {
622  return EvaluateRegression(methodTag, aux).at(tgtNumber);
623  }
624  catch (std::out_of_range &) {
625  Log() << kWARNING << "Regression could not be evaluated for target-number " << tgtNumber << Endl;
626  return 0;
627  }
628 }
629 
630 
631 
632 ////////////////////////////////////////////////////////////////////////////////
633 /// evaluates MVA for given set of input variables
634 
635 const std::vector< Float_t >& TMVA::Reader::EvaluateMulticlass( const TString& methodTag, Double_t aux )
636 {
637  IMethod* method = 0;
638 
639  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
640  if (it == fMethodMap.end()) {
641  Log() << kINFO << "<EvaluateMVA> unknown method in map; "
642  << "you looked for \"" << methodTag << "\" within available methods: " << Endl;
643  for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "--> " << it->first << Endl;
644  Log() << "Check calling string" << kFATAL << Endl;
645  }
646  else method = it->second;
647 
648  MethodBase * kl = dynamic_cast<TMVA::MethodBase*>(method);
649 
650  if(kl==0)
651  Log() << kFATAL << methodTag << " is not a method" << Endl;
652  // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
653  // it is not again checked in each of these subsequent calls..
654 
655  const Event* ev = kl->GetEvent();
656  for (UInt_t i=0; i<ev->GetNVariables(); i++){
657  if (TMath::IsNaN(ev->GetValue(i))) {
658  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;
659  }
660  }
661 
662  return this->EvaluateMulticlass( kl, aux );
663 }
664 
665 ////////////////////////////////////////////////////////////////////////////////
666 /// evaluates the multiclass MVA
667 /// check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
668 /// it is not again checked in each of these subsequent calls..
669 
670 const std::vector< Float_t >& TMVA::Reader::EvaluateMulticlass( MethodBase* method, Double_t /*aux*/ )
671 {
672  const Event* ev = method->GetEvent();
673  for (UInt_t i=0; i<ev->GetNVariables(); i++){
674  if (TMath::IsNaN(ev->GetValue(i))) {
675  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;
676  }
677  }
678  return method->GetMulticlassValues();
679 }
680 
681 
682 ////////////////////////////////////////////////////////////////////////////////
683 /// evaluates the multiclass MVA
684 
685 Float_t TMVA::Reader::EvaluateMulticlass( UInt_t clsNumber, const TString& methodTag, Double_t aux )
686 {
687  try {
688  return EvaluateMulticlass(methodTag, aux).at(clsNumber);
689  }
690  catch (std::out_of_range &) {
691  Log() << kWARNING << "Multiclass could not be evaluated for class-number " << clsNumber << Endl;
692  return 0;
693  }
694 }
695 
696 
697 ////////////////////////////////////////////////////////////////////////////////
698 /// return pointer to method with tag "methodTag"
699 
700 TMVA::IMethod* TMVA::Reader::FindMVA( const TString& methodTag )
701 {
702  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
703  if (it != fMethodMap.end()) return it->second;
704  Log() << kERROR << "Method " << methodTag << " not found!" << Endl;
705  return 0;
706 }
707 
708 ////////////////////////////////////////////////////////////////////////////////
709 /// special function for Cuts to avoid dynamic_casts in ROOT macros,
710 /// which are not properly handled by CINT
711 
712 TMVA::MethodCuts* TMVA::Reader::FindCutsMVA( const TString& methodTag )
713 {
714  return dynamic_cast<MethodCuts*>(FindMVA(methodTag));
715 }
716 
717 ////////////////////////////////////////////////////////////////////////////////
718 /// evaluates probability of MVA for given set of input variables
719 
720 Double_t TMVA::Reader::GetProba( const TString& methodTag, Double_t ap_sig, Double_t mvaVal )
721 {
722  IMethod* method = 0;
723  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
724  if (it == fMethodMap.end()) {
725  for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "M" << it->first << Endl;
726  Log() << kFATAL << "<EvaluateMVA> unknown classifier in map: " << method << "; "
727  << "you looked for " << methodTag<< " while the available methods are : " << Endl;
728  }
729  else method = it->second;
730 
731  MethodBase* kl = dynamic_cast<MethodBase*>(method);
732  if(kl==0) return -1;
733  // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
734  // it is not again checked in each of these subsequent calls..
735  const Event* ev = kl->GetEvent();
736  for (UInt_t i=0; i<ev->GetNVariables(); i++){
737  if (TMath::IsNaN(ev->GetValue(i))) {
738  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;
739  return -999;
740  }
741  }
742 
743  if (mvaVal == -9999999) mvaVal = kl->GetMvaValue();
744 
745  return kl->GetProba( mvaVal, ap_sig );
746 }
747 
748 ////////////////////////////////////////////////////////////////////////////////
749 /// evaluates the MVA's rarity
750 
751 Double_t TMVA::Reader::GetRarity( const TString& methodTag, Double_t mvaVal )
752 {
753  IMethod* method = 0;
754  std::map<TString, IMethod*>::iterator it = fMethodMap.find( methodTag );
755  if (it == fMethodMap.end()) {
756  for (it = fMethodMap.begin(); it!=fMethodMap.end(); ++it) Log() << "M" << it->first << Endl;
757  Log() << kFATAL << "<EvaluateMVA> unknown classifier in map: \"" << method << "\"; "
758  << "you looked for \"" << methodTag<< "\" while the available methods are : " << Endl;
759  }
760  else method = it->second;
761 
762  MethodBase* kl = dynamic_cast<MethodBase*>(method);
763  if(kl==0) return -1;
764  // check for NaN in event data: (note: in the factory, this check was done already at the creation of the datasets, hence
765  // it is not again checked in each of these subsequent calls..
766  const Event* ev = kl->GetEvent();
767  for (UInt_t i=0; i<ev->GetNVariables(); i++){
768  if (TMath::IsNaN(ev->GetValue(i))) {
769  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;
770  return -999;
771  }
772  }
773 
774  if (mvaVal == -9999999) mvaVal = kl->GetMvaValue();
775 
776  return kl->GetRarity( mvaVal );
777 }
778 
779 // ---------------------------------------------------------------------------------------
780 // ----- methods related to the decoding of the input variable names ---------------------
781 // ---------------------------------------------------------------------------------------
782 
783 ////////////////////////////////////////////////////////////////////////////////
784 /// decodes "name1:name2:..." form
785 
786 void TMVA::Reader::DecodeVarNames( const std::string& varNames )
787 {
788  size_t ipos = 0, f = 0;
789  while (f != varNames.length()) {
790  f = varNames.find( ':', ipos );
791  if (f > varNames.length()) f = varNames.length();
792  std::string subs = varNames.substr( ipos, f-ipos ); ipos = f+1;
793  DataInfo().AddVariable( subs.c_str() );
794  }
795 }
796 
797 ////////////////////////////////////////////////////////////////////////////////
798 /// decodes "name1:name2:..." form
799 
800 void TMVA::Reader::DecodeVarNames( const TString& varNames )
801 {
802  TString format;
803  Int_t n = varNames.Length();
804  TString format_obj;
805 
806  for (int i=0; i< n+1 ; i++) {
807  format.Append(varNames(i));
808  if (varNames(i) == ':' || i == n) {
809  format.Chop();
810  format_obj = format;
811  format_obj.ReplaceAll("@","");
812  DataInfo().AddVariable( format_obj );
813  format.Resize(0);
814  }
815  }
816 }