File indexing completed on 2024-05-15 22:43:57
0001
0002 #include "FWCore/Framework/interface/Frameworkfwd.h"
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/Framework/interface/one/EDFilter.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010 #include "FWCore/Utilities/interface/ESGetToken.h"
0011 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0012 #include "DataFormats/HcalDigi/interface/QIE11DataFrame.h"
0013 #include "DataFormats/HcalDetId/interface/HcalGenericDetId.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include <fstream>
0017
0018 class HBHEstuckADCfilter : public edm::one::EDFilter<> {
0019 public:
0020 explicit HBHEstuckADCfilter(const edm::ParameterSet&);
0021 ~HBHEstuckADCfilter() override;
0022
0023 static void fillDescriptions(edm::ConfigurationDescriptions&);
0024
0025 private:
0026 bool filter(edm::Event&, const edm::EventSetup&) override;
0027 void endJob() override;
0028
0029 edm::EDGetTokenT<QIE11DigiCollection> tok_qie11_;
0030 int thresholdADC_;
0031 bool writeList_;
0032 std::ofstream outfile_;
0033 };
0034
0035 HBHEstuckADCfilter::HBHEstuckADCfilter(const edm::ParameterSet& conf)
0036 : tok_qie11_(consumes<QIE11DigiCollection>(conf.getParameter<edm::InputTag>("digiLabel"))),
0037 thresholdADC_(conf.getParameter<int>("thresholdADC")),
0038 writeList_(conf.getParameter<bool>("writeList")) {
0039 if (writeList_)
0040 outfile_.open("events_list_stuckADC.txt");
0041 }
0042
0043 HBHEstuckADCfilter::~HBHEstuckADCfilter() {}
0044
0045 bool HBHEstuckADCfilter::filter(edm::Event& ev, const edm::EventSetup& set) {
0046 edm::Handle<QIE11DigiCollection> theDigis;
0047 ev.getByToken(tok_qie11_, theDigis);
0048
0049 bool result = true;
0050 for (QIE11DigiCollection::const_iterator itr = theDigis->begin(); itr != theDigis->end(); itr++) {
0051 int tsize = (*itr).size();
0052 const QIE11DataFrame frame = *itr;
0053
0054 bool flag = true;
0055 int adc0 = (frame[0]).adc();
0056 if (adc0 < thresholdADC_)
0057 flag = false;
0058 else {
0059 for (int i = 1; i < tsize; i++) {
0060 if ((frame[i]).adc() != adc0) {
0061 flag = false;
0062 break;
0063 }
0064 }
0065 }
0066
0067
0068 if (flag) {
0069 const HcalDetId cell(itr->id());
0070 edm::LogWarning("HBHEstuckADCfilter") << "stuck ADC = " << adc0 << " in " << cell << std::endl;
0071 result = false;
0072 }
0073 }
0074 if (!result && writeList_)
0075 outfile_ << ev.id().run() << ":" << ev.luminosityBlock() << ":" << ev.id().event() << std::endl;
0076
0077 return result;
0078 }
0079
0080 void HBHEstuckADCfilter::endJob() {
0081 if (writeList_)
0082 outfile_.close();
0083 }
0084
0085 void HBHEstuckADCfilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0086 edm::ParameterSetDescription desc;
0087 desc.add<edm::InputTag>("digiLabel", edm::InputTag("hcalDigis"));
0088 desc.add<int>("thresholdADC", 100);
0089 desc.add<bool>("writeList", true);
0090 descriptions.add("hbhestuckADCfilter", desc);
0091 }
0092
0093
0094 DEFINE_FWK_MODULE(HBHEstuckADCfilter);