File indexing completed on 2023-03-17 10:59:52
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/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/one/EDFilter.h"
0025
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030
0031 #include "FWCore/Framework/interface/ESHandle.h"
0032
0033 #include <cmath>
0034 #include <iostream>
0035 #include <sstream>
0036 #include <fstream>
0037
0038 #include "FWCore/Framework/interface/Run.h"
0039 #include "FWCore/Framework/interface/LuminosityBlock.h"
0040
0041 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0042
0043 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0044 #include "CalibFormats/HcalObjects/interface/HcalCoderDb.h"
0045
0046
0047
0048
0049
0050 class HcalLaserHBHEFilter2012 : public edm::one::EDFilter<> {
0051 public:
0052 explicit HcalLaserHBHEFilter2012(const edm::ParameterSet&);
0053 ~HcalLaserHBHEFilter2012() override;
0054
0055 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0056
0057 private:
0058 bool filter(edm::Event&, const edm::EventSetup&) override;
0059 void endJob() override;
0060
0061
0062
0063 bool verbose_;
0064
0065 std::string prefix_;
0066
0067 int minCalibChannelsHBHELaser_;
0068
0069
0070
0071
0072 double minFracDiffHBHELaser_;
0073 edm::InputTag digiLabel_;
0074 edm::EDGetTokenT<HcalCalibDigiCollection> tok_calib_;
0075 edm::EDGetTokenT<HBHEDigiCollection> tok_hbhe_;
0076
0077 double HBHEcalibThreshold_;
0078 std::vector<int> CalibTS_;
0079 bool WriteBadToFile_;
0080 bool forceFilterTrue_;
0081 std::ofstream outfile_;
0082 };
0083
0084
0085
0086
0087
0088
0089
0090
0091
0092
0093
0094
0095 HcalLaserHBHEFilter2012::HcalLaserHBHEFilter2012(const edm::ParameterSet& ps) {
0096
0097 verbose_ = ps.getUntrackedParameter<bool>("verbose", false);
0098 prefix_ = ps.getUntrackedParameter<std::string>("prefix", "");
0099 minCalibChannelsHBHELaser_ = ps.getUntrackedParameter<int>("minCalibChannelsHBHELaser", 20);
0100 minFracDiffHBHELaser_ = ps.getUntrackedParameter<double>("minFracDiffHBHELaser", 0.3);
0101 edm::InputTag digi_default("hcalDigis");
0102 digiLabel_ = ps.getUntrackedParameter<edm::InputTag>("digiLabel", digi_default);
0103
0104 tok_calib_ = consumes<HcalCalibDigiCollection>(digiLabel_);
0105 tok_hbhe_ = consumes<HBHEDigiCollection>(digiLabel_);
0106
0107 HBHEcalibThreshold_ = ps.getUntrackedParameter<double>("HBHEcalibThreshold", 15);
0108 std::vector<int> dummyTS;
0109 for (int i = 3; i <= 6; ++i)
0110 dummyTS.push_back(i);
0111 CalibTS_ = ps.getUntrackedParameter<std::vector<int> >("CalibTS", dummyTS);
0112 WriteBadToFile_ = ps.getUntrackedParameter<bool>("WriteBadToFile", false);
0113 if (WriteBadToFile_)
0114 outfile_.open("badHcalLaserList_hbhefilter.txt");
0115 forceFilterTrue_ = ps.getUntrackedParameter<bool>("forceFilterTrue", false);
0116 }
0117
0118 HcalLaserHBHEFilter2012::~HcalLaserHBHEFilter2012() {
0119
0120
0121 }
0122
0123
0124
0125
0126
0127
0128 bool HcalLaserHBHEFilter2012::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0129
0130
0131 edm::Handle<HcalCalibDigiCollection> calib_digi;
0132 if (!(iEvent.getByToken(tok_calib_, calib_digi))) {
0133 edm::LogWarning("HcalLaserHBHEFilter2012") << digiLabel_ << " calib_digi not available";
0134 return true;
0135 }
0136
0137 if (!(calib_digi.isValid())) {
0138 edm::LogWarning("HcalLaserHBHEFilter2012") << digiLabel_ << " calib_digi is not valid";
0139 return true;
0140 }
0141
0142 edm::Handle<HBHEDigiCollection> hbhe_digi;
0143 if (!(iEvent.getByToken(tok_hbhe_, hbhe_digi))) {
0144 edm::LogWarning("HcalLaserHBHEFilter2012") << digiLabel_ << " hbhe_digi not available";
0145 return true;
0146 }
0147
0148 if (!(hbhe_digi.isValid())) {
0149 edm::LogWarning("HcalLaserHBHEFilter2012") << digiLabel_ << " hbhe_digi is not valid";
0150 return true;
0151 }
0152
0153
0154 int ncalibHBHE = 0;
0155
0156 for (HcalCalibDigiCollection::const_iterator Calibiter = calib_digi->begin(); Calibiter != calib_digi->end();
0157 ++Calibiter) {
0158 const HcalCalibDataFrame digi = (const HcalCalibDataFrame)(*Calibiter);
0159 if (digi.zsMarkAndPass())
0160 continue;
0161 HcalCalibDetId myid = (HcalCalibDetId)digi.id();
0162 if (myid.hcalSubdet() != HcalBarrel && myid.hcalSubdet() != HcalEndcap)
0163 continue;
0164 if (myid.calibFlavor() == HcalCalibDetId::HOCrosstalk)
0165 continue;
0166
0167
0168 double thischarge = 0;
0169 for (unsigned int i = 0; i < CalibTS_.size(); ++i) {
0170 thischarge += digi[CalibTS_[i]].nominal_fC();
0171 if (thischarge > HBHEcalibThreshold_) {
0172 ++ncalibHBHE;
0173 break;
0174 }
0175 }
0176
0177 if (ncalibHBHE >= minCalibChannelsHBHELaser_) {
0178 if (verbose_)
0179 std::cout << prefix_ << iEvent.id().run() << ":" << iEvent.luminosityBlock() << ":" << iEvent.id().event()
0180 << std::endl;
0181 if (WriteBadToFile_)
0182 outfile_ << iEvent.id().run() << ":" << iEvent.luminosityBlock() << ":" << iEvent.id().event() << std::endl;
0183 if (forceFilterTrue_)
0184 return true;
0185 else
0186 return false;
0187 }
0188 }
0189
0190
0191
0192
0193 double badrbxfrac = 0.;
0194 double goodrbxfrac = 0.;
0195 int Nbad = 72 * 3;
0196 int Ngood = 2592 * 2 - Nbad;
0197 for (HBHEDigiCollection::const_iterator hbhe = hbhe_digi->begin(); hbhe != hbhe_digi->end(); ++hbhe) {
0198 const HBHEDataFrame digi = (const HBHEDataFrame)(*hbhe);
0199 HcalDetId myid = (HcalDetId)digi.id();
0200 bool isbad = false;
0201 if (myid.subdet() == HcalBarrel && myid.ieta() < 0) {
0202 if (myid.iphi() >= 15 && myid.iphi() <= 18)
0203 isbad = true;
0204 else if (myid.iphi() >= 27 && myid.iphi() <= 34)
0205 isbad = true;
0206 }
0207 if (isbad == true)
0208 badrbxfrac += 1.;
0209 else
0210 goodrbxfrac += 1.;
0211 }
0212 goodrbxfrac /= Ngood;
0213 badrbxfrac /= Nbad;
0214 if (goodrbxfrac - badrbxfrac > minFracDiffHBHELaser_) {
0215 if (verbose_)
0216 std::cout << prefix_ << iEvent.id().run() << ":" << iEvent.luminosityBlock() << ":" << iEvent.id().event()
0217 << std::endl;
0218 if (WriteBadToFile_)
0219 outfile_ << iEvent.id().run() << ":" << iEvent.luminosityBlock() << ":" << iEvent.id().event() << std::endl;
0220
0221 if (forceFilterTrue_)
0222 return true;
0223 else
0224 return false;
0225 }
0226
0227 return true;
0228
0229 }
0230
0231
0232 void HcalLaserHBHEFilter2012::endJob() {
0233 if (WriteBadToFile_)
0234 outfile_.close();
0235 }
0236
0237
0238 void HcalLaserHBHEFilter2012::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0239
0240
0241 edm::ParameterSetDescription desc;
0242 desc.setUnknown();
0243 descriptions.addDefault(desc);
0244 }
0245
0246 DEFINE_FWK_MODULE(HcalLaserHBHEFilter2012);