File indexing completed on 2021-02-14 14:24:22
0001
0002
0003
0004
0005
0006
0007
0008 #include "RecoLocalCalo/EcalRecProducers/plugins/EcalTPSkimmer.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "FWCore/Framework/interface/ESHandle.h"
0012
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014
0015 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0016 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0017 #include "DataFormats/EcalDetId/interface/EcalTrigTowerDetId.h"
0018
0019 EcalTPSkimmer::EcalTPSkimmer(const edm::ParameterSet& ps) {
0020 skipModule_ = ps.getParameter<bool>("skipModule");
0021
0022 doBarrel_ = ps.getParameter<bool>("doBarrel");
0023 doEndcap_ = ps.getParameter<bool>("doEndcap");
0024
0025 chStatusToSelectTP_ = ps.getParameter<std::vector<uint32_t> >("chStatusToSelectTP");
0026
0027 tpOutputCollection_ = ps.getParameter<std::string>("tpOutputCollection");
0028 tpInputToken_ = consumes<EcalTrigPrimDigiCollection>(ps.getParameter<edm::InputTag>("tpInputCollection"));
0029 ttMapToken_ = esConsumes<EcalTrigTowerConstituentsMap, IdealGeometryRecord>();
0030 if (not skipModule_) {
0031 chStatusToken_ = esConsumes<EcalChannelStatus, EcalChannelStatusRcd>();
0032 }
0033 produces<EcalTrigPrimDigiCollection>(tpOutputCollection_);
0034 }
0035
0036 EcalTPSkimmer::~EcalTPSkimmer() {}
0037
0038 void EcalTPSkimmer::produce(edm::Event& evt, const edm::EventSetup& es) {
0039 insertedTP_.clear();
0040
0041 using namespace edm;
0042
0043 ttMap_ = es.getHandle(ttMapToken_);
0044
0045
0046 auto tpOut = std::make_unique<EcalTrigPrimDigiCollection>();
0047
0048 if (skipModule_) {
0049 evt.put(std::move(tpOut), tpOutputCollection_);
0050 return;
0051 }
0052
0053 edm::ESHandle<EcalChannelStatus> chStatus = es.getHandle(chStatusToken_);
0054
0055 edm::Handle<EcalTrigPrimDigiCollection> tpIn;
0056 evt.getByToken(tpInputToken_, tpIn);
0057
0058 if (doBarrel_) {
0059 EcalChannelStatusMap::const_iterator chit;
0060 uint16_t code = 0;
0061 for (int i = 0; i < EBDetId::kSizeForDenseIndexing; ++i) {
0062 if (!EBDetId::validDenseIndex(i))
0063 continue;
0064 EBDetId id = EBDetId::detIdFromDenseIndex(i);
0065 chit = chStatus->find(id);
0066
0067 if (chit != chStatus->end()) {
0068 code = (*chit).getStatusCode();
0069 if (std::find(chStatusToSelectTP_.begin(), chStatusToSelectTP_.end(), code) != chStatusToSelectTP_.end()) {
0070
0071 EcalTrigTowerDetId ttDetId(((EBDetId)id).tower());
0072
0073 if (!alreadyInserted(ttDetId))
0074 insertTP(ttDetId, tpIn, *tpOut);
0075 }
0076 } else {
0077 edm::LogError("EcalDetIdToBeRecoveredProducer") << "No channel status found for xtal " << id.rawId()
0078 << "! something wrong with EcalChannelStatus in your DB? ";
0079 }
0080 }
0081 }
0082
0083 if (doEndcap_) {
0084 EcalChannelStatusMap::const_iterator chit;
0085 uint16_t code = 0;
0086 for (int i = 0; i < EEDetId::kSizeForDenseIndexing; ++i) {
0087 if (!EEDetId::validDenseIndex(i))
0088 continue;
0089 EEDetId id = EEDetId::detIdFromDenseIndex(i);
0090 chit = chStatus->find(id);
0091
0092 if (chit != chStatus->end()) {
0093 code = (*chit).getStatusCode();
0094 if (std::find(chStatusToSelectTP_.begin(), chStatusToSelectTP_.end(), code) != chStatusToSelectTP_.end()) {
0095
0096 EcalTrigTowerDetId ttDetId = ttMap_->towerOf(id);
0097
0098 if (!alreadyInserted(ttDetId))
0099 insertTP(ttDetId, tpIn, *tpOut);
0100 }
0101 } else {
0102 edm::LogError("EcalDetIdToBeRecoveredProducer") << "No channel status found for xtal " << id.rawId()
0103 << "! something wrong with EcalChannelStatus in your DB? ";
0104 }
0105 }
0106 }
0107
0108
0109 LogInfo("EcalTPSkimmer") << "total # of TP inserted: " << tpOut->size();
0110
0111 evt.put(std::move(tpOut), tpOutputCollection_);
0112 }
0113
0114 bool EcalTPSkimmer::alreadyInserted(EcalTrigTowerDetId ttId) { return (insertedTP_.find(ttId) != insertedTP_.end()); }
0115
0116 void EcalTPSkimmer::insertTP(EcalTrigTowerDetId ttId,
0117 edm::Handle<EcalTrigPrimDigiCollection>& tpIn,
0118 EcalTrigPrimDigiCollection& tpOut) {
0119 EcalTrigPrimDigiCollection::const_iterator tpIt = tpIn->find(ttId);
0120 if (tpIt != tpIn->end()) {
0121 tpOut.push_back(*tpIt);
0122 insertedTP_.insert(ttId);
0123 }
0124 }
0125
0126 #include "FWCore/Framework/interface/MakerMacros.h"
0127 DEFINE_FWK_MODULE(EcalTPSkimmer);