16 typename propagator_options_t>
25 auto res = getDistanceAndMomentum(gctx, trkParams, vtxPos, deltaR, momDir);
36 typename propagator_options_t>
46 auto res = getDistanceAndMomentum(gctx, trkParams, vtxPos, deltaR, momDir);
56 Vector3D corrDeltaR = deltaR - (deltaR.dot(momDir)) * momDir;
59 Vector3D perpDir = momDir.cross(corrDeltaR);
63 thePlane.matrix().block(0, 0, 3, 1) = corrDeltaR;
64 thePlane.matrix().block(0, 1, 3, 1) = perpDir;
65 thePlane.matrix().block(0, 2, 3, 1) = momDir;
67 thePlane.matrix().block(0, 3, 3, 1) = vtxPos;
69 std::shared_ptr<PlaneSurface> planeSurface =
70 Surface::makeShared<PlaneSurface>(
71 std::make_shared<Transform3D>(thePlane));
74 propagator_options_t pOptions(gctx, mctx);
78 auto result = m_cfg.propagator->propagate(trkParams, *planeSurface, pOptions);
80 return std::move((*result).endParameters);
82 return result.error();
87 typename propagator_options_t>
93 if (trkParams ==
nullptr) {
94 return VertexingError::EmptyInput;
105 Vector3D xDirPlane = myRotation.col(0);
106 Vector3D yDirPlane = myRotation.col(1);
109 Vector3D vertexLocPlane = vertexPos - myTranslation;
112 Vector2D vertexLocXY{vertexLocPlane.dot(xDirPlane),
113 vertexLocPlane.dot(yDirPlane)};
125 return myXYpos.dot(myWeightXY * myXYpos);
128 template <
typename input_track_t,
typename propagator_t,
129 typename propagator_options_t>
132 propagator_options_t>::performNewtonApproximation(
const Vector3D& trkPos,
136 double sinNewPhi = -std::sin(phi);
137 double cosNewPhi = std::cos(phi);
140 bool hasConverged =
false;
142 double cotTheta = 1. / std::tan(theta);
145 while (!hasConverged && nIter < m_cfg.maxIterations) {
146 double x0 = trkPos.x();
147 double y0 = trkPos.y();
148 double z0 = trkPos.z();
150 double xc = vtxPos.x();
151 double yc = vtxPos.y();
152 double zc = vtxPos.z();
154 double derivative = (x0 - xc) * (-r * cosNewPhi) +
155 (y0 - yc) * r * sinNewPhi +
156 (z0 - zc - r * phi * cotTheta) * (-r * cotTheta);
157 double secDerivative = r * (-(x0 - xc) * sinNewPhi - (y0 - yc) * cosNewPhi +
158 r * cotTheta * cotTheta);
160 if (secDerivative < 0.) {
161 return VertexingError::NumericFailure;
164 double deltaPhi = -derivative / secDerivative;
167 sinNewPhi = -std::sin(phi);
168 cosNewPhi = std::cos(phi);
172 if (
std::abs(deltaPhi) < m_cfg.precision) {
179 return VertexingError::NotConverged;
184 template <
typename input_track_t,
typename propagator_t,
185 typename propagator_options_t>
200 double cotTheta = 1. / std::tan(theta);
203 double bZ = m_cfg.bField.getField(trkSurfaceCenter)[
eZ];
208 if (bZ == 0. ||
std::abs(qOvP) < m_cfg.minQoP) {
212 r = std::sin(theta) * (1. / qOvP) / bZ;
216 (d0 - r) * std::cos(phi),
217 z0 + r * phi * cotTheta);
221 auto res = performNewtonApproximation(vec0, vtxPos, phi, theta, r);
229 momDir =
Vector3D(std::sin(theta) * std::cos(phi),
230 std::sin(phi) * std::sin(theta), std::cos(theta));
234 vec0 + r *
Vector3D(-std::sin(phi), std::cos(phi), -cotTheta * phi);
237 deltaR = pointCA3d - vtxPos;
242 template <
typename input_track_t,
typename propagator_t,
243 typename propagator_options_t>
254 const std::shared_ptr<PerigeeSurface> perigeeSurface =
255 Surface::makeShared<PerigeeSurface>(vtx.
position());
258 propagator_options_t pOptions(gctx, mctx);
262 auto result = m_cfg.propagator->propagate(track, *perigeeSurface, pOptions);
265 return result.error();
268 const auto& propRes = *result;
269 const auto& params = propRes.endParameters->parameters();
278 const auto& perigeeCov = *(propRes.endParameters->covariance());
284 newIPandSigma.
IPd0 = d0;
285 double d0_PVcontrib = d0JacXY.transpose() * (vrtXYCov * d0JacXY);
286 if (d0_PVcontrib >= 0) {
287 newIPandSigma.
sigmad0 = std::sqrt(
289 newIPandSigma.
PVsigmad0 = std::sqrt(d0_PVcontrib);
304 ActsVectorD<2> z0JacZ0Theta(std::sin(theta), z0 * std::cos(theta));
309 z0JacZ0Theta.transpose() * (covPerigeeZ0Theta * z0JacZ0Theta) +
310 std::sin(theta) * vtxZZCov * std::sin(theta));
313 std::sqrt(std::sin(theta) * vtxZZCov * std::sin(theta));
317 newIPandSigma.
PVsigmaz0 = std::sqrt(vtxZZCov);
321 double sigma2z0sinTheta =
322 (z0JacZ0Theta.transpose() * (covPerigeeZ0Theta * z0JacZ0Theta));
332 return newIPandSigma;