File indexing completed on 2024-04-06 12:25:40
0001 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0002 #include "FWCore/Framework/interface/MakerMacros.h"
0003 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
0004 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
0005 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include <FWCore/Framework/interface/Event.h>
0010
0011 class EcalSeverityLevelAlgoAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0012 public:
0013 explicit EcalSeverityLevelAlgoAnalyzer(const edm::ParameterSet& ps);
0014 ~EcalSeverityLevelAlgoAnalyzer() override = default;
0015
0016 void analyze(edm::Event const& iEvent, const edm::EventSetup& iSetup) override;
0017
0018 private:
0019 const edm::EDGetTokenT<EcalRecHitCollection> recHitsToken_;
0020 const edm::ESGetToken<EcalSeverityLevelAlgo, EcalSeverityLevelAlgoRcd> severityLevelAlgoToken_;
0021
0022 TH1F* shisto_;
0023 TH1F* fhisto_;
0024 };
0025
0026 DEFINE_FWK_MODULE(EcalSeverityLevelAlgoAnalyzer);
0027
0028 EcalSeverityLevelAlgoAnalyzer::EcalSeverityLevelAlgoAnalyzer(const edm::ParameterSet& ps)
0029 : recHitsToken_(consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit", "EcalRecHitsEB"))),
0030 severityLevelAlgoToken_(esConsumes()) {
0031 usesResource(TFileService::kSharedResource);
0032
0033 edm::Service<TFileService> fs;
0034 shisto_ = fs->make<TH1F>("SeverityLevel", "SeverityLevel", 6, 0, 6);
0035 fhisto_ = fs->make<TH1F>("Flags", "Flags", 32, 0, 32);
0036 }
0037
0038 void EcalSeverityLevelAlgoAnalyzer::analyze(edm::Event const& iEvent, const edm::EventSetup& iSetup) {
0039 edm::Handle<EcalRecHitCollection> rechits;
0040 iEvent.getByToken(recHitsToken_, rechits);
0041
0042 const auto& sevlv = iSetup.getData(severityLevelAlgoToken_);
0043
0044 EcalRecHitCollection::const_iterator rechit = rechits->begin();
0045 for (; rechit != rechits->end(); ++rechit) {
0046 for (int flag = 0; flag < 32; ++flag) {
0047 if (rechit->checkFlag(flag))
0048 fhisto_->Fill(flag);
0049 }
0050
0051 int severity = sevlv.severityLevel(rechit->id(), *rechits);
0052 shisto_->Fill(severity);
0053 }
0054 }