File indexing completed on 2023-10-25 10:00:08
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022 #include "FWCore/Framework/interface/Frameworkfwd.h"
0023 #include "FWCore/Framework/interface/stream/EDProducer.h"
0024 #include "FWCore/Utilities/interface/ESGetToken.h"
0025 #include "FWCore/Framework/interface/Event.h"
0026 #include "FWCore/Framework/interface/MakerMacros.h"
0027 #include "FWCore/Utilities/interface/Exception.h"
0028
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 #include "FWCore/Utilities/interface/StreamID.h"
0031
0032 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0033
0034 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
0035 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
0036 #include "CondFormats/HcalObjects/interface/HFPhase1PMTParams.h"
0037 #include "CondFormats/DataRecord/interface/HFPhase1PMTParamsRcd.h"
0038
0039 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputer.h"
0040 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputerRcd.h"
0041
0042
0043 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalHF_PETalgorithm.h"
0044 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalHF_S9S1algorithm.h"
0045 #include "RecoLocalCalo/HcalRecAlgos/interface/HFStripFilter.h"
0046
0047
0048 #include "RecoLocalCalo/HcalRecAlgos/interface/parseHFPhase1AlgoDescription.h"
0049
0050
0051
0052
0053 class HFPhase1Reconstructor : public edm::stream::EDProducer<> {
0054 public:
0055 explicit HFPhase1Reconstructor(const edm::ParameterSet&);
0056 ~HFPhase1Reconstructor() override;
0057
0058 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0059
0060 private:
0061 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0062 void produce(edm::Event&, const edm::EventSetup&) override;
0063
0064
0065 std::string algoConfigClass_;
0066 bool setNoiseFlags_;
0067 bool runHFStripFilter_;
0068 bool useChannelQualityFromDB_;
0069 bool checkChannelQualityForDepth3and4_;
0070
0071
0072 edm::EDGetTokenT<HFPreRecHitCollection> tok_PreRecHit_;
0073 std::unique_ptr<AbsHFPhase1Algo> reco_;
0074 std::unique_ptr<AbsHcalAlgoData> recoConfig_;
0075 edm::ESGetToken<HFPhase1PMTParams, HFPhase1PMTParamsRcd> recoConfigToken_;
0076
0077
0078 std::unique_ptr<HcalHF_S9S1algorithm> hfS9S1_;
0079 std::unique_ptr<HcalHF_S9S1algorithm> hfS8S1_;
0080 std::unique_ptr<HcalHF_PETalgorithm> hfPET_;
0081 std::unique_ptr<HFStripFilter> hfStripFilter_;
0082
0083
0084 edm::ESGetToken<HcalDbService, HcalDbRecord> conditionsToken_;
0085 edm::ESGetToken<HcalChannelQuality, HcalChannelQualityRcd> qualToken_;
0086 edm::ESGetToken<HcalSeverityLevelComputer, HcalSeverityLevelComputerRcd> sevToken_;
0087 };
0088
0089
0090
0091
0092 HFPhase1Reconstructor::HFPhase1Reconstructor(const edm::ParameterSet& conf)
0093 : algoConfigClass_(conf.getParameter<std::string>("algoConfigClass")),
0094 setNoiseFlags_(conf.getParameter<bool>("setNoiseFlags")),
0095 runHFStripFilter_(conf.getParameter<bool>("runHFStripFilter")),
0096 useChannelQualityFromDB_(conf.getParameter<bool>("useChannelQualityFromDB")),
0097 checkChannelQualityForDepth3and4_(conf.getParameter<bool>("checkChannelQualityForDepth3and4")),
0098 reco_(parseHFPhase1AlgoDescription(conf.getParameter<edm::ParameterSet>("algorithm"))) {
0099
0100 if (!reco_.get())
0101 throw cms::Exception("HFPhase1BadConfig") << "Invalid HFPhase1Reconstructor algorithm configuration" << std::endl;
0102
0103 if (reco_->isConfigurable()) {
0104 if ("HFPhase1PMTParams" != algoConfigClass_) {
0105 throw cms::Exception("HBHEPhase1BadConfig")
0106 << "Invalid HBHEPhase1Reconstructor \"algoConfigClass\" parameter value \"" << algoConfigClass_ << '"'
0107 << std::endl;
0108 }
0109 recoConfigToken_ = esConsumes<edm::Transition::BeginRun>();
0110 }
0111
0112
0113 if (setNoiseFlags_) {
0114 const edm::ParameterSet& psS9S1 = conf.getParameter<edm::ParameterSet>("S9S1stat");
0115 hfS9S1_ = std::make_unique<HcalHF_S9S1algorithm>(psS9S1.getParameter<std::vector<double> >("short_optimumSlope"),
0116 psS9S1.getParameter<std::vector<double> >("shortEnergyParams"),
0117 psS9S1.getParameter<std::vector<double> >("shortETParams"),
0118 psS9S1.getParameter<std::vector<double> >("long_optimumSlope"),
0119 psS9S1.getParameter<std::vector<double> >("longEnergyParams"),
0120 psS9S1.getParameter<std::vector<double> >("longETParams"),
0121 psS9S1.getParameter<int>("HcalAcceptSeverityLevel"),
0122 psS9S1.getParameter<bool>("isS8S1"));
0123
0124 const edm::ParameterSet& psS8S1 = conf.getParameter<edm::ParameterSet>("S8S1stat");
0125 hfS8S1_ = std::make_unique<HcalHF_S9S1algorithm>(psS8S1.getParameter<std::vector<double> >("short_optimumSlope"),
0126 psS8S1.getParameter<std::vector<double> >("shortEnergyParams"),
0127 psS8S1.getParameter<std::vector<double> >("shortETParams"),
0128 psS8S1.getParameter<std::vector<double> >("long_optimumSlope"),
0129 psS8S1.getParameter<std::vector<double> >("longEnergyParams"),
0130 psS8S1.getParameter<std::vector<double> >("longETParams"),
0131 psS8S1.getParameter<int>("HcalAcceptSeverityLevel"),
0132 psS8S1.getParameter<bool>("isS8S1"));
0133
0134 const edm::ParameterSet& psPET = conf.getParameter<edm::ParameterSet>("PETstat");
0135 hfPET_ = std::make_unique<HcalHF_PETalgorithm>(psPET.getParameter<std::vector<double> >("short_R"),
0136 psPET.getParameter<std::vector<double> >("shortEnergyParams"),
0137 psPET.getParameter<std::vector<double> >("shortETParams"),
0138 psPET.getParameter<std::vector<double> >("long_R"),
0139 psPET.getParameter<std::vector<double> >("longEnergyParams"),
0140 psPET.getParameter<std::vector<double> >("longETParams"),
0141 psPET.getParameter<int>("HcalAcceptSeverityLevel"),
0142 psPET.getParameter<std::vector<double> >("short_R_29"),
0143 psPET.getParameter<std::vector<double> >("long_R_29"));
0144
0145
0146 if (runHFStripFilter_) {
0147 const edm::ParameterSet& psStripFilter = conf.getParameter<edm::ParameterSet>("HFStripFilter");
0148 hfStripFilter_ = HFStripFilter::parseParameterSet(psStripFilter);
0149 }
0150 }
0151
0152
0153 tok_PreRecHit_ = consumes<HFPreRecHitCollection>(conf.getParameter<edm::InputTag>("inputLabel"));
0154
0155
0156 produces<HFRecHitCollection>();
0157
0158
0159 conditionsToken_ = esConsumes<HcalDbService, HcalDbRecord>();
0160 qualToken_ = esConsumes<HcalChannelQuality, HcalChannelQualityRcd>(edm::ESInputTag("", "withTopo"));
0161 sevToken_ = esConsumes<HcalSeverityLevelComputer, HcalSeverityLevelComputerRcd>();
0162 }
0163
0164 HFPhase1Reconstructor::~HFPhase1Reconstructor() {
0165
0166
0167 }
0168
0169 void HFPhase1Reconstructor::beginRun(const edm::Run& r, const edm::EventSetup& es) {
0170 if (reco_->isConfigurable()) {
0171 recoConfig_ = std::make_unique<HFPhase1PMTParams>(es.getData(recoConfigToken_));
0172 if (!reco_->configure(recoConfig_.get()))
0173 throw cms::Exception("HFPhase1BadConfig")
0174 << "Failed to configure HFPhase1Reconstructor algorithm from EventSetup" << std::endl;
0175 }
0176 }
0177
0178
0179 void HFPhase1Reconstructor::produce(edm::Event& e, const edm::EventSetup& eventSetup) {
0180 using namespace edm;
0181
0182
0183 const HcalDbService* conditions = &eventSetup.getData(conditionsToken_);
0184 const HcalChannelQuality* myqual = &eventSetup.getData(qualToken_);
0185 const HcalSeverityLevelComputer* mySeverity = &eventSetup.getData(sevToken_);
0186
0187
0188 Handle<HFPreRecHitCollection> preRecHits;
0189 e.getByToken(tok_PreRecHit_, preRecHits);
0190
0191
0192 std::unique_ptr<HFRecHitCollection> rec(std::make_unique<HFRecHitCollection>());
0193 rec->reserve(preRecHits->size());
0194
0195
0196 for (HFPreRecHitCollection::const_iterator it = preRecHits->begin(); it != preRecHits->end(); ++it) {
0197
0198
0199
0200 const bool thisIsSingleAnodePMT = false;
0201
0202
0203 bool taggedBadByDb[2] = {false, false};
0204 if (useChannelQualityFromDB_) {
0205 if (checkChannelQualityForDepth3and4_ && !thisIsSingleAnodePMT) {
0206 HcalDetId anodeIds[2];
0207 anodeIds[0] = it->id();
0208 anodeIds[1] = anodeIds[0].secondAnodeId();
0209 for (unsigned i = 0; i < 2; ++i) {
0210 const HcalChannelStatus* mydigistatus = myqual->getValues(anodeIds[i].rawId());
0211 taggedBadByDb[i] = mySeverity->dropChannel(mydigistatus->getValue());
0212 }
0213 } else {
0214 const HcalChannelStatus* mydigistatus = myqual->getValues(it->id().rawId());
0215 const bool b = mySeverity->dropChannel(mydigistatus->getValue());
0216 taggedBadByDb[0] = b;
0217 taggedBadByDb[1] = b;
0218 }
0219 }
0220
0221
0222 const HFRecHit& rh =
0223 reco_->reconstruct(*it, conditions->getHcalCalibrations(it->id()), taggedBadByDb, thisIsSingleAnodePMT);
0224
0225
0226
0227 if (rh.id().rawId())
0228 rec->push_back(rh);
0229 }
0230
0231
0232
0233
0234
0235
0236
0237 if (setNoiseFlags_) {
0238
0239 for (HFRecHitCollection::iterator i = rec->begin(); i != rec->end(); ++i) {
0240 int depth = i->id().depth();
0241 int ieta = i->id().ieta();
0242
0243 if (depth == 2 || abs(ieta) == 29)
0244 hfPET_->HFSetFlagFromPET(*i, *rec, myqual, mySeverity);
0245 }
0246
0247
0248 for (HFRecHitCollection::iterator i = rec->begin(); i != rec->end(); ++i) {
0249 int depth = i->id().depth();
0250 int ieta = i->id().ieta();
0251
0252 if (depth == 2 || abs(ieta) == 29)
0253 hfS8S1_->HFSetFlagFromS9S1(*i, *rec, myqual, mySeverity);
0254 }
0255
0256
0257 for (HFRecHitCollection::iterator i = rec->begin(); i != rec->end(); ++i) {
0258 int depth = i->id().depth();
0259 int ieta = i->id().ieta();
0260
0261 if (depth == 1 && abs(ieta) != 29)
0262 hfS9S1_->HFSetFlagFromS9S1(*i, *rec, myqual, mySeverity);
0263 }
0264
0265
0266 if (runHFStripFilter_)
0267 hfStripFilter_->runFilter(*rec, myqual);
0268 }
0269
0270
0271 e.put(std::move(rec));
0272 }
0273
0274 #define add_param_set(name) \
0275 edm::ParameterSetDescription name; \
0276 name.setAllowAnything(); \
0277 desc.add<edm::ParameterSetDescription>(#name, name)
0278
0279
0280 void HFPhase1Reconstructor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0281 edm::ParameterSetDescription desc;
0282
0283 desc.add<edm::InputTag>("inputLabel");
0284 desc.add<std::string>("algoConfigClass");
0285 desc.add<bool>("setNoiseFlags");
0286 desc.add<bool>("runHFStripFilter", false);
0287 desc.add<bool>("useChannelQualityFromDB");
0288 desc.add<bool>("checkChannelQualityForDepth3and4");
0289 desc.add<edm::ParameterSetDescription>("algorithm", fillDescriptionForParseHFPhase1AlgoDescription());
0290 desc.add<edm::ParameterSetDescription>("HFStripFilter", HFStripFilter::fillDescription());
0291
0292 add_param_set(S9S1stat);
0293 add_param_set(S8S1stat);
0294 add_param_set(PETstat);
0295
0296 descriptions.addDefault(desc);
0297 }
0298
0299
0300 DEFINE_FWK_MODULE(HFPhase1Reconstructor);