217 CADEBUG(int32_t ii; printf(
"\nRefitting track\n"));
218 typename internal::refitTrackTypes<S>::propagator prop;
220 float TrackParCovChi2 = 0.f;
221 convertTrack<S, T, typename internal::refitTrackTypes<S>::propagator>(trk,
trkX, prop, &TrackParCovChi2);
224 if constexpr (std::is_same_v<T, GPUTPCGMMergedTrack>) {
227 int32_t
leg = mPtrackHits[
trkX.FirstClusterRef() +
trkX.NClusters() - 1].leg;
228 for (int32_t
i =
trkX.NClusters() - 2;
i > 0;
i--) {
229 if (mPtrackHits[
trkX.FirstClusterRef() +
i].leg !=
leg) {
235 tOffset =
trkX.GetParam().GetTZOffset();
236 }
else if constexpr (std::is_same_v<T, TrackTPC>) {
238 tOffset =
trkX.getTime0();
239 }
else if constexpr (std::is_same_v<T, TrackParCovWithArgs>) {
241 tOffset =
trkX.time0;
243 static_assert(
"Invalid template");
245 if constexpr (std::is_same_v<S, GPUTPCGMTrackParam>) {
246 CADEBUG(printf(
"\t%21sInit Alpha %8.3f , X %8.3f - Y %8.3f, Z %8.3f - QPt %7.2f (%7.2f), SP %5.2f (%5.2f) --- Cov sY %8.3f sZ %8.3f sSP %8.3f sPt %8.3f\n",
"", prop.GetAlpha(), trk.GetX(), trk.Par()[0], trk.Par()[1], trk.Par()[4], prop.GetQPt0(), trk.Par()[2], prop.GetSinPhi0(), sqrtf(trk.Cov()[0]), sqrtf(trk.Cov()[2]), sqrtf(trk.Cov()[5]), sqrtf(trk.Cov()[14])));
249 int32_t direction = outward ? -1 : 1;
253 uint8_t sector = 255,
row = 255;
254 int32_t lastSector = -1, currentSector = -1, currentRow = -1;
255 int16_t clusterState = 0, nextState = 0;
257 float sumInvSqrtCharge = 0.f;
258 int32_t nAvgCharge = 0;
260 for (int32_t
i =
start;
i != stop;
i += cl ? 0 : direction) {
262 float time = 0.f, invCharge = 0.f, invSqrtCharge = 0.f;
267 if constexpr (std::is_same_v<T, GPUTPCGMMergedTrack>) {
268 const auto& hit = mPtrackHits[
trkX.FirstClusterRef() +
i];
269 cl = &mPclusterNative->clustersLinear[hit.num];
272 if (
i + direction != stop) {
280 nextState = mPclusterState[hit.num];
281 }
else if constexpr (std::is_same_v<T, TrackTPC>) {
282 cl = &
trkX.getCluster(mPtrackHitReferences,
i, *mPclusterNative, sector,
row);
283 nextState = mPclusterState[cl - mPclusterNative->clustersLinear];
284 }
else if constexpr (std::is_same_v<T, TrackParCovWithArgs>) {
285 cl = &TrackTPC::getCluster(mPtrackHitReferences,
i, *mPclusterNative, sector,
row,
trkX.clusRef);
286 nextState = mPclusterState[cl - mPclusterNative->clustersLinear];
288 static_assert(
"Invalid template");
291 if (
clusters == 0 || (
row == currentRow && sector == currentSector)) {
298 mPfastTransformHelper->Transform(sector,
row, cl->getPad(), cl->getTime(),
x,
y,
z, tOffset);
299 CADEBUG(printf(
"\tHit %3d/%3d Row %3d: Cluster Alpha %8.3f %3d, X %8.3f - Y %8.3f, Z %8.3f - State %d\n", ii,
count,
row, mPparam->Alpha(sector), (int32_t)sector,
x,
y,
z, (int32_t)nextState));
301 currentSector = sector;
303 clusterState = nextState;
304 time = cl->getTime();
305 invSqrtCharge = CAMath::InvSqrt(cl->
qMax);
306 invCharge = (1.f / cl->
qMax);
309 mPfastTransformHelper->Transform(sector,
row, cl->getPad(), cl->getTime(), xx, yy, zz, tOffset);
310 CADEBUG(printf(
"\tHit %3d/%3d Row %3d: Cluster Alpha %8.3f %3d, X %8.3f - Y %8.3f, Z %8.3f - State %d\n", ii,
count,
row, mPparam->Alpha(sector), (int32_t)sector, xx, yy, zz, (int32_t)nextState));
315 clusterState |= nextState;
319 if (
i + direction != stop) {
332 CADEBUG(printf(
"\tMerged Hit Row %3d: Cluster Alpha %8.3f %3d, X %8.3f - Y %8.3f, Z %8.3f - State %d\n",
row, mPparam->Alpha(sector), (int32_t)sector,
x,
y,
z, (int32_t)clusterState));
335 float invAvgCharge = (sumInvSqrtCharge += invSqrtCharge) / ++nAvgCharge;
336 invAvgCharge *= invAvgCharge;
338 if constexpr (std::is_same_v<S, GPUTPCGMTrackParam>) {
339 if (prop.PropagateToXAlpha(
x, mPparam->Alpha(currentSector), !outward)) {
344 trk.ResetCovariance();
345 }
else if (lastSector != -1 && (lastSector < 18) != (sector < 18)) {
346 if (mPparam->rec.tpc.addErrorsCECrossing) {
347 if (mPparam->rec.tpc.addErrorsCECrossing >= 2) {
348 trk.AddCovDiagErrorsWithCorrelations(mPparam->rec.tpc.errorsCECrossing);
350 trk.AddCovDiagErrors(mPparam->rec.tpc.errorsCECrossing);
352 }
else if (trk.Cov()[2] < 0.5f) {
356 CADEBUG(printf(
"\t%21sPropaga Alpha %8.3f , X %8.3f - Y %8.3f, Z %8.3f - QPt %7.2f (%7.2f), SP %5.2f (%5.2f) --- Res %8.3f %8.3f --- Cov sY %8.3f sZ %8.3f sSP %8.3f sPt %8.3f - YPt %8.3f\n",
"", prop.GetAlpha(),
x, trk.Par()[0], trk.Par()[1], trk.Par()[4], prop.GetQPt0(), trk.Par()[2], prop.GetSinPhi0(), trk.Par()[0] -
y, trk.Par()[1] -
z, sqrtf(trk.Cov()[0]), sqrtf(trk.Cov()[2]), sqrtf(trk.Cov()[5]), sqrtf(trk.Cov()[14]), trk.Cov()[10]));
358 if (prop.Update(
y,
z,
row, *mPparam, clusterState, 0,
nullptr,
true, sector,
time, invAvgCharge, invCharge)) {
362 trk.ConstrainSinPhi();
363 CADEBUG(printf(
"\t%21sFit Alpha %8.3f , X %8.3f - Y %8.3f, Z %8.3f - QPt %7.2f (%7.2f), SP %5.2f (%5.2f), DzDs %5.2f %16s --- Cov sY %8.3f sZ %8.3f sSP %8.3f sPt %8.3f - YPt %8.3f\n",
"", prop.GetAlpha(),
x, trk.Par()[0], trk.Par()[1], trk.Par()[4], prop.GetQPt0(), trk.Par()[2], prop.GetSinPhi0(), trk.Par()[3],
"", sqrtf(trk.Cov()[0]), sqrtf(trk.Cov()[2]), sqrtf(trk.Cov()[5]), sqrtf(trk.Cov()[14]), trk.Cov()[10]));
364 }
else if constexpr (std::is_same_v<S, TrackParCov>) {
365 if (!trk.rotate(mPparam->Alpha(currentSector))) {
373 if (lastSector != -1 && (lastSector < 18) != (sector < 18)) {
374 if (mPparam->rec.tpc.addErrorsCECrossing) {
375 trk.updateCov(mPparam->rec.tpc.errorsCECrossing, mPparam->rec.tpc.addErrorsCECrossing >= 2);
376 }
else if (trk.getCov()[2] < 0.5f) {
381 trk.resetCovariance();
382 float bzkG = prop->getNominalBz(), qptB5Scale = CAMath::Abs(bzkG) > 0.1f ? CAMath::Abs(bzkG) / 5.006680f : 1.f;
383 float q2pt2 = trk.getQ2Pt() * trk.getQ2Pt(), q2pt2Wgh = q2pt2 * qptB5Scale * qptB5Scale;
384 float err2 = (100.f + q2pt2Wgh) / (1.f + q2pt2Wgh) * q2pt2;
385 trk.setCov(err2, 14);
386 TrackParCovChi2 = 0.f;
388 CADEBUG(printf(
"\t%21sPropaga Alpha %8.3f , X %8.3f - Y %8.3f, Z %8.3f - QPt %7.2f (%7.2f), SP %5.2f (%5.2f) --- Res %8.3f %8.3f --- Cov sY %8.3f sZ %8.3f sSP %8.3f sPt %8.3f - YPt %8.3f\n",
"", trk.getAlpha(),
x, trk.getParams()[0], trk.getParams()[1], trk.getParams()[4], trk.getParams()[4], trk.getParams()[2], trk.getParams()[2], trk.getParams()[0] -
y, trk.getParams()[1] -
z, sqrtf(trk.getCov()[0]), sqrtf(trk.getCov()[2]), sqrtf(trk.getCov()[5]), sqrtf(trk.getCov()[14]), trk.getCov()[10]));
389 std::array<float, 2> p = {
y,
z};
390 std::array<float, 3>
c = {0, 0, 0};
391 GPUTPCGMPropagator::GetErr2(
c[0],
c[2], *mPparam, getPar(trk)[2], getPar(trk)[3],
z,
x,
y, currentRow, clusterState, sector,
time, invAvgCharge, invCharge,
false);
392 TrackParCovChi2 += trk.getPredictedChi2(p,
c);
393 if (!trk.update(p,
c)) {
397 CADEBUG(printf(
"\t%21sFit Alpha %8.3f , X %8.3f - Y %8.3f, Z %8.3f - QPt %7.2f (%7.2f), SP %5.2f (%5.2f), DzDs %5.2f %16s --- Cov sY %8.3f sZ %8.3f sSP %8.3f sPt %8.3f - YPt %8.3f\n",
"", trk.getAlpha(),
x, trk.getParams()[0], trk.getParams()[1], trk.getParams()[4], trk.getParams()[4], trk.getParams()[2], trk.getParams()[2], trk.getParams()[3],
"", sqrtf(trk.getCov()[0]), sqrtf(trk.getCov()[2]), sqrtf(trk.getCov()[5]), sqrtf(trk.getCov()[14]), trk.getCov()[10]));
399 static_assert(
"Invalid template");
404 if constexpr (std::is_same_v<S, GPUTPCGMTrackParam>) {
405 float alpha = prop.GetAlpha();
406 trk.MoveToReference(prop, *mPparam,
alpha);
407 trk.NormalizeAlpha(
alpha);
408 prop.SetAlpha(
alpha);
409 }
else if constexpr (std::is_same_v<S, TrackParCov>) {
410 static constexpr float kDeg2Rad = M_PI / 180.f;
411 static constexpr float kSectAngle = 2 * M_PI / 18.f;
412 if (mPparam->rec.tpc.trackReferenceX <= 500) {
413 if (prop->PropagateToXBxByBz(trk, mPparam->rec.tpc.trackReferenceX)) {
414 if (CAMath::Abs(trk.getY()) > trk.getX() * CAMath::Tan(kSectAngle / 2.f)) {
415 float newAlpha = trk.getAlpha() + CAMath::Round(CAMath::ATan2(trk.getY(), trk.getX()) / kDeg2Rad / 20.f) * kSectAngle;
416 GPUTPCGMTrackParam::NormalizeAlpha(newAlpha);
417 trk.rotate(newAlpha) && prop->PropagateToXBxByBz(trk, mPparam->rec.tpc.trackReferenceX);
422 static_assert(
"Invalid template");
425 convertTrack<T, S, typename internal::refitTrackTypes<S>::propagator>(
trkX, trk, prop, &TrackParCovChi2);