CMS 3D CMS Logo

EcalDeadChannelRecoveryBDTG.cc
Go to the documentation of this file.
1 // Original Authors: S. Taroni, N. Marinelli
2 // University of Notre Dame - US
3 // Created:
4 //
5 //
6 //
7 
9 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h" // can I use a egammatools here?
11 #include <iostream>
12 
13 #include <iostream>
14 #include <memory>
15 
16 #include <fstream>
17 #include <ostream>
18 #include <string>
19 
20 template <>
22  for (int i = 0; i < 9; ++i) {
23  if (i == 4)
24  continue;
25  reader->AddVariable("E" + std::to_string(i + 1) + "/(E1+E2+E3+E4+E6+E7+E8+E9)", &(mx_.rEn[i]));
26  }
27  reader->AddVariable("E1+E2+E3+E4+E6+E7+E8+E9", &(mx_.sumE8));
28  for (int i = 0; i < 9; ++i)
29  reader->AddVariable("iEta" + std::to_string(i + 1), &(mx_.ieta[i]));
30  for (int i = 0; i < 9; ++i)
31  reader->AddVariable("iPhi" + std::to_string(i + 1), &(mx_.iphi[i]));
32 }
33 template <>
35  readerNoCrack = std::make_unique<TMVA::Reader>("!Color:!Silent");
36  readerCrack = std::make_unique<TMVA::Reader>("!Color:!Silent");
37 
40 
43 }
44 
45 template <typename T>
47 
48 template <typename T>
50 
51 template <>
53  bdtWeightFileNoCracks_ = ps.getParameter<edm::FileInPath>("bdtWeightFileNoCracks");
54  bdtWeightFileCracks_ = ps.getParameter<edm::FileInPath>("bdtWeightFileCracks");
55 
56  loadFile();
57 }
58 
59 template <>
61 
62 template <>
64  const EEDetId id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag) {
65  return 0;
66 }
67 
68 template <>
70  const EBDetId id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag) {
71  bool isCrack = false;
72  int cellIndex = 0.;
73  double neighTotEn = 0.;
74  float val = 0.;
75 
76  //find the matrix around id
77  std::vector<DetId> m3x3aroundDC = EcalClusterTools::matrixDetId(topology_, id, 1);
78  if (m3x3aroundDC.size() < 9) {
79  *acceptFlag = false;
80  return 0;
81  }
82 
83  // Loop over all cells in the vector "NxNaroundDC", and for each cell find it's energy
84  // (from the EcalRecHits collection).
85  for (auto const &theCells : m3x3aroundDC) {
86  EBDetId cell = EBDetId(theCells);
87  if (cell == id) {
88  int iEtaCentral = std::abs(cell.ieta());
89  int iPhiCentral = cell.iphi();
90 
91  if (iEtaCentral < 2 || std::abs(iEtaCentral - 25) < 2 || std::abs(iEtaCentral - 45) < 2 ||
92  std::abs(iEtaCentral - 65) < 2 || iEtaCentral > 83 || (int(iPhiCentral + 0.5) % 20 == 0))
93  isCrack = true;
94  }
95  if (!cell.null()) {
96  EcalRecHitCollection::const_iterator goS_it = hit_collection.find(cell);
97  if (goS_it != hit_collection.end() && cell != id) {
98  if (goS_it->energy() < single8Cut) {
99  *acceptFlag = false;
100  return 0.;
101  } else {
102  neighTotEn += goS_it->energy();
103  mx_.rEn[cellIndex] = goS_it->energy();
104  mx_.iphi[cellIndex] = cell.iphi();
105  mx_.ieta[cellIndex] = cell.ieta();
106  cellIndex++;
107  }
108  } else if (cell == id) { // the cell is the central one
109  mx_.rEn[cellIndex] = 0;
110  cellIndex++;
111  } else { //goS_it is not in the rechitcollection
112  *acceptFlag = false;
113  return 0.;
114  }
115  } else { //cell is null
116  *acceptFlag = false;
117  return 0.;
118  }
119  }
120  if (cellIndex > 0 && neighTotEn >= single8Cut * 8. && neighTotEn >= sum8Cut) {
121  bool allneighs = true;
122  mx_.sumE8 = neighTotEn;
123  for (unsigned int icell = 0; icell < 9; icell++) {
124  if (mx_.rEn[icell] < single8Cut && icell != 4) {
125  allneighs = false;
126  }
127  mx_.rEn[icell] = mx_.rEn[icell] / neighTotEn;
128  }
129  if (allneighs == true) {
130  // evaluate the regression
131  if (isCrack) {
132  val = exp((readerCrack->EvaluateRegression("BDTG"))[0]);
133  } else {
134  val = exp((readerNoCrack->EvaluateRegression("BDTG"))[0]);
135  }
136 
137  *acceptFlag = true;
138  //return the estimated energy
139  return val;
140  } else {
141  *acceptFlag = false;
142  return 0;
143  }
144  } else {
145  *acceptFlag = false;
146  return 0;
147  }
148 }
149 
151 template class EcalDeadChannelRecoveryBDTG<EEDetId>; //not used.
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
std::string fullPath() const
Definition: FileInPath.cc:161
void setParameters(const edm::ParameterSet &ps)
int iphi() const
get the crystal iphi
Definition: EBDetId.h:51
std::vector< EcalRecHit >::const_iterator const_iterator
std::string to_string(const V &value)
Definition: OMSAccess.h:71
reader
Definition: DQM.py:105
int ieta() const
get the crystal ieta
Definition: EBDetId.h:49
constexpr bool null() const
is this a null id ?
Definition: DetId.h:59
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::unique_ptr< TMVA::Reader > readerCrack
const_iterator end() const
TMVA::IMethod * loadTMVAWeights(TMVA::Reader *reader, const std::string &method, const std::string &weightFile, bool verbose=false)
double recover(const DetIdT id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag)
iterator find(key_type k)
std::unique_ptr< TMVA::Reader > readerNoCrack
void addVariables(TMVA::Reader *reader)