CMS 3D CMS Logo

HLTEgammaL1MatchFilterPairs.cc
Go to the documentation of this file.
1 
9 
10 //#include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h"
11 
13 
15 
17 
20 
25 
26 #include <vector>
27 #include <cmath>
28 
29 #define TWOPI 2*M_PI
30 //
31 // constructors and destructor
32 //
34 {
35  candIsolatedTag_ = iConfig.getParameter< edm::InputTag > ("candIsolatedTag");
36  l1IsolatedTag_ = iConfig.getParameter< edm::InputTag > ("l1IsolatedTag");
37  candNonIsolatedTag_ = iConfig.getParameter< edm::InputTag > ("candNonIsolatedTag");
38  l1NonIsolatedTag_ = iConfig.getParameter< edm::InputTag > ("l1NonIsolatedTag");
39  L1SeedFilterTag_ = iConfig.getParameter< edm::InputTag > ("L1SeedFilterTag");
40 
41  AlsoNonIsolatedFirst_ = iConfig.getParameter<bool>("AlsoNonIsolatedFirst");
42  AlsoNonIsolatedSecond_ = iConfig.getParameter<bool>("AlsoNonIsolatedSecond");
43 
44  region_eta_size_ = iConfig.getParameter<double> ("region_eta_size");
45  region_eta_size_ecap_ = iConfig.getParameter<double> ("region_eta_size_ecap");
46  region_phi_size_ = iConfig.getParameter<double> ("region_phi_size");
47  barrel_end_ = iConfig.getParameter<double> ("barrel_end");
48  endcap_end_ = iConfig.getParameter<double> ("endcap_end");
49 
50  candIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candIsolatedTag_);
51  candNonIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candNonIsolatedTag_);
52  L1SeedFilterToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(L1SeedFilterTag_);
53 }
54 
56 
57 void
61  desc.add<edm::InputTag>("candIsolatedTag",edm::InputTag("hltRecoIsolatedEcalCandidate"));
62  desc.add<edm::InputTag>("l1IsolatedTag",edm::InputTag("l1extraParticles","Isolated"));
63  desc.add<edm::InputTag>("candNonIsolatedTag",edm::InputTag("hltRecoNonIsolatedEcalCandidate"));
64  desc.add<edm::InputTag>("l1NonIsolatedTag",edm::InputTag("l1extraParticles","NonIsolated"));
65  desc.add<edm::InputTag>("L1SeedFilterTag",edm::InputTag("theL1SeedFilter"));
66  desc.add<bool>("AlsoNonIsolatedFirst",false);
67  desc.add<bool>("AlsoNonIsolatedSecond",false);
68  desc.add<double>("region_eta_size",0.522);
69  desc.add<double>("region_eta_size_ecap",1.0);
70  desc.add<double>("region_phi_size",1.044);
71  desc.add<double>("barrel_end",1.4791);
72  desc.add<double>("endcap_end",2.65);
73  descriptions.add("hltEgammaL1MatchFilterPairs",desc);
74 }
75 
76 
77 
78 // ------------ method called to produce the data ------------
79 bool
81 {
82 
83  using namespace trigger;
84  using namespace l1extra;
85  std::vector < std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > > thePairs;
86 
88  iEvent.getByToken(candIsolatedToken_,recoIsolecalcands);
90  iEvent.getByToken(candNonIsolatedToken_,recoNonIsolecalcands);
91 
92  // create pairs <L1Iso,L1Iso> and optionally <L1Iso, L1NonIso>
93  for (auto recoecalcand1= recoIsolecalcands->begin(); recoecalcand1!=recoIsolecalcands->end(); recoecalcand1++) {
94  edm::Ref<reco::RecoEcalCandidateCollection> ref1 = edm::Ref<reco::RecoEcalCandidateCollection>(recoIsolecalcands, distance(recoIsolecalcands->begin(),recoecalcand1) );
95  for (auto recoecalcand2= recoIsolecalcands->begin(); recoecalcand2!=recoIsolecalcands->end(); recoecalcand2++) {
96  edm::Ref<reco::RecoEcalCandidateCollection> ref2 = edm::Ref<reco::RecoEcalCandidateCollection>(recoIsolecalcands, distance(recoIsolecalcands->begin(),recoecalcand2) );
97  if( &(*ref1) != &(*ref2) ) {thePairs.push_back(std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > (ref1,ref2) );}
98  }
100  for (auto recoecalcand2= recoNonIsolecalcands->begin(); recoecalcand2!=recoNonIsolecalcands->end(); recoecalcand2++) {
101  edm::Ref<reco::RecoEcalCandidateCollection> ref2 = edm::Ref<reco::RecoEcalCandidateCollection>(recoNonIsolecalcands, distance(recoNonIsolecalcands->begin(),recoecalcand2) );
102  if( &(*ref1) != &(*ref2) ) {thePairs.push_back(std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > (ref1,ref2) );}
103  }
104  }
105  }
106 
107 
108  // create pairs <L1NonIso,L1Iso> and optionally <L1NonIso, L1NonIso>
110  for (auto recoecalcand1= recoNonIsolecalcands->begin(); recoecalcand1!=recoNonIsolecalcands->end(); recoecalcand1++) {
111  edm::Ref<reco::RecoEcalCandidateCollection> ref1 = edm::Ref<reco::RecoEcalCandidateCollection>(recoNonIsolecalcands, distance(recoNonIsolecalcands->begin(),recoecalcand1) );
112  for (auto recoecalcand2= recoIsolecalcands->begin(); recoecalcand2!=recoIsolecalcands->end(); recoecalcand2++) {
113  edm::Ref<reco::RecoEcalCandidateCollection> ref2 = edm::Ref<reco::RecoEcalCandidateCollection>(recoIsolecalcands, distance(recoIsolecalcands->begin(),recoecalcand2) );
114  if( &(*ref1) != &(*ref2) ) {thePairs.push_back(std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > (ref1,ref2) );}
115  }
117  for (auto recoecalcand2= recoNonIsolecalcands->begin(); recoecalcand2!=recoNonIsolecalcands->end(); recoecalcand2++) {
118  edm::Ref<reco::RecoEcalCandidateCollection> ref2 = edm::Ref<reco::RecoEcalCandidateCollection>(recoNonIsolecalcands, distance(recoNonIsolecalcands->begin(),recoecalcand2) );
119  if( &(*ref1) != &(*ref2) ) {thePairs.push_back(std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > (ref1,ref2) );}
120  }
121  }
122  }
123  }
124 
125 
126  // Get the CaloGeometry
127  edm::ESHandle<L1CaloGeometry> l1CaloGeom ;
128  iSetup.get<L1CaloGeometryRecord>().get(l1CaloGeom) ;
129 
130 
131  // look at all candidates, check cuts and add to filter object
132  int n(0);
133 
135  iEvent.getByToken (L1SeedFilterToken_,L1SeedOutput);
136 
137  std::vector<l1extra::L1EmParticleRef > l1EGIso;
138  L1SeedOutput->getObjects(TriggerL1IsoEG, l1EGIso);
139  std::vector<l1extra::L1EmParticleRef > l1EGNonIso;
140  L1SeedOutput->getObjects(TriggerL1NoIsoEG, l1EGNonIso);
141 
142 // std::cout<<"L1EGIso size: "<<l1EGIso.size()<<std::endl;
143 // for (unsigned int i=0; i<l1EGIso.size(); i++){std::cout<<"L1EGIso Et Eta phi: "<<l1EGIso[i]->et()<<" "<<l1EGIso[i]->eta()<<" "<<l1EGIso[i]->phi()<<std::endl;}
144 // std::cout<<"L1EGNonIso size: "<<l1EGNonIso.size()<<std::endl;
145 // for (unsigned int i=0; i<l1EGNonIso.size(); i++){std::cout<<"L1EGNonIso Et Eta phi: "<<l1EGNonIso[i]->et()<<" "<<l1EGNonIso[i]->eta()<<" "<<l1EGNonIso[i]->phi()<<std::endl;}
146 // std::cout<<"Lpair vector size: "<<thePairs.size()<<std::endl;
147 
148  std::vector < std::pair< edm::Ref<reco::RecoEcalCandidateCollection>, edm::Ref<reco::RecoEcalCandidateCollection> > >::iterator pairsIt;
149  for (pairsIt = thePairs.begin(); pairsIt != thePairs.end(); pairsIt++){
150 // edm::Ref<reco::RecoEcalCandidateCollection> r1 = pairsIt->first;
151 // edm::Ref<reco::RecoEcalCandidateCollection> r2 = pairsIt->second;
152 // std::cout<<"1) Et Eta phi: "<<r1->et()<<" "<<r1->eta()<<" "<<r1->phi()<<" 2) Et eta phi: "<<r2->et()<<" "<<r2->eta()<<" "<<r2->phi()<<std::endl;
153 
154  if ( CheckL1Matching(pairsIt->first,l1EGIso,l1EGNonIso) && CheckL1Matching(pairsIt->second,l1EGIso,l1EGNonIso) ){
155  filterproduct.addObject(TriggerCluster, pairsIt->first);
156  filterproduct.addObject(TriggerCluster, pairsIt->second);
157  n++;
158  }
159  }
160 
161 
162  // std::cout<<"#####################################################"<<std::endl;
163  // filter decision
164  bool accept(n>=1);
165 
166  return accept;
167 }
168 
169 bool HLTEgammaL1MatchFilterPairs::CheckL1Matching(edm::Ref<reco::RecoEcalCandidateCollection> ref, std::vector<l1extra::L1EmParticleRef >& l1EGIso, std::vector<l1extra::L1EmParticleRef >& l1EGNonIso) const {
170 
171  for (auto & i : l1EGIso) {
172  //ORCA matching method
173  double etaBinLow = 0.;
174  double etaBinHigh = 0.;
175  if(fabs(ref->eta()) < barrel_end_){
176  etaBinLow = i->eta() - region_eta_size_/2.;
177  etaBinHigh = etaBinLow + region_eta_size_;
178  }
179  else{
180  etaBinLow = i->eta() - region_eta_size_ecap_/2.;
181  etaBinHigh = etaBinLow + region_eta_size_ecap_;
182  }
183 
184  float deltaphi=fabs(ref->phi() -i->phi());
185  if(deltaphi>TWOPI) deltaphi-=TWOPI;
186  if(deltaphi>M_PI) deltaphi=TWOPI-deltaphi;
187 
188  if(ref->eta() < etaBinHigh && ref->eta() > etaBinLow && deltaphi <region_phi_size_/2. ) {return true;}
189 
190  }
191 
192  for (auto & i : l1EGNonIso) {
193  //ORCA matching method
194  double etaBinLow = 0.;
195  double etaBinHigh = 0.;
196  if(fabs(ref->eta()) < barrel_end_){
197  etaBinLow = i->eta() - region_eta_size_/2.;
198  etaBinHigh = etaBinLow + region_eta_size_;
199  }
200  else{
201  etaBinLow = i->eta() - region_eta_size_ecap_/2.;
202  etaBinHigh = etaBinLow + region_eta_size_ecap_;
203  }
204 
205  float deltaphi=fabs(ref->phi() - i->phi());
206  if(deltaphi>TWOPI) deltaphi-=TWOPI;
207  if(deltaphi>M_PI) deltaphi=TWOPI-deltaphi;
208 
209  if(ref->eta() < etaBinHigh && ref->eta() > etaBinLow && deltaphi <region_phi_size_/2. ) {return true;}
210 
211  }
212 
213  return false;
214 }
215 
216 // 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
edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > L1SeedFilterToken_
bool accept(const edm::Event &event, const edm::TriggerResults &triggerTable, const std::string &triggerPath)
Definition: TopDQMHelpers.h:30
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
ParameterDescriptionBase * add(U const &iLabel, T const &value)
#define M_PI
~HLTEgammaL1MatchFilterPairs() override
edm::EDGetTokenT< reco::RecoEcalCandidateCollection > candNonIsolatedToken_
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void add(std::string const &label, ParameterSetDescription const &psetDescription)
edm::EDGetTokenT< reco::RecoEcalCandidateCollection > candIsolatedToken_
T get() const
Definition: EventSetup.h:71
bool CheckL1Matching(edm::Ref< reco::RecoEcalCandidateCollection >ref, std::vector< l1extra::L1EmParticleRef > &l1EGIso, std::vector< l1extra::L1EmParticleRef > &l1EGNonIso) const
HLTEgammaL1MatchFilterPairs(const edm::ParameterSet &)