CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
KfTrackProducerBase.cc
Go to the documentation of this file.
2 // system include files
3 #include <memory>
4 // user include files
11 
13 
15 
17 #include "TrajectoryToResiduals.h"
18 
20  const Propagator* prop,
21  const MeasurementTracker* measTk,
22  std::auto_ptr<TrackingRecHitCollection>& selHits,
23  std::auto_ptr<reco::TrackCollection>& selTracks,
24  std::auto_ptr<reco::TrackExtraCollection>& selTrackExtras,
25  std::auto_ptr<std::vector<Trajectory> >& selTrajectories,
26  AlgoProductCollection& algoResults)
27 {
28 
31 
35  edm::Ref< std::vector<Trajectory> >::key_type iTjRef = 0;
36  std::map<unsigned int, unsigned int> tjTkMap;
37 
38  for(AlgoProductCollection::iterator i=algoResults.begin(); i!=algoResults.end();i++){
39  Trajectory * theTraj = (*i).first;
40  if(trajectoryInEvent_) {
41  selTrajectories->push_back(*theTraj);
42  iTjRef++;
43  }
44 
45  // const TrajectoryFitter::RecHitContainer& transHits = theTraj->recHits(useSplitting); // NO: the return type in Trajectory is by VALUE
47 
48  reco::Track * theTrack = (*i).second.first;
49 
50  // Hits are going to be re-sorted along momentum few lines later.
51  // Therefore the direction stored in the TrackExtra
52  // has to be "alongMomentum" as well. Anyway, this direction can be differnt from the one of the orignal
53  // seed! The name seedDirection() for the Track's method (and the corresponding data member) is
54  // misleading and should be changed into something like "hitsDirection()". TO BE FIXED!
56 
57  LogDebug("TrackProducer") << "In KfTrackProducerBase::putInEvt - seedDir=" << seedDir;
58 
59  reco::Track t = * theTrack;
60  selTracks->push_back( t );
61  iTkRef++;
62 
63  // Store indices in local map (starts at 0)
64  if(trajectoryInEvent_) tjTkMap[iTjRef-1] = iTkRef-1;
65 
66  //sets the outermost and innermost TSOSs
67  TrajectoryStateOnSurface outertsos;
68  TrajectoryStateOnSurface innertsos;
69  unsigned int innerId, outerId;
70 
71  // --- NOTA BENE: the convention is to sort hits and measurements "along the momentum".
72  // This is consistent with innermost and outermost labels only for tracks from LHC collision
73  if (theTraj->direction() == alongMomentum) {
74  outertsos = theTraj->lastMeasurement().updatedState();
75  innertsos = theTraj->firstMeasurement().updatedState();
76  outerId = theTraj->lastMeasurement().recHit()->geographicalId().rawId();
77  innerId = theTraj->firstMeasurement().recHit()->geographicalId().rawId();
78  } else {
79  outertsos = theTraj->firstMeasurement().updatedState();
80  innertsos = theTraj->lastMeasurement().updatedState();
81  outerId = theTraj->firstMeasurement().recHit()->geographicalId().rawId();
82  innerId = theTraj->lastMeasurement().recHit()->geographicalId().rawId();
83  }
84  // ---
85  //build the TrackExtra
86  GlobalPoint v = outertsos.globalParameters().position();
87  GlobalVector p = outertsos.globalParameters().momentum();
88  math::XYZVector outmom( p.x(), p.y(), p.z() );
89  math::XYZPoint outpos( v.x(), v.y(), v.z() );
90  v = innertsos.globalParameters().position();
91  p = innertsos.globalParameters().momentum();
92  math::XYZVector inmom( p.x(), p.y(), p.z() );
93  math::XYZPoint inpos( v.x(), v.y(), v.z() );
94 
95  reco::TrackExtraRef teref= reco::TrackExtraRef ( rTrackExtras, idx ++ );
96  reco::Track & track = selTracks->back();
97  track.setExtra( teref );
98 
99  //======= I want to set the second hitPattern here =============
100  if (theSchool.isValid())
101  {
102  NavigationSetter setter( *theSchool );
103  setSecondHitPattern(theTraj,track,prop,measTk);
104  }
105  //==============================================================
106 
107  selTrackExtras->push_back( reco::TrackExtra (outpos, outmom, true, inpos, inmom, true,
108  outertsos.curvilinearError(), outerId,
109  innertsos.curvilinearError(), innerId,
110  seedDir, theTraj->seedRef()));
111 
112 
113  reco::TrackExtra & tx = selTrackExtras->back();
114 
115  // --- NOTA BENE: the convention is to sort hits and measurements "along the momentum".
116  // This is consistent with innermost and outermost labels only for tracks from LHC collisions
117  size_t ih = 0;
118  if (theTraj->direction() == alongMomentum) {
119  for( TrajectoryFitter::RecHitContainer::const_iterator j = transHits.begin();
120  j != transHits.end(); j ++ ) {
121  if ((**j).hit()!=0){
122  TrackingRecHit * hit = (**j).hit()->clone();
123  track.setHitPattern( * hit, ih ++ );
124  selHits->push_back( hit );
125  tx.add( TrackingRecHitRef( rHits, hidx ++ ) );
126  }
127  }
128  }else{
129  for( TrajectoryFitter::RecHitContainer::const_iterator j = transHits.end()-1;
130  j != transHits.begin()-1; --j ) {
131  if ((**j).hit()!=0){
132  TrackingRecHit * hit = (**j).hit()->clone();
133  track.setHitPattern( * hit, ih ++ );
134  selHits->push_back( hit );
135  tx.add( TrackingRecHitRef( rHits, hidx ++ ) );
136  }
137  }
138  }
139  // ----
140  tx.setResiduals(trajectoryToResiduals(*theTraj));
141 
142  delete theTrack;
143  delete theTraj;
144  }
145 
146  // Now we can re-set refs to hits, as they have already been cloned
147  if (rekeyClusterRefs_) {
149  for (TrackingRecHitCollection::iterator it = selHits->begin(), ed = selHits->end(); it != ed; ++it) {
150  refSetter.reKey(&*it);
151  }
152  }
153 
154  LogTrace("TrackingRegressionTest") << "========== TrackProducer Info ===================";
155  LogTrace("TrackingRegressionTest") << "number of finalTracks: " << selTracks->size();
156  for (reco::TrackCollection::const_iterator it = selTracks->begin(); it != selTracks->end(); it++) {
157  LogTrace("TrackingRegressionTest") << "track's n valid and invalid hit, chi2, pt, eta : "
158  << it->found() << " , "
159  << it->lost() <<" , "
160  << it->normalizedChi2() << " , "
161  << it->pt() << " , "
162  << it->eta() ;
163  }
164  LogTrace("TrackingRegressionTest") << "=================================================";
165 
166 
167  rTracks_ = evt.put( selTracks );
168  evt.put( selTrackExtras );
169  evt.put( selHits );
170 
171  if(trajectoryInEvent_) {
172  edm::OrphanHandle<std::vector<Trajectory> > rTrajs = evt.put(selTrajectories);
173 
174  // Now Create traj<->tracks association map
175  std::auto_ptr<TrajTrackAssociationCollection> trajTrackMap( new TrajTrackAssociationCollection() );
176  for ( std::map<unsigned int, unsigned int>::iterator i = tjTkMap.begin();
177  i != tjTkMap.end(); i++ ) {
178  edm::Ref<std::vector<Trajectory> > trajRef( rTrajs, (*i).first );
179  edm::Ref<reco::TrackCollection> tkRef( rTracks_, (*i).second );
180  trajTrackMap->insert( edm::Ref<std::vector<Trajectory> >( rTrajs, (*i).first ),
181  edm::Ref<reco::TrackCollection>( rTracks_, (*i).second ) );
182  }
183  evt.put( trajTrackMap );
184  }
185 }
186 
#define LogDebug(id)
int i
Definition: DBlmapReader.cc:9
ConstRecHitPointer const & recHit() const
void setSecondHitPattern(Trajectory *traj, reco::Track &track, const Propagator *prop, const MeasurementTracker *measTk)
const CurvilinearTrajectoryError & curvilinearError() const
void reKey(TrackingRecHit *hit) const
T y() const
Definition: PV3DBase.h:63
edm::Ref< TrackExtraCollection > TrackExtraRef
persistent reference to a TrackExtra
Definition: TrackExtraFwd.h:13
PropagationDirection
ConstRecHitContainer recHits(bool splitting=false) const
Definition: Trajectory.cc:67
PropagationDirection const & direction() const
Definition: Trajectory.cc:196
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:94
edm::Ref< TrackingRecHitCollection > TrackingRecHitRef
persistent reference to a TrackingRecHit
TrajectoryMeasurement const & lastMeasurement() const
Definition: Trajectory.h:193
T z() const
Definition: PV3DBase.h:64
int j
Definition: DBlmapReader.cc:9
edm::RefToBase< TrajectorySeed > seedRef(void) const
Definition: Trajectory.h:308
std::vector< AlgoProduct > AlgoProductCollection
edm::AssociationMap< edm::OneToOne< std::vector< Trajectory >, reco::TrackCollection, unsigned short > > TrajTrackAssociationCollection
#define LogTrace(id)
RefProd< PROD > getRefBeforePut()
Definition: Event.h:106
void setHitPattern(const C &c)
set hit patterns from vector of hit references
Definition: TrackBase.h:246
std::vector< TrackExtra > TrackExtraCollection
collection of TrackExtra objects
Definition: TrackExtraFwd.h:9
void setExtra(const TrackExtraRef &ref)
set reference to &quot;extra&quot; object
Definition: Track.h:95
TrajectoryMeasurement const & firstMeasurement() const
Definition: Trajectory.h:206
const GlobalTrajectoryParameters & globalParameters() const
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:13
tuple idx
DEBUGGING if hasattr(process,&quot;trackMonIterativeTracking2012&quot;): print &quot;trackMonIterativeTracking2012 D...
void add(const TrackingRecHitRef &r)
add a reference to a RecHit
virtual void putInEvt(edm::Event &, const Propagator *prop, const MeasurementTracker *measTk, std::auto_ptr< TrackingRecHitCollection > &, std::auto_ptr< reco::TrackCollection > &, std::auto_ptr< reco::TrackExtraCollection > &, std::auto_ptr< std::vector< Trajectory > > &, AlgoProductCollection &)
Put produced collections in the event.
edm::ESHandle< NavigationSchool > theSchool
TrajectoryStateOnSurface const & updatedState() const
reco::TrackResiduals trajectoryToResiduals(const Trajectory &trajectory, enum reco::TrackResiduals::ResidualType type)
boost::remove_cv< typename boost::remove_reference< argument_type >::type >::type key_type
Definition: Ref.h:170
bool isValid() const
Definition: ESHandle.h:37
edm::OrphanHandle< TrackCollection > rTracks_
T x() const
Definition: PV3DBase.h:62
Trajectory::RecHitContainer RecHitContainer
void setResiduals(const TrackResiduals &r)
set the residuals
Definition: TrackExtra.h:131