CMS 3D CMS Logo

PhotonMVAEstimator.cc
Go to the documentation of this file.
4 
7  , mvaVarMngr_ (conf.getParameter<std::string>("variableDefinition"))
8 {
9 
10  //
11  // Construct the MVA estimators
12  //
13  if (getTag() == "Run2Spring16NonTrigV1") {
14  effectiveAreas_ = std::make_unique<EffectiveAreas>((conf.getParameter<edm::FileInPath>("effAreasConfigFile")).fullPath());
15  phoIsoPtScalingCoeff_ = conf.getParameter<std::vector<double >>("phoIsoPtScalingCoeff");
16  phoIsoCutoff_ = conf.getParameter<double>("phoIsoCutoff");
17  }
18 
19  const auto weightFileNames = conf.getParameter<std::vector<std::string> >("weightFileNames");
20  const auto categoryCutStrings = conf.getParameter<std::vector<std::string> >("categoryCuts");
21 
22  if( (int)(categoryCutStrings.size()) != getNCategories() )
23  throw cms::Exception("MVA config failure: ")
24  << "wrong number of category cuts in PhotonMVAEstimator" << getTag() << std::endl;
25 
26  for (auto const& cut : categoryCutStrings) categoryFunctions_.emplace_back(cut);
27 
28  // Initialize GBRForests
29  if( static_cast<int>(weightFileNames.size()) != getNCategories() )
30  throw cms::Exception("MVA config failure: ")
31  << "wrong number of weightfiles in PhotonMVAEstimator" << getTag() << std::endl;
32 
33  gbrForests_.clear();
34  // Create a TMVA reader object for each category
35  for(int i=0; i<getNCategories(); i++){
36 
37  std::vector<int> variablesInCategory;
38 
39  std::vector<std::string> variableNamesInCategory;
40  gbrForests_.push_back(createGBRForest(weightFileNames[i], variableNamesInCategory));
41 
42  nVariables_.push_back(variableNamesInCategory.size());
43 
44  variables_.push_back(variablesInCategory);
45 
46  for (int j=0; j<nVariables_[i];++j) {
47  int index = mvaVarMngr_.getVarIndex(variableNamesInCategory[j]);
48  if(index == -1) {
49  throw cms::Exception("MVA config failure: ")
50  << "Concerning PhotonMVAEstimator" << getTag() << std::endl
51  << "Variable " << variableNamesInCategory[j]
52  << " not found in variable definition file!" << std::endl;
53  }
54  variables_[i].push_back(index);
55 
56  }
57  }
58 }
59 
61 mvaValue(const reco::Candidate* candPtr, std::vector<float> const& auxVars, int &iCategory) const {
62 
63  const reco::Photon* phoPtr = dynamic_cast<const reco::Photon*>(candPtr);
64  if( phoPtr == nullptr) {
65  throw cms::Exception("MVA failure: ")
66  << " given particle is expected to be reco::Photon or pat::Photon," << std::endl
67  << " but appears to be neither" << std::endl;
68  }
69 
70  iCategory = findCategory( phoPtr );
71 
72  std::vector<float> vars;
73 
74  for (int i = 0; i < nVariables_[iCategory]; ++i) {
75  vars.push_back(mvaVarMngr_.getValue(variables_[iCategory][i], *phoPtr, auxVars));
76  }
77 
78  // Special case for Spring16!
79  if (getTag() == "Run2Spring16NonTrigV1" and iCategory == 1) { // Endcap category
80  // Raw value for EB only, because of loss of transparency in EE
81  // for endcap MVA only in 2016
82  double eA = effectiveAreas_->getEffectiveArea( std::abs(phoPtr->superCluster()->eta()) );
83  double phoIsoCorr = vars[10] - eA*(double)vars[9] - phoIsoPtScalingCoeff_.at(1) * phoPtr->pt();
84  vars[10] = TMath::Max( phoIsoCorr, phoIsoCutoff_);
85  }
86 
87  if(isDebug()) {
88  std::cout << " *** Inside PhotonMVAEstimator" << getTag() << std::endl;
89  std::cout << " category " << iCategory << std::endl;
90  for (int i = 0; i < nVariables_[iCategory]; ++i) {
91  std::cout << " " << mvaVarMngr_.getName(variables_[iCategory][i]) << " " << vars[i] << std::endl;
92  }
93  }
94 
95  const float response = gbrForests_.at(iCategory)->GetResponse(vars.data());
96 
97  if(isDebug()) {
98  std::cout << " ### MVA " << response << std::endl << std::endl;
99  }
100 
101  return response;
102 }
103 
105 
106  const reco::Photon* phoPtr = dynamic_cast<const reco::Photon*>(candPtr);
107  if( phoPtr == nullptr ) {
108  throw cms::Exception("MVA failure: ")
109  << " given particle is expected to be reco::Photon or pat::Photon," << std::endl
110  << " but appears to be neither" << std::endl;
111  }
112 
113  return findCategory(*phoPtr);
114 
115 }
116 
118 
119  for (int i = 0; i < getNCategories(); ++i) {
120  if (categoryFunctions_[i](photon)) return i;
121  }
122 
123  edm::LogWarning ("MVA warning") <<
124  "category not defined for particle with pt " << photon.pt() << " GeV, eta " <<
125  photon.superCluster()->eta() << " in PhotonMVAEstimator" << getTag();
126 
127  return -1;
128 
129 }
130 
133  "PhotonMVAEstimator");
T getParameter(std::string const &) const
const std::string & getName(int index) const
float mvaValue(const reco::Candidate *candPtr, std::vector< float > const &auxVars, int &iCategory) const override
PhotonMVAEstimator(const edm::ParameterSet &conf)
MVAVariableManager< reco::Photon > mvaVarMngr_
double pt() const final
transverse momentum
float getValue(int index, const ParticleType &particle, const std::vector< float > &auxVariables) const
reco::SuperClusterRef superCluster() const override
Ref to SuperCluster.
const std::string & getTag() const
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::vector< std::unique_ptr< const GBRForest > > gbrForests_
T Max(T a, T b)
Definition: MathUtil.h:44
std::vector< int > nVariables_
std::vector< double > phoIsoPtScalingCoeff_
int getVarIndex(const std::string &name)
std::unique_ptr< EffectiveAreas > effectiveAreas_
#define DEFINE_EDM_PLUGIN(factory, type, name)
std::vector< std::vector< int > > variables_
vars
Definition: DeepTauId.cc:77
std::vector< ThreadSafeStringCut< StringCutObjectSelector< reco::Photon >, reco::Photon > > categoryFunctions_
int findCategory(const reco::Candidate *candPtr) const override
std::unique_ptr< const GBRForest > createGBRForest(const std::string &weightsFile)