HcalRawToDigi/plugins/HcalLaserHBHEFilter2012.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    HcalLaserHBHEFilter2012
0004 // Class:      HcalLaserHBHEFilter2012
0005 //
0006 /**\class HcalLaserHBHEFilter2012 HcalLaserHBHEFilter2012.cc UserCode/HcalLaserHBHEFilter2012/src/HcalLaserHBHEFilter2012.cc
0007
0008  Description: [filters out HBHE laser events based on number of HBHE calib digis above threshold and the distribution of HBHE calib digis]
0009
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:
0015 //         Created:  Fri Oct 19 13:15:44 EDT 2012
0016 //
0017 //
0018
0019 // system include files
0020 #include <memory>
0021
0022 // user include files
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 // class declaration
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   // ----------member data ---------------------------
0062   // if set to true, then the run:LS:event for any event failing the cut will be printed out
0063   bool verbose_;
0064   // prefix will be printed before any event if verbose mode is true, in order to make searching for events easier
0065   std::string prefix_;
0066   // set minimum number of HBHE Calib events that causes an event to be considered a bad (i.e., HBHE laser) event
0067   int minCalibChannelsHBHELaser_;
0068   // minimum difference in fractional occupancies between 'good' and 'bad' HBHE regions (i.e., regions whose
0069   // RBXes receive laser signals and those whose RBXes see no laser) necessary to declare an event as a laser
0070   // event.  In laser events, good fractional occupancy is generally near 1, while bad fractional occupancy
0071   // is considerably less
0072   double minFracDiffHBHELaser_;
0073   edm::InputTag digiLabel_;
0074   edm::EDGetTokenT<HcalCalibDigiCollection> tok_calib_;
0075   edm::EDGetTokenT<HBHEDigiCollection> tok_hbhe_;
0076
0077   double HBHEcalibThreshold_;  // minimum integrated charge needed for a hit to count as an occupied calib channel
0078   std::vector<int> CalibTS_;   // time slices used when integrating calib charges
0079   bool WriteBadToFile_;
0080   bool forceFilterTrue_;
0081   std::ofstream outfile_;
0082 };
0083
0084 //
0085 // constants, enums and typedefs
0086 //
0087
0088 //
0089 // static data member definitions
0090 //
0091
0092 //
0093 // constructors and destructor
0094 //
0095 HcalLaserHBHEFilter2012::HcalLaserHBHEFilter2012(const edm::ParameterSet& ps) {
0096   //now do what ever initialization is needed
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 }  // HcalLaserHBHEFilter2012::HcalLaserHBHEFilter2012  constructor
0117
0118 HcalLaserHBHEFilter2012::~HcalLaserHBHEFilter2012() {
0119   // do anything here that needs to be done at destruction time
0120   // (e.g. close files, deallocate resources etc.)
0121 }
0122
0123 //
0124 // member functions
0125 //
0126
0127 // ------------ method called on each new Event  ------------
0128 bool HcalLaserHBHEFilter2012::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0129   // Step 1:: try to get calib digi and HBHE collections.
0130   // Return true if collection not found?  Or false?  What should default behavior be?
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   // Step 2:  Count HBHE digi calib channels
0154   int ncalibHBHE = 0;  // this will track number of HBHE digi channels
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;  // skip digis labeled as "mark and pass" in NZS events
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     // Compute charge in current channel (for relevant TS only)
0167     // If total charge in channel exceeds threshold, increment count of calib channels
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;  // if special input boolean set, always return true, regardless of filter decision
0185       else
0186         return false;
0187     }
0188   }
0189
0190   // Step 3:  Look at distribution of HBHE hits
0191
0192   // Count digis in good, bad RBXes.  ('bad' RBXes see no laser signal)
0193   double badrbxfrac = 0.;
0194   double goodrbxfrac = 0.;
0195   int Nbad = 72 * 3;            // 3 bad RBXes, 72 channels each
0196   int Ngood = 2592 * 2 - Nbad;  // remaining HBHE channels are 'good'
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;  // assume channel is not bad
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;  // if special input boolean set, always return true, regardless of filter decision
0223     else
0224       return false;
0225   }
0226   // Step 4:  HBHE laser tests passed.  return true
0227   return true;
0228
0229 }  // HcalLaserHBHEFilter2012::filter
0230
0231 // ------------ method called once each job just after ending the event loop  ------------
0232 void HcalLaserHBHEFilter2012::endJob() {
0233   if (WriteBadToFile_)
0234     outfile_.close();
0235 }
0236
0237 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0238 void HcalLaserHBHEFilter2012::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0239   //The following says we do not know what parameters are allowed so do no validation
0240   // Please change this to state exactly what you do use, even if it is no parameters
0241   edm::ParameterSetDescription desc;
0242   desc.setUnknown();
0243   descriptions.addDefault(desc);
0244 }
0245 //define this as a plug-in
0246 DEFINE_FWK_MODULE(HcalLaserHBHEFilter2012);