CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ElectronRegressionValueMapProducer.cc
Go to the documentation of this file.
3 
6 
8 
11 
14 
18 
19 #include "TVector2.h"
20 
21 #include <memory>
22 #include <vector>
23 
24 namespace {
25  constexpr char sigmaIEtaIPhi_[] = "sigmaIEtaIPhi";
26  constexpr char eMax_[] = "eMax";
27  constexpr char e2nd_[] = "e2nd";
28  constexpr char eTop_[] = "eTop";
29  constexpr char eBottom_[] = "eBottom";
30  constexpr char eLeft_[] = "eLeft";
31  constexpr char eRight_[] = "eRight";
32  constexpr char clusterMaxDR_[] = "clusterMaxDR";
33  constexpr char clusterMaxDRDPhi_[] = "clusterMaxDRDPhi";
34  constexpr char clusterMaxDRDEta_[] = "clusterMaxDRDEta";
35  constexpr char clusterMaxDRRawEnergy_[] = "clusterMaxDRRawEnergy";
36  constexpr char clusterRawEnergy0_[] = "clusterRawEnergy0";
37  constexpr char clusterRawEnergy1_[] = "clusterRawEnergy1";
38  constexpr char clusterRawEnergy2_[] = "clusterRawEnergy2";
39  constexpr char clusterDPhiToSeed0_[] = "clusterDPhiToSeed0";
40  constexpr char clusterDPhiToSeed1_[] = "clusterDPhiToSeed1";
41  constexpr char clusterDPhiToSeed2_[] = "clusterDPhiToSeed2";
42  constexpr char clusterDEtaToSeed0_[] = "clusterDEtaToSeed0";
43  constexpr char clusterDEtaToSeed1_[] = "clusterDEtaToSeed1";
44  constexpr char clusterDEtaToSeed2_[] = "clusterDEtaToSeed2";
45  constexpr char eleIPhi_[] = "iPhi";
46  constexpr char eleIEta_[] = "iEta";
47  constexpr char eleCryPhi_[] = "cryPhi";
48  constexpr char eleCryEta_[] = "cryEta";
49 }
50 
52 
53  public:
54 
57 
58  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
59 
60  private:
61 
62  virtual void produce(edm::Event&, const edm::EventSetup&) override;
63 
66  const std::vector<float> & values,
67  const std::string & label) const ;
68 
69  void writeValueMap(edm::Event &iEvent,
71  const std::vector<int> & values,
72  const std::string & label) const ;
73 
74  std::unique_ptr<EcalClusterLazyToolsBase> lazyTools;
75 
76  // for AOD case
81 
82  // for miniAOD case
87 
88  const bool use_full5x5_;
89 };
90 
92  use_full5x5_(iConfig.getParameter<bool>("useFull5x5")) {
93 
94  //
95  // Declare consummables, handle both AOD and miniAOD case
96  //
97  ebReducedRecHitCollection_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
98  ("ebReducedRecHitCollection"));
99  ebReducedRecHitCollectionMiniAOD_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
100  ("ebReducedRecHitCollectionMiniAOD"));
101 
102  eeReducedRecHitCollection_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
103  ("eeReducedRecHitCollection"));
104  eeReducedRecHitCollectionMiniAOD_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
105  ("eeReducedRecHitCollectionMiniAOD"));
106 
107  esReducedRecHitCollection_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
108  ("esReducedRecHitCollection"));
109  esReducedRecHitCollectionMiniAOD_ = mayConsume<EcalRecHitCollection>(iConfig.getParameter<edm::InputTag>
110  ("esReducedRecHitCollectionMiniAOD"));
111 
112  src_ = mayConsume<edm::View<reco::GsfElectron> >(iConfig.getParameter<edm::InputTag>("src"));
113  srcMiniAOD_ = mayConsume<edm::View<reco::GsfElectron> >(iConfig.getParameter<edm::InputTag>("srcMiniAOD"));
114 
115  produces<edm::ValueMap<float> >(sigmaIEtaIPhi_);
116  produces<edm::ValueMap<float> >(eMax_);
117  produces<edm::ValueMap<float> >(e2nd_);
118  produces<edm::ValueMap<float> >(eTop_);
119  produces<edm::ValueMap<float> >(eBottom_);
120  produces<edm::ValueMap<float> >(eLeft_);
121  produces<edm::ValueMap<float> >(eRight_);
122  produces<edm::ValueMap<float> >(clusterMaxDR_);
123  produces<edm::ValueMap<float> >(clusterMaxDRDPhi_);
124  produces<edm::ValueMap<float> >(clusterMaxDRDEta_);
125  produces<edm::ValueMap<float> >(clusterMaxDRRawEnergy_);
126  produces<edm::ValueMap<float> >(clusterRawEnergy0_);
127  produces<edm::ValueMap<float> >(clusterRawEnergy1_);
128  produces<edm::ValueMap<float> >(clusterRawEnergy2_);
129  produces<edm::ValueMap<float> >(clusterDPhiToSeed0_);
130  produces<edm::ValueMap<float> >(clusterDPhiToSeed1_);
131  produces<edm::ValueMap<float> >(clusterDPhiToSeed2_);
132  produces<edm::ValueMap<float> >(clusterDEtaToSeed0_);
133  produces<edm::ValueMap<float> >(clusterDEtaToSeed1_);
134  produces<edm::ValueMap<float> >(clusterDEtaToSeed2_);
135  produces<edm::ValueMap<int> >(eleIPhi_);
136  produces<edm::ValueMap<int> >(eleIEta_);
137  produces<edm::ValueMap<float> >(eleCryPhi_);
138  produces<edm::ValueMap<float> >(eleCryEta_);
139 }
140 
142 }
143 
144 template<typename LazyTools>
145 inline void calculateValues(EcalClusterLazyToolsBase* tools_tocast,
146  const edm::Ptr<reco::GsfElectron>& iEle,
147  const edm::EventSetup& iSetup,
148  std::vector<float>& vsigmaIEtaIPhi,
149  std::vector<float>& veMax,
150  std::vector<float>& ve2nd,
151  std::vector<float>& veTop,
152  std::vector<float>& veBottom,
153  std::vector<float>& veLeft,
154  std::vector<float>& veRight,
155  std::vector<float>& vclusterMaxDR,
156  std::vector<float>& vclusterMaxDRDPhi,
157  std::vector<float>& vclusterMaxDRDEta,
158  std::vector<float>& vclusterMaxDRRawEnergy,
159  std::vector<float>& vclusterRawEnergy0,
160  std::vector<float>& vclusterRawEnergy1,
161  std::vector<float>& vclusterRawEnergy2,
162  std::vector<float>& vclusterDPhiToSeed0,
163  std::vector<float>& vclusterDPhiToSeed1,
164  std::vector<float>& vclusterDPhiToSeed2,
165  std::vector<float>& vclusterDEtaToSeed0,
166  std::vector<float>& vclusterDEtaToSeed1,
167  std::vector<float>& vclusterDEtaToSeed2,
168  std::vector<int>& veleIPhi,
169  std::vector<int>& veleIEta,
170  std::vector<float>& veleCryPhi,
171  std::vector<float>& veleCryEta) {
172  LazyTools* tools = static_cast<LazyTools*>(tools_tocast);
173 
174  const auto& the_sc = iEle->superCluster();
175  const auto& theseed = the_sc->seed();
176 
177  std::vector<float> vCov = tools->localCovariances( *theseed );
178 
179  const float eMax = tools->eMax( *theseed );
180  const float e2nd = tools->e2nd( *theseed );
181  const float eTop = tools->eTop( *theseed );
182  const float eLeft = tools->eLeft( *theseed );
183  const float eRight = tools->eRight( *theseed );
184  const float eBottom = tools->eBottom( *theseed );
185 
186  float dummy;
187  int iPhi;
188  int iEta;
189  float cryPhi;
190  float cryEta;
191  EcalClusterLocal _ecalLocal;
192  if (iEle->isEB())
193  _ecalLocal.localCoordsEB(*theseed, iSetup, cryEta, cryPhi, iEta, iPhi, dummy, dummy);
194  else
195  _ecalLocal.localCoordsEE(*theseed, iSetup, cryEta, cryPhi, iEta, iPhi, dummy, dummy);
196 
197  double see = (isnan(vCov[0]) ? 0. : sqrt(vCov[0]));
198  double spp = (isnan(vCov[2]) ? 0. : sqrt(vCov[2]));
199  double sep;
200  if (see*spp > 0)
201  sep = vCov[1] / (see * spp);
202  else if (vCov[1] > 0)
203  sep = 1.0;
204  else
205  sep = -1.0;
206 
207  vsigmaIEtaIPhi.push_back(sep);
208  veMax.push_back(eMax);
209  ve2nd.push_back(e2nd);
210  veTop.push_back(eTop);
211  veBottom.push_back(eBottom);
212  veLeft.push_back(eLeft);
213  veRight.push_back(eRight);
214  veleIPhi.push_back(iPhi);
215  veleIEta.push_back(iEta);
216  veleCryPhi.push_back(cryPhi);
217  veleCryEta.push_back(cryEta);
218 
219  // loop over all clusters that aren't the seed
220  auto clusend = the_sc->clustersEnd();
221  int numberOfClusters = the_sc->clusters().size();
222 
223  std::vector<float> _clusterRawEnergy;
224  _clusterRawEnergy.resize(std::max(3, numberOfClusters), 0);
225  std::vector<float> _clusterDEtaToSeed;
226  _clusterDEtaToSeed.resize(std::max(3, numberOfClusters), 0);
227  std::vector<float> _clusterDPhiToSeed;
228  _clusterDPhiToSeed.resize(std::max(3, numberOfClusters), 0);
229  float _clusterMaxDR = 999.;
230  float _clusterMaxDRDPhi = 999.;
231  float _clusterMaxDRDEta = 999.;
232  float _clusterMaxDRRawEnergy = 0.;
233 
234  size_t iclus = 0;
235  float maxDR = 0;
237  for( auto clus = the_sc->clustersBegin(); clus != clusend; ++clus ) {
238  pclus = *clus;
239 
240  if(theseed == pclus )
241  continue;
242  _clusterRawEnergy.push_back(pclus->energy());
243  _clusterDPhiToSeed.push_back(reco::deltaPhi(pclus->phi(),theseed->phi()));
244  _clusterDEtaToSeed.push_back(pclus->eta() - theseed->eta());
245 
246  // find cluster with max dR
247  if(reco::deltaR(*pclus, *theseed) > maxDR) {
248  maxDR = reco::deltaR(*pclus, *theseed);
249  _clusterMaxDR = maxDR;
250  _clusterMaxDRDPhi = _clusterDPhiToSeed[iclus];
251  _clusterMaxDRDEta = _clusterDEtaToSeed[iclus];
252  _clusterMaxDRRawEnergy = _clusterRawEnergy[iclus];
253  }
254 
255  ++iclus;
256  }
257 
258  vclusterMaxDR.push_back(_clusterMaxDR);
259  vclusterMaxDRDPhi.push_back(_clusterMaxDRDPhi);
260  vclusterMaxDRDEta.push_back(_clusterMaxDRDEta);
261  vclusterMaxDRRawEnergy.push_back(_clusterMaxDRRawEnergy);
262  vclusterRawEnergy0.push_back(_clusterRawEnergy[0]);
263  vclusterRawEnergy1.push_back(_clusterRawEnergy[1]);
264  vclusterRawEnergy2.push_back(_clusterRawEnergy[2]);
265  vclusterDPhiToSeed0.push_back(_clusterDPhiToSeed[0]);
266  vclusterDPhiToSeed1.push_back(_clusterDPhiToSeed[1]);
267  vclusterDPhiToSeed2.push_back(_clusterDPhiToSeed[2]);
268  vclusterDEtaToSeed0.push_back(_clusterDEtaToSeed[0]);
269  vclusterDEtaToSeed1.push_back(_clusterDEtaToSeed[1]);
270  vclusterDEtaToSeed2.push_back(_clusterDEtaToSeed[2]);
271 }
272 
274 
275  using namespace edm;
276 
278 
279  // Retrieve the collection of electrons from the event.
280  // If we fail to retrieve the collection with the standard AOD
281  // name, we next look for the one with the stndard miniAOD name.
282  bool isAOD = true;
283  iEvent.getByToken(src_, src);
284 
285  if( !src.isValid() ) {
286  isAOD = false;
287  iEvent.getByToken(srcMiniAOD_,src);
288  }
289 
290  // configure lazy tools
292 
293  if( isAOD ) {
297  } else {
301  }
302 
303  if( use_full5x5_ ) {
304  lazyTools = std::make_unique<noZS::EcalClusterLazyTools>(iEvent, iSetup,
305  ebrh, eerh, esrh );
306  } else {
307  lazyTools = std::make_unique<EcalClusterLazyTools>(iEvent, iSetup,
308  ebrh, eerh, esrh );
309  }
310 
311  std::vector<float> sigmaIEtaIPhi;
312  std::vector<float> eMax;
313  std::vector<float> e2nd;
314  std::vector<float> eTop;
315  std::vector<float> eBottom;
316  std::vector<float> eLeft;
317  std::vector<float> eRight;
318  std::vector<float> clusterMaxDR;
319  std::vector<float> clusterMaxDRDPhi;
320  std::vector<float> clusterMaxDRDEta;
321  std::vector<float> clusterMaxDRRawEnergy;
322  std::vector<float> clusterRawEnergy0;
323  std::vector<float> clusterRawEnergy1;
324  std::vector<float> clusterRawEnergy2;
325  std::vector<float> clusterDPhiToSeed0;
326  std::vector<float> clusterDPhiToSeed1;
327  std::vector<float> clusterDPhiToSeed2;
328  std::vector<float> clusterDEtaToSeed0;
329  std::vector<float> clusterDEtaToSeed1;
330  std::vector<float> clusterDEtaToSeed2;
331  std::vector<int> eleIPhi;
332  std::vector<int> eleIEta;
333  std::vector<float> eleCryPhi;
334  std::vector<float> eleCryEta;
335 
336  // reco::GsfElectron::superCluster() is virtual so we can exploit polymorphism
337  for (size_t i = 0; i < src->size(); ++i){
338  auto iEle = src->ptrAt(i);
339 
340  if( use_full5x5_ ) {
341  calculateValues<noZS::EcalClusterLazyTools>(lazyTools.get(),
342  iEle,
343  iSetup,
344  sigmaIEtaIPhi,
345  eMax,
346  e2nd,
347  eTop,
348  eBottom,
349  eLeft,
350  eRight,
351  clusterMaxDR,
352  clusterMaxDRDPhi,
353  clusterMaxDRDEta,
354  clusterMaxDRRawEnergy,
355  clusterRawEnergy0,
356  clusterRawEnergy1,
357  clusterRawEnergy2,
358  clusterDPhiToSeed0,
359  clusterDPhiToSeed1,
360  clusterDPhiToSeed2,
361  clusterDEtaToSeed0,
362  clusterDEtaToSeed1,
363  clusterDEtaToSeed2,
364  eleIPhi,
365  eleIEta,
366  eleCryPhi,
367  eleCryEta);
368  } else {
369  calculateValues<EcalClusterLazyTools>(lazyTools.get(),
370  iEle,
371  iSetup,
372  sigmaIEtaIPhi,
373  eMax,
374  e2nd,
375  eTop,
376  eBottom,
377  eLeft,
378  eRight,
379  clusterMaxDR,
380  clusterMaxDRDPhi,
381  clusterMaxDRDEta,
382  clusterMaxDRRawEnergy,
383  clusterRawEnergy0,
384  clusterRawEnergy1,
385  clusterRawEnergy2,
386  clusterDPhiToSeed0,
387  clusterDPhiToSeed1,
388  clusterDPhiToSeed2,
389  clusterDEtaToSeed0,
390  clusterDEtaToSeed1,
391  clusterDEtaToSeed2,
392  eleIPhi,
393  eleIEta,
394  eleCryPhi,
395  eleCryEta);
396  }
397  }
398 
399  writeValueMap(iEvent, src, sigmaIEtaIPhi, sigmaIEtaIPhi_);
400  writeValueMap(iEvent, src, eMax ,eMax_);
401  writeValueMap(iEvent, src, e2nd ,e2nd_);
402  writeValueMap(iEvent, src, eTop ,eTop_);
403  writeValueMap(iEvent, src, eBottom ,eBottom_);
404  writeValueMap(iEvent, src, eLeft ,eLeft_);
405  writeValueMap(iEvent, src, eRight ,eRight_);
406  writeValueMap(iEvent, src, clusterMaxDR, clusterMaxDR_);
407  writeValueMap(iEvent, src, clusterMaxDRDPhi, clusterMaxDRDPhi_);
408  writeValueMap(iEvent, src, clusterMaxDRDEta, clusterMaxDRDEta_);
409  writeValueMap(iEvent, src, clusterMaxDRRawEnergy,clusterMaxDRRawEnergy_);
410  writeValueMap(iEvent, src, clusterRawEnergy0, clusterRawEnergy0_);
411  writeValueMap(iEvent, src, clusterRawEnergy1, clusterRawEnergy1_);
412  writeValueMap(iEvent, src, clusterRawEnergy2, clusterRawEnergy2_);
413  writeValueMap(iEvent, src, clusterDPhiToSeed0, clusterDPhiToSeed0_);
414  writeValueMap(iEvent, src, clusterDPhiToSeed1, clusterDPhiToSeed1_);
415  writeValueMap(iEvent, src, clusterDPhiToSeed2, clusterDPhiToSeed2_);
416  writeValueMap(iEvent, src, clusterDEtaToSeed0, clusterDEtaToSeed0_);
417  writeValueMap(iEvent, src, clusterDEtaToSeed1, clusterDEtaToSeed1_);
418  writeValueMap(iEvent, src, clusterDEtaToSeed2, clusterDEtaToSeed2_);
419  writeValueMap(iEvent, src, eleIPhi ,eleIPhi_);
420  writeValueMap(iEvent, src, eleIEta ,eleIEta_);
421  writeValueMap(iEvent, src, eleCryPhi ,eleCryPhi_);
422  writeValueMap(iEvent, src, eleCryEta ,eleCryEta_);
423  lazyTools.reset();
424 }
425 
428  const std::vector<float> & values,
429  const std::string & label) const
430 {
431  using namespace edm;
432  using namespace std;
433  auto_ptr<ValueMap<float> > valMap(new ValueMap<float>());
434  edm::ValueMap<float>::Filler filler(*valMap);
435  filler.insert(handle, values.begin(), values.end());
436  filler.fill();
437  iEvent.put(valMap, label);
438 }
439 
442  const std::vector<int> & values,
443  const std::string & label) const
444 {
445  using namespace edm;
446  using namespace std;
447  auto_ptr<ValueMap<int> > valMap(new ValueMap<int>());
448  edm::ValueMap<int>::Filler filler(*valMap);
449  filler.insert(handle, values.begin(), values.end());
450  filler.fill();
451  iEvent.put(valMap, label);
452 }
453 
455  //The following says we do not know what parameters are allowed so do no validation
456  // Please change this to state exactly what you do use, even if it is no parameters
458  desc.setUnknown();
459  descriptions.addDefault(desc);
460 }
461 
T getParameter(std::string const &) const
edm::EDGetTokenT< EcalRecHitCollection > eeReducedRecHitCollection_
int i
Definition: DBlmapReader.cc:9
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
edm::EDGetTokenT< EcalRecHitCollection > esReducedRecHitCollection_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:457
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
void insert(const H &h, I begin, I end)
Definition: ValueMap.h:52
double deltaR(const T1 &t1, const T2 &t2)
Definition: deltaR.h:48
#define constexpr
edm::EDGetTokenT< EcalRecHitCollection > esReducedRecHitCollectionMiniAOD_
int iEvent
Definition: GenABIO.cc:230
void addDefault(ParameterSetDescription const &psetDescription)
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:115
bool isnan(float x)
Definition: math.h:13
edm::EDGetTokenT< EcalRecHitCollection > ebReducedRecHitCollection_
T sqrt(T t)
Definition: SSEVec.h:48
edm::EDGetTokenT< EcalRecHitCollection > eeReducedRecHitCollectionMiniAOD_
tuple handle
Definition: patZpeak.py:22
void writeValueMap(edm::Event &iEvent, const edm::Handle< edm::View< reco::GsfElectron > > &handle, const std::vector< float > &values, const std::string &label) const
std::unique_ptr< EcalClusterLazyToolsBase > lazyTools
double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:12
void calculateValues(EcalClusterLazyToolsBase *tools_tocast, const edm::Ptr< reco::GsfElectron > &iEle, const edm::EventSetup &iSetup, std::vector< float > &vsigmaIEtaIPhi, std::vector< float > &veMax, std::vector< float > &ve2nd, std::vector< float > &veTop, std::vector< float > &veBottom, std::vector< float > &veLeft, std::vector< float > &veRight, std::vector< float > &vclusterMaxDR, std::vector< float > &vclusterMaxDRDPhi, std::vector< float > &vclusterMaxDRDEta, std::vector< float > &vclusterMaxDRRawEnergy, std::vector< float > &vclusterRawEnergy0, std::vector< float > &vclusterRawEnergy1, std::vector< float > &vclusterRawEnergy2, std::vector< float > &vclusterDPhiToSeed0, std::vector< float > &vclusterDPhiToSeed1, std::vector< float > &vclusterDPhiToSeed2, std::vector< float > &vclusterDEtaToSeed0, std::vector< float > &vclusterDEtaToSeed1, std::vector< float > &vclusterDEtaToSeed2, std::vector< int > &veleIPhi, std::vector< int > &veleIEta, std::vector< float > &veleCryPhi, std::vector< float > &veleCryEta)
edm::EDGetTokenT< EcalRecHitCollection > ebReducedRecHitCollectionMiniAOD_
virtual void produce(edm::Event &, const edm::EventSetup &) override
void localCoordsEB(const reco::CaloCluster &bclus, const edm::EventSetup &es, float &etacry, float &phicry, int &ieta, int &iphi, float &thetatilt, float &phitilt) const
void localCoordsEE(const reco::CaloCluster &bclus, const edm::EventSetup &es, float &xcry, float &ycry, int &ix, int &iy, float &thetatilt, float &phitilt) const