CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
V0Fitter.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: V0Producer
4 // Class: V0Fitter
5 //
13 //
14 // Original Author: Brian Drell
15 // Created: Fri May 18 22:57:40 CEST 2007
16 //
17 //
18 
19 #include "V0Fitter.h"
20 
28 #include <Math/Functions.h>
29 #include <Math/SVector.h>
30 #include <Math/SMatrix.h>
31 #include <typeinfo>
32 #include <memory>
35 
36 // pdg mass constants
37 namespace {
38  const double piMass = 0.13957018;
39  const double piMassSquared = piMass*piMass;
40  const double protonMass = 0.938272046;
41  const double protonMassSquared = protonMass*protonMass;
42  const double kShortMass = 0.497614;
43  const double lambdaMass = 1.115683;
44 }
45 
46 typedef ROOT::Math::SMatrix<double, 3, 3, ROOT::Math::MatRepSym<double, 3> > SMatrixSym3D;
47 typedef ROOT::Math::SVector<double, 3> SVector3;
48 
50 {
51  token_beamSpot = iC.consumes<reco::BeamSpot>(theParameters.getParameter<edm::InputTag>("beamSpot"));
52  useVertex_ = theParameters.getParameter<bool>("useVertex");
53  token_vertices = iC.consumes<std::vector<reco::Vertex>>(theParameters.getParameter<edm::InputTag>("vertices"));
54 
55  token_tracks = iC.consumes<reco::TrackCollection>(theParameters.getParameter<edm::InputTag>("trackRecoAlgorithm"));
56  vertexFitter_ = theParameters.getParameter<bool>("vertexFitter");
57  useRefTracks_ = theParameters.getParameter<bool>("useRefTracks");
58 
59  // whether to reconstruct KShorts
60  doKShorts_ = theParameters.getParameter<bool>("doKShorts");
61  // whether to reconstruct Lambdas
62  doLambdas_ = theParameters.getParameter<bool>("doLambdas");
63 
64  // cuts on initial track selection
65  tkChi2Cut_ = theParameters.getParameter<double>("tkChi2Cut");
66  tkNHitsCut_ = theParameters.getParameter<int>("tkNHitsCut");
67  tkPtCut_ = theParameters.getParameter<double>("tkPtCut");
68  tkIPSigXYCut_ = theParameters.getParameter<double>("tkIPSigXYCut");
69  tkIPSigZCut_ = theParameters.getParameter<double>("tkIPSigZCut");
70 
71  // cuts on vertex
72  vtxChi2Cut_ = theParameters.getParameter<double>("vtxChi2Cut");
73  vtxDecaySigXYZCut_ = theParameters.getParameter<double>("vtxDecaySigXYZCut");
74  vtxDecaySigXYCut_ = theParameters.getParameter<double>("vtxDecaySigXYCut");
75  // miscellaneous cuts
76  tkDCACut_ = theParameters.getParameter<double>("tkDCACut");
77  mPiPiCut_ = theParameters.getParameter<double>("mPiPiCut");
78  innerHitPosCut_ = theParameters.getParameter<double>("innerHitPosCut");
79  cosThetaXYCut_ = theParameters.getParameter<double>("cosThetaXYCut");
80  cosThetaXYZCut_ = theParameters.getParameter<double>("cosThetaXYZCut");
81  // cuts on the V0 candidate mass
82  kShortMassCut_ = theParameters.getParameter<double>("kShortMassCut");
83  lambdaMassCut_ = theParameters.getParameter<double>("lambdaMassCut");
84 }
85 
86 // method containing the algorithm for vertex reconstruction
89 {
90  using std::vector;
91 
93  iEvent.getByToken(token_tracks, theTrackHandle);
94  if (!theTrackHandle->size()) return;
95  const reco::TrackCollection* theTrackCollection = theTrackHandle.product();
96 
97  edm::Handle<reco::BeamSpot> theBeamSpotHandle;
98  iEvent.getByToken(token_beamSpot, theBeamSpotHandle);
99  const reco::BeamSpot* theBeamSpot = theBeamSpotHandle.product();
100  math::XYZPoint referencePos(theBeamSpot->position());
101 
102  if (useVertex_) {
104  iEvent.getByToken(token_vertices, vertices);
105  const reco::Vertex & vertex = (*vertices)[0];
106  referencePos = vertex.position();
107  }
108 
109  edm::ESHandle<MagneticField> theMagneticFieldHandle;
110  iSetup.get<IdealMagneticFieldRecord>().get(theMagneticFieldHandle);
111  const MagneticField* theMagneticField = theMagneticFieldHandle.product();
112 
113  std::vector<reco::TrackRef> theTrackRefs;
114  std::vector<reco::TransientTrack> theTransTracks;
115 
116  // fill vectors of TransientTracks and TrackRefs after applying preselection cuts
117  for (reco::TrackCollection::const_iterator iTk = theTrackCollection->begin(); iTk != theTrackCollection->end(); ++iTk) {
118  const reco::Track* tmpTrack = &(*iTk);
119  double ipsigXY;
120  if (useVertex_) {
121  ipsigXY = std::abs(tmpTrack->dxy(referencePos)/tmpTrack->dxyError());
122  } else {
123  ipsigXY = std::abs(tmpTrack->dxy(*theBeamSpot)/tmpTrack->dxyError());
124  }
125  double ipsigZ = std::abs(tmpTrack->dz(referencePos)/tmpTrack->dzError());
126  if (tmpTrack->normalizedChi2() < tkChi2Cut_ && tmpTrack->numberOfValidHits() >= tkNHitsCut_ &&
127  tmpTrack->pt() > tkPtCut_ && ipsigXY > tkIPSigXYCut_ && ipsigZ > tkIPSigZCut_) {
128  reco::TrackRef tmpRef(theTrackHandle, std::distance(theTrackCollection->begin(), iTk));
129  theTrackRefs.push_back(std::move(tmpRef));
130  reco::TransientTrack tmpTransient(*tmpRef, theMagneticField);
131  theTransTracks.push_back(std::move(tmpTransient));
132  }
133  }
134  // good tracks have now been selected for vertexing
135 
136  // loop over tracks and vertex good charged track pairs
137  for (unsigned int trdx1 = 0; trdx1 < theTrackRefs.size(); ++trdx1) {
138  for (unsigned int trdx2 = trdx1 + 1; trdx2 < theTrackRefs.size(); ++trdx2) {
139 
140  reco::TrackRef positiveTrackRef;
141  reco::TrackRef negativeTrackRef;
142  reco::TransientTrack* posTransTkPtr = nullptr;
143  reco::TransientTrack* negTransTkPtr = nullptr;
144 
145  if (theTrackRefs[trdx1]->charge() < 0. && theTrackRefs[trdx2]->charge() > 0.) {
146  negativeTrackRef = theTrackRefs[trdx1];
147  positiveTrackRef = theTrackRefs[trdx2];
148  negTransTkPtr = &theTransTracks[trdx1];
149  posTransTkPtr = &theTransTracks[trdx2];
150  } else if (theTrackRefs[trdx1]->charge() > 0. && theTrackRefs[trdx2]->charge() < 0.) {
151  negativeTrackRef = theTrackRefs[trdx2];
152  positiveTrackRef = theTrackRefs[trdx1];
153  negTransTkPtr = &theTransTracks[trdx2];
154  posTransTkPtr = &theTransTracks[trdx1];
155  } else {
156  continue;
157  }
158 
159  // measure distance between tracks at their closest approach
160  if (!posTransTkPtr->impactPointTSCP().isValid() || !negTransTkPtr->impactPointTSCP().isValid()) continue;
161  FreeTrajectoryState const & posState = posTransTkPtr->impactPointTSCP().theState();
162  FreeTrajectoryState const & negState = negTransTkPtr->impactPointTSCP().theState();
164  cApp.calculate(posState, negState);
165  if (!cApp.status()) continue;
166  float dca = std::abs(cApp.distance());
167  if (dca > tkDCACut_) continue;
168 
169  // the POCA should at least be in the sensitive volume
170  GlobalPoint cxPt = cApp.crossingPoint();
171  if (sqrt(cxPt.x()*cxPt.x() + cxPt.y()*cxPt.y()) > 120. || std::abs(cxPt.z()) > 300.) continue;
172 
173  // the tracks should at least point in the same quadrant
174  TrajectoryStateClosestToPoint const & posTSCP = posTransTkPtr->trajectoryStateClosestToPoint(cxPt);
175  TrajectoryStateClosestToPoint const & negTSCP = negTransTkPtr->trajectoryStateClosestToPoint(cxPt);
176  if (!posTSCP.isValid() || !negTSCP.isValid()) continue;
177  if (posTSCP.momentum().dot(negTSCP.momentum()) < 0) continue;
178 
179  // calculate mPiPi
180  double totalE = sqrt(posTSCP.momentum().mag2() + piMassSquared) + sqrt(negTSCP.momentum().mag2() + piMassSquared);
181  double totalESq = totalE*totalE;
182  double totalPSq = (posTSCP.momentum() + negTSCP.momentum()).mag2();
183  double mass = sqrt(totalESq - totalPSq);
184  if (mass > mPiPiCut_) continue;
185 
186  // Fill the vector of TransientTracks to send to KVF
187  std::vector<reco::TransientTrack> transTracks;
188  transTracks.reserve(2);
189  transTracks.push_back(*posTransTkPtr);
190  transTracks.push_back(*negTransTkPtr);
191 
192  // create the vertex fitter object and vertex the tracks
193  TransientVertex theRecoVertex;
194  if (vertexFitter_) {
195  KalmanVertexFitter theKalmanFitter(useRefTracks_ == 0 ? false : true);
196  theRecoVertex = theKalmanFitter.vertex(transTracks);
197  } else if (!vertexFitter_) {
198  useRefTracks_ = false;
199  AdaptiveVertexFitter theAdaptiveFitter;
200  theRecoVertex = theAdaptiveFitter.vertex(transTracks);
201  }
202  if (!theRecoVertex.isValid()) continue;
203 
204  reco::Vertex theVtx = theRecoVertex;
205  if (theVtx.normalizedChi2() > vtxChi2Cut_) continue;
206  GlobalPoint vtxPos(theVtx.x(), theVtx.y(), theVtx.z());
207 
208  // 2D decay significance
209  const SMatrixSym3D totalCov = theBeamSpot->rotatedCovariance3D() + theVtx.covariance();
210  SVector3 distVecXY(vtxPos.x()-referencePos.x(), vtxPos.y()-referencePos.y(), 0.);
211  double distMagXY = ROOT::Math::Mag(distVecXY);
212  double sigmaDistMagXY = sqrt(ROOT::Math::Similarity(totalCov, distVecXY)) / distMagXY;
213  if (distMagXY/sigmaDistMagXY < vtxDecaySigXYCut_) continue;
214 
215  // 3D decay significance
216  SVector3 distVecXYZ(vtxPos.x()-referencePos.x(), vtxPos.y()-referencePos.y(), vtxPos.z()-referencePos.z());
217  double distMagXYZ = ROOT::Math::Mag(distVecXYZ);
218  double sigmaDistMagXYZ = sqrt(ROOT::Math::Similarity(totalCov, distVecXYZ)) / distMagXYZ;
219  if (distMagXYZ/sigmaDistMagXYZ < vtxDecaySigXYZCut_) continue;
220 
221  // make sure the vertex radius is within the inner track hit radius
222  if (innerHitPosCut_ > 0. && positiveTrackRef->innerOk()) {
223  reco::Vertex::Point posTkHitPos = positiveTrackRef->innerPosition();
224  double posTkHitPosD2 = (posTkHitPos.x()-referencePos.x())*(posTkHitPos.x()-referencePos.x()) +
225  (posTkHitPos.y()-referencePos.y())*(posTkHitPos.y()-referencePos.y());
226  if (sqrt(posTkHitPosD2) < (distMagXY - sigmaDistMagXY*innerHitPosCut_)) continue;
227  }
228  if (innerHitPosCut_ > 0. && negativeTrackRef->innerOk()) {
229  reco::Vertex::Point negTkHitPos = negativeTrackRef->innerPosition();
230  double negTkHitPosD2 = (negTkHitPos.x()-referencePos.x())*(negTkHitPos.x()-referencePos.x()) +
231  (negTkHitPos.y()-referencePos.y())*(negTkHitPos.y()-referencePos.y());
232  if (sqrt(negTkHitPosD2) < (distMagXY - sigmaDistMagXY*innerHitPosCut_)) continue;
233  }
234 
235  std::auto_ptr<TrajectoryStateClosestToPoint> trajPlus;
236  std::auto_ptr<TrajectoryStateClosestToPoint> trajMins;
237  std::vector<reco::TransientTrack> theRefTracks;
238  if (theRecoVertex.hasRefittedTracks()) {
239  theRefTracks = theRecoVertex.refittedTracks();
240  }
241 
242  if (useRefTracks_ && theRefTracks.size() > 1) {
243  reco::TransientTrack* thePositiveRefTrack = 0;
244  reco::TransientTrack* theNegativeRefTrack = 0;
245  for (std::vector<reco::TransientTrack>::iterator iTrack = theRefTracks.begin(); iTrack != theRefTracks.end(); ++iTrack) {
246  if (iTrack->track().charge() > 0.) {
247  thePositiveRefTrack = &*iTrack;
248  } else if (iTrack->track().charge() < 0.) {
249  theNegativeRefTrack = &*iTrack;
250  }
251  }
252  if (thePositiveRefTrack == 0 || theNegativeRefTrack == 0) continue;
253  trajPlus.reset(new TrajectoryStateClosestToPoint(thePositiveRefTrack->trajectoryStateClosestToPoint(vtxPos)));
254  trajMins.reset(new TrajectoryStateClosestToPoint(theNegativeRefTrack->trajectoryStateClosestToPoint(vtxPos)));
255  } else {
256  trajPlus.reset(new TrajectoryStateClosestToPoint(posTransTkPtr->trajectoryStateClosestToPoint(vtxPos)));
257  trajMins.reset(new TrajectoryStateClosestToPoint(negTransTkPtr->trajectoryStateClosestToPoint(vtxPos)));
258  }
259 
260  if (trajPlus.get() == 0 || trajMins.get() == 0 || !trajPlus->isValid() || !trajMins->isValid()) continue;
261 
262  GlobalVector positiveP(trajPlus->momentum());
263  GlobalVector negativeP(trajMins->momentum());
264  GlobalVector totalP(positiveP + negativeP);
265 
266  // 2D pointing angle
267  double dx = theVtx.x()-referencePos.x();
268  double dy = theVtx.y()-referencePos.y();
269  double px = totalP.x();
270  double py = totalP.y();
271  double angleXY = (dx*px+dy*py)/(sqrt(dx*dx+dy*dy)*sqrt(px*px+py*py));
272  if (angleXY < cosThetaXYCut_) continue;
273 
274  // 3D pointing angle
275  double dz = theVtx.z()-referencePos.z();
276  double pz = totalP.z();
277  double angleXYZ = (dx*px+dy*py+dz*pz)/(sqrt(dx*dx+dy*dy+dz*dz)*sqrt(px*px+py*py+pz*pz));
278  if (angleXYZ < cosThetaXYZCut_) continue;
279 
280  // calculate total energy of V0 3 ways: assume it's a kShort, a Lambda, or a LambdaBar.
281  double piPlusE = sqrt(positiveP.mag2() + piMassSquared);
282  double piMinusE = sqrt(negativeP.mag2() + piMassSquared);
283  double protonE = sqrt(positiveP.mag2() + protonMassSquared);
284  double antiProtonE = sqrt(negativeP.mag2() + protonMassSquared);
285  double kShortETot = piPlusE + piMinusE;
286  double lambdaEtot = protonE + piMinusE;
287  double lambdaBarEtot = antiProtonE + piPlusE;
288 
289  // Create momentum 4-vectors for the 3 candidate types
290  const reco::Particle::LorentzVector kShortP4(totalP.x(), totalP.y(), totalP.z(), kShortETot);
291  const reco::Particle::LorentzVector lambdaP4(totalP.x(), totalP.y(), totalP.z(), lambdaEtot);
292  const reco::Particle::LorentzVector lambdaBarP4(totalP.x(), totalP.y(), totalP.z(), lambdaBarEtot);
293 
294  reco::Particle::Point vtx(theVtx.x(), theVtx.y(), theVtx.z());
295  const reco::Vertex::CovarianceMatrix vtxCov(theVtx.covariance());
296  double vtxChi2(theVtx.chi2());
297  double vtxNdof(theVtx.ndof());
298 
299  // Create the VertexCompositeCandidate object that will be stored in the Event
300  reco::VertexCompositeCandidate* theKshort = nullptr;
301  reco::VertexCompositeCandidate* theLambda = nullptr;
302  reco::VertexCompositeCandidate* theLambdaBar = nullptr;
303 
304  if (doKShorts_) {
305  theKshort = new reco::VertexCompositeCandidate(0, kShortP4, vtx, vtxCov, vtxChi2, vtxNdof);
306  }
307  if (doLambdas_) {
308  if (positiveP.mag2() > negativeP.mag2()) {
309  theLambda = new reco::VertexCompositeCandidate(0, lambdaP4, vtx, vtxCov, vtxChi2, vtxNdof);
310  } else {
311  theLambdaBar = new reco::VertexCompositeCandidate(0, lambdaBarP4, vtx, vtxCov, vtxChi2, vtxNdof);
312  }
313  }
314 
315  // Create daughter candidates for the VertexCompositeCandidates
316  reco::RecoChargedCandidate thePiPlusCand(
317  1, reco::Particle::LorentzVector(positiveP.x(), positiveP.y(), positiveP.z(), piPlusE), vtx);
318  thePiPlusCand.setTrack(positiveTrackRef);
319 
320  reco::RecoChargedCandidate thePiMinusCand(
321  -1, reco::Particle::LorentzVector(negativeP.x(), negativeP.y(), negativeP.z(), piMinusE), vtx);
322  thePiMinusCand.setTrack(negativeTrackRef);
323 
324  reco::RecoChargedCandidate theProtonCand(
325  1, reco::Particle::LorentzVector(positiveP.x(), positiveP.y(), positiveP.z(), protonE), vtx);
326  theProtonCand.setTrack(positiveTrackRef);
327 
328  reco::RecoChargedCandidate theAntiProtonCand(
329  -1, reco::Particle::LorentzVector(negativeP.x(), negativeP.y(), negativeP.z(), antiProtonE), vtx);
330  theAntiProtonCand.setTrack(negativeTrackRef);
331 
332  AddFourMomenta addp4;
333  // Store the daughter Candidates in the VertexCompositeCandidates if they pass mass cuts
334  if (doKShorts_) {
335  theKshort->addDaughter(thePiPlusCand);
336  theKshort->addDaughter(thePiMinusCand);
337  theKshort->setPdgId(310);
338  addp4.set(*theKshort);
339  if (theKshort->mass() < kShortMass + kShortMassCut_ && theKshort->mass() > kShortMass - kShortMassCut_) {
340  theKshorts.push_back(std::move(*theKshort));
341  }
342  }
343  if (doLambdas_ && theLambda) {
344  theLambda->addDaughter(theProtonCand);
345  theLambda->addDaughter(thePiMinusCand);
346  theLambda->setPdgId(3122);
347  addp4.set( *theLambda );
348  if (theLambda->mass() < lambdaMass + lambdaMassCut_ && theLambda->mass() > lambdaMass - lambdaMassCut_) {
349  theLambdas.push_back(std::move(*theLambda));
350  }
351  } else if (doLambdas_ && theLambdaBar) {
352  theLambdaBar->addDaughter(theAntiProtonCand);
353  theLambdaBar->addDaughter(thePiPlusCand);
354  theLambdaBar->setPdgId(-3122);
355  addp4.set(*theLambdaBar);
356  if (theLambdaBar->mass() < lambdaMass + lambdaMassCut_ && theLambdaBar->mass() > lambdaMass - lambdaMassCut_) {
357  theLambdas.push_back(std::move(*theLambdaBar));
358  }
359  }
360 
361  delete theKshort;
362  delete theLambda;
363  delete theLambdaBar;
364  theKshort = theLambda = theLambdaBar = nullptr;
365 
366  }
367  }
368 }
369 
bool doLambdas_
Definition: V0Fitter.h:58
double tkChi2Cut_
Definition: V0Fitter.h:61
T getParameter(std::string const &) const
T mag2() const
Definition: PV3DBase.h:66
std::vector< VertexCompositeCandidate > VertexCompositeCandidateCollection
collection of Candidate objects
const double protonMass
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > SMatrixSym3D
Definition: V0Fitter.cc:46
double cosThetaXYZCut_
Definition: V0Fitter.h:75
double lambdaMassCut_
Definition: V0Fitter.h:78
double tkIPSigXYCut_
Definition: V0Fitter.h:64
bool vertexFitter_
Definition: V0Fitter.h:55
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:524
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:457
TrajectoryStateClosestToPoint impactPointTSCP() const
const FreeTrajectoryState & theState() const
double dxyError() const
error on dxy
Definition: TrackBase.h:759
edm::EDGetTokenT< std::vector< reco::Vertex > > token_vertices
Definition: V0Fitter.h:83
double y() const
y coordinate
Definition: Vertex.h:110
virtual CachingVertex< 5 > vertex(const std::vector< reco::TransientTrack > &tracks) const
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:13
T y() const
Definition: PV3DBase.h:63
PreciseFloatType< T, U >::Type dot(const Vector3DBase< U, FrameTag > &v) const
Definition: Vector3DBase.h:107
double covariance(int i, int j) const
(i, j)-th element of error matrix, i, j = 0, ... 2
Definition: Vertex.h:123
bool hasRefittedTracks() const
double vtxDecaySigXYZCut_
Definition: V0Fitter.h:69
const Point & position() const
position
Definition: Vertex.h:106
math::Error< dimension >::type CovarianceMatrix
covariance error matrix (3x3)
Definition: Vertex.h:45
virtual GlobalPoint crossingPoint() const
double tkDCACut_
Definition: V0Fitter.h:71
int tkNHitsCut_
Definition: V0Fitter.h:62
virtual double mass() const
mass
int iEvent
Definition: GenABIO.cc:230
bool useRefTracks_
Definition: V0Fitter.h:56
edm::EDGetTokenT< reco::TrackCollection > token_tracks
Definition: V0Fitter.h:80
T sqrt(T t)
Definition: SSEVec.h:48
double pt() const
track transverse momentum
Definition: TrackBase.h:584
T z() const
Definition: PV3DBase.h:64
T mag2() const
The vector magnitude squared. Equivalent to vec.dot(vec)
virtual CachingVertex< 5 > vertex(const std::vector< reco::TransientTrack > &) const
math::XYZPoint Point
point in the space
Definition: Particle.h:25
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double chi2() const
chi-squares
Definition: Vertex.h:95
double z() const
y coordinate
Definition: Vertex.h:112
void fitAll(const edm::Event &iEvent, const edm::EventSetup &iSetup, reco::VertexCompositeCandidateCollection &k, reco::VertexCompositeCandidateCollection &l)
Definition: V0Fitter.cc:87
bool useVertex_
Definition: V0Fitter.h:82
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
unsigned short numberOfValidHits() const
number of valid hits found
Definition: TrackBase.h:783
double mPiPiCut_
Definition: V0Fitter.h:72
double kShortMassCut_
Definition: V0Fitter.h:77
V0Fitter(const edm::ParameterSet &theParams, edm::ConsumesCollector &&iC)
Definition: V0Fitter.cc:49
double tkIPSigZCut_
Definition: V0Fitter.h:65
void addDaughter(const Candidate &, const std::string &s="")
add a clone of the passed candidate as daughter
double ndof() const
Definition: Vertex.h:102
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:572
double dzError() const
error on dz
Definition: TrackBase.h:777
double cosThetaXYCut_
Definition: V0Fitter.h:74
double x() const
x coordinate
Definition: Vertex.h:108
double innerHitPosCut_
Definition: V0Fitter.h:73
double vtxDecaySigXYCut_
Definition: V0Fitter.h:68
virtual bool calculate(const TrajectoryStateOnSurface &sta, const TrajectoryStateOnSurface &stb)
T const * product() const
Definition: Handle.h:81
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
TrajectoryStateClosestToPoint trajectoryStateClosestToPoint(const GlobalPoint &point) const
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:86
double tkPtCut_
Definition: V0Fitter.h:63
double vtxChi2Cut_
Definition: V0Fitter.h:67
virtual void setPdgId(int pdgId)
const Point & position() const
position
Definition: BeamSpot.h:62
edm::EDGetTokenT< reco::BeamSpot > token_beamSpot
Definition: V0Fitter.h:81
std::vector< reco::TransientTrack > const & refittedTracks() const
void set(reco::Candidate &c) const
set up a candidate
double normalizedChi2() const
chi-squared divided by n.d.o.f.
Definition: Vertex.h:104
Covariance3DMatrix rotatedCovariance3D() const
Definition: BeamSpot.cc:78
ROOT::Math::SVector< double, 3 > SVector3
Definition: V0Fitter.cc:47
void setTrack(const reco::TrackRef &r)
set reference to track
virtual float distance() const
double dxy() const
dxy parameter. (This is the transverse impact parameter w.r.t. to (0,0,0) ONLY if refPoint is close t...
Definition: TrackBase.h:554
T x() const
Definition: PV3DBase.h:62
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:21
bool isValid() const
virtual bool status() const
bool doKShorts_
Definition: V0Fitter.h:57