METFilters/plugins/HFNoisyHitsFilter.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    RecoMET/METFilters
0004 // Class:      HFNoisyHitsFilter
0005 //
0006 /**\class HFNoisyHitsFilter HFNoisyHitsFilter.cc RecoMET/METFilters/plugins/HFNoisyHitsFilter.cc
0007
0008  Description: [one line class summary]
0009
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Laurent Thomas
0015 //         Created:  Tue, 01 Sep 2020 11:24:33 GMT
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/global/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/Utilities/interface/Exception.h"
0032
0033 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0034 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0035
0036 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0037
0038 #include "DataFormats/HcalRecHit/interface/HFRecHit.h"
0039 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0040 #include "DataFormats/METReco/interface/HcalPhase1FlagLabels.h"
0041 //
0042 // class declaration
0043 //
0044
0045 class HFNoisyHitsFilter : public edm::global::EDFilter<> {
0046 public:
0047   explicit HFNoisyHitsFilter(const edm::ParameterSet&);
0048   ~HFNoisyHitsFilter() override {}
0049   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0050
0051 private:
0052   bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0053   std::vector<HcalPhase1FlagLabels::HFStatusFlag> getNoiseBits() const;
0054   const edm::EDGetTokenT<HFRecHitCollection> hfhits_token_;
0055   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geom_token_;
0056   const double rechitPtThreshold_;
0057   const std::vector<std::string> listOfNoises_;
0058   const bool taggingMode_;
0059   const bool debug_;
0060   std::vector<HcalPhase1FlagLabels::HFStatusFlag> noiseBits_;
0061 };
0062
0063 //
0064 // constructors and destructor
0065 //
0066 HFNoisyHitsFilter::HFNoisyHitsFilter(const edm::ParameterSet& iConfig)
0067     : hfhits_token_(consumes<HFRecHitCollection>(iConfig.getParameter<edm::InputTag>("hfrechits"))),
0068       geom_token_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0069       rechitPtThreshold_(iConfig.getParameter<double>("rechitPtThreshold")),
0070       listOfNoises_(iConfig.getParameter<std::vector<std::string>>("listOfNoises")),
0071       taggingMode_(iConfig.getParameter<bool>("taggingMode")),
0072       debug_(iConfig.getParameter<bool>("debug")) {
0073   noiseBits_ = getNoiseBits();
0074   produces<bool>();
0075 }
0076
0077 //
0078 // member functions
0079 //
0080
0081 // ------------ method called on each new Event  ------------
0082 bool HFNoisyHitsFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0083   using namespace edm;
0084   bool pass = true;
0085
0086   // Calo Geometry - needed for computing E_t
0087   const CaloGeometry& geo = iSetup.getData(geom_token_);
0088
0089   auto const& hfHits = iEvent.get(hfhits_token_);
0090
0091   //Loop over the HF rechits. If one of them has Et>X and fires one the noise bits, declare the event as bad
0092   for (auto const& hfhit : hfHits) {
0093     float ene = hfhit.energy();
0094     float et = 0;
0095     // compute transverse energy
0096     const GlobalPoint& poshf = geo.getPosition(hfhit.detid());
0097     float pf = poshf.perp() / poshf.mag();
0098     et = ene * pf;
0099     if (et < rechitPtThreshold_)
0100       continue;
0101     int hitFlags = hfhit.flags();
0102     for (auto noiseBit : noiseBits_) {
0103       if ((hitFlags >> noiseBit) & 1) {
0104         pass = false;
0105         break;
0106       }
0107     }
0108     if (!pass)
0109       break;
0110   }
0111   iEvent.put(std::make_unique<bool>(pass));
0112   if (debug_)
0113     LogDebug("HFNoisyHitsFilter") << "Passing filter? " << pass;
0114   return taggingMode_ || pass;
0115 }
0116
0117 std::vector<HcalPhase1FlagLabels::HFStatusFlag> HFNoisyHitsFilter::getNoiseBits() const {
0118   std::vector<HcalPhase1FlagLabels::HFStatusFlag> result;
0119   for (auto const& noise : listOfNoises_) {
0120     if (noise == "HFLongShort")
0121       result.push_back(HcalPhase1FlagLabels::HFLongShort);
0122     else if (noise == "HFS8S1Ratio")
0123       result.push_back(HcalPhase1FlagLabels::HFS8S1Ratio);
0124     else if (noise == "HFPET")
0125       result.push_back(HcalPhase1FlagLabels::HFPET);
0126     else if (noise == "HFSignalAsymmetry")
0127       result.push_back(HcalPhase1FlagLabels::HFSignalAsymmetry);
0128     else if (noise == "HFAnomalousHit")
0129       result.push_back(HcalPhase1FlagLabels::HFAnomalousHit);
0130     else
0131       throw cms::Exception("Error") << "Couldn't find the bit index associated to this string: " << noise;
0132   }
0133
0134   return result;
0135 }
0136
0137 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0138 void HFNoisyHitsFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0139   edm::ParameterSetDescription desc;
0140   desc.add<edm::InputTag>("hfrechits", {"reducedHcalRecHits:hfreco"});
0141   desc.add<double>("rechitPtThreshold", 20.);
0142   desc.add<std::vector<std::string>>("listOfNoises", {"HFLongShort", "HFS8S1Ratio", "HFPET", "HFSignalAsymmetry"});
0143   desc.add<bool>("taggingMode", false);
0144   desc.add<bool>("debug", false);
0145   descriptions.add("hfNoisyHitsFilter", desc);
0146 }
0147 //define this as a plug-in
0148 DEFINE_FWK_MODULE(HFNoisyHitsFilter);