CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TrackProducerAlgorithm.cc
Go to the documentation of this file.
2 
5 
8 
12 
21 
24 
27 // #include "TrackingTools/MaterialEffects/interface/PropagatorWithMaterial.h"
31 
32 template <> bool
34  const Propagator * thePropagator,
35  AlgoProductCollection& algoResults,
37  TrajectoryStateOnSurface& theTSOS,
38  const TrajectorySeed& seed,
39  float ndof,
40  const reco::BeamSpot& bs,
41  SeedRef seedRef,
42  int qualityMask,signed char nLoops)
43 {
44  //variable declarations
45  std::vector<Trajectory> trajVec;
46  reco::Track * theTrack;
47  Trajectory * theTraj;
48  PropagationDirection seedDir = seed.direction();
49 
50  //perform the fit: the result's size is 1 if it succeded, 0 if fails
51  if(nLoops>0)
52  trajVec = theFitter->fit(seed, hits, theTSOS,TrajectoryFitter::looper);
53  else
54  trajVec = theFitter->fit(seed, hits, theTSOS,TrajectoryFitter::standard);
55 
56  LogDebug("TrackProducer") <<" FITTER FOUND "<< trajVec.size() << " TRAJECTORIES" <<"\n";
57  TrajectoryStateOnSurface innertsos;
58 
59  if (trajVec.size() != 0){
60  theTraj = new Trajectory( trajVec.front() );
61  theTraj->setSeedRef(seedRef);
62 
63  if (theTraj->direction() == alongMomentum) {
64  innertsos = theTraj->firstMeasurement().updatedState();
65  } else {
66  innertsos = theTraj->lastMeasurement().updatedState();
67  }
68 
69  ndof = 0;
71  theTraj->validRecHits(validHits);
72  for (TransientTrackingRecHit::RecHitContainer::iterator h=validHits.begin();h!=validHits.end();++h)
73  ndof = ndof + ((*h)->dimension())*((*h)->weight());
74  if (theTSOS.magneticField()->inTesla(GlobalPoint(0,0,0)).mag2()<DBL_MIN) ndof = ndof - 4;
75  else ndof = ndof - 5;
76 
77  TSCBLBuilderNoMaterial tscblBuilder;
78  // const FreeTrajectoryState & stateForProjectionToBeamLine=*innertsos.freeState();
79  const TrajectoryStateOnSurface & stateForProjectionToBeamLineOnSurface = theTraj->closestMeasurement(GlobalPoint(bs.x0(),bs.y0(),bs.z0())).updatedState();
80  if (!stateForProjectionToBeamLineOnSurface.isValid()){
81  edm::LogError("CannotPropagateToBeamLine")<<"the state on the closest measurement isnot valid. skipping track.";
82  delete theTraj;
83  return false;
84  }
85  const FreeTrajectoryState & stateForProjectionToBeamLine=*stateForProjectionToBeamLineOnSurface.freeState();
86 
87  LogDebug("TrackProducer") << "stateForProjectionToBeamLine=" << stateForProjectionToBeamLine;
88 
89  TrajectoryStateClosestToBeamLine tscbl = tscblBuilder(stateForProjectionToBeamLine,bs);
90 
91  if (tscbl.isValid()==false) {
92  delete theTraj;
93  return false;
94  }
95 
97  math::XYZPoint pos( v.x(), v.y(), v.z() );
99  math::XYZVector mom( p.x(), p.y(), p.z() );
100 
101  LogDebug("TrackProducer") << "pos=" << v << " mom=" << p << " pt=" << p.perp() << " mag=" << p.mag();
102 
103  theTrack = new reco::Track(theTraj->chiSquared(),
104  int(ndof),//FIXME fix weight() in TrackingRecHit
105  pos, mom, tscbl.trackStateAtPCA().charge(),
107  algo_);
108 
109  theTrack->setQualityMask(qualityMask);
110  theTrack->setNLoops(nLoops);
111 
112  LogDebug("TrackProducer") << "theTrack->pt()=" << theTrack->pt();
113 
114  LogDebug("TrackProducer") <<"track done\n";
115 
116  AlgoProduct aProduct(theTraj,std::make_pair(theTrack,seedDir));
117  algoResults.push_back(aProduct);
118 
119  return true;
120  }
121  else return false;
122 }
123 
124 template <> bool
126  const Propagator * thePropagator,
127  AlgoProductCollection& algoResults,
129  TrajectoryStateOnSurface& theTSOS,
130  const TrajectorySeed& seed,
131  float ndof,
132  const reco::BeamSpot& bs,
133  SeedRef seedRef,
134  int qualityMask,signed char nLoops)
135 {
136  //variable declarations
137  std::vector<Trajectory> trajVec;
138  reco::GsfTrack * theTrack;
139  Trajectory * theTraj;
140  PropagationDirection seedDir = seed.direction();
141 
142  //perform the fit: the result's size is 1 if it succeded, 0 if fails
143  if(nLoops>0)
144  trajVec = theFitter->fit(seed, hits, theTSOS,TrajectoryFitter::looper);
145  else
146  trajVec = theFitter->fit(seed, hits, theTSOS,TrajectoryFitter::standard);
147 
148 
149 
150  LogDebug("GsfTrackProducer") <<" FITTER FOUND "<< trajVec.size() << " TRAJECTORIES" <<"\n";
151 
152  TrajectoryStateOnSurface innertsos;
153  TrajectoryStateOnSurface outertsos;
154 
155  if (trajVec.size() != 0){
156 
157  theTraj = new Trajectory( trajVec.front() );
158  theTraj->setSeedRef(seedRef);
159 
160  if (theTraj->direction() == alongMomentum) {
161  innertsos = theTraj->firstMeasurement().updatedState();
162  outertsos = theTraj->lastMeasurement().updatedState();
163  } else {
164  innertsos = theTraj->lastMeasurement().updatedState();
165  outertsos = theTraj->firstMeasurement().updatedState();
166  }
167 // std::cout
168 // << "Nr. of first / last states = "
169 // << innertsos.components().size() << " "
170 // << outertsos.components().size() << std::endl;
171 // std::vector<TrajectoryStateOnSurface> components =
172 // innertsos.components();
173 // double sinTheta =
174 // sin(innertsos.globalMomentum().theta());
175 // for ( std::vector<TrajectoryStateOnSurface>::const_iterator ic=components.begin();
176 // ic!=components.end(); ic++ ) {
177 // std::cout << " comp " << ic-components.begin() << " "
178 // << (*ic).weight() << " "
179 // << (*ic).localParameters().vector()[0]/sinTheta << " "
180 // << sqrt((*ic).localError().matrix()[0][0])/sinTheta << std::endl;
181 // }
182 
183  ndof = 0;
185  theTraj->validRecHits(validHits);
186  for (TransientTrackingRecHit::RecHitContainer::iterator h=validHits.begin();h!=validHits.end();++h)
187  ndof = ndof + ((*h)->dimension())*((*h)->weight());
188  if (theTSOS.magneticField()->inTesla(GlobalPoint(0,0,0)).mag2()<DBL_MIN) ndof = ndof - 4;
189  else ndof = ndof - 5;
190 
191  TSCBLBuilderNoMaterial tscblBuilder;
192  // const FreeTrajectoryState & stateForProjectionToBeamLine=*innertsos.freeState();
193  const TrajectoryStateOnSurface & stateForProjectionToBeamLineOnSurface = theTraj->closestMeasurement(GlobalPoint(bs.x0(),bs.y0(),bs.z0())).updatedState();
194  if (!stateForProjectionToBeamLineOnSurface.isValid()){
195  edm::LogError("CannotPropagateToBeamLine")<<"the state on the closest measurement isnot valid. skipping track.";
196  delete theTraj;
197  return false;
198  } const FreeTrajectoryState & stateForProjectionToBeamLine=*theTraj->closestMeasurement(GlobalPoint(bs.x0(),bs.y0(),bs.z0())).updatedState().freeState();
199 
200  LogDebug("TrackProducer") << "stateForProjectionToBeamLine=" << stateForProjectionToBeamLine;
201 
202  TrajectoryStateClosestToBeamLine tscbl = tscblBuilder(stateForProjectionToBeamLine,bs);
203 
204  if (tscbl.isValid()==false) {
205  delete theTraj;
206  return false;
207  }
208 
210  math::XYZPoint pos( v.x(), v.y(), v.z() );
212  math::XYZVector mom( p.x(), p.y(), p.z() );
213 
214  LogDebug("TrackProducer") << "pos=" << v << " mom=" << p << " pt=" << p.perp() << " mag=" << p.mag();
215 
216  theTrack = new reco::GsfTrack(theTraj->chiSquared(),
217  int(ndof),//FIXME fix weight() in TrackingRecHit
218  // theTraj->foundHits(),//FIXME to be fixed in Trajectory.h
219  // 0, //FIXME no corresponding method in trajectory.h
220  // theTraj->lostHits(),//FIXME to be fixed in Trajectory.h
221  pos, mom, tscbl.trackStateAtPCA().charge(), tscbl.trackStateAtPCA().curvilinearError());
222  theTrack->setAlgorithm(algo_);
223 
224  LogDebug("GsfTrackProducer") <<"track done\n";
225 
226  AlgoProduct aProduct(theTraj,std::make_pair(theTrack,seedDir));
227  LogDebug("GsfTrackProducer") <<"track done1\n";
228  algoResults.push_back(aProduct);
229  LogDebug("GsfTrackProducer") <<"track done2\n";
230 
231  return true;
232  }
233  else return false;
234 }
#define LogDebug(id)
PropagationDirection direction() const
void setQualityMask(int qualMask)
Definition: TrackBase.h:287
double z0() const
z coordinate
Definition: BeamSpot.h:69
T mag2() const
Definition: PV3DBase.h:65
T perp() const
Definition: PV3DBase.h:71
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
T y() const
Definition: PV3DBase.h:62
PropagationDirection
TrackCharge charge() const
std::vector< AlgoProduct > AlgoProductCollection
TrajectoryMeasurement const & closestMeasurement(GlobalPoint) const
Definition: Trajectory.cc:269
const MagneticField * magneticField() const
const CurvilinearTrajectoryError & curvilinearError() const
std::vector< ConstRecHitPointer > RecHitContainer
bool buildTrack(const TrajectoryFitter *, const Propagator *, AlgoProductCollection &, TransientTrackingRecHit::RecHitContainer &, TrajectoryStateOnSurface &, const TrajectorySeed &, float, const reco::BeamSpot &, SeedRef seedRef=SeedRef(), int qualityMask=0, signed char nLoops=0)
Construct Tracks to be put in the event.
PropagationDirection const & direction() const
Definition: Trajectory.cc:195
void setNLoops(signed char value)
Definition: TrackBase.h:290
std::pair< Trajectory *, std::pair< reco::Track *, PropagationDirection > > AlgoProduct
FreeTrajectoryState * freeState(bool withErrors=true) const
double pt() const
track transverse momentum
Definition: TrackBase.h:131
TrajectoryMeasurement const & lastMeasurement() const
Definition: Trajectory.h:181
T z() const
Definition: PV3DBase.h:63
void validRecHits(ConstRecHitContainer &cont) const
Definition: Trajectory.cc:187
TrajectoryStateOnSurface updatedState() const
GlobalVector momentum() const
GlobalPoint position() const
TrajectoryMeasurement const & firstMeasurement() const
Definition: Trajectory.h:194
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
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
reco::TrackBase::TrackAlgorithm algo_
void setSeedRef(const edm::RefToBase< TrajectorySeed > &seedRef)
Definition: Trajectory.h:298
virtual std::vector< Trajectory > fit(const Trajectory &) const =0
void setAlgorithm(const TrackAlgorithm a, bool set=true)
position index
Definition: TrackBase.h:273
double y0() const
y coordinate
Definition: BeamSpot.h:67
T x() const
Definition: PV3DBase.h:61
mathSSE::Vec4< T > v
double chiSquared() const
Definition: Trajectory.h:242
double x0() const
x coordinate
Definition: BeamSpot.h:65