test
CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ElectronSeedGenerator.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: EgammaElectronAlgos
4 // Class: ElectronSeedGenerator.
5 //
13 //
14 // Original Author: Ursula Berthon, Claude Charlot
15 // Created: Mon Mar 27 13:22:06 CEST 2006
16 //
17 
20 
26 
27 //#include "DataFormats/EgammaReco/interface/EcalCluster.h"
29 
31 //#include "DataFormats/BeamSpot/interface/BeamSpot.h"
32 //#include "DataFormats/Common/interface/Handle.h"
33 
37 
41 
43 
45 
46 
47 #include <vector>
48 #include <utility>
49 
52  : dynamicphiroad_(pset.getParameter<bool>("dynamicPhiRoad")),
53  fromTrackerSeeds_(pset.getParameter<bool>("fromTrackerSeeds")),
54  useRecoVertex_(false),
55  verticesTag_(ts.token_vtx),
56  beamSpotTag_(ts.token_bs),
57  lowPtThreshold_(pset.getParameter<double>("LowPtThreshold")),
58  highPtThreshold_(pset.getParameter<double>("HighPtThreshold")),
59  nSigmasDeltaZ1_(pset.getParameter<double>("nSigmasDeltaZ1")),
60  deltaZ1WithVertex_(0.5),
61  sizeWindowENeg_(pset.getParameter<double>("SizeWindowENeg")),
62  deltaPhi1Low_(pset.getParameter<double>("DeltaPhi1Low")),
63  deltaPhi1High_(pset.getParameter<double>("DeltaPhi1High")),
64  deltaPhi1Coef1_(0.), deltaPhi1Coef2_(0.),
65  myMatchEle(0), myMatchPos(0),
66  thePropagator(0),
67  theMeasurementTracker(0),
68  theMeasurementTrackerEventTag(ts.token_measTrkEvt),
69  theSetup(0),
70  cacheIDMagField_(0),/*cacheIDGeom_(0),*/cacheIDNavSchool_(0),cacheIDCkfComp_(0),cacheIDTrkGeom_(0)
71 {
72  // so that deltaPhi1 = deltaPhi1Coef1_ + deltaPhi1Coef2_/clusterEnergyT
73  if (dynamicphiroad_)
74  {
77  }
78 
79  theMeasurementTrackerName = pset.getParameter<std::string>("measurementTrackerName");
80 
81  // use of reco vertex
82  useRecoVertex_ = pset.getParameter<bool>("useRecoVertex");
83  deltaZ1WithVertex_ = pset.getParameter<double>("deltaZ1WithVertex");
84 
85  // new B/F configurables
86  deltaPhi2B_ = pset.getParameter<double>("DeltaPhi2B") ;
87  deltaPhi2F_ = pset.getParameter<double>("DeltaPhi2F") ;
88 
89  phiMin2B_ = pset.getParameter<double>("PhiMin2B") ;
90  phiMin2F_ = pset.getParameter<double>("PhiMin2F") ;
91 
92  phiMax2B_ = pset.getParameter<double>("PhiMax2B") ;
93  phiMax2F_ = pset.getParameter<double>("PhiMax2F") ;
94 
95  // Instantiate the pixel hit matchers
97  ( pset.getParameter<double>("ePhiMin1"),
98  pset.getParameter<double>("ePhiMax1"),
99  phiMin2B_, phiMax2B_, phiMin2F_, phiMax2F_,
100  pset.getParameter<double>("z2MinB"),
101  pset.getParameter<double>("z2MaxB"),
102  pset.getParameter<double>("r2MinF"),
103  pset.getParameter<double>("r2MaxF"),
104  pset.getParameter<double>("rMinI"),
105  pset.getParameter<double>("rMaxI"),
106  pset.getParameter<bool>("searchInTIDTEC") ) ;
107 
109  ( pset.getParameter<double>("pPhiMin1"),
110  pset.getParameter<double>("pPhiMax1"),
111  phiMin2B_, phiMax2B_, phiMin2F_, phiMax2F_,
112  pset.getParameter<double>("z2MinB"),
113  pset.getParameter<double>("z2MaxB"),
114  pset.getParameter<double>("r2MinF"),
115  pset.getParameter<double>("r2MaxF"),
116  pset.getParameter<double>("rMinI"),
117  pset.getParameter<double>("rMaxI"),
118  pset.getParameter<bool>("searchInTIDTEC") ) ;
119 
120  theUpdator = new KFUpdator() ;
121 }
122 
124  {
125  delete myMatchEle ;
126  delete myMatchPos ;
127  delete thePropagator ;
128  delete theUpdator ;
129  }
130 
132 
133  // get records if necessary (called once per event)
134  bool tochange=false;
135 
138  cacheIDMagField_=setup.get<IdealMagneticFieldRecord>().cacheIdentifier();
139  if (thePropagator) delete thePropagator;
141  tochange=true;
142  }
143 
145  edm::ESHandle<MeasurementTracker> measurementTrackerHandle;
146  setup.get<CkfComponentsRecord>().get(theMeasurementTrackerName,measurementTrackerHandle);
147  cacheIDCkfComp_=setup.get<CkfComponentsRecord>().cacheIdentifier();
148  theMeasurementTracker = measurementTrackerHandle.product();
149  tochange=true;
150  }
151 
152  //edm::ESHandle<TrackerGeometry> trackerGeometryHandle;
154  cacheIDTrkGeom_=setup.get<TrackerDigiGeometryRecord>().cacheIdentifier();
156  tochange=true; //FIXME
157  }
158 
159  if (tochange) {
162  }
163 
166  setup.get<NavigationSchoolRecord>().get("SimpleNavigationSchool", nav);
167  cacheIDNavSchool_=setup.get<NavigationSchoolRecord>().cacheIdentifier();
169  }
170 
171 // if (cacheIDGeom_!=setup.get<TrackerRecoGeometryRecord>().cacheIdentifier()) {
172 // setup.get<TrackerRecoGeometryRecord>().get( theGeomSearchTracker );
173 // cacheIDGeom_=setup.get<TrackerRecoGeometryRecord>().cacheIdentifier();
174 // }
175 
176 }
177 
178 void display_seed( const std::string & title1, const std::string & title2, const reco::ElectronSeed & seed, edm::ESHandle<TrackerGeometry> trackerGeometry )
179  {
180  const PTrajectoryStateOnDet & startingState = seed.startingState() ;
181  const LocalTrajectoryParameters & parameters = startingState.parameters() ;
182  std::cout<<title1
183  <<" ("<<seed.subDet2()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<")"
184  <<" ("<<seed.direction()<<"/"<<startingState.detId()<<"/"<<startingState.surfaceSide()<<"/"<<parameters.charge()<<"/"<<parameters.position()<<"/"<<parameters.momentum()<<")"
185  <<std::endl ;
186  }
187 
188 bool equivalent( const TrajectorySeed & s1, const TrajectorySeed & s2 )
189  {
190  if (s1.nHits()!=s2.nHits()) return false ;
191 
192  unsigned int nHits ;
193  TrajectorySeed::range r1 = s1.recHits(), r2 = s2.recHits() ;
195  for ( i1=r1.first, i2=r2.first, nHits=0 ; i1!=r1.second ; ++i1, ++i2, ++nHits )
196  {
197  if ( !i1->isValid() || !i2->isValid() ) return false ;
198  if ( i1->geographicalId()!=i2->geographicalId() ) return false ;
199  if ( ! ( i1->localPosition()==i2->localPosition() ) ) return false ;
200  }
201 
202  return true ;
203  }
204 
207  const reco::SuperClusterRefVector & sclRefs, const std::vector<float> & hoe1s, const std::vector<float> & hoe2s,
209  {
210  theInitialSeedColl = seeds ;
211 // bool duplicateTrajectorySeeds =false ;
212 // unsigned int i,j ;
213 // for (i=0;i<seeds->size();++i)
214 // for (j=i+1;j<seeds->size();++j)
215 // {
216 // const TrajectorySeed & s1 =(*seeds)[i] ;
217 // const TrajectorySeed & s2 =(*seeds)[j] ;
218 // if ( equivalent(s1,s2) )
219 // {
220 // const PTrajectoryStateOnDet & ss1 = s1.startingState() ;
221 // const LocalTrajectoryParameters & p1 = ss1.parameters() ;
222 // const PTrajectoryStateOnDet & ss2 = s2.startingState() ;
223 // const LocalTrajectoryParameters & p2 = ss2.parameters() ;
224 // duplicateTrajectorySeeds = true ;
225 // std::cout<<"Same hits for "
226 // <<"\n s["<<i<<"] ("<<s1.direction()<<"/"<<ss1.detId()<<"/"<<ss1.surfaceSide()<<"/"<<p1.charge()<<"/"<<p1.position()<<"/"<<p1.momentum()<<")"
227 // <<"\n s["<<j<<"] ("<<s2.direction()<<"/"<<ss2.detId()<<"/"<<ss2.surfaceSide()<<"/"<<p2.charge()<<"/"<<p2.position()<<"/"<<p2.momentum()<<")"
228 // <<std::endl ;
229 // }
230 // }
231 // if (duplicateTrajectorySeeds)
232 // { edm::LogWarning("ElectronSeedGenerator|DuplicateTrajectorySeeds")<<"We see several identical trajectory seeds." ; }
233 
234  //Retrieve tracker topology from geometry
236  setup.get<TrackerTopologyRcd>().get(tTopoHand);
237  const TrackerTopology *tTopo=tTopoHand.product();
238 
239  theSetup= &setup;
240 
241 
242  // Step A: set Event for the TrajectoryBuilder
244  e.getByToken(theMeasurementTrackerEventTag, data);
245  myMatchEle->setEvent(*data);
246  myMatchPos->setEvent(*data);
247 
248  // get initial TrajectorySeeds if necessary
249  // if (fromTrackerSeeds_) e.getByToken(initialSeeds_, theInitialSeedColl);
250 
251  // get the beamspot from the Event:
252  //e.getByType(theBeamSpot);
253  e.getByToken(beamSpotTag_,theBeamSpot);
254 
255  // if required get the vertices
256  if (useRecoVertex_) e.getByToken(verticesTag_,theVertices);
257 
258  if (!fromTrackerSeeds_)
259  { throw cms::Exception("NotSupported") << "Here in ElectronSeedGenerator " << __FILE__ << ":" << __LINE__ << " I would like to do theMeasurementTracker->update(e); but that no longer makes sense.\n";
260  }
261 
262  for (unsigned int i=0;i<sclRefs.size();++i) {
263  // Find the seeds
264  recHits_.clear();
265 
266  LogDebug ("run") << "new cluster, calling seedsFromThisCluster";
267  seedsFromThisCluster(sclRefs[i],hoe1s[i],hoe2s[i],out,tTopo);
268  }
269 
270  LogDebug ("run") << ": For event "<<e.id();
271  LogDebug ("run") <<"Nr of superclusters after filter: "<<sclRefs.size()
272  <<", no. of ElectronSeeds found = " << out.size();
273 }
274 
277  float hoe1, float hoe2,
278  reco::ElectronSeedCollection & out, const TrackerTopology *tTopo )
279 {
280  float clusterEnergy = seedCluster->energy() ;
281  GlobalPoint clusterPos
282  ( seedCluster->position().x(),
283  seedCluster->position().y(),
284  seedCluster->position().z() ) ;
285  reco::ElectronSeed::CaloClusterRef caloCluster(seedCluster) ;
286 
287  //LogDebug("") << "[ElectronSeedGenerator::seedsFromThisCluster] new supercluster with energy: " << clusterEnergy ;
288  //LogDebug("") << "[ElectronSeedGenerator::seedsFromThisCluster] and position: " << clusterPos ;
289 
290  if (dynamicphiroad_)
291  {
292  float clusterEnergyT = clusterEnergy / cosh( EleRelPoint(clusterPos,theBeamSpot->position()).eta() ) ;
293 
294  float deltaPhi1 ;
295  if (clusterEnergyT < lowPtThreshold_)
296  { deltaPhi1= deltaPhi1Low_ ; }
297  else if (clusterEnergyT > highPtThreshold_)
298  { deltaPhi1= deltaPhi1High_ ; }
299  else
300  { deltaPhi1 = deltaPhi1Coef1_ + deltaPhi1Coef2_/clusterEnergyT ; }
301 
302  float ephimin1 = -deltaPhi1*sizeWindowENeg_ ;
303  float ephimax1 = deltaPhi1*(1.-sizeWindowENeg_);
304  float pphimin1 = -deltaPhi1*(1.-sizeWindowENeg_);
305  float pphimax1 = deltaPhi1*sizeWindowENeg_;
306 
307  float phimin2B = -deltaPhi2B_/2. ;
308  float phimax2B = deltaPhi2B_/2. ;
309  float phimin2F = -deltaPhi2F_/2. ;
310  float phimax2F = deltaPhi2F_/2. ;
311 
312 
313  myMatchEle->set1stLayer(ephimin1,ephimax1);
314  myMatchPos->set1stLayer(pphimin1,pphimax1);
315  myMatchEle->set2ndLayer(phimin2B,phimax2B, phimin2F,phimax2F);
316  myMatchPos->set2ndLayer(phimin2B,phimax2B, phimin2F,phimax2F);
317  }
318 
320 
321  if (!useRecoVertex_) // here use the beam spot position
322  {
323  double sigmaZ=theBeamSpot->sigmaZ();
324  double sigmaZ0Error=theBeamSpot->sigmaZ0Error();
325  double sq=sqrt(sigmaZ*sigmaZ+sigmaZ0Error*sigmaZ0Error);
326  double myZmin1=theBeamSpot->position().z()-nSigmasDeltaZ1_*sq;
327  double myZmax1=theBeamSpot->position().z()+nSigmasDeltaZ1_*sq;
328 
329  GlobalPoint vertexPos ;
330  ele_convert(theBeamSpot->position(),vertexPos) ;
331 
332  myMatchEle->set1stLayerZRange(myZmin1,myZmax1);
333  myMatchPos->set1stLayerZRange(myZmin1,myZmax1);
334 
335  if (!fromTrackerSeeds_)
336  {
337  // try electron
338  std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> > elePixelHits
339  = myMatchEle->compatibleHits(clusterPos,vertexPos,
340  clusterEnergy,-1., tTopo, *theNavigationSchool) ;
341  GlobalPoint eleVertex(theBeamSpot->position().x(),theBeamSpot->position().y(),myMatchEle->getVertex()) ;
342  seedsFromRecHits(elePixelHits,dir,eleVertex,caloCluster,out,false) ;
343  // try positron
344  std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> > posPixelHits
345  = myMatchPos->compatibleHits(clusterPos,vertexPos,clusterEnergy,1.,tTopo, *theNavigationSchool) ;
346  GlobalPoint posVertex(theBeamSpot->position().x(),theBeamSpot->position().y(),myMatchPos->getVertex()) ;
347  seedsFromRecHits(posPixelHits,dir,posVertex,caloCluster,out,true) ;
348  }
349  else
350  {
351  // try electron
352  std::vector<SeedWithInfo> elePixelSeeds
353  = myMatchEle->compatibleSeeds(theInitialSeedColl,clusterPos,vertexPos,clusterEnergy,-1.) ;
354  seedsFromTrajectorySeeds(elePixelSeeds,caloCluster,hoe1,hoe2,out,false) ;
355  // try positron
356  std::vector<SeedWithInfo> posPixelSeeds
357  = myMatchPos->compatibleSeeds(theInitialSeedColl,clusterPos,vertexPos,clusterEnergy,1.) ;
358  seedsFromTrajectorySeeds(posPixelSeeds,caloCluster,hoe1,hoe2,out,true) ;
359  }
360 
361  }
362  else // here we use the reco vertices
363  {
364 
365  myMatchEle->setUseRecoVertex(true) ; //Hit matchers need to know that the vertex is known
366  myMatchPos->setUseRecoVertex(true) ;
367 
368  const std::vector<reco::Vertex> * vtxCollection = theVertices.product() ;
369  std::vector<reco::Vertex>::const_iterator vtxIter ;
370  for (vtxIter = vtxCollection->begin(); vtxIter != vtxCollection->end() ; vtxIter++)
371  {
372 
373  GlobalPoint vertexPos(vtxIter->position().x(),vtxIter->position().y(),vtxIter->position().z());
374  double myZmin1, myZmax1 ;
375  if (vertexPos.z()==theBeamSpot->position().z())
376  { // in case vetex not found
377  double sigmaZ=theBeamSpot->sigmaZ();
378  double sigmaZ0Error=theBeamSpot->sigmaZ0Error();
379  double sq=sqrt(sigmaZ*sigmaZ+sigmaZ0Error*sigmaZ0Error);
380  myZmin1=theBeamSpot->position().z()-nSigmasDeltaZ1_*sq;
381  myZmax1=theBeamSpot->position().z()+nSigmasDeltaZ1_*sq;
382  }
383  else
384  { // a vertex has been recoed
385  myZmin1=vtxIter->position().z()-deltaZ1WithVertex_;
386  myZmax1=vtxIter->position().z()+deltaZ1WithVertex_;
387  }
388 
389  myMatchEle->set1stLayerZRange(myZmin1,myZmax1);
390  myMatchPos->set1stLayerZRange(myZmin1,myZmax1);
391 
392  if (!fromTrackerSeeds_)
393  {
394  // try electron
395  std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> > elePixelHits
396  = myMatchEle->compatibleHits(clusterPos,vertexPos,clusterEnergy,-1.,tTopo, *theNavigationSchool) ;
397  seedsFromRecHits(elePixelHits,dir,vertexPos,caloCluster,out,false) ;
398  // try positron
399  std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> > posPixelHits
400  = myMatchPos->compatibleHits(clusterPos,vertexPos,clusterEnergy,1.,tTopo, *theNavigationSchool) ;
401  seedsFromRecHits(posPixelHits,dir,vertexPos,caloCluster,out,true) ;
402  }
403  else
404  {
405  // try electron
406  std::vector<SeedWithInfo> elePixelSeeds
407  = myMatchEle->compatibleSeeds(theInitialSeedColl,clusterPos,vertexPos,clusterEnergy,-1.) ;
408  seedsFromTrajectorySeeds(elePixelSeeds,caloCluster,hoe1,hoe2,out,false) ;
409  // try positron
410  std::vector<SeedWithInfo> posPixelSeeds
411  = myMatchPos->compatibleSeeds(theInitialSeedColl,clusterPos,vertexPos,clusterEnergy,1.) ;
412  seedsFromTrajectorySeeds(posPixelSeeds,caloCluster,hoe1,hoe2,out,true) ;
413  }
414  }
415  }
416 
417  return ;
418  }
419 
421  ( std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> > & pixelHits,
423  const GlobalPoint & vertexPos, const reco::ElectronSeed::CaloClusterRef & cluster,
425  bool positron )
426  {
427  if (!pixelHits.empty())
428  { LogDebug("ElectronSeedGenerator") << "Compatible "<<(positron?"positron":"electron")<<" hits found." ; }
429 
430  std::vector<std::pair<RecHitWithDist,ConstRecHitPointer> >::iterator v ;
431  for ( v = pixelHits.begin() ; v != pixelHits.end() ; v++ )
432  {
433  if (!positron)
434  { (*v).first.invert() ; }
435  if (!prepareElTrackSeed((*v).first.recHit(),(*v).second,vertexPos))
436  { continue ; }
437  reco::ElectronSeed seed(pts_,recHits_,dir) ;
438  seed.setCaloCluster(cluster) ;
439  addSeed(seed,0,positron,out) ;
440  }
441  }
442 
444  ( const std::vector<SeedWithInfo> & pixelSeeds,
445  const reco::ElectronSeed::CaloClusterRef & cluster,
446  float hoe1, float hoe2,
448  bool positron )
449  {
450  if (!pixelSeeds.empty())
451  { LogDebug("ElectronSeedGenerator") << "Compatible "<<(positron?"positron":"electron")<<" seeds found." ; }
452 
453  std::vector<SeedWithInfo>::const_iterator s;
454  for ( s = pixelSeeds.begin() ; s != pixelSeeds.end() ; s++ )
455  {
456  reco::ElectronSeed seed(s->seed()) ;
457  seed.setCaloCluster(cluster,s->hitsMask(),s->subDet2(),s->subDet1(),hoe1,hoe2) ;
458  addSeed(seed,&*s,positron,out) ;
459  }
460  }
461 
464  const SeedWithInfo * info,
465  bool positron,
467  {
468  if (!info)
469  { out.push_back(seed) ; return ; }
470 
471  if (positron)
472  { seed.setPosAttributes(info->dRz2(),info->dPhi2(),info->dRz1(),info->dPhi1()) ; }
473  else
474  { seed.setNegAttributes(info->dRz2(),info->dPhi2(),info->dRz1(),info->dPhi1()) ; }
475  reco::ElectronSeedCollection::iterator resItr ;
476  for ( resItr=out.begin() ; resItr!=out.end() ; ++resItr )
477  {
478  if ( (seed.caloCluster()==resItr->caloCluster()) &&
479  (seed.hitsMask()==resItr->hitsMask()) &&
480  equivalent(seed,*resItr) )
481  {
482  if (positron)
483  {
484  if ( resItr->dRz2Pos()==std::numeric_limits<float>::infinity() &&
485  resItr->dRz2()!=std::numeric_limits<float>::infinity() )
486  {
487  resItr->setPosAttributes(info->dRz2(),info->dPhi2(),info->dRz1(),info->dPhi1()) ;
488  seed.setNegAttributes(resItr->dRz2(),resItr->dPhi2(),resItr->dRz1(),resItr->dPhi1()) ;
489  break ;
490  }
491  else
492  {
493  if ( resItr->dRz2Pos()!=std::numeric_limits<float>::infinity() )
494  {
495  if ( resItr->dRz2Pos()!=seed.dRz2Pos() )
496  {
497  edm::LogWarning("ElectronSeedGenerator|BadValue")
498  <<"this similar old seed already has another dRz2Pos"
499  <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
500  <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
501  }
502 // else
503 // {
504 // edm::LogWarning("ElectronSeedGenerator|UnexpectedValue")
505 // <<"this old seed already knows its dRz2Pos, we suspect duplicates in input trajectry seeds"
506 // <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
507 // <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
508 // }
509  }
510 // if (resItr->dRz2()==std::numeric_limits<float>::infinity())
511 // {
512 // edm::LogWarning("ElectronSeedGenerator|BadValue")
513 // <<"this old seed has no dRz2, we suspect duplicates in input trajectry seeds"
514 // <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
515 // <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
516 // }
517  }
518  }
519  else
520  {
521  if ( resItr->dRz2()==std::numeric_limits<float>::infinity()
522  && resItr->dRz2Pos()!=std::numeric_limits<float>::infinity() )
523  {
524  resItr->setNegAttributes(info->dRz2(),info->dPhi2(),info->dRz1(),info->dPhi1()) ;
525  seed.setPosAttributes(resItr->dRz2Pos(),resItr->dPhi2Pos(),resItr->dRz1Pos(),resItr->dPhi1Pos()) ;
526  break ;
527  }
528  else
529  {
530  if ( resItr->dRz2()!=std::numeric_limits<float>::infinity() )
531  {
532  if (resItr->dRz2()!=seed.dRz2())
533  {
534  edm::LogWarning("ElectronSeedGenerator|BadValue")
535  <<"this old seed already has another dRz2"
536  <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
537  <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
538  }
539  // else
540  // {
541  // edm::LogWarning("ElectronSeedGenerator|UnexpectedValue")
542  // <<"this old seed already knows its dRz2, we suspect duplicates in input trajectry seeds"
543  // <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
544  // <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
545  // seed.setPosAttributes(resItr->dRz2Pos(),resItr->dPhi2Pos(),resItr->dRz1Pos(),resItr->dPhi1Pos()) ;
546  // }
547  }
548 // if (resItr->dRz2Pos()==std::numeric_limits<float>::infinity())
549 // {
550 // edm::LogWarning("ElectronSeedGenerator|BadValue")
551 // <<"this old seed has no dRz2Pos"
552 // <<"\nold seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)resItr->hitsMask()<<"/"<<resItr->dRz2()<<"/"<<resItr->dPhi2()<<"/"<<resItr->dRz2Pos()<<"/"<<resItr->dPhi2Pos()
553 // <<"\nnew seed mask/dRz2/dPhi2/dRz2Pos/dPhi2Pos: "<<(unsigned int)seed.hitsMask()<<"/"<<seed.dRz2()<<"/"<<seed.dPhi2()<<"/"<<seed.dRz2Pos()<<"/"<<seed.dPhi2Pos() ;
554 // }
555  }
556  }
557  }
558  }
559 
560  out.push_back(seed) ;
561  }
562 
564  ( ConstRecHitPointer innerhit,
565  ConstRecHitPointer outerhit,
566  const GlobalPoint& vertexPos )
567  {
568 
569  // debug prints
570  LogDebug("") <<"[ElectronSeedGenerator::prepareElTrackSeed] inner PixelHit x,y,z "<<innerhit->globalPosition();
571  LogDebug("") <<"[ElectronSeedGenerator::prepareElTrackSeed] outer PixelHit x,y,z "<<outerhit->globalPosition();
572 
573  recHits_.clear();
574 
575  SiPixelRecHit *pixhit=0;
576  SiStripMatchedRecHit2D *striphit=0;
577  const SiPixelRecHit* constpixhit = dynamic_cast <const SiPixelRecHit*> (innerhit->hit());
578  if (constpixhit) {
579  pixhit=new SiPixelRecHit(*constpixhit);
580  recHits_.push_back(pixhit);
581  } else return false;
582  constpixhit = dynamic_cast <const SiPixelRecHit *> (outerhit->hit());
583  if (constpixhit) {
584  pixhit=new SiPixelRecHit(*constpixhit);
585  recHits_.push_back(pixhit);
586  } else {
587  const SiStripMatchedRecHit2D * conststriphit=dynamic_cast <const SiStripMatchedRecHit2D *> (outerhit->hit());
588  if (conststriphit) {
589  striphit = new SiStripMatchedRecHit2D(*conststriphit);
590  recHits_.push_back(striphit);
591  } else return false;
592  }
593 
595 
596  // FIXME to be optimized outside the loop
598  theSetup->get<IdealMagneticFieldRecord>().get(bfield);
599  float nomField = bfield->nominalValue();
600 
601  // make a spiral
602  FastHelix helix(outerhit->globalPosition(),innerhit->globalPosition(),vertexPos,nomField,&*bfield);
603  if ( !helix.isValid()) {
604  return false;
605  }
606  FreeTrajectoryState fts(helix.stateAtVertex());
607  TSOS propagatedState = thePropagator->propagate(fts,innerhit->det()->surface()) ;
608  if (!propagatedState.isValid())
609  return false;
610  TSOS updatedState = theUpdator->update(propagatedState, *innerhit);
611 
612  TSOS propagatedState_out = thePropagator->propagate(updatedState,outerhit->det()->surface()) ;
613  if (!propagatedState_out.isValid())
614  return false;
615  TSOS updatedState_out = theUpdator->update(propagatedState_out, *outerhit);
616  // debug prints
617  LogDebug("") <<"[ElectronSeedGenerator::prepareElTrackSeed] final TSOS, position: "<<updatedState_out.globalPosition()<<" momentum: "<<updatedState_out.globalMomentum();
618  LogDebug("") <<"[ElectronSeedGenerator::prepareElTrackSeed] final TSOS Pt: "<<updatedState_out.globalMomentum().perp();
619  pts_ = trajectoryStateTransform::persistentState(updatedState_out, outerhit->geographicalId().rawId());
620 
621  return true;
622  }
#define LogDebug(id)
PropagationDirection direction() const
void setCaloCluster(const CaloClusterRef &, unsigned char hitsMask=0, int subDet2=0, int subDet1=0, float hoe1=std::numeric_limits< float >::infinity(), float hoe2=std::numeric_limits< float >::infinity())
Definition: ElectronSeed.cc:71
T getParameter(std::string const &) const
unsigned long long cacheIdentifier() const
unsigned long long cacheIDNavSchool_
int i
Definition: DBlmapReader.cc:9
static const TGPicture * info(bool iBackgroundIsBlack)
T perp() const
Definition: PV3DBase.h:72
PixelHitMatcher * myMatchEle
void setupES(const edm::EventSetup &setup)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:462
LocalPoint position() const
Local x and y position coordinates.
TransientTrackingRecHit::ConstRecHitPointer ConstRecHitPointer
void setES(const MagneticField *, const MeasurementTracker *theMeasurementTracker, const TrackerGeometry *trackerGeometry)
const NavigationSchool * theNavigationSchool
float dRz1() const
PTrajectoryStateOnDet persistentState(const TrajectoryStateOnSurface &ts, unsigned int detid)
float dRz2() const
GlobalPoint globalPosition() const
unsigned long long cacheIDMagField_
PropagationDirection
void run(edm::Event &, const edm::EventSetup &setup, const reco::SuperClusterRefVector &, const std::vector< float > &hoe1s, const std::vector< float > &hoe2s, TrajectorySeedCollection *seeds, reco::ElectronSeedCollection &)
unsigned long long cacheIDTrkGeom_
tuple s2
Definition: indexGen.py:106
CaloClusterRef caloCluster() const
Definition: ElectronSeed.h:70
float dPhi2() const
Definition: ElectronSeed.h:74
float dRz2Pos() const
Definition: ElectronSeed.h:75
const SurfaceType & surface() const
edm::ESHandle< TrackerGeometry > theTrackerGeometry
void seedsFromThisCluster(edm::Ref< reco::SuperClusterCollection > seedCluster, float hoe1, float hoe2, reco::ElectronSeedCollection &out, const TrackerTopology *tTopo)
std::vector< TrajectorySeed > TrajectorySeedCollection
recHitContainer::const_iterator const_iterator
T sqrt(T t)
Definition: SSEVec.h:18
const double infinity
std::pair< const_iterator, const_iterator > range
void setNegAttributes(float dRz2=std::numeric_limits< float >::infinity(), float dPhi2=std::numeric_limits< float >::infinity(), float dRz1=std::numeric_limits< float >::infinity(), float dPhi1=std::numeric_limits< float >::infinity())
Definition: ElectronSeed.cc:86
unsigned int detId() const
PropagatorWithMaterial * thePropagator
LocalVector momentum() const
Momentum vector in the local frame.
std::vector< ElectronSeed > ElectronSeedCollection
collection of ElectronSeed objects
PixelHitMatcher * myMatchPos
float dRz2() const
Definition: ElectronSeed.h:73
void seedsFromRecHits(std::vector< std::pair< RecHitWithDist, ConstRecHitPointer > > &elePixelHits, PropagationDirection &dir, const GlobalPoint &vertexPos, const reco::ElectronSeed::CaloClusterRef &cluster, reco::ElectronSeedCollection &out, bool positron)
float dPhi1() const
float dPhi2() const
bool prepareElTrackSeed(ConstRecHitPointer outerhit, ConstRecHitPointer innerhit, const GlobalPoint &vertexPos)
TrajectoryStateOnSurface TSOS
Definition: TestHits.cc:19
PTrajectoryStateOnDet const & startingState() const
float dPhi2Pos() const
Definition: ElectronSeed.h:76
void addSeed(reco::ElectronSeed &seed, const SeedWithInfo *info, bool positron, reco::ElectronSeedCollection &out)
const T & get() const
Definition: EventSetup.h:56
T const * product() const
Definition: ESHandle.h:86
range recHits() const
ElectronSeedGenerator(const edm::ParameterSet &, const Tokens &)
void ele_convert(const Type1 &obj1, Type2 &obj2)
TrackCharge charge() const
Charge (-1, 0 or 1)
edm::ESHandle< MagneticField > theMagField
return(e1-e2)*(e1-e2)+dp *dp
edm::EventID id() const
Definition: EventBase.h:59
bool equivalent(const TrajectorySeed &s1, const TrajectorySeed &s2)
unsigned int nHits() const
GlobalVector globalMomentum() const
void display_seed(const std::string &title1, const std::string &title2, const reco::ElectronSeed &seed, edm::ESHandle< TrackerGeometry > trackerGeometry)
const MeasurementTracker * theMeasurementTracker
size_type size() const
Size of the RefVector.
Definition: RefVector.h:107
int subDet2() const
Definition: ElectronSeed.h:72
tuple cout
Definition: gather_cfg.py:145
dbl *** dir
Definition: mlp_gen.cc:35
volatile std::atomic< bool > shutdown_flag false
void setPosAttributes(float dRz2=std::numeric_limits< float >::infinity(), float dPhi2=std::numeric_limits< float >::infinity(), float dRz1=std::numeric_limits< float >::infinity(), float dPhi1=std::numeric_limits< float >::infinity())
Definition: ElectronSeed.cc:95
void seedsFromTrajectorySeeds(const std::vector< SeedWithInfo > &elePixelSeeds, const reco::ElectronSeed::CaloClusterRef &cluster, float hoe1, float hoe2, reco::ElectronSeedCollection &out, bool positron)
const LocalTrajectoryParameters & parameters() const
Our base class.
Definition: SiPixelRecHit.h:23
unsigned char hitsMask() const
Definition: ElectronSeed.h:71
unsigned long long cacheIDCkfComp_