CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
RecoTauMVAHelper.cc
Go to the documentation of this file.
3 
4 #include <boost/foreach.hpp>
5 #include <boost/bind.hpp>
6 #include <sstream>
7 
11 
13 
14 namespace reco { namespace tau {
15 
17  const std::string& eslabel,
18  const edm::ParameterSet& pluginOptions):
19  name_(name), eslabel_(eslabel), pluginOptions_(pluginOptions) {}
20 
22  const edm::EventSetup &es) {
23  // Update our MVA from the DB
26  if (eslabel_.size()) {
27  es.get<TauTagMVAComputerRcd>().get(eslabel_, handle);
28  } else {
29  es.get<TauTagMVAComputerRcd>().get(handle);
30  }
31  const MVAComputerContainer *container = handle.product();
32  // Load our MVA
33  bool reload = computer_.update(container, name_.c_str());
34  // If the MVA changed, update our list of discriminant plugins
35  if (reload && computer_.get())
36  loadDiscriminantPlugins(container->find(name_));
37  // Update the event info for all of our discriminators
38  BOOST_FOREACH(PluginMap::value_type plugin, plugins_) {
39  plugin.second->setup(evt, es);
40  }
41 }
42 
45  typedef std::vector<PhysicsTools::Calibration::Variable> VarList;
46  // List of input variables for this MVA.
47  const VarList &vars = comp.inputSet;
48  // Load the plugin for each of the Var if needed
49  BOOST_FOREACH(const VarList::value_type& var, vars) {
50  // Check to make sure it isn't a magic variable
51  if (std::strncmp(var.name.c_str(), "__", 2) != 0) {
52  // If we haven't added yet, build it.
53  PhysicsTools::AtomicId varId(var.name);
54  if (!plugins_.count(varId)) {
56  if (pluginOptions_.exists(var.name)) {
57  options = pluginOptions_.getParameter<edm::ParameterSet>(var.name);
58  };
59  // Make sure it has a name (required by base class)
60  if (!options.exists("name"))
61  options.addParameter("name", "MVA_" + var.name);
62  // Check if we want to specify the plugin name manually. This is
63  // required for things like the discriminant from discriminators, which
64  // take an InputTag. If we want to have more than one, we have to be
65  // able take the MVA name (like FlightPathSig) and map it to
66  // RecoTauDiscriminantFromDiscriminator[disc input tag = flight path sig]
67  // Otherwise we just keep our regular plugin mapping.
68  std::string pluginName = reco::tau::discPluginName(var.name);
69  if (options.exists("plugin")) {
70  pluginName = options.getParameter<std::string>("plugin");
71  }
72  plugins_.insert(
74  pluginName, options));
75  }
76  }
77  }
78 }
79 
81  // Loop over the relevant discriminators and the output
82  for (PluginMap::const_iterator plugin = plugins_.begin();
83  plugin != plugins_.end(); ++plugin) {
84  PhysicsTools::AtomicId id = plugin->first;
85  std::vector<double> pluginOutput = (plugin->second)->operator()(tau);
86  // Check for nans
87  for(size_t instance = 0; instance < pluginOutput.size(); ++instance) {
88  if (std::isnan(pluginOutput[instance])) {
89  std::ostringstream error;
90  error << "A nan was detected in"
91  << " the tau MVA variable " << id << " returning zero instead!"
92  << " The PFTau: " << *tau << std::endl;
93  tau->dump(error);
94  edm::LogError("CorruptedMVAInput") << error.str();
95  pluginOutput[instance] = 0.0;
96  }
97  }
98  //std::cout << "id: " << id << " first: " << pluginOutput[0] << std::endl;
99  // Build values and copy into values vector
100  std::for_each(pluginOutput.begin(), pluginOutput.end(),
102  boost::ref(values_), id, _1));
103  }
104 }
105 
106 // Get values
109  values_.clear();
110  fillValues(tau);
111  return values_;
112 }
113 
114 // Apply the MVA to a given tau
116  // Clear output
117  values_.clear();
118  // Build the values
119  fillValues(tau);
120  // Call the MVA
121  return computer_->eval(values_);
122 }
123 
125  double weight) const {
126  static const PhysicsTools::AtomicId kTargetId("__TARGET__");
127  static const PhysicsTools::AtomicId kWeightId("__WEIGHT__");
128  if (!computer_)
129  return;
130  values_.clear();
131  values_.add(kTargetId, target);
132  values_.add(kWeightId, weight);
133  // Build the discriminant values
134  fillValues(tau);
136 }
137 
138 }} // end namespace reco::tau
T getParameter(std::string const &) const
void setEvent(const edm::Event &evt, const edm::EventSetup &es)
double eval(Iterator_t first, Iterator_t last) const
evaluate variables given by a range of iterators given by first and last
void fillValues(const reco::PFTauRef &tau) const
bool exists(std::string const &parameterName) const
checks if a parameter exists
Cheap generic unique keyword identifier class.
Definition: AtomicId.h:32
std::string discPluginName(const std::string &mvaName)
PhysicsTools::Variable::ValueList values_
bool isnan(float x)
Definition: math.h:13
void addParameter(std::string const &name, T const &value)
Definition: ParameterSet.h:139
tuple handle
Definition: patZpeak.py:22
double operator()(const PFTauRef &tau) const
edm::ParameterSet pluginOptions_
const PhysicsTools::Variable::ValueList & discriminants(const PFTauRef &tau) const
Container::value_type value_type
Helper class that can contain an list of identifier-value pairs.
Definition: Variable.h:82
RecoTauMVAHelper(const std::string &name, const std::string &eslabel, const edm::ParameterSet &pluginOptions)
bool update(const Calibration::MVAComputer *computer)
void train(const PFTauRef &tau, bool target, double weight=1.0) const
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:62
std::vector< Variable > inputSet
Definition: MVAComputer.h:177
PhysicsTools::MVAComputerCache computer_
void loadDiscriminantPlugins(const PhysicsTools::Calibration::MVAComputer &computer)
SurfaceDeformation * create(int type, const std::vector< double > &params)
T get(const Candidate &c)
Definition: component.h:56
void add(AtomicId id, double value)
Definition: Variable.h:107