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 <string>
0021 #include <iostream>
0022 #include <memory>
0023
0024
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026
0027 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/MakerMacros.h"
0030 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032
0033 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0034 #include "DataFormats/HcalDigi/interface/HcalCalibrationEventTypes.h"
0035 #include "EventFilter/HcalRawToDigi/interface/HcalDCCHeader.h"
0036 #include "EventFilter/HcalRawToDigi/interface/HcalHTRData.h"
0037
0038 #include "HLTHcalNZSFilter.h"
0039
0040
0041
0042
0043 HLTHcalNZSFilter::HLTHcalNZSFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0044
0045
0046 dataInputTag_ = iConfig.getParameter<edm::InputTag>("InputTag");
0047 dataInputToken_ = consumes<FEDRawDataCollection>(dataInputTag_);
0048 }
0049
0050 HLTHcalNZSFilter::~HLTHcalNZSFilter() {
0051
0052
0053 }
0054
0055 void HLTHcalNZSFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0056 edm::ParameterSetDescription desc;
0057 makeHLTFilterDescription(desc);
0058 desc.add<edm::InputTag>("InputTag", edm::InputTag("source"));
0059 descriptions.add("hltHcalNZSFilter", desc);
0060 }
0061
0062
0063
0064
0065
0066
0067 bool HLTHcalNZSFilter::hltFilter(edm::Event& iEvent,
0068 const edm::EventSetup& iSetup,
0069 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0070 using namespace edm;
0071
0072
0073 if (!iEvent.isRealData())
0074 return false;
0075
0076 edm::Handle<FEDRawDataCollection> rawdata;
0077 iEvent.getByToken(dataInputToken_, rawdata);
0078
0079 int nFEDs = 0;
0080 int nNZSfed = 0;
0081 int nZSfed = 0;
0082 for (int i = FEDNumbering::MINHCALFEDID; i <= FEDNumbering::MAXHCALFEDID; i++) {
0083 const FEDRawData& fedData = rawdata->FEDData(i);
0084 if (fedData.size() < 24)
0085 continue;
0086 nFEDs++;
0087
0088
0089 HcalHTRData htr;
0090 const HcalDCCHeader* dccHeader = (const HcalDCCHeader*)(fedData.data());
0091 int nZS = 0;
0092 int nUS = 0;
0093 int nSpigot = 0;
0094 for (int spigot = 0; spigot < HcalDCCHeader::SPIGOT_COUNT; spigot++) {
0095 if (!dccHeader->getSpigotPresent(spigot))
0096 continue;
0097
0098
0099 dccHeader->getSpigotData(spigot, htr, fedData.size());
0100 if ((htr.getFirmwareFlavor() & 0xE0) == 0x80)
0101 continue;
0102
0103 nSpigot++;
0104
0105 if (htr.isUnsuppressed())
0106 nUS++;
0107 else
0108 nZS++;
0109 }
0110
0111 if (nUS == nSpigot)
0112 nNZSfed++;
0113 else {
0114 nZSfed++;
0115 if (nUS > 0)
0116 LogWarning("HLTHcalNZSFilter") << "Mixture of ZS(" << nZS << ") and NZS(" << nUS << ") spigots in FED " << i;
0117 }
0118 }
0119
0120 if ((nNZSfed == nFEDs) && (nFEDs > 0)) {
0121 return true;
0122 } else {
0123 if (nNZSfed > 0)
0124 LogWarning("HLTHcalNZSFilter") << "Mixture of ZS(" << nZSfed << ") and NZS(" << nNZSfed << ") FEDs in this event";
0125 return false;
0126 }
0127 }
0128
0129
0130 #include "FWCore/Framework/interface/MakerMacros.h"
0131 DEFINE_FWK_MODULE(HLTHcalNZSFilter);