CMS 3D CMS Logo

HFNoisyHitsFilter.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: RecoMET/METFilters
4 // Class: HFNoisyHitsFilter
5 //
13 //
14 // Original Author: Laurent Thomas
15 // Created: Tue, 01 Sep 2020 11:24:33 GMT
16 //
17 //
18 
19 // system include files
20 #include <memory>
21 
22 // user include files
25 
28 
30 
32 
35 
37 
41 //
42 // class declaration
43 //
44 
46 public:
47  explicit HFNoisyHitsFilter(const edm::ParameterSet&);
48  ~HFNoisyHitsFilter() override {}
49  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
50 
51 private:
52  bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
53  std::vector<HcalPhase1FlagLabels::HFStatusFlag> getNoiseBits() const;
56  const double rechitPtThreshold_;
57  const std::vector<std::string> listOfNoises_;
58  const bool taggingMode_;
59  const bool debug_;
60  std::vector<HcalPhase1FlagLabels::HFStatusFlag> noiseBits_;
61 };
62 
63 //
64 // constructors and destructor
65 //
67  : hfhits_token_(consumes<HFRecHitCollection>(iConfig.getParameter<edm::InputTag>("hfrechits"))),
69  rechitPtThreshold_(iConfig.getParameter<double>("rechitPtThreshold")),
70  listOfNoises_(iConfig.getParameter<std::vector<std::string>>("listOfNoises")),
71  taggingMode_(iConfig.getParameter<bool>("taggingMode")),
72  debug_(iConfig.getParameter<bool>("debug")) {
74  produces<bool>();
75 }
76 
77 //
78 // member functions
79 //
80 
81 // ------------ method called on each new Event ------------
83  using namespace edm;
84  bool pass = true;
85 
86  // Calo Geometry - needed for computing E_t
87  const CaloGeometry& geo = iSetup.getData(geom_token_);
88 
89  auto const& hfHits = iEvent.get(hfhits_token_);
90 
91  //Loop over the HF rechits. If one of them has Et>X and fires one the noise bits, declare the event as bad
92  for (auto const& hfhit : hfHits) {
93  float ene = hfhit.energy();
94  float et = 0;
95  // compute transverse energy
96  const GlobalPoint& poshf = geo.getPosition(hfhit.detid());
97  float pf = poshf.perp() / poshf.mag();
98  et = ene * pf;
99  if (et < rechitPtThreshold_)
100  continue;
101  int hitFlags = hfhit.flags();
102  for (auto noiseBit : noiseBits_) {
103  if ((hitFlags >> noiseBit) & 1) {
104  pass = false;
105  break;
106  }
107  }
108  if (!pass)
109  break;
110  }
111  iEvent.put(std::make_unique<bool>(pass));
112  if (debug_)
113  LogDebug("HFNoisyHitsFilter") << "Passing filter? " << pass;
114  return taggingMode_ || pass;
115 }
116 
117 std::vector<HcalPhase1FlagLabels::HFStatusFlag> HFNoisyHitsFilter::getNoiseBits() const {
118  std::vector<HcalPhase1FlagLabels::HFStatusFlag> result;
119  for (auto const& noise : listOfNoises_) {
120  if (noise == "HFLongShort")
121  result.push_back(HcalPhase1FlagLabels::HFLongShort);
122  else if (noise == "HFS8S1Ratio")
123  result.push_back(HcalPhase1FlagLabels::HFS8S1Ratio);
124  else if (noise == "HFPET")
125  result.push_back(HcalPhase1FlagLabels::HFPET);
126  else if (noise == "HFSignalAsymmetry")
127  result.push_back(HcalPhase1FlagLabels::HFSignalAsymmetry);
128  else if (noise == "HFAnomalousHit")
129  result.push_back(HcalPhase1FlagLabels::HFAnomalousHit);
130  else
131  throw cms::Exception("Error") << "Couldn't find the bit index associated to this string: " << noise;
132  }
133 
134  return result;
135 }
136 
137 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
140  desc.add<edm::InputTag>("hfrechits", {"reducedHcalRecHits:hfreco"});
141  desc.add<double>("rechitPtThreshold", 20.);
142  desc.add<std::vector<std::string>>("listOfNoises", {"HFLongShort", "HFS8S1Ratio", "HFPET", "HFSignalAsymmetry"});
143  desc.add<bool>("taggingMode", false);
144  desc.add<bool>("debug", false);
145  descriptions.add("hfNoisyHitsFilter", desc);
146 }
147 //define this as a plug-in
#define LogDebug(id)
HFNoisyHitsFilter(const edm::ParameterSet &)
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:125
const edm::ESGetToken< CaloGeometry, CaloGeometryRecord > geom_token_
T perp() const
Definition: PV3DBase.h:72
~HFNoisyHitsFilter() override
const edm::EDGetTokenT< HFRecHitCollection > hfhits_token_
const std::vector< std::string > listOfNoises_
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
std::vector< HcalPhase1FlagLabels::HFStatusFlag > noiseBits_
bool getData(T &iHolder) const
Definition: EventSetup.h:111
const double rechitPtThreshold_
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
T mag() const
Definition: PV3DBase.h:67
bool filter(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
GlobalPoint getPosition(const DetId &id) const
Get the position of a given detector id.
Definition: CaloGeometry.cc:74
bool get(ProductID const &oid, Handle< PROD > &result) const
Definition: Event.h:326
std::vector< HcalPhase1FlagLabels::HFStatusFlag > getNoiseBits() const
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void add(std::string const &label, ParameterSetDescription const &psetDescription)
et
define resolution functions of each parameter
HLT enums.