CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
EcalListOfFEDSProducer.cc
Go to the documentation of this file.
1 
2 //#include <FWCore/Framework/interface/Handle.h>
6 
9 
11 
12 // Ecal Mapping
16 
17 
18 // Level 1 Trigger
23 
24 // EgammaCoreTools
27 
28 // Muon stuff
32 
33 // Jets stuff
36 #include <vector>
37 
38 using namespace l1extra;
39 
40 
42 
43  debug_ = pset.getUntrackedParameter<bool>("debug");
44 
45  Pi0ListToIgnore_ = pset.getParameter<edm::InputTag>("Pi0ListToIgnore");
46 
47  EGamma_ = pset.getUntrackedParameter<bool>("EGamma",false);
48  Muon_ = pset.getUntrackedParameter<bool>("Muon",false);
49  Jets_ = pset.getUntrackedParameter<bool>("Jets",false);
50 
51  if (EGamma_ && Muon_) {
52  throw cms::Exception("EcalListOfFEDSProducer") <<
53  " Wrong configuration : EGamma and Muon should not be true at the same time." ;
54  }
55 
56 
57  if (EGamma_) {
58 
59  EMl1TagIsolated_ = pset.getUntrackedParameter<edm::InputTag>("EM_l1TagIsolated");
60  EMl1TagNonIsolated_ = pset.getUntrackedParameter<edm::InputTag>("EM_l1TagNonIsolated");
61  EMdoIsolated_ = pset.getUntrackedParameter<bool>("EM_doIsolated",true);
62  EMdoNonIsolated_ = pset.getUntrackedParameter<bool>("EM_doNonIsolated",true);
63  EMregionEtaMargin_ = pset.getUntrackedParameter<double>("EM_regionEtaMargin",0.25);
64  EMregionPhiMargin_ = pset.getUntrackedParameter<double>("EM_regionPhiMargin",0.40);
65  Ptmin_iso_ = pset.getUntrackedParameter<double>("Ptmin_iso",0.);
66  Ptmin_noniso_ = pset.getUntrackedParameter<double>("Ptmin_noniso",0.);
67  }
68 
69  if (Muon_) {
70  MUregionEtaMargin_ = pset.getUntrackedParameter<double>("MU_regionEtaMargin",1.0);
71  MUregionPhiMargin_ = pset.getUntrackedParameter<double>("MU_regionPhiMargin",1.0);
72  Ptmin_muon_ = pset.getUntrackedParameter<double>("Ptmin_muon",0.);
73  MuonSource_ = pset.getUntrackedParameter<edm::InputTag>("MuonSource");
74  }
75 
76  if (Jets_) {
77  JETSregionEtaMargin_ = pset.getUntrackedParameter<double>("JETS_regionEtaMargin",1.0);
78  JETSregionPhiMargin_ = pset.getUntrackedParameter<double>("JETS_regionPhiMargin",1.0);
79  Ptmin_jets_ = pset.getUntrackedParameter<double>("Ptmin_jets",0.);
80  CentralSource_ = pset.getUntrackedParameter<edm::InputTag>("CentralSource");
81  ForwardSource_ = pset.getUntrackedParameter<edm::InputTag>("ForwardSource");
82  TauSource_ = pset.getUntrackedParameter<edm::InputTag>("TauSource");
83  JETSdoCentral_ = pset.getUntrackedParameter<bool>("JETS_doCentral",true);
84  JETSdoForward_ = pset.getUntrackedParameter<bool>("JETS_doForward",true);
85  JETSdoTau_ = pset.getUntrackedParameter<bool>("JETS_doTau",true);
86  }
87 
88 
89  OutputLabel_ = pset.getUntrackedParameter<std::string>("OutputLabel");
90 
91  TheMapping = new EcalElectronicsMapping();
92  first_ = true;
93 
94  produces<EcalListOfFEDS>(OutputLabel_);
95 }
96 
97 
98 
100  delete TheMapping;
101 }
102 
103 
105 }
106 
108 }
109 
111 
112  if (first_) {
114  iSetup.get< EcalMappingRcd >().get(ecalmapping);
115  const EcalElectronicsMapping* TheMapping_ = ecalmapping.product();
116  *TheMapping = *TheMapping_;
117  first_ = false;
118  }
119 
120  std::auto_ptr<EcalListOfFEDS> productAddress(new EcalListOfFEDS);
121 
122  std::vector<int> feds; // the list of FEDS produced by this module
123 
124 
125  // ---- First, get the list of ECAL FEDs which have already been unpacked
126 
127  std::vector< edm::Handle<EcalListOfFEDS> > FEDs_Done;
128  std::vector<int> Done;
129  e.getManyByType(FEDs_Done);
130  unsigned int nDone = FEDs_Done.size();
131  if (debug_) std::cout << " ECAL unpacking module has already run " << nDone << " times. " << std::endl;
132  for (unsigned int id=0; id < nDone; id++) {
133  // ignore the FEDs coming from unpacking in pi0 paths
134  if( Pi0ListToIgnore_.label() == FEDs_Done[id].provenance()->moduleLabel() ){continue;}
135  std::vector<int> done = FEDs_Done[id]-> GetList();
136  for (int jd=0; jd < (int)done.size(); jd++) {
137  Done.push_back(done[jd] - FEDNumbering::MINECALFEDID);
138  }
139  }
140  if (debug_) std::cout << " For this event, " << Done.size() << " ECAL FEDs have already been unpacked." << std::endl;
141 
142 
143  if (EGamma_) {
144  // feds = Egamma(e, iSetup, Done);
145  Egamma(e, iSetup, Done, feds);
146  }
147 
148  if (Muon_) {
149  // feds = Muon(e, iSetup, Done);
150  Muon(e, iSetup, Done, feds);
151  }
152 
153  if (Jets_) {
154  // feds = Jets(e, iSetup, Done);
155  Jets(e, iSetup, Done, feds);
156  }
157 
158  if ( !EGamma_ && !Muon_ && ! Jets_) {
159  for (int i=1; i <= 54; i++) {
160  if ( std::find(Done.begin(), Done.end(), i) == Done.end())
161  feds.push_back(i);
162  }
163  }
164 
165  int nf = (int)feds.size();
166  for (int i=0; i <nf; i++) {
167  feds[i] += FEDNumbering::MINECALFEDID;
168  if (debug_) std::cout << "Will unpack FED " << feds[i] << std::endl;
169  }
170 
171  if (debug_ && nf < 1 )
172  std::cout << " Warning : no ECAL FED to unpack for Run " << e.id().run() << " Event " << e.id().event() << std::endl;
173 
174 
175  productAddress.get() -> SetList(feds);
176  e.put(productAddress,OutputLabel_);
177 }
178 
179 
180 // std::vector<int> EcalListOfFEDSProducer::Egamma(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done) {
181 void EcalListOfFEDSProducer::Egamma(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done, std::vector<int>& FEDs ) {
182 
183  // std::vector<int> FEDs;
184 
185  if (debug_) std::cout << std::endl << std::endl << " enter in EcalListOfFEDSProducer::Egamma" << std::endl;
186 
187  //Get the L1 EM Particle Collection
188  //Get the L1 EM Particle Collection
190  if(EMdoIsolated_)
191  e.getByLabel(EMl1TagIsolated_, emIsolColl);
192  //Get the L1 EM Particle Collection
194  if (EMdoNonIsolated_)
195  e.getByLabel(EMl1TagNonIsolated_, emNonIsolColl);
196 
197  // Get the CaloGeometry
198  edm::ESHandle<L1CaloGeometry> l1CaloGeom ;
199  es.get<L1CaloGeometryRecord>().get(l1CaloGeom) ;
200 
201  if(EMdoIsolated_) {
202 
203  for( l1extra::L1EmParticleCollection::const_iterator emItr = emIsolColl->begin();
204  emItr != emIsolColl->end() ;++emItr ){
205 
206  float pt = emItr -> pt();
207  if (pt < Ptmin_iso_ ) continue;
208  if (debug_) std::cout << " Here is an L1 isoEM candidate of pt " << pt << std::endl;
209  // Access the GCT hardware object corresponding to the L1Extra EM object.
210  int etaIndex = emItr->gctEmCand()->etaIndex() ;
211  int phiIndex = emItr->gctEmCand()->phiIndex() ;
212  // Use the L1CaloGeometry to find the eta, phi bin boundaries.
213  double etaLow = l1CaloGeom->etaBinLowEdge( etaIndex ) ;
214  double etaHigh = l1CaloGeom->etaBinHighEdge( etaIndex ) ;
215  double phiLow = l1CaloGeom->emJetPhiBinLowEdge( phiIndex ) ;
216  double phiHigh = l1CaloGeom->emJetPhiBinHighEdge( phiIndex ) ;
217 
218  std::vector<int> feds = ListOfFEDS(etaLow, etaHigh, phiLow, phiHigh, EMregionEtaMargin_, EMregionPhiMargin_);
219  for (int i=0; i < (int)feds.size(); i++) {
220  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
221  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
222  }
223 
224  } // end loop on L1EmParticleCollection
225 
226  } // endif doIsolated_
227 
228 
229  if (EMdoNonIsolated_) {
230 
231  for( l1extra::L1EmParticleCollection::const_iterator emItr = emNonIsolColl->begin();
232  emItr != emNonIsolColl->end() ;++emItr ){
233 
234  float pt = emItr -> pt();
235  if (debug_) std::cout << " Here is an L1 nonisoEM candidate of pt " << pt << std::endl;
236  if (pt < Ptmin_noniso_ ) continue;
237  // Access the GCT hardware object corresponding to the L1Extra EM object.
238  int etaIndex = emItr->gctEmCand()->etaIndex() ;
239  int phiIndex = emItr->gctEmCand()->phiIndex() ;
240  // std::cout << " etaIndex phiIndex " << etaIndex << " " << phiIndex << std::endl;
241  // Use the L1CaloGeometry to find the eta, phi bin boundaries.
242  double etaLow = l1CaloGeom->etaBinLowEdge( etaIndex ) ;
243  double etaHigh = l1CaloGeom->etaBinHighEdge( etaIndex ) ;
244  double phiLow = l1CaloGeom->emJetPhiBinLowEdge( phiIndex ) ;
245  double phiHigh = l1CaloGeom->emJetPhiBinHighEdge( phiIndex ) ;
246 
247  std::vector<int> feds = ListOfFEDS(etaLow, etaHigh, phiLow, phiHigh, EMregionEtaMargin_, EMregionPhiMargin_);
248  for (int i=0; i < (int)feds.size(); i++) {
249  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
250  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
251 
252  }
253 
254  } // end loop on L1EmParticleCollection
255  }
256 
257  // std::cout << "end of get list of feds " << std::endl;
258 
259  if (debug_) {
260  std::cout << std::endl;
261  for (int i=0; i < (int)FEDs.size(); i++) {
262  std::cout << "Egamma: unpack FED " << FEDs[i] << std::endl;
263  }
264  std::cout << "Number of FEDS is " << FEDs.size() << std::endl;
265  }
266 
267  // return FEDs;
268 
269 }
270 
271 
272 // std::vector<int> EcalListOfFEDSProducer::Muon(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done) {
273 void EcalListOfFEDSProducer::Muon(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done, std::vector<int>& FEDs) {
274 
275  // std::vector<int> FEDs;
276 
277  if (debug_) std::cout << std::endl << std::endl << " enter in EcalListOfFEDSProducer::Muon" << std::endl;
278 
280  e.getByLabel(MuonSource_, muColl);
281 
282 
283  double epsilon = 0.01;
284 
285  for (L1MuonParticleCollection::const_iterator it=muColl->begin(); it != muColl->end(); it++) {
286 
287  const L1MuGMTExtendedCand muonCand = (*it).gmtMuonCand();
288  double pt = (*it).pt();
289  double eta = (*it).eta();
290  double phi = (*it).phi();
291 
292  if (debug_) std::cout << " here is a L1 muon Seed with (eta,phi) = " <<
293  eta << " " << phi << " and pt " << pt << std::endl;
294  if (pt < Ptmin_muon_ ) continue;
295 
296  std::vector<int> feds = ListOfFEDS(eta, eta, phi-epsilon, phi+epsilon, MUregionEtaMargin_, MUregionPhiMargin_);
297 
298  for (int i=0; i < (int)feds.size(); i++) {
299  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
300  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
301  }
302  }
303 
304  if (debug_) {
305  std::cout << std::endl;
306  for (int i=0; i < (int)FEDs.size(); i++) {
307  std::cout << "Muons: unpack FED " << FEDs[i] << std::endl;
308  }
309  std::cout << "Number of FEDS is " << FEDs.size() << std::endl;
310  }
311 
312 
313  // return FEDs;
314 
315 }
316 
317 
318 
319 // std::vector<int> EcalListOfFEDSProducer::Jets(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done) {
320 void EcalListOfFEDSProducer::Jets(edm::Event& e, const edm::EventSetup& es, std::vector<int>& done, std::vector<int>& FEDs) {
321 
322  // std::vector<int> FEDs;
323 
324  if (debug_) std::cout << std::endl << std::endl << " enter in EcalListOfFEDSProducer::Jets" << std::endl;
325  double epsilon = 0.01;
326 
327  // Get the CaloGeometry
328  // edm::ESHandle<L1CaloGeometry> l1CaloGeom ;
329  // es.get<L1CaloGeometryRecord>().get(l1CaloGeom) ;
330 
331  if (JETSdoCentral_) {
332 
334  e.getByLabel(CentralSource_,jetColl);
335 
336  for (L1JetParticleCollection::const_iterator it=jetColl->begin(); it != jetColl->end(); it++) {
337 
338  double pt = it -> pt();
339  double eta = it -> eta();
340  double phi = it -> phi();
341 
342  if (debug_) std::cout << " here is a L1 CentralJet Seed with (eta,phi) = " <<
343  eta << " " << phi << " and pt " << pt << std::endl;
344  // int etaIndex = it->gctJetCand()->etaIndex() ;
345  // int phiIndex = it->gctJetCand()->phiIndex() ;
346  // double etaLow = l1CaloGeom->etaBinLowEdge( etaIndex ) ;
347  // double etaHigh = l1CaloGeom->etaBinHighEdge( etaIndex ) ;
348  // double phiLow = l1CaloGeom->emJetPhiBinLowEdge( phiIndex ) ;
349  // double phiHigh = l1CaloGeom->emJetPhiBinHighEdge( phiIndex ) ;
350  // if (debug_) std::cout << " etaLow .. " << etaLow << " " << etaHigh << " " << phiLow << " " << phiHigh << std::endl;
351 
352  if (pt < Ptmin_jets_ ) continue;
353 
354  std::vector<int> feds = ListOfFEDS(eta, eta, phi-epsilon, phi+epsilon, JETSregionEtaMargin_, JETSregionPhiMargin_);
355 
356  for (int i=0; i < (int)feds.size(); i++) {
357  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
358  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
359  }
360  }
361  }
362 
363  if (JETSdoForward_) {
364 
366  e.getByLabel(ForwardSource_,jetColl);
367 
368  for (L1JetParticleCollection::const_iterator it=jetColl->begin(); it != jetColl->end(); it++) {
369 
370  double pt = it -> pt();
371  double eta = it -> eta();
372  double phi = it -> phi();
373 
374  if (debug_) std::cout << " here is a L1 ForwardJet Seed with (eta,phi) = " <<
375  eta << " " << phi << " and pt " << pt << std::endl;
376  if (pt < Ptmin_jets_ ) continue;
377 
378  std::vector<int> feds = ListOfFEDS(eta, eta, phi-epsilon, phi+epsilon, JETSregionEtaMargin_, JETSregionPhiMargin_);
379 
380  for (int i=0; i < (int)feds.size(); i++) {
381  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
382  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
383  }
384  }
385  }
386 
387  if (JETSdoTau_) {
388 
390  e.getByLabel(TauSource_,jetColl);
391 
392  for (L1JetParticleCollection::const_iterator it=jetColl->begin(); it != jetColl->end(); it++) {
393 
394  double pt = it -> pt();
395  double eta = it -> eta();
396  double phi = it -> phi();
397 
398  if (debug_) std::cout << " here is a L1 TauJet Seed with (eta,phi) = " <<
399  eta << " " << phi << " and pt " << pt << std::endl;
400  if (pt < Ptmin_jets_ ) continue;
401 
402  std::vector<int> feds = ListOfFEDS(eta, eta, phi-epsilon, phi+epsilon, JETSregionEtaMargin_, JETSregionPhiMargin_);
403 
404  for (int i=0; i < (int)feds.size(); i++) {
405  if ( std::find(FEDs.begin(), FEDs.end(), feds[i]) == FEDs.end() &&
406  std::find(done.begin(), done.end(), feds[i]) == done.end() ) FEDs.push_back(feds[i]);
407  }
408  }
409  }
410 
411 
412 
413  if (debug_) {
414  std::cout << std::endl;
415  for (int i=0; i < (int)FEDs.size(); i++) {
416  std::cout << "Jets: unpack FED " << FEDs[i] << std::endl;
417  }
418  std::cout << "Number of FEDS is " << FEDs.size() << std::endl;
419  }
420 
421 
422  // return FEDs;
423 
424 }
425 
426 
427 std::vector<int> EcalListOfFEDSProducer::ListOfFEDS(double etaLow, double etaHigh, double phiLow,
428  double phiHigh, double etamargin, double phimargin)
429 {
430 
431  std::vector<int> FEDs;
432 
433  if (phimargin > Geom::pi()) phimargin = Geom::pi() ;
434 
435 
436  if (debug_) std::cout << " etaLow etaHigh phiLow phiHigh " << etaLow << " " <<
437  etaHigh << " " << phiLow << " " << phiHigh << std::endl;
438 
439  etaLow -= etamargin;
440  etaHigh += etamargin;
441  double phiMinus = phiLow - phimargin;
442  double phiPlus = phiHigh + phimargin;
443 
444  bool all = false;
445  double dd = fabs(phiPlus-phiMinus);
446  if (debug_) std::cout << " dd = " << dd << std::endl;
447  if (dd > 2.*Geom::pi() ) all = true;
448 
449  while (phiPlus > Geom::pi()) { phiPlus -= 2.*Geom::pi() ; }
450  while (phiMinus < 0) { phiMinus += 2.*Geom::pi() ; }
451  if ( phiMinus > Geom::pi()) phiMinus -= 2.*Geom::pi() ;
452 
453  double dphi = phiPlus - phiMinus;
454  if (dphi < 0) dphi += 2.*Geom::pi() ;
455  if (debug_) std::cout << "dphi = " << dphi << std::endl;
456  if (dphi > Geom::pi()) {
457  int fed_low1 = TheMapping -> GetFED(etaLow,phiMinus*180./Geom::pi());
458  int fed_low2 = TheMapping -> GetFED(etaLow,phiPlus*180./Geom::pi());
459  if (debug_) std::cout << "fed_low1 fed_low2 " << fed_low1 << " " << fed_low2 << std::endl;
460  if (fed_low1 == fed_low2) all = true;
461  int fed_hi1 = TheMapping -> GetFED(etaHigh,phiMinus*180./Geom::pi());
462  int fed_hi2 = TheMapping -> GetFED(etaHigh,phiPlus*180./Geom::pi());
463  if (debug_) std::cout << "fed_hi1 fed_hi2 " << fed_hi1 << " " << fed_hi2 << std::endl;
464  if (fed_hi1 == fed_hi2) all = true;
465  }
466 
467 
468  if (all) {
469  if (debug_) std::cout << " unpack everything in phi ! " << std::endl;
470  phiMinus = -20 * Geom::pi() / 180.; // -20 deg
471  phiPlus = -40 * Geom::pi() / 180.; // -20 deg
472  }
473 
474  if (debug_) std::cout << " with margins : " << etaLow << " " << etaHigh << " " <<
475  phiMinus << " " << phiPlus << std::endl;
476 
477 
478  const EcalEtaPhiRegion ecalregion(etaLow,etaHigh,phiMinus,phiPlus);
479 
480  FEDs = TheMapping -> GetListofFEDs(ecalregion);
481 
482 /*
483  if (debug_) {
484  int nn = (int)FEDs.size();
485  for (int ii=0; ii < nn; ii++) {
486  std::cout << "unpack fed " << FEDs[ii] << std::endl;
487  }
488  }
489 */
490 
491  return FEDs;
492 
493 }
494 
495 
496 
497 /*
498 bool EcalListOfFEDSProducer::PhiInbetween(double phiPlus,double phiMinus,double phiLow) {
499 
500  bool all = false;
501  // first, everything back between 0 and 2pi
502 
503  if (phiPlus < 0) phiPlus += 2.*Geom::pi();
504  if (phiMinus < 0) phiMinus += 2.*Geom::pi();
505  if (phiLow < 0) phiLow += 2.*Geom::pi();
506 
507  if (phiMinus <= phiPlus && phiPlus <= phiLow) all = true;
508  return all;
509 
510 }
511 */
512 
513 
514 
515 
516 
517 
518 
RunNumber_t run() const
Definition: EventID.h:42
void getManyByType(std::vector< Handle< PROD > > &results) const
Definition: Event.h:408
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:44
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
void produce(edm::Event &e, const edm::EventSetup &c)
void Muon(edm::Event &e, const edm::EventSetup &es, std::vector< int > &done, std::vector< int > &FEDs)
T eta() const
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
void Egamma(edm::Event &e, const edm::EventSetup &es, std::vector< int > &done, std::vector< int > &FEDs)
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:85
tuple Jets
Definition: METSkim_cff.py:17
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
void Jets(edm::Event &e, const edm::EventSetup &es, std::vector< int > &done, std::vector< int > &FEDs)
EcalListOfFEDSProducer(const edm::ParameterSet &pset)
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:62
std::vector< int > ListOfFEDS(double etaLow, double etaHigh, double phiLow, double phiHigh, double etamargin, double phimargin)
edm::EventID id() const
Definition: EventBase.h:56
double pi()
Definition: Pi.h:31
tuple cout
Definition: gather_cfg.py:121
const double epsilon
Definition: DDAxes.h:10