84 tbb::blocked_range<int>(0, mTrkParams[iteration].TrackletsPerRoad()),
85 [&](
const tbb::blocked_range<int>& Layers) {
86 for (
int iLayer = Layers.begin(); iLayer < Layers.end(); ++iLayer) {
87 gsl::span<const Cluster> layer0 = mTimeFrame->getClustersOnLayer(rof0, iLayer);
91 float meanDeltaR{mTrkParams[iteration].LayerRadii[iLayer + 1] - mTrkParams[iteration].LayerRadii[iLayer]};
93 const int currentLayerClustersNum{
static_cast<int>(layer0.size())};
94 for (
int iCluster{0}; iCluster < currentLayerClustersNum; ++iCluster) {
95 const Cluster& currentCluster{layer0[iCluster]};
96 const int currentSortedIndex{mTimeFrame->getSortedIndex(rof0, iLayer, iCluster)};
98 if (mTimeFrame->isClusterUsed(iLayer, currentCluster.clusterId)) {
101 const float inverseR0{1.f / currentCluster.
radius};
103 for (
int iV{startVtx}; iV < endVtx; ++iV) {
104 auto& primaryVertex{primaryVertices[iV]};
105 if (primaryVertex.isFlagSet(2) && iteration != 3) {
108 const float resolution = o2::gpu::CAMath::Sqrt(Sq(mTrkParams[iteration].PVres) / primaryVertex.getNContributors() + Sq(mTimeFrame->getPositionResolution(iLayer)));
110 const float tanLambda{(currentCluster.zCoordinate - primaryVertex.getZ()) * inverseR0};
112 const float zAtRmin{tanLambda * (mTimeFrame->getMinR(iLayer + 1) - currentCluster.radius) + currentCluster.zCoordinate};
113 const float zAtRmax{tanLambda * (mTimeFrame->getMaxR(iLayer + 1) - currentCluster.radius) + currentCluster.zCoordinate};
115 const float sqInverseDeltaZ0{1.f / (Sq(currentCluster.zCoordinate - primaryVertex.getZ()) + 2.e-8f)};
116 const float sigmaZ{o2::gpu::CAMath::Sqrt(Sq(resolution) * Sq(tanLambda) * ((Sq(inverseR0) + sqInverseDeltaZ0) * Sq(meanDeltaR) + 1.f) + Sq(meanDeltaR * mTimeFrame->getMSangle(iLayer)))};
118 const int4 selectedBinsRect{getBinsRect(currentCluster, iLayer + 1, zAtRmin, zAtRmax, sigmaZ * mTrkParams[iteration].NSigmaCut, mTimeFrame->getPhiCut(iLayer))};
119 if (selectedBinsRect.x == 0 && selectedBinsRect.y == 0 && selectedBinsRect.z == 0 && selectedBinsRect.w == 0) {
123 int phiBinsNum{selectedBinsRect.
w - selectedBinsRect.y + 1};
125 if (phiBinsNum < 0) {
126 phiBinsNum += mTrkParams[iteration].PhiBins;
129 for (
int rof1{minRof}; rof1 <= maxRof; ++rof1) {
130 auto layer1 = mTimeFrame->getClustersOnLayer(rof1, iLayer + 1);
131 if (layer1.empty()) {
134 for (
int iPhiCount{0}; iPhiCount < phiBinsNum; iPhiCount++) {
135 int iPhiBin = (selectedBinsRect.y + iPhiCount) % mTrkParams[iteration].PhiBins;
136 const int firstBinIndex{mTimeFrame->mIndexTableUtils.getBinIndex(selectedBinsRect.x, iPhiBin)};
137 const int maxBinIndex{firstBinIndex + selectedBinsRect.z - selectedBinsRect.x + 1};
139 if (firstBinIndex < 0 || firstBinIndex > mTimeFrame->getIndexTable(rof1, iLayer + 1).size() ||
140 maxBinIndex < 0 || maxBinIndex > mTimeFrame->getIndexTable(rof1, iLayer + 1).size()) {
141 std::cout << iLayer <<
"\t" << iCluster <<
"\t" << zAtRmin <<
"\t" << zAtRmax <<
"\t" << sigmaZ * mTrkParams[iteration].NSigmaCut <<
"\t" << mTimeFrame->getPhiCut(iLayer) << std::endl;
142 std::cout << currentCluster.zCoordinate <<
"\t" << primaryVertex.getZ() <<
"\t" << currentCluster.radius << std::endl;
143 std::cout << mTimeFrame->getMinR(iLayer + 1) <<
"\t" << currentCluster.radius <<
"\t" << currentCluster.zCoordinate << std::endl;
144 std::cout <<
"Illegal access to IndexTable " << firstBinIndex <<
"\t" << maxBinIndex <<
"\t" << selectedBinsRect.z <<
"\t" << selectedBinsRect.x << std::endl;
148 const int firstRowClusterIndex = mTimeFrame->getIndexTable(rof1, iLayer + 1)[firstBinIndex];
149 const int maxRowClusterIndex = mTimeFrame->getIndexTable(rof1, iLayer + 1)[maxBinIndex];
150 for (
int iNextCluster{firstRowClusterIndex}; iNextCluster < maxRowClusterIndex; ++iNextCluster) {
151 if (iNextCluster >= (
int)layer1.size()) {
155 const Cluster& nextCluster{layer1[iNextCluster]};
156 if (mTimeFrame->isClusterUsed(iLayer + 1, nextCluster.clusterId)) {
160 const float deltaPhi{o2::gpu::GPUCommonMath::Abs(currentCluster.phi - nextCluster.phi)};
161 const float deltaZ{o2::gpu::GPUCommonMath::Abs(tanLambda * (nextCluster.radius - currentCluster.radius) +
162 currentCluster.zCoordinate - nextCluster.zCoordinate)};
164#ifdef OPTIMISATION_OUTPUT
166 int currentId{currentCluster.clusterId};
167 int nextId{nextCluster.clusterId};
168 for (
auto& lab1 : mTimeFrame->getClusterLabels(iLayer, currentId)) {
169 for (
auto& lab2 : mTimeFrame->getClusterLabels(iLayer + 1, nextId)) {
170 if (lab1 == lab2 && lab1.isValid()) {
175 if (
label.isValid()) {
179 off << std::format(
"{}\t{:d}\t{}\t{}\t{}\t{}", iLayer,
label.isValid(), (tanLambda * (nextCluster.radius - currentCluster.radius) + currentCluster.zCoordinate - nextCluster.zCoordinate) / sigmaZ, tanLambda, resolution, sigmaZ) << std::endl;
182 if (deltaZ / sigmaZ < mTrkParams[iteration].NSigmaCut &&
183 (deltaPhi < mTimeFrame->getPhiCut(iLayer) ||
186 mTimeFrame->getTrackletsLookupTable()[iLayer - 1][currentSortedIndex]++;
188 const float phi{o2::gpu::GPUCommonMath::ATan2(currentCluster.yCoordinate - nextCluster.yCoordinate,
189 currentCluster.xCoordinate - nextCluster.xCoordinate)};
190 const float tanL{(currentCluster.zCoordinate - nextCluster.zCoordinate) /
191 (currentCluster.radius - nextCluster.radius)};
192 mTimeFrame->getTracklets()[iLayer].emplace_back(currentSortedIndex, mTimeFrame->getSortedIndex(rof1, iLayer + 1, iNextCluster), tanL, phi, rof0, rof1);
204 auto sortTracklets = [](
const Tracklet&
a,
const Tracklet&
b) ->
bool {
205 return a.firstClusterIndex <
b.firstClusterIndex || (
a.firstClusterIndex ==
b.firstClusterIndex &&
a.secondClusterIndex <
b.secondClusterIndex);
207 auto equalTracklets = [](
const Tracklet&
a,
const Tracklet&
b) ->
bool {
208 return a.firstClusterIndex ==
b.firstClusterIndex &&
a.secondClusterIndex ==
b.secondClusterIndex;
211 mTaskArena.execute([&] {
213 tbb::blocked_range<int>(0, mTrkParams[iteration].
CellsPerRoad()),
214 [&](
const tbb::blocked_range<int>& Layers) {
215 for (
int iLayer = Layers.begin(); iLayer < Layers.end(); ++iLayer) {
217 auto& trkl{mTimeFrame->getTracklets()[iLayer + 1]};
218 tbb::parallel_sort(trkl.begin(), trkl.end(), sortTracklets);
220 trkl.erase(std::unique(trkl.begin(), trkl.end(), equalTracklets), trkl.end());
221 trkl.shrink_to_fit();
223 auto& lut{mTimeFrame->getTrackletsLookupTable()[iLayer]};
224 std::fill(lut.begin(), lut.end(), 0);
228 for (
const auto& tkl : trkl) {
229 lut[tkl.firstClusterIndex]++;
231 std::exclusive_scan(lut.begin(), lut.end(), lut.begin(), 0);
232 lut.push_back(trkl.size());
239 auto& trklt0 = mTimeFrame->getTracklets()[0];
240 mTaskArena.execute([&] { tbb::parallel_sort(trklt0.begin(), trklt0.end(), sortTracklets); });
241 trklt0.erase(std::unique(trklt0.begin(), trklt0.end(), equalTracklets), trklt0.end());
242 trklt0.shrink_to_fit();
245 if (mTimeFrame->hasMCinformation()) {
246 for (
int iLayer{0}; iLayer < mTrkParams[iteration].TrackletsPerRoad(); ++iLayer) {
247 for (
auto& trk : mTimeFrame->getTracklets()[iLayer]) {
249 int currentId{mTimeFrame->getClusters()[iLayer][trk.firstClusterIndex].clusterId};
250 int nextId{mTimeFrame->getClusters()[iLayer + 1][trk.secondClusterIndex].clusterId};
251 for (
auto& lab1 : mTimeFrame->getClusterLabels(iLayer, currentId)) {
252 for (
auto& lab2 : mTimeFrame->getClusterLabels(iLayer + 1, nextId)) {
253 if (lab1 == lab2 && lab1.isValid()) {
258 if (
label.isValid()) {
262 mTimeFrame->getTrackletsLabel(iLayer).emplace_back(
label);
268template <
int nLayers>
271#ifdef OPTIMISATION_OUTPUT
273 std::ofstream off(std::format(
"cells{}.txt", iter++));
276 constexpr float radl = 9.36f;
277 constexpr float rho = 2.33f;
279 for (
int iLayer = 0; iLayer < mTrkParams[iteration].CellsPerRoad(); ++iLayer) {
284 if (mTimeFrame->hasMCinformation()) {
289 mTaskArena.execute([&] {
291 tbb::blocked_range<int>(0, mTrkParams[iteration].CellsPerRoad()),
292 [&](
const tbb::blocked_range<int>& Layers) {
293 for (
int iLayer = Layers.begin(); iLayer < Layers.end(); ++iLayer) {
295 if (mTimeFrame->getTracklets()[iLayer + 1].empty() ||
296 mTimeFrame->getTracklets()[iLayer].empty()) {
300#ifdef OPTIMISATION_OUTPUT
301 float resolution{o2::gpu::CAMath::Sqrt(0.5f * (mTrkParams[iteration].SystErrorZ2[iLayer] + mTrkParams[iteration].SystErrorZ2[iLayer + 1] + mTrkParams[iteration].SystErrorZ2[iLayer + 2] + mTrkParams[iteration].SystErrorY2[iLayer] + mTrkParams[iteration].SystErrorY2[iLayer + 1] + mTrkParams[iteration].SystErrorY2[iLayer + 2])) / mTrkParams[iteration].LayerResolution[iLayer]};
302 resolution = resolution > 1.e-12 ? resolution : 1.f;
306 const int currentLayerTrackletsNum{
static_cast<int>(mTimeFrame->getTracklets()[iLayer].size())};
309 tbb::blocked_range<int>(0, currentLayerTrackletsNum),
310 [&](
const tbb::blocked_range<int>& Tracklets) {
311 for (
int iTracklet = Tracklets.begin(); iTracklet < Tracklets.end(); ++iTracklet) {
312 const Tracklet& currentTracklet{mTimeFrame->getTracklets()[iLayer][iTracklet]};
314 const int nextLayerFirstTrackletIndex{
315 mTimeFrame->getTrackletsLookupTable()[iLayer][nextLayerClusterIndex]};
316 const int nextLayerLastTrackletIndex{
317 mTimeFrame->getTrackletsLookupTable()[iLayer][nextLayerClusterIndex + 1]};
319 if (nextLayerFirstTrackletIndex == nextLayerLastTrackletIndex) {
324 for (
int iNextTracklet{nextLayerFirstTrackletIndex}; iNextTracklet < nextLayerLastTrackletIndex; ++iNextTracklet) {
325 if (mTimeFrame->getTracklets()[iLayer + 1][iNextTracklet].firstClusterIndex != nextLayerClusterIndex) {
328 const Tracklet& nextTracklet{mTimeFrame->getTracklets()[iLayer + 1][iNextTracklet]};
329 const float deltaTanLambda{std::abs(currentTracklet.tanLambda - nextTracklet.tanLambda)};
331#ifdef OPTIMISATION_OUTPUT
332 bool good{mTimeFrame->getTrackletsLabel(iLayer)[iTracklet] == mTimeFrame->getTrackletsLabel(iLayer + 1)[iNextTracklet]};
333 float signedDelta{currentTracklet.
tanLambda - nextTracklet.tanLambda};
334 off << std::format(
"{}\t{:d}\t{}\t{}\t{}\t{}", iLayer, good, signedDelta, signedDelta / (mTrkParams[iteration].CellDeltaTanLambdaSigma), tanLambda, resolution) << std::endl;
337 if (deltaTanLambda / mTrkParams[iteration].CellDeltaTanLambdaSigma < mTrkParams[iteration].NSigmaCut) {
341 mTimeFrame->getClusters()[iLayer][currentTracklet.firstClusterIndex].clusterId,
342 mTimeFrame->getClusters()[iLayer + 1][nextTracklet.firstClusterIndex].clusterId,
343 mTimeFrame->getClusters()[iLayer + 2][nextTracklet.secondClusterIndex].clusterId};
344 const auto& cluster1_glo = mTimeFrame->getUnsortedClusters()[iLayer][clusId[0]];
345 const auto& cluster2_glo = mTimeFrame->getUnsortedClusters()[iLayer + 1][clusId[1]];
346 const auto& cluster3_tf = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer + 2)[clusId[2]];
347 auto track{buildTrackSeed(cluster1_glo, cluster2_glo, cluster3_tf)};
351 for (
int iC{2}; iC--;) {
352 const TrackingFrameInfo& trackingHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer + iC)[clusId[iC]];
362 if (!track.correctForMaterial(mTrkParams[0].LayerxX0[iLayer + iC], mTrkParams[0].LayerxX0[iLayer] * radl * rho,
true)) {
367 if (!iC && predChi2 > mTrkParams[iteration].MaxChi2ClusterAttachment) {
383 perTrackletCount[iTracklet] = foundCells;
388 std::exclusive_scan(perTrackletCount.begin(), perTrackletCount.end(), perTrackletCount.begin(), 0);
389 auto totalCells{perTrackletCount.back()};
390 if (totalCells == 0) {
393 auto& layerCells = mTimeFrame->getCells()[iLayer];
394 layerCells.resize(totalCells);
397 tbb::blocked_range<int>(0, currentLayerTrackletsNum),
398 [&](
const tbb::blocked_range<int>& Tracklets) {
399 for (
int iTracklet = Tracklets.begin(); iTracklet < Tracklets.end(); ++iTracklet) {
400 if (perTrackletCount[iTracklet] == perTrackletCount[iTracklet + 1]) {
404 const Tracklet& currentTracklet{mTimeFrame->getTracklets()[iLayer][iTracklet]};
406 const int nextLayerFirstTrackletIndex{
407 mTimeFrame->getTrackletsLookupTable()[iLayer][nextLayerClusterIndex]};
408 const int nextLayerLastTrackletIndex{
409 mTimeFrame->getTrackletsLookupTable()[iLayer][nextLayerClusterIndex + 1]};
411 int position = perTrackletCount[iTracklet];
412 for (
int iNextTracklet{nextLayerFirstTrackletIndex}; iNextTracklet < nextLayerLastTrackletIndex; ++iNextTracklet) {
413 if (mTimeFrame->getTracklets()[iLayer + 1][iNextTracklet].firstClusterIndex != nextLayerClusterIndex) {
416 const Tracklet& nextTracklet{mTimeFrame->getTracklets()[iLayer + 1][iNextTracklet]};
417 const float deltaTanLambda{std::abs(currentTracklet.tanLambda - nextTracklet.tanLambda)};
419 if (deltaTanLambda / mTrkParams[iteration].CellDeltaTanLambdaSigma < mTrkParams[iteration].NSigmaCut) {
423 mTimeFrame->getClusters()[iLayer][currentTracklet.firstClusterIndex].clusterId,
424 mTimeFrame->getClusters()[iLayer + 1][nextTracklet.firstClusterIndex].clusterId,
425 mTimeFrame->getClusters()[iLayer + 2][nextTracklet.secondClusterIndex].clusterId};
426 const auto& cluster1_glo = mTimeFrame->getUnsortedClusters()[iLayer][clusId[0]];
427 const auto& cluster2_glo = mTimeFrame->getUnsortedClusters()[iLayer + 1][clusId[1]];
428 const auto& cluster3_tf = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer + 2)[clusId[2]];
429 auto track{buildTrackSeed(cluster1_glo, cluster2_glo, cluster3_tf)};
433 for (
int iC{2}; iC--;) {
434 const TrackingFrameInfo& trackingHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer + iC)[clusId[iC]];
444 if (!track.correctForMaterial(mTrkParams[0].LayerxX0[iLayer + iC], mTrkParams[0].LayerxX0[iLayer] * radl * rho,
true)) {
449 if (!iC && predChi2 > mTrkParams[iteration].MaxChi2ClusterAttachment) {
461 layerCells[position++] =
CellSeed(iLayer, clusId[0], clusId[1], clusId[2], iTracklet, iNextTracklet, track, chi2);
469 auto& lut = mTimeFrame->getCellsLookupTable()[iLayer - 1];
470 lut.resize(currentLayerTrackletsNum + 1);
471 std::copy_n(perTrackletCount.begin(), currentLayerTrackletsNum + 1, lut.begin());
478 if (mTimeFrame->hasMCinformation()) {
479 for (
int iLayer{0}; iLayer < mTrkParams[iteration].CellsPerRoad(); ++iLayer) {
480 for (
auto& cell : mTimeFrame->getCells()[iLayer]) {
481 MCCompLabel currentLab{mTimeFrame->getTrackletsLabel(iLayer)[cell.getFirstTrackletIndex()]};
482 MCCompLabel nextLab{mTimeFrame->getTrackletsLabel(iLayer + 1)[cell.getSecondTrackletIndex()]};
483 mTimeFrame->getCellsLabel(iLayer).emplace_back(currentLab == nextLab ? currentLab :
MCCompLabel());
489 for (
int iLayer{0}; iLayer < mTrkParams[iteration].CellsPerRoad(); ++iLayer) {
490 std::cout <<
"Cells on layer " << iLayer <<
" " << mTimeFrame->getCells()[iLayer].size() << std::endl;
495template <
int nLayers>
498#ifdef OPTIMISATION_OUTPUT
499 std::ofstream off(std::format(
"cellneighs{}.txt", iteration));
501 for (
int iLayer{0}; iLayer < mTrkParams[iteration].CellsPerRoad() - 1; ++iLayer) {
502 const int nextLayerCellsNum{
static_cast<int>(mTimeFrame->getCells()[iLayer + 1].size())};
505 if (mTimeFrame->getCells()[iLayer + 1].empty() ||
506 mTimeFrame->getCellsLookupTable()[iLayer].empty()) {
510 mTaskArena.execute([&] {
511 int layerCellsNum{
static_cast<int>(mTimeFrame->getCells()[iLayer].
size())};
515 tbb::blocked_range<int>(0, layerCellsNum),
516 [&](
const tbb::blocked_range<int>& Cells) {
517 for (
int iCell = Cells.begin(); iCell < Cells.end(); ++iCell) {
518 const auto& currentCellSeed{mTimeFrame->getCells()[iLayer][iCell]};
519 const int nextLayerTrackletIndex{currentCellSeed.getSecondTrackletIndex()};
520 const int nextLayerFirstCellIndex{mTimeFrame->getCellsLookupTable()[iLayer][nextLayerTrackletIndex]};
521 const int nextLayerLastCellIndex{mTimeFrame->getCellsLookupTable()[iLayer][nextLayerTrackletIndex + 1]};
523 int foundNextCells{0};
524 for (
int iNextCell{nextLayerFirstCellIndex}; iNextCell < nextLayerLastCellIndex; ++iNextCell) {
525 auto nextCellSeed{mTimeFrame->getCells()[iLayer + 1][iNextCell]};
526 if (nextCellSeed.getFirstTrackletIndex() != nextLayerTrackletIndex) {
530 if (!nextCellSeed.rotate(currentCellSeed.getAlpha()) ||
531 !nextCellSeed.propagateTo(currentCellSeed.getX(), getBz())) {
534 float chi2 = currentCellSeed.getPredictedChi2(nextCellSeed);
536#ifdef OPTIMISATION_OUTPUT
537 bool good{mTimeFrame->getCellsLabel(iLayer)[iCell] == mTimeFrame->getCellsLabel(iLayer + 1)[iNextCell]};
538 off << std::format(
"{}\t{:d}\t{}", iLayer, good, chi2) << std::endl;
541 if (chi2 > mTrkParams[0].MaxChi2ClusterAttachment) {
546 perCellCount[iCell] = foundNextCells;
550 std::exclusive_scan(perCellCount.begin(), perCellCount.end(), perCellCount.begin(), 0);
551 int totalCellNeighbours = perCellCount.back();
552 if (totalCellNeighbours == 0) {
558 int cell{-1}, nextCell{-1},
level{-1};
561 cellsNeighbours.resize(totalCellNeighbours);
564 tbb::blocked_range<int>(0, layerCellsNum),
565 [&](
const tbb::blocked_range<int>& Cells) {
566 for (
int iCell = Cells.begin(); iCell < Cells.end(); ++iCell) {
567 if (perCellCount[iCell] == perCellCount[iCell + 1]) {
570 const auto& currentCellSeed{mTimeFrame->getCells()[iLayer][iCell]};
571 const int nextLayerTrackletIndex{currentCellSeed.getSecondTrackletIndex()};
572 const int nextLayerFirstCellIndex{mTimeFrame->getCellsLookupTable()[iLayer][nextLayerTrackletIndex]};
573 const int nextLayerLastCellIndex{mTimeFrame->getCellsLookupTable()[iLayer][nextLayerTrackletIndex + 1]};
575 int position = perCellCount[iCell];
576 for (
int iNextCell{nextLayerFirstCellIndex}; iNextCell < nextLayerLastCellIndex; ++iNextCell) {
577 auto nextCellSeed{mTimeFrame->getCells()[iLayer + 1][iNextCell]};
578 if (nextCellSeed.getFirstTrackletIndex() != nextLayerTrackletIndex) {
582 if (!nextCellSeed.rotate(currentCellSeed.getAlpha()) ||
583 !nextCellSeed.propagateTo(currentCellSeed.getX(), getBz())) {
587 float chi2 = currentCellSeed.getPredictedChi2(nextCellSeed);
588 if (chi2 > mTrkParams[0].MaxChi2ClusterAttachment) {
592 cellsNeighbours[position++] = {iCell, iNextCell, currentCellSeed.getLevel() + 1};
597 tbb::parallel_sort(cellsNeighbours.begin(), cellsNeighbours.end(), [](
const auto&
a,
const auto&
b) {
598 return a.nextCell < b.nextCell;
601 auto& cellsNeighbourLUT = mTimeFrame->getCellsNeighboursLUT()[iLayer];
602 cellsNeighbourLUT.assign(nextLayerCellsNum, 0);
603 for (
const auto& neigh : cellsNeighbours) {
604 ++cellsNeighbourLUT[neigh.nextCell];
606 std::inclusive_scan(cellsNeighbourLUT.begin(), cellsNeighbourLUT.end(), cellsNeighbourLUT.begin());
608 mTimeFrame->getCellsNeighbours()[iLayer].reserve(totalCellNeighbours);
609 std::ranges::transform(cellsNeighbours, std::back_inserter(mTimeFrame->getCellsNeighbours()[iLayer]), [](
const auto& neigh) { return neigh.cell; });
611 auto it = cellsNeighbours.begin();
612 while (it != cellsNeighbours.end()) {
613 const int current_nextCell = it->nextCell;
614 auto group_end = std::find_if_not(it, cellsNeighbours.end(),
615 [current_nextCell](
const auto& nb) { return nb.nextCell == current_nextCell; });
616 const auto max_level_it = std::max_element(it, group_end,
617 [](
const auto&
a,
const auto&
b) {
return a.level <
b.level; });
618 mTimeFrame->getCells()[iLayer + 1][current_nextCell].setLevel(max_level_it->level);
625template <
int nLayers>
628 CA_DEBUGGER(std::cout <<
"Processing neighbours layer " << iLayer <<
" level " << iLevel <<
", size of the cell seeds: " << currentCellSeed.size() << std::endl);
632 int failed[5]{0, 0, 0, 0, 0}, attempts{0}, failedByMismatch{0};
635 mTaskArena.execute([&] {
638 tbb::blocked_range<int>(0, (
int)currentCellSeed.size()),
639 [&](
const tbb::blocked_range<int>& Cells) {
640 for (int iCell = Cells.begin(); iCell < Cells.end(); ++iCell) {
641 const CellSeed& currentCell{currentCellSeed[iCell]};
643 if (currentCell.getLevel() != iLevel) {
646 if (currentCellId.empty() && (mTimeFrame->isClusterUsed(iLayer, currentCell.getFirstClusterIndex()) ||
647 mTimeFrame->isClusterUsed(iLayer + 1, currentCell.getSecondClusterIndex()) ||
648 mTimeFrame->isClusterUsed(iLayer + 2, currentCell.getThirdClusterIndex()))) {
651 const int cellId = currentCellId.empty() ? iCell : currentCellId[iCell];
652 const int startNeighbourId{cellId ? mTimeFrame->getCellsNeighboursLUT()[iLayer - 1][cellId - 1] : 0};
653 const int endNeighbourId{mTimeFrame->getCellsNeighboursLUT()[iLayer - 1][cellId]};
655 for (int iNeighbourCell{startNeighbourId}; iNeighbourCell < endNeighbourId; ++iNeighbourCell) {
656 CA_DEBUGGER(attempts++);
657 const int neighbourCellId = mTimeFrame->getCellsNeighbours()[iLayer - 1][iNeighbourCell];
658 const CellSeed& neighbourCell = mTimeFrame->getCells()[iLayer - 1][neighbourCellId];
659 if (neighbourCell.getSecondTrackletIndex() != currentCell.getFirstTrackletIndex()) {
660 CA_DEBUGGER(failedByMismatch++);
663 if (mTimeFrame->isClusterUsed(iLayer - 1, neighbourCell.getFirstClusterIndex())) {
666 if (currentCell.getLevel() - 1 != neighbourCell.getLevel()) {
667 CA_DEBUGGER(failed[0]++);
671 CellSeed seed{currentCell};
672 auto& trHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer - 1)[neighbourCell.getFirstClusterIndex()];
674 if (!seed.rotate(trHit.alphaTrackingFrame)) {
675 CA_DEBUGGER(failed[1]++);
679 if (!propagator->propagateToX(seed, trHit.xTrackingFrame, getBz(), o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, mCorrType)) {
680 CA_DEBUGGER(failed[2]++);
684 if (mCorrType == o2::base::PropagatorF::MatCorrType::USEMatCorrNONE) {
687 if (!seed.correctForMaterial(mTrkParams[0].LayerxX0[iLayer - 1], mTrkParams[0].LayerxX0[iLayer - 1] * radl * rho, true)) {
692 auto predChi2{seed.getPredictedChi2Quiet(trHit.positionTrackingFrame, trHit.covarianceTrackingFrame)};
693 if ((predChi2 > mTrkParams[0].MaxChi2ClusterAttachment) || predChi2 < 0.f) {
694 CA_DEBUGGER(failed[3]++);
697 seed.setChi2(seed.getChi2() + predChi2);
698 if (!seed.o2::track::TrackParCov::update(trHit.positionTrackingFrame, trHit.covarianceTrackingFrame)) {
699 CA_DEBUGGER(failed[4]++);
704 perCellCount[iCell] = foundSeeds;
708 std::exclusive_scan(perCellCount.begin(), perCellCount.end(), perCellCount.begin(), 0);
709 auto totalNeighbours{perCellCount.back()};
710 if (totalNeighbours == 0) {
713 updatedCellSeeds.resize(totalNeighbours);
714 updatedCellsIds.resize(totalNeighbours);
717 tbb::blocked_range<int>(0, (
int)currentCellSeed.size()),
718 [&](
const tbb::blocked_range<int>& Cells) {
719 for (int iCell = Cells.begin(); iCell < Cells.end(); ++iCell) {
720 if (perCellCount[iCell] == perCellCount[iCell + 1]) {
725 const CellSeed& currentCell{currentCellSeed[iCell]};
726 const int cellId = currentCellId.empty() ? iCell : currentCellId[iCell];
727 const int startNeighbourId{cellId ? mTimeFrame->getCellsNeighboursLUT()[iLayer - 1][cellId - 1] : 0};
728 const int endNeighbourId{mTimeFrame->getCellsNeighboursLUT()[iLayer - 1][cellId]};
730 int offset = perCellCount[iCell];
731 for (int iNeighbourCell{startNeighbourId}; iNeighbourCell < endNeighbourId; ++iNeighbourCell) {
732 const int neighbourCellId = mTimeFrame->getCellsNeighbours()[iLayer - 1][iNeighbourCell];
733 const CellSeed& neighbourCell = mTimeFrame->getCells()[iLayer - 1][neighbourCellId];
734 if (neighbourCell.getSecondTrackletIndex() != currentCell.getFirstTrackletIndex() ||
735 mTimeFrame->isClusterUsed(iLayer - 1, neighbourCell.getFirstClusterIndex()) ||
736 currentCell.getLevel() - 1 != neighbourCell.getLevel()) {
740 auto seed = currentCell;
742 const auto& trHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer - 1)[neighbourCell.getFirstClusterIndex()];
743 if (!seed.rotate(trHit.alphaTrackingFrame) || !propagator->propagateToX(seed, trHit.xTrackingFrame, getBz(), o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, mCorrType)) {
747 if (mCorrType == o2::base::PropagatorF::MatCorrType::USEMatCorrNONE) {
750 if (!seed.correctForMaterial(mTrkParams[0].LayerxX0[iLayer - 1], mTrkParams[0].LayerxX0[iLayer - 1] * radl * rho, true)) {
755 auto predChi2{seed.getPredictedChi2Quiet(trHit.positionTrackingFrame, trHit.covarianceTrackingFrame)};
756 if ((predChi2 > mTrkParams[0].MaxChi2ClusterAttachment) || predChi2 < 0.f) {
759 seed.setChi2(seed.getChi2() + predChi2);
760 if (!seed.o2::track::TrackParCov::update(trHit.positionTrackingFrame, trHit.covarianceTrackingFrame)) {
764 seed.getClusters()[iLayer - 1] = neighbourCell.getFirstClusterIndex();
765 seed.setLevel(neighbourCell.getLevel());
766 seed.setFirstTrackletIndex(neighbourCell.getFirstTrackletIndex());
767 seed.setSecondTrackletIndex(neighbourCell.getSecondTrackletIndex());
769 updatedCellSeeds[offset] = seed;
770 updatedCellsIds[offset++] = neighbourCellId;
777 std::cout <<
"\t\t- Found " << updatedCellSeeds.size() <<
" cell seeds out of " << attempts <<
" attempts" << std::endl;
778 std::cout <<
"\t\t\t> " <<
failed[0] <<
" failed because of level" << std::endl;
779 std::cout <<
"\t\t\t> " <<
failed[1] <<
" failed because of rotation" << std::endl;
780 std::cout <<
"\t\t\t> " <<
failed[2] <<
" failed because of propagation" << std::endl;
781 std::cout <<
"\t\t\t> " <<
failed[3] <<
" failed because of chi2 cut" << std::endl;
782 std::cout <<
"\t\t\t> " <<
failed[4] <<
" failed because of update" << std::endl;
783 std::cout <<
"\t\t\t> " << failedByMismatch <<
" failed because of mismatch" << std::endl;
787template <
int nLayers>
790 CA_DEBUGGER(std::cout <<
"Finding roads, iteration " << iteration << std::endl);
792 for (
int startLevel{mTrkParams[iteration].CellsPerRoad()}; startLevel >= mTrkParams[iteration].CellMinimumLevel(); --startLevel) {
793 CA_DEBUGGER(std::cout <<
"\t > Processing level " << startLevel << std::endl);
794 auto seedFilter = [&](
const CellSeed& seed) {
795 return seed.getQ2Pt() <= 1.e3 && seed.getChi2() <= mTrkParams[0].MaxChi2NDF * ((startLevel + 2) * 2 - 5);
798 for (
int startLayer{mTrkParams[iteration].CellsPerRoad() - 1}; startLayer >= startLevel - 1; --startLayer) {
799 if ((mTrkParams[iteration].StartLayerMask & (1 << (startLayer + 2))) == 0) {
802 CA_DEBUGGER(std::cout <<
"\t\t > Starting processing layer " << startLayer << std::endl);
806 processNeighbours(startLayer, startLevel, mTimeFrame->getCells()[startLayer], lastCellId, updatedCellSeed, updatedCellId);
808 int level = startLevel;
809 for (
int iLayer{startLayer - 1}; iLayer > 0 &&
level > 2; --iLayer) {
810 lastCellSeed.swap(updatedCellSeed);
811 lastCellId.swap(updatedCellId);
814 processNeighbours(iLayer, --
level, lastCellSeed, lastCellId, updatedCellSeed, updatedCellId);
819 if (!updatedCellSeed.empty()) {
820 trackSeeds.reserve(trackSeeds.size() + std::count_if(updatedCellSeed.begin(), updatedCellSeed.end(), seedFilter));
821 std::copy_if(updatedCellSeed.begin(), updatedCellSeed.end(), std::back_inserter(trackSeeds), seedFilter);
825 if (trackSeeds.empty()) {
830 mTaskArena.execute([&] {
833 tbb::blocked_range<size_t>(
size_t(0), trackSeeds.size()),
834 [&](
const tbb::blocked_range<size_t>& Seeds) {
835 for (int iSeed = Seeds.begin(); iSeed < Seeds.end(); ++iSeed) {
836 const CellSeed& seed{trackSeeds[iSeed]};
837 TrackITSExt temporaryTrack{seed};
838 temporaryTrack.resetCovariance();
839 temporaryTrack.setChi2(0);
840 for (int iL{0}; iL < 7; ++iL) {
841 temporaryTrack.setExternalClusterIndex(iL, seed.getCluster(iL), seed.getCluster(iL) != constants::its::UnusedIndex);
844 bool fitSuccess = fitTrack(temporaryTrack, 0, mTrkParams[0].NLayers, 1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF);
848 temporaryTrack.getParamOut() = temporaryTrack.getParamIn();
849 temporaryTrack.resetCovariance();
850 temporaryTrack.setChi2(0);
851 fitSuccess = fitTrack(temporaryTrack, mTrkParams[0].NLayers - 1, -1, -1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF, 50.f);
852 if (!fitSuccess || temporaryTrack.getPt() < mTrkParams[iteration].MinPt[mTrkParams[iteration].NLayers - temporaryTrack.getNClusters()]) {
855 ++perSeedCount[iSeed];
858 std::exclusive_scan(perSeedCount.begin(), perSeedCount.end(), perSeedCount.begin(), 0);
859 auto totalTracks{perSeedCount.back()};
860 if (totalTracks == 0) {
863 tracks.resize(totalTracks);
866 tbb::blocked_range<int>(0, (
int)trackSeeds.size()),
867 [&](
const tbb::blocked_range<int>& Seeds) {
868 for (int iSeed = Seeds.begin(); iSeed < Seeds.end(); ++iSeed) {
869 if (perSeedCount[iSeed] == perSeedCount[iSeed + 1]) {
872 const CellSeed& seed{trackSeeds[iSeed]};
873 auto& trk = tracks[perSeedCount[iSeed]] = TrackITSExt(seed);
874 trk.resetCovariance();
876 for (int iL{0}; iL < 7; ++iL) {
877 trk.setExternalClusterIndex(iL, seed.getCluster(iL), seed.getCluster(iL) != constants::its::UnusedIndex);
880 bool fitSuccess = fitTrack(trk, 0, mTrkParams[0].NLayers, 1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF);
884 trk.getParamOut() = trk.getParamIn();
885 trk.resetCovariance();
887 fitTrack(trk, mTrkParams[0].NLayers - 1, -1, -1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF, 50.f);
892 tbb::parallel_sort(tracks.begin(), tracks.end(), [](
const auto&
a,
const auto&
b) {
893 return a.getChi2() < b.getChi2();
897 for (
auto& track : tracks) {
899 bool isFirstShared{
false};
900 for (
int iLayer{0}; iLayer < mTrkParams[0].NLayers; ++iLayer) {
904 nShared +=
int(mTimeFrame->isClusterUsed(iLayer, track.getClusterIndex(iLayer)));
905 isFirstShared |= !iLayer && mTimeFrame->isClusterUsed(iLayer, track.getClusterIndex(iLayer));
908 if (nShared > mTrkParams[0].ClusterSharing) {
912 std::array<int, 3> rofs{INT_MAX, INT_MAX, INT_MAX};
913 for (
int iLayer{0}; iLayer < mTrkParams[0].NLayers; ++iLayer) {
917 mTimeFrame->markUsedCluster(iLayer, track.getClusterIndex(iLayer));
918 int currentROF = mTimeFrame->getClusterROF(iLayer, track.getClusterIndex(iLayer));
919 for (
int iR{0}; iR < 3; ++iR) {
920 if (rofs[iR] == INT_MAX) {
921 rofs[iR] = currentROF;
923 if (rofs[iR] == currentROF) {
928 if (rofs[2] != INT_MAX) {
931 track.setUserField(0);
932 track.getParamOut().setUserField(0);
933 if (rofs[1] != INT_MAX) {
934 track.setNextROFbit();
936 mTimeFrame->getTracks(o2::gpu::CAMath::Min(rofs[0], rofs[1])).emplace_back(track);
941template <
int nLayers>
944 for (
int rof{0}; rof < mTimeFrame->getNrof(); ++rof) {
945 for (
auto& track : mTimeFrame->getTracks(rof)) {
949 if ((mTrkParams[iteration].UseTrackFollowerMix || mTrkParams[iteration].UseTrackFollowerTop) && track.getLastClusterLayer() != mTrkParams[iteration].NLayers - 1) {
950 success = success || trackFollowing(&track, rof,
true, iteration);
952 if ((mTrkParams[iteration].UseTrackFollowerMix || (mTrkParams[iteration].UseTrackFollowerBot && !success)) && track.getFirstClusterLayer() != 0) {
953 success = success || trackFollowing(&track, rof,
false, iteration);
957 track.resetCovariance();
959 bool fitSuccess = fitTrack(track, 0, mTrkParams[iteration].NLayers, 1, mTrkParams[iteration].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF);
964 track.getParamOut() = track;
965 track.resetCovariance();
967 fitSuccess = fitTrack(track, mTrkParams[iteration].NLayers - 1, -1, -1, mTrkParams[iteration].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF, 50.);
972 mTimeFrame->mNExtendedTracks++;
973 mTimeFrame->mNExtendedUsedClusters += track.getNClusters() - backup.getNClusters();
974 auto pattern = track.getPattern();
975 auto diff = (
pattern & ~backup.getPattern()) & 0xff;
979 for (
int iLayer{0}; iLayer < mTrkParams[iteration].NLayers; ++iLayer) {
980 if (track.getClusterIndex(iLayer) == constants::its::UnusedIndex) {
983 mTimeFrame->markUsedCluster(iLayer, track.getClusterIndex(iLayer));
990template <
int nLayers>
994 mTimeFrame->fillPrimaryVerticesXandAlpha();
996 for (
auto& cell : mTimeFrame->getCells()[0]) {
997 auto& cluster3_glo = mTimeFrame->getClusters()[2][cell.getThirdClusterIndex()];
998 auto& cluster2_glo = mTimeFrame->getClusters()[1][cell.getSecondClusterIndex()];
999 auto& cluster1_glo = mTimeFrame->getClusters()[0][cell.getFirstClusterIndex()];
1000 if (mTimeFrame->isClusterUsed(2, cluster1_glo.clusterId) ||
1001 mTimeFrame->isClusterUsed(1, cluster2_glo.clusterId) ||
1002 mTimeFrame->isClusterUsed(0, cluster3_glo.clusterId)) {
1006 std::array<int, 3> rofs{
1007 mTimeFrame->getClusterROF(2, cluster3_glo.clusterId),
1008 mTimeFrame->getClusterROF(1, cluster2_glo.clusterId),
1009 mTimeFrame->getClusterROF(0, cluster1_glo.clusterId)};
1010 if (rofs[0] != rofs[1] && rofs[1] != rofs[2] && rofs[0] != rofs[2]) {
1015 if (rofs[1] == rofs[2]) {
1019 auto pvs{mTimeFrame->getPrimaryVertices(rof)};
1020 auto pvsXAlpha{mTimeFrame->getPrimaryVerticesXAlpha(rof)};
1022 const auto& cluster3_tf = mTimeFrame->getTrackingFrameInfoOnLayer(2)[cluster3_glo.clusterId];
1023 TrackITSExt temporaryTrack{buildTrackSeed(cluster1_glo, cluster2_glo, cluster3_tf)};
1024 temporaryTrack.setExternalClusterIndex(0, cluster1_glo.clusterId,
true);
1025 temporaryTrack.setExternalClusterIndex(1, cluster2_glo.clusterId,
true);
1026 temporaryTrack.setExternalClusterIndex(2, cluster3_glo.clusterId,
true);
1029 bool fitSuccess = fitTrack(temporaryTrack, 1, -1, -1);
1035 TrackITSExt bestTrack{temporaryTrack}, backup{temporaryTrack};
1036 float bestChi2{std::numeric_limits<float>::max()};
1037 for (
int iV{0}; iV < (
int)pvs.size(); ++iV) {
1038 temporaryTrack = backup;
1039 if (!temporaryTrack.rotate(pvsXAlpha[iV][1])) {
1042 if (!propagator->propagateTo(temporaryTrack, pvsXAlpha[iV][0],
true)) {
1046 float pvRes{mTrkParams[0].PVres / o2::gpu::CAMath::Sqrt(
float(pvs[iV].getNContributors()))};
1047 const float posVtx[2]{0.f, pvs[iV].getZ()};
1048 const float covVtx[3]{pvRes, 0.f, pvRes};
1049 float chi2 = temporaryTrack.getPredictedChi2Quiet(posVtx, covVtx);
1050 if (chi2 < bestChi2) {
1051 if (!temporaryTrack.track::TrackParCov::update(posVtx, covVtx)) {
1054 bestTrack = temporaryTrack;
1059 bestTrack.resetCovariance();
1060 bestTrack.setChi2(0.f);
1061 fitSuccess = fitTrack(bestTrack, 0, mTrkParams[0].NLayers, 1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF);
1065 bestTrack.getParamOut() = bestTrack;
1066 bestTrack.resetCovariance();
1067 bestTrack.setChi2(0.f);
1068 fitSuccess = fitTrack(bestTrack, mTrkParams[0].NLayers - 1, -1, -1, mTrkParams[0].MaxChi2ClusterAttachment, mTrkParams[0].MaxChi2NDF, 50.);
1072 mTimeFrame->markUsedCluster(0, bestTrack.getClusterIndex(0));
1073 mTimeFrame->markUsedCluster(1, bestTrack.getClusterIndex(1));
1074 mTimeFrame->markUsedCluster(2, bestTrack.getClusterIndex(2));
1075 mTimeFrame->getTracks(rof).emplace_back(bestTrack);
1079template <
int nLayers>
1084 for (
int iLayer{
start}; iLayer !=
end; iLayer += step) {
1085 if (track.getClusterIndex(iLayer) == constants::its::UnusedIndex) {
1088 const TrackingFrameInfo& trackingHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer)[track.getClusterIndex(iLayer)];
1090 if (!track.rotate(trackingHit.alphaTrackingFrame)) {
1098 if (mCorrType == o2::base::PropagatorF::MatCorrType::USEMatCorrNONE) {
1099 constexpr float radl = 9.36f;
1100 constexpr float rho = 2.33f;
1101 if (!track.correctForMaterial(mTrkParams[0].LayerxX0[iLayer], mTrkParams[0].LayerxX0[iLayer] * radl * rho,
true)) {
1106 auto predChi2{track.getPredictedChi2Quiet(trackingHit.positionTrackingFrame, trackingHit.covarianceTrackingFrame)};
1107 if ((nCl >= 3 && predChi2 > chi2clcut) || predChi2 < 0.f) {
1110 track.setChi2(track.getChi2() + predChi2);
1111 if (!track.o2::track::TrackParCov::update(trackingHit.positionTrackingFrame, trackingHit.covarianceTrackingFrame)) {
1116 return std::abs(track.getQ2Pt()) < maxQoverPt && track.getChi2() < chi2ndfcut * (nCl * 2 - 5);
1119template <
int nLayers>
1123 const int step = -1 + outward * 2;
1124 const int end = outward ? mTrkParams[iteration].NLayers - 1 : 0;
1126 for (
size_t iHypo{0}; iHypo < hypotheses.size(); ++iHypo) {
1127 auto hypo{hypotheses[iHypo]};
1128 int iLayer =
static_cast<int>(outward ? hypo.getLastClusterLayer() : hypo.getFirstClusterLayer());
1130 while (iLayer !=
end) {
1132 const float r = mTrkParams[iteration].LayerRadii[iLayer];
1139 auto& hypoParam{outward ? hypo.getParamOut() : hypo.getParamIn()};
1140 if (!propInstance->propagateToX(hypoParam,
x, mTimeFrame->getBz(), PropagatorF::MAX_SIN_PHI,
1141 PropagatorF::MAX_STEP, mTrkParams[iteration].CorrType)) {
1145 if (mTrkParams[iteration].CorrType == PropagatorF::MatCorrType::USEMatCorrNONE) {
1146 constexpr float radl = 9.36f;
1147 constexpr float rho = 2.33f;
1148 if (!hypoParam.correctForMaterial(mTrkParams[iteration].LayerxX0[iLayer], mTrkParams[iteration].LayerxX0[iLayer] * radl * rho,
true)) {
1154 const float phi{hypoParam.getPhi()};
1155 const float ePhi{o2::gpu::CAMath::Sqrt(hypoParam.getSigmaSnp2() / hypoParam.getCsp2())};
1156 const float z{hypoParam.getZ()};
1157 const float eZ{o2::gpu::CAMath::Sqrt(hypoParam.getSigmaZ2())};
1158 const int4 selectedBinsRect{getBinsRect(iLayer, phi, mTrkParams[iteration].TrackFollowerNSigmaCutPhi * ePhi,
z, mTrkParams[iteration].TrackFollowerNSigmaCutZ * eZ)};
1159 if (selectedBinsRect.x == 0 && selectedBinsRect.y == 0 && selectedBinsRect.z == 0 && selectedBinsRect.w == 0) {
1163 int phiBinsNum{selectedBinsRect.
w - selectedBinsRect.y + 1};
1165 if (phiBinsNum < 0) {
1166 phiBinsNum += mTrkParams[iteration].PhiBins;
1169 gsl::span<const Cluster> layer1 = mTimeFrame->getClustersOnLayer(rof, iLayer);
1170 if (layer1.empty()) {
1175 for (
int iPhiCount = 0; iPhiCount < phiBinsNum; iPhiCount++) {
1176 int iPhiBin = (selectedBinsRect.y + iPhiCount) % mTrkParams[iteration].PhiBins;
1177 const int firstBinIndex{mTimeFrame->mIndexTableUtils.getBinIndex(selectedBinsRect.x, iPhiBin)};
1178 const int maxBinIndex{firstBinIndex + selectedBinsRect.z - selectedBinsRect.x + 1};
1179 const int firstRowClusterIndex = mTimeFrame->getIndexTable(rof, iLayer)[firstBinIndex];
1180 const int maxRowClusterIndex = mTimeFrame->getIndexTable(rof, iLayer)[maxBinIndex];
1182 for (
int iNextCluster{firstRowClusterIndex}; iNextCluster < maxRowClusterIndex; ++iNextCluster) {
1183 if (iNextCluster >= (
int)layer1.size()) {
1186 const Cluster& nextCluster{layer1[iNextCluster]};
1188 if (mTimeFrame->isClusterUsed(iLayer, nextCluster.clusterId)) {
1192 const TrackingFrameInfo& trackingHit = mTimeFrame->getTrackingFrameInfoOnLayer(iLayer)[nextCluster.clusterId];
1194 auto tbupdated{hypo};
1195 auto& tbuParams = outward ? tbupdated.getParamOut() : tbupdated.getParamIn();
1200 if (!propInstance->propagateToX(tbuParams, trackingHit.
xTrackingFrame, mTimeFrame->getBz(),
1201 PropagatorF::MAX_SIN_PHI, PropagatorF::MAX_STEP, PropagatorF::MatCorrType::USEMatCorrNONE)) {
1206 if (predChi2 >= track->getChi2() * mTrkParams[iteration].NSigmaCut) {
1213 tbupdated.setChi2(tbupdated.getChi2() + predChi2);
1214 tbupdated.setExternalClusterIndex(iLayer, nextCluster.clusterId,
true);
1215 hypotheses.emplace_back(tbupdated);
1222 bool swapped{
false};
1223 for (
auto& hypo : hypotheses) {
1224 if (hypo.isBetter(*bestHypo, track->getChi2() * mTrkParams[iteration].NSigmaCut)) {
1235template <
int nLayers>
1238 float ca{-999.f}, sa{-999.f};
1249 float tgp{1.f}, crv{1.f}, snp{-999.f}, tgl12{-999.f}, tgl23{-999.f}, q2pt{1.f / track::kMostProbablePt}, q2pt2{1.f}, sg2q2pt{-999.f};
1251 tgp = o2::gpu::CAMath::ATan2(y3 -
y1, x3 -
x1);
1252 snp = tgp / o2::gpu::CAMath::Sqrt(1.f + tgp * tgp);
1254 crv = math_utils::computeCurvature(x3, y3, x2, y2,
x1,
y1);
1255 snp = crv * (
x3 - math_utils::computeCurvatureCentreX(x3, y3, x2, y2,
x1,
y1));
1259 tgl12 = math_utils::computeTanDipAngle(
x1,
y1, x2, y2, z1, z2);
1260 tgl23 = math_utils::computeTanDipAngle(x2, y2, x3, y3, z2, z3);
1261 sg2q2pt = track::kC1Pt2max * (q2pt2 > 0.0005f ? (q2pt2 < 1.f ? q2pt2 : 1.f) : 0.0005f);
1262 return {tf3.
xTrackingFrame, tf3.
alphaTrackingFrame, {
y3, z3, snp, 0.5f * (tgl12 + tgl23), q2pt}, {tf3.
covarianceTrackingFrame[0], tf3.
covarianceTrackingFrame[1], tf3.
covarianceTrackingFrame[2], 0.f, 0.f, track::kCSnp2max, 0.f, 0.f, 0.f, track::kCTgl2max, 0.f, 0.f, 0.f, 0.f, sg2q2pt}};
1265template <
int nLayers>
1269 mIsZeroField = std::abs(mBz) < 0.01;
1270 mTimeFrame->
setBz(bz);
1273template <
int nLayers>
1279template <
int nLayers>
1282 if (mNThreads ==
n && mTaskArena.is_active()) {
1285 mNThreads =
n > 0 ?
n : 1;
1286#if defined(OPTIMISATION_OUTPUT) || defined(CA_DEBUG)
1289 mTaskArena.initialize(mNThreads);
1290 LOGP(info,
"Setting tracker with {} threads.", mNThreads);