CMS 3D CMS Logo

HLTEcalTowerFilter.cc
Go to the documentation of this file.
1 
13 
14 #include <memory>
15 
24 
25 
26 //
27 // class declaration
28 //
29 
30 class HLTEcalTowerFilter : public HLTFilter {
31 public:
32  explicit HLTEcalTowerFilter(const edm::ParameterSet &);
33  ~HLTEcalTowerFilter() override;
34  static void fillDescriptions(edm::ConfigurationDescriptions & descriptions);
35 
36 private:
37  bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs & filterproduct) const override;
38 
40  edm::InputTag inputTag_; // input tag identifying product
41  double min_E_; // energy threshold in GeV
42  double max_Eta_; // maximum eta
43  int min_N_; // minimum number
44 
45 };
46 
47 //
48 // constructors and destructor
49 //
51  inputTag_ (config.getParameter<edm::InputTag>("inputTag")),
52  min_E_ (config.getParameter<double> ("MinE" )),
53  max_Eta_ (config.getParameter<double> ("MaxEta" )),
54  min_N_ (config.getParameter<int> ("MinN" ))
55 {
56  inputToken_ = consumes<CaloTowerCollection>(inputTag_);
57  LogDebug("") << "Input/ecut/etacut/ncut : "
58  << inputTag_.encode() << " "
59  << min_E_ << " "
60  << max_Eta_ << " "
61  << min_N_ ;
62 }
63 
65 
66 void
70  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltTowerMakerForEcal"));
71  desc.add<double>("MinE",10.);
72  desc.add<double>("MaxEta",3.);
73  desc.add<int>("MinN",1);
74  descriptions.add("hltEcalTowerFilter",desc);
75 }
76 
77 //
78 // member functions
79 //
80 
81 // ------------ method called to produce the data ------------
82  bool
84 {
85  using namespace std;
86  using namespace edm;
87  using namespace reco;
88  using namespace trigger;
89 
90  // All HLT filters must create and fill an HLT filter object,
91  // recording any reconstructed physics objects satisfying (or not)
92  // this HLT filter, and place it in the Event.
93 
94  // The filter object
95  if (saveTags()) filterproduct.addCollectionTag(inputTag_);
96 
97  // get hold of collection of objects
99  event.getByToken(inputToken_, towers);
100 
101  LogDebug("HLTEcalTowerFilter") << "Number of towers: " << towers->size();
102 
103  // look at all objects, check cuts and add to filter object
104  int n = 0;
105  for (auto const & i : *towers) {
106  if (i.emEnergy() >= min_E_ and fabs(i.eta()) <= max_Eta_) {
107  ++n;
108  //edm::Ref<CaloTowerCollection> ref(towers, std::distance(towers->begin(), i));
109  //filterproduct.addObject(TriggerJet, ref);
110  }
111  }
112 
113  LogDebug("HLTEcalTowerFilter") << "Number of towers with eta < " << max_Eta_ << " and energy > " << min_E_ << ": " << n;
114 
115  // filter decision
116  bool accept(n>=min_N_);
117 
118  return accept;
119 }
120 
121 // define as a framework module
#define LogDebug(id)
HLTEcalTowerFilter(const edm::ParameterSet &)
~HLTEcalTowerFilter() override
def setup(process, global_tag, zero_tesla=False)
Definition: GeneralSetup.py:2
Definition: config.py:1
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
std::string encode() const
Definition: InputTag.cc:159
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
edm::EDGetTokenT< CaloTowerCollection > inputToken_
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
size_type size() const
Definition: event.py:1