File indexing completed on 2023-03-17 11:09:45
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/MessageLogger/interface/MessageLogger.h"
0024 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0026
0027 #include "DataFormats/Common/interface/Handle.h"
0028 #include "DataFormats/Common/interface/RefToBase.h"
0029 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0030 #include "DataFormats/HcalDigi/interface/HcalCalibrationEventTypes.h"
0031 #include "EventFilter/HcalRawToDigi/interface/HcalDCCHeader.h"
0032 #include "HLTHcalLaserMisfireFilter.h"
0033
0034
0035
0036
0037 HLTHcalLaserMisfireFilter::HLTHcalLaserMisfireFilter(const edm::ParameterSet& config) {
0038 inputHBHE_ = config.getParameter<edm::InputTag>("InputHBHE");
0039 inputHF_ = config.getParameter<edm::InputTag>("InputHF");
0040 minFracDiffHBHELaser_ = config.getParameter<double>("minFracDiffHBHELaser");
0041 minFracHFLaser_ = config.getParameter<double>("minFracHFLaser");
0042 minADCHBHE_ = config.getParameter<int>("minADCHBHE");
0043 minADCHF_ = config.getParameter<int>("minADCHF"), testMode_ = config.getUntrackedParameter<bool>("testMode", false);
0044
0045 inputTokenHBHE_ = consumes<HBHEDigiCollection>(inputHBHE_);
0046 inputTokenHF_ = consumes<QIE10DigiCollection>(inputHF_);
0047 }
0048
0049 HLTHcalLaserMisfireFilter::~HLTHcalLaserMisfireFilter() {}
0050
0051 void HLTHcalLaserMisfireFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0052 edm::ParameterSetDescription desc;
0053 desc.add<edm::InputTag>("InputHBHE", edm::InputTag("source"));
0054 desc.add<edm::InputTag>("InputHF", edm::InputTag("source"));
0055 desc.add<int>("minADCHBHE", 10);
0056 desc.add<int>("minADCHF", 10);
0057 desc.add<double>("minFracDiffHBHELaser", 0.3);
0058 desc.add<double>("minFracHFLaser", 0.3);
0059 desc.addUntracked<bool>("testMode", false);
0060 descriptions.add("hltHcalLaserMisfireFilter", desc);
0061 }
0062
0063 bool HLTHcalLaserMisfireFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0064 edm::Handle<HBHEDigiCollection> hbhe_digi;
0065 iEvent.getByToken(inputTokenHBHE_, hbhe_digi);
0066
0067 edm::Handle<QIE10DigiCollection> hf_digi;
0068 iEvent.getByToken(inputTokenHF_, hf_digi);
0069
0070
0071 double badrbxfracHBHE(0), goodrbxfracHBHE(0), allrbxfracHF(0);
0072 int NbadHBHE = 72 * 3;
0073 int NgoodHBHE = 2592 * 2 - NbadHBHE;
0074 int NallHF = 864 * 4;
0075
0076 for (auto hbhe = hbhe_digi->begin(); hbhe != hbhe_digi->end(); ++hbhe) {
0077 const HBHEDataFrame digi = (const HBHEDataFrame)(*hbhe);
0078 HcalDetId myid = (HcalDetId)digi.id();
0079 bool isbad(false);
0080
0081 bool passCut(false);
0082 int maxdigiHB(0);
0083 for (int i = 0; i < digi.size(); i++)
0084 if (digi.sample(i).adc() > maxdigiHB)
0085 maxdigiHB = digi.sample(i).adc();
0086 if (maxdigiHB > minADCHBHE_)
0087 passCut = true;
0088
0089
0090
0091
0092 if (myid.subdet() == HcalBarrel && myid.ieta() < 0) {
0093 if (myid.iphi() >= 15 && myid.iphi() <= 18)
0094 isbad = true;
0095 else if (myid.iphi() >= 27 && myid.iphi() <= 34)
0096 isbad = true;
0097 }
0098
0099 if (passCut) {
0100 if (isbad) {
0101 badrbxfracHBHE += 1.;
0102 } else
0103 goodrbxfracHBHE += 1.;
0104 }
0105 }
0106 goodrbxfracHBHE /= NgoodHBHE;
0107 badrbxfracHBHE /= NbadHBHE;
0108
0109 for (auto hf = hf_digi->begin(); hf != hf_digi->end(); ++hf) {
0110 const QIE10DataFrame digi = (const QIE10DataFrame)(*hf);
0111 bool passCut(false);
0112 int maxdigiHF(0);
0113 for (int i = 0; i < digi.samples(); i++)
0114 if (digi[i].adc() > maxdigiHF)
0115 maxdigiHF = digi[i].adc();
0116 if (maxdigiHF > minADCHF_)
0117 passCut = true;
0118
0119 if (passCut) {
0120 allrbxfracHF += 1.;
0121 }
0122 }
0123 allrbxfracHF /= NallHF;
0124
0125 if (testMode_)
0126 edm::LogVerbatim("Report") << "******************************************************************\n"
0127 << "goodrbxfracHBHE: " << goodrbxfracHBHE << " badrbxfracHBHE: " << badrbxfracHBHE
0128 << " Size " << hbhe_digi->size() << "\n"
0129 << "allrbxfracHF: " << allrbxfracHF << " Size " << hf_digi->size()
0130 << "\n******************************************************************";
0131
0132 if (((goodrbxfracHBHE - badrbxfracHBHE) < minFracDiffHBHELaser_) || (allrbxfracHF < minFracHFLaser_))
0133 return false;
0134
0135 return true;
0136 }
0137
0138
0139 #include "FWCore/Framework/interface/MakerMacros.h"
0140 DEFINE_FWK_MODULE(HLTHcalLaserMisfireFilter);