File indexing completed on 2024-04-06 12:25:52
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021
0022
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/stream/EDProducer.h"
0025
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030
0031 #include "FWCore/Utilities/interface/ESGetToken.h"
0032 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0033 #include "Geometry/Records/interface/HcalRecNumberingRecord.h"
0034 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0035 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0036 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0037 #include "DataFormats/DetId/interface/DetIdCollection.h"
0038
0039 #include "CondFormats/HcalObjects/interface/HcalChannelQuality.h"
0040 #include "CondFormats/DataRecord/interface/HcalChannelQualityRcd.h"
0041
0042 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputer.h"
0043 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputerRcd.h"
0044
0045
0046
0047
0048
0049 class HcalHitSelection : public edm::stream::EDProducer<> {
0050 public:
0051 explicit HcalHitSelection(const edm::ParameterSet&);
0052 ~HcalHitSelection() override;
0053
0054 private:
0055 void produce(edm::Event&, const edm::EventSetup&) override;
0056
0057 edm::InputTag hbheTag, hoTag, hfTag;
0058 edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0059 edm::EDGetTokenT<HORecHitCollection> tok_ho_;
0060 edm::EDGetTokenT<HFRecHitCollection> tok_hf_;
0061 std::vector<edm::EDGetTokenT<DetIdCollection> > toks_did_;
0062 int hoSeverityLevel;
0063 std::vector<edm::InputTag> interestingDetIdCollections;
0064 const HcalTopology* theHcalTopology_;
0065
0066
0067 const HcalChannelQuality* theHcalChStatus;
0068 const HcalSeverityLevelComputer* theHcalSevLvlComputer;
0069 std::set<DetId> toBeKept;
0070 template <typename CollectionType>
0071 void skim(const edm::Handle<CollectionType>& input, CollectionType& output, int severityThreshold = 0) const;
0072
0073
0074 edm::ESGetToken<HcalTopology, HcalRecNumberingRecord> htopoToken_;
0075 edm::ESGetToken<HcalChannelQuality, HcalChannelQualityRcd> qualToken_;
0076 edm::ESGetToken<HcalSeverityLevelComputer, HcalSeverityLevelComputerRcd> sevToken_;
0077 };
0078
0079 template <class CollectionType>
0080 void HcalHitSelection::skim(const edm::Handle<CollectionType>& input,
0081 CollectionType& output,
0082 int severityThreshold) const {
0083 output.reserve(input->size());
0084 typename CollectionType::const_iterator begin = input->begin();
0085 typename CollectionType::const_iterator end = input->end();
0086 typename CollectionType::const_iterator hit = begin;
0087
0088 for (; hit != end; ++hit) {
0089
0090 HcalDetId id = hit->detid();
0091 if (theHcalTopology_->getMergePositionFlag() && id.subdet() == HcalEndcap) {
0092 id = theHcalTopology_->idFront(id);
0093 }
0094 const uint32_t& recHitFlag = hit->flags();
0095
0096
0097 const uint32_t& dbStatusFlag = theHcalChStatus->getValues(id)->getValue();
0098 int severityLevel = theHcalSevLvlComputer->getSeverityLevel(id, recHitFlag, dbStatusFlag);
0099
0100 if (severityLevel > severityThreshold) {
0101 output.push_back(*hit);
0102 } else {
0103
0104 if (toBeKept.find(id) != toBeKept.end())
0105 output.push_back(*hit);
0106 }
0107 }
0108 }
0109
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119
0120
0121 HcalHitSelection::HcalHitSelection(const edm::ParameterSet& iConfig)
0122 : hbheTag(iConfig.getParameter<edm::InputTag>("hbheTag")),
0123 hoTag(iConfig.getParameter<edm::InputTag>("hoTag")),
0124 hfTag(iConfig.getParameter<edm::InputTag>("hfTag")),
0125 theHcalTopology_(nullptr),
0126 theHcalChStatus(nullptr),
0127 theHcalSevLvlComputer(nullptr) {
0128
0129 tok_hbhe_ = consumes<HBHERecHitCollection>(hbheTag);
0130 tok_hf_ = consumes<HFRecHitCollection>(hfTag);
0131 tok_ho_ = consumes<HORecHitCollection>(hoTag);
0132
0133 interestingDetIdCollections = iConfig.getParameter<std::vector<edm::InputTag> >("interestingDetIds");
0134
0135 const unsigned nLabels = interestingDetIdCollections.size();
0136 for (unsigned i = 0; i != nLabels; i++)
0137 toks_did_.push_back(consumes<DetIdCollection>(interestingDetIdCollections[i]));
0138
0139 hoSeverityLevel = iConfig.getParameter<int>("hoSeverityLevel");
0140
0141 produces<HBHERecHitCollection>(hbheTag.label());
0142 produces<HFRecHitCollection>(hfTag.label());
0143 produces<HORecHitCollection>(hoTag.label());
0144
0145
0146 htopoToken_ = esConsumes<HcalTopology, HcalRecNumberingRecord>();
0147 qualToken_ = esConsumes<HcalChannelQuality, HcalChannelQualityRcd>(edm::ESInputTag("", "withTopo"));
0148 sevToken_ = esConsumes<HcalSeverityLevelComputer, HcalSeverityLevelComputerRcd>();
0149 }
0150
0151 HcalHitSelection::~HcalHitSelection() {
0152
0153
0154 }
0155
0156
0157
0158
0159
0160
0161 void HcalHitSelection::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0162 theHcalChStatus = &iSetup.getData(qualToken_);
0163 theHcalSevLvlComputer = &iSetup.getData(sevToken_);
0164 theHcalTopology_ = &iSetup.getData(htopoToken_);
0165
0166 edm::Handle<HBHERecHitCollection> hbhe;
0167 edm::Handle<HFRecHitCollection> hf;
0168 edm::Handle<HORecHitCollection> ho;
0169
0170 iEvent.getByToken(tok_hbhe_, hbhe);
0171 iEvent.getByToken(tok_hf_, hf);
0172 iEvent.getByToken(tok_ho_, ho);
0173
0174 toBeKept.clear();
0175 edm::Handle<DetIdCollection> detId;
0176 for (unsigned int t = 0; t < toks_did_.size(); ++t) {
0177 iEvent.getByToken(toks_did_[t], detId);
0178 if (!detId.isValid()) {
0179 edm::LogError("MissingInput") << "the collection of interesting detIds:" << interestingDetIdCollections[t]
0180 << " is not found.";
0181 continue;
0182 }
0183 toBeKept.insert(detId->begin(), detId->end());
0184 }
0185
0186 auto hbhe_out = std::make_unique<HBHERecHitCollection>();
0187 skim(hbhe, *hbhe_out);
0188 iEvent.put(std::move(hbhe_out), hbheTag.label());
0189
0190 auto hf_out = std::make_unique<HFRecHitCollection>();
0191 skim(hf, *hf_out);
0192 iEvent.put(std::move(hf_out), hfTag.label());
0193
0194 auto ho_out = std::make_unique<HORecHitCollection>();
0195 skim(ho, *ho_out, hoSeverityLevel);
0196 iEvent.put(std::move(ho_out), hoTag.label());
0197 }
0198
0199
0200 DEFINE_FWK_MODULE(HcalHitSelection);