CMS 3D CMS Logo

HLTEgammaDoubleEtDeltaPhiFilter.cc
Go to the documentation of this file.
1 
10 
12 
16 
18 
19 //
20 // constructors and destructor
21 //
23 {
24  inputTag_ = iConfig.getParameter< edm::InputTag > ("inputTag");
25  etcut_ = iConfig.getParameter<double> ("etcut");
26  minDeltaPhi_ = iConfig.getParameter<double> ("minDeltaPhi");
27  l1EGTag_ = iConfig.getParameter< edm::InputTag > ("l1EGCand");
28  inputToken_ = consumes<trigger::TriggerFilterObjectWithRefs> (inputTag_);
29 }
30 
32 
33 void
37  desc.add<edm::InputTag>("inputTag", edm::InputTag("hltDoublePhotonEt5L1MatchFilterRegional"));
38  desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
39  desc.add<double>("etcut", 5.0);
40  desc.add<double>("minDeltaPhi", 2.5);
41  descriptions.add("hltEgammaDoubleEtDeltaPhiFilter", desc);
42 }
43 
44 // ------------ method called to produce the data ------------
45 bool
47 {
48  using namespace trigger;
49  // The filter object
50  if (saveTags()) {
51  filterproduct.addCollectionTag(l1EGTag_);
52  }
53 
54  // get hold of filtered candidates
56  iEvent.getByToken (inputToken_,PrevFilterOutput);
57 
58  std::vector<edm::Ref<reco::RecoEcalCandidateCollection> > recoecalcands;
59  PrevFilterOutput->getObjects(TriggerCluster, recoecalcands);
60  if(recoecalcands.empty()) PrevFilterOutput->getObjects(TriggerPhoton,recoecalcands); //we dont know if its type trigger cluster or trigger photon
61 
62  // Refs to the two Candidate objects used to calculate deltaPhi
64 
65  // look at all candidates, check cuts
66  int n(0);
67  for(auto & ref : recoecalcands) {
68  if( ref->et() >= etcut_) {
69  ++n;
70  if(n==1) ref1 = ref;
71  if(n==2) ref2 = ref;
72  }
73  }
74 
75  // if there are only two Candidate objects, calculate deltaPhi
76  double deltaPhi(0.0);
77  if(n==2) {
78  deltaPhi = fabs(ref1->phi()-ref2->phi());
79  if(deltaPhi>M_PI) deltaPhi = 2*M_PI - deltaPhi;
80 
81  filterproduct.addObject(TriggerCluster, ref1);
82  filterproduct.addObject(TriggerCluster, ref2);
83  }
84 
85  // filter decision
86  bool accept(n==2 && deltaPhi>minDeltaPhi_);
87 
88  return accept;
89 }
90 
91 // define as a framework module
92 // #include "FWCore/Framework/interface/MakerMacros.h"
93 // DEFINE_FWK_MODULE(HLTEgammaDoubleEtDeltaPhiFilter);
94 
95 
96 // declare this class as a framework plugin
T getParameter(std::string const &) const
void getObjects(Vids &ids, VRphoton &refs) const
various physics-level getters:
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
HLTEgammaDoubleEtDeltaPhiFilter(const edm::ParameterSet &)
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
int iEvent
Definition: GenABIO.cc:224
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > inputToken_
ParameterDescriptionBase * add(U const &iLabel, T const &value)
#define M_PI
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)
bool saveTags() const
Definition: HLTFilter.h:45