Skip to content

Commit

Permalink
Cherry pick files from PR27175 for Backport
Browse files Browse the repository at this point in the history
  • Loading branch information
nancymarinelli authored and Rajdeep Mohan Chatterjee committed Jan 30, 2020
1 parent 781d8a3 commit 09ae771
Show file tree
Hide file tree
Showing 10 changed files with 986 additions and 786 deletions.
1 change: 1 addition & 0 deletions RecoLocalCalo/EcalDeadChannelRecoveryAlgos/BuildFile.xml
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
<use name="Geometry/CaloGeometry"/>
<use name="Geometry/CaloTopology"/>
<use name="Geometry/Records"/>
<use name="CommonTools/MVAUtils"/>
<use name="root"/>
<use name="rootmath"/>
<use name="rootcore"/>
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -7,16 +7,25 @@
#include "DataFormats/EcalDetId/interface/EBDetId.h"
#include <string>

#include "RecoLocalCalo/EcalDeadChannelRecoveryAlgos/interface/EcalDeadChannelRecoveryNN.h"
#include "FWCore/ParameterSet/interface/ParameterSetDescription.h"

template <typename DetIdT> class EcalDeadChannelRecoveryAlgos {
public:
#include "RecoLocalCalo/EcalDeadChannelRecoveryAlgos/interface/EcalDeadChannelRecoveryBDTG.h"

#include "FWCore/ParameterSet/interface/ParameterSet.h"

template <typename DetIdT>
class EcalDeadChannelRecoveryAlgos {
public:
void setParameters(const edm::ParameterSet &ps);
void setCaloTopology(const CaloTopology *topology);
EcalRecHit correct(const DetIdT id,
const EcalRecHitCollection &hit_collection,
std::string algo, double Sum8Cut, bool *AccFlag);
float correct(const DetIdT id,
const EcalRecHitCollection &hit_collection,
std::string algo,
double single8Cut,
double sum8Cut,
bool *accFlag);

private:
EcalDeadChannelRecoveryNN<DetIdT> nn;
private:
EcalDeadChannelRecoveryBDTG<DetIdT> bdtg_;
};
#endif // RecoLocalCalo_EcalDeadChannelRecoveryAlgos_EcalDeadChannelRecoveryAlgos_HH
#endif
Original file line number Diff line number Diff line change
@@ -0,0 +1,52 @@
#ifndef RecoLocalCalo_EcalDeadChannelRecoveryAlgos_EcalDeadChannelRecoveryBDTG_H
#define RecoLocalCalo_EcalDeadChannelRecoveryAlgos_EcalDeadChannelRecoveryBDTG_H

#include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"

#include "DataFormats/EcalDetId/interface/EBDetId.h"
#include "DataFormats/EcalDetId/interface/EEDetId.h"

#include "Geometry/CaloTopology/interface/CaloTopology.h"
#include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"

#include "CommonTools/MVAUtils/interface/TMVAZipReader.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "TMVA/Reader.h"

#include <string>
#include <memory>

template <typename DetIdT>
class EcalDeadChannelRecoveryBDTG {
public:
EcalDeadChannelRecoveryBDTG();
~EcalDeadChannelRecoveryBDTG();

void setParameters(const edm::ParameterSet &ps);
void setCaloTopology(const CaloTopology *topo) { topology_ = topo; };

double recover(
const DetIdT id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag);

void loadFile();
void addVariables(TMVA::Reader *reader);

private:
const CaloTopology *topology_;
struct XtalMatrix {
std::array<float, 9> rEn, ieta, iphi;
float sumE8;
};

XtalMatrix mx_;

edm::FileInPath bdtWeightFileNoCracks_;
edm::FileInPath bdtWeightFileCracks_;

TMVA::Reader *readerNoCrack;
TMVA::Reader *readerCrack;
};

#endif
Original file line number Diff line number Diff line change
Expand Up @@ -8,33 +8,41 @@
// Feb 14 2013: Implementation of the criterion to select the "correct"
// max. cont. crystal.
//
//modified by S.Taroni, N. Marinelli 11 June 2019

#include "RecoLocalCalo/EcalDeadChannelRecoveryAlgos/interface/EcalDeadChannelRecoveryAlgos.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

template <typename T>
void EcalDeadChannelRecoveryAlgos<T>::setCaloTopology(
const CaloTopology *topo) {
nn.setCaloTopology(topo);
void EcalDeadChannelRecoveryAlgos<T>::setParameters(const edm::ParameterSet &ps) {
bdtg_.setParameters(ps);
}

template <typename T>
EcalRecHit EcalDeadChannelRecoveryAlgos<T>::correct(
const T id, const EcalRecHitCollection &hit_collection, std::string algo,
double Sum8Cut, bool *AcceptFlag) {
// recover as single dead channel
double NewEnergy = 0.0;
void EcalDeadChannelRecoveryAlgos<T>::setCaloTopology(const CaloTopology *topo) {
bdtg_.setCaloTopology(topo);
}

if (algo == "NeuralNetworks") {
NewEnergy = this->nn.recover(id, hit_collection, Sum8Cut, AcceptFlag);
template <typename T>
float EcalDeadChannelRecoveryAlgos<T>::correct(const T id,
const EcalRecHitCollection &hit_collection,
std::string algo,
double single8Cut,
double sum8Cut,
bool *acceptFlag) {
// recover as single dead channel
double newEnergy = 0.0;
if (algo == "BDTG") {
*acceptFlag = false;
newEnergy = this->bdtg_.recover(id, hit_collection, single8Cut, sum8Cut, acceptFlag); //ADD here
if (newEnergy > 0.)
*acceptFlag = true; //bdtg set to 0 if there is more than one channel in the matrix that is not reponding
} else {
edm::LogError("EcalDeadChannelRecoveryAlgos")
<< "Invalid algorithm for dead channel recovery.";
*AcceptFlag = false;
edm::LogError("EcalDeadChannelRecoveryAlgos") << "Invalid algorithm for dead channel recovery.";
*acceptFlag = false;
}

uint32_t flag = 0;
return EcalRecHit(id, NewEnergy, 0, flag);
return newEnergy;
}

template class EcalDeadChannelRecoveryAlgos<EBDetId>;
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,150 @@
// Original Authors: S. Taroni, N. Marinelli
// University of Notre Dame - US
// Created:
//
//
//

#include "RecoLocalCalo/EcalDeadChannelRecoveryAlgos/interface/EcalDeadChannelRecoveryBDTG.h"
#include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h" // can I use a egammatools here?
#include "FWCore/ParameterSet/interface/FileInPath.h"
#include <iostream>

#include <iostream>
#include <ostream>
#include <string>
#include <fstream>

template <>
void EcalDeadChannelRecoveryBDTG<EBDetId>::addVariables(TMVA::Reader *reader) {
for (int i = 0; i < 9; ++i) {
if (i == 4)
continue;
reader->AddVariable("E" + std::to_string(i + 1) + "/(E1+E2+E3+E4+E6+E7+E8+E9)", &(mx_.rEn[i]));
}
reader->AddVariable("E1+E2+E3+E4+E6+E7+E8+E9", &(mx_.sumE8));
for (int i = 0; i < 9; ++i)
reader->AddVariable("iEta" + std::to_string(i + 1), &(mx_.ieta[i]));
for (int i = 0; i < 9; ++i)
reader->AddVariable("iPhi" + std::to_string(i + 1), &(mx_.iphi[i]));
}
template <>
void EcalDeadChannelRecoveryBDTG<EBDetId>::loadFile() {
readerNoCrack = new TMVA::Reader("!Color:!Silent");
readerCrack = new TMVA::Reader("!Color:!Silent");

this->addVariables(readerNoCrack);
this->addVariables(readerCrack);

reco::details::loadTMVAWeights(readerNoCrack, "BDTG", bdtWeightFileNoCracks_.fullPath());
reco::details::loadTMVAWeights(readerCrack, "BDTG", bdtWeightFileCracks_.fullPath());
}

template <typename T>
EcalDeadChannelRecoveryBDTG<T>::EcalDeadChannelRecoveryBDTG() {}

template <typename T>
EcalDeadChannelRecoveryBDTG<T>::~EcalDeadChannelRecoveryBDTG() {}

template <>
void EcalDeadChannelRecoveryBDTG<EBDetId>::setParameters(const edm::ParameterSet &ps) {
bdtWeightFileNoCracks_ = ps.getParameter<edm::FileInPath>("bdtWeightFileNoCracks");
bdtWeightFileCracks_ = ps.getParameter<edm::FileInPath>("bdtWeightFileCracks");

this->loadFile();
}

template <>
void EcalDeadChannelRecoveryBDTG<EEDetId>::setParameters(const edm::ParameterSet &ps) {}

template <>
double EcalDeadChannelRecoveryBDTG<EEDetId>::recover(
const EEDetId id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag) {
return 0;
}

template <>
double EcalDeadChannelRecoveryBDTG<EBDetId>::recover(
const EBDetId id, const EcalRecHitCollection &hit_collection, double single8Cut, double sum8Cut, bool *acceptFlag) {
bool isCrack = false;
int cellIndex = 0.;
double neighTotEn = 0.;
float val = 0.;

//find the matrix around id
std::vector<DetId> m3x3aroundDC = EcalClusterTools::matrixDetId(topology_, id, 1);
if (m3x3aroundDC.size() < 9) {
*acceptFlag = false;
return 0;
}

// Loop over all cells in the vector "NxNaroundDC", and for each cell find it's energy
// (from the EcalRecHits collection).
for (auto const &theCells : m3x3aroundDC) {
EBDetId cell = EBDetId(theCells);
if (cell == id) {

int iEtaCentral = std::abs(cell.ieta());
int iPhiCentral = cell.iphi();

if (iEtaCentral < 2 || std::abs(iEtaCentral - 25) < 2 || std::abs(iEtaCentral - 45) < 2 ||
std::abs(iEtaCentral - 65) < 2 || iEtaCentral > 83 || (int(iPhiCentral + 0.5) % 20 == 0))
isCrack = true;
}
if (!cell.null()) {
EcalRecHitCollection::const_iterator goS_it = hit_collection.find(cell);
if (goS_it != hit_collection.end() && cell!=id) {
if (goS_it->energy() < single8Cut) {
*acceptFlag = false;
return 0.;
} else {
neighTotEn += goS_it->energy();
mx_.rEn[cellIndex] = goS_it->energy();
mx_.iphi[cellIndex] = cell.iphi();
mx_.ieta[cellIndex] = cell.ieta();
cellIndex++;
}
} else if (cell==id) { // the cell is the central one
mx_.rEn[cellIndex] = 0;
cellIndex++;
}else { //goS_it is not in the rechitcollection
*acceptFlag = false;
return 0.;
}
} else { //cell is null
*acceptFlag = false;
return 0.;
}
}
if (cellIndex > 0 && neighTotEn >= single8Cut * 8. && neighTotEn >= sum8Cut) {
bool allneighs = true;
mx_.sumE8 = neighTotEn;
for (unsigned int icell = 0; icell < 9; icell++) {
if (mx_.rEn[icell] < single8Cut && icell != 4) {
allneighs = false;
}
mx_.rEn[icell] = mx_.rEn[icell] / neighTotEn;
}
if (allneighs == true) {
// evaluate the regression
if (isCrack) {
val = exp((readerCrack->EvaluateRegression("BDTG"))[0]);
} else {
val = exp((readerNoCrack->EvaluateRegression("BDTG"))[0]);
}

*acceptFlag = true;
//return the estimated energy
return val;
} else {
*acceptFlag = false;
return 0;
}
} else {
*acceptFlag = false;
return 0;
}
}

template class EcalDeadChannelRecoveryBDTG<EBDetId>;
template class EcalDeadChannelRecoveryBDTG<EEDetId>; //not used.
Loading

0 comments on commit 09ae771

Please sign in to comment.