From 194fa51ba4f3354207a7cf58be2dcc6c0c1a2c31 Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Tue, 23 Apr 2024 19:26:39 +0200 Subject: [PATCH 1/6] Add multi rof vertexer idea --- .../include/ITStracking/Configuration.h | 2 +- .../tracking/include/ITStracking/TimeFrame.h | 6 + .../include/ITStracking/TrackingConfigParam.h | 2 + .../ITSMFT/ITS/tracking/src/Vertexer.cxx | 1 + .../ITS/tracking/src/VertexerTraits.cxx | 133 ++++++++++-------- 5 files changed, 88 insertions(+), 56 deletions(-) diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h index a700dc1e806c0..8453dce2a11d5 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h @@ -111,7 +111,7 @@ struct VertexingParameters { std::vector LayerRadii = {2.33959f, 3.14076f, 3.91924f, 19.6213f, 24.5597f, 34.388f, 39.3329f}; int ZBins{1}; int PhiBins{128}; - + int deltaRof = 0; float zCut = 0.002f; float phiCut = 0.005f; float pairCut = 0.04f; diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h index 3d99192d8f478..be07f3cd8d351 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h @@ -151,6 +151,7 @@ class TimeFrame bool isClusterUsed(int layer, int clusterId) const; void markUsedCluster(int layer, int clusterId); + gsl::span getUsedClusters(const int layer); std::vector>& getTracklets(); std::vector>& getTrackletsLookupTable(); @@ -519,6 +520,11 @@ inline bool TimeFrame::isClusterUsed(int layer, int clusterId) const return mUsedClusters[layer][clusterId]; } +inline gsl::span TimeFrame::getUsedClusters(const int layer) +{ + return {&mUsedClusters[layer][0], static_cast::size_type>(mUsedClusters[layer].size())}; +} + inline void TimeFrame::markUsedCluster(int layer, int clusterId) { mUsedClusters[layer][clusterId] = true; } inline std::vector>& TimeFrame::getTracklets() diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h index 36a5fd63b12d1..5ccb74babea6d 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h @@ -23,6 +23,8 @@ namespace its struct VertexerParamConfig : public o2::conf::ConfigurableParamHelper { bool allowSingleContribClusters = false; + // Number of ROFs to be considered for the vertexing + int deltaRof = 0; // geometrical cuts float zCut = 0.002f; diff --git a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx index acc59a66a7c2c..8a4a2cf7664a2 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx @@ -71,6 +71,7 @@ void Vertexer::getGlobalConfiguration() auto& grc = o2::its::GpuRecoParamConfig::Instance(); VertexingParameters verPar; + verPar.deltaRof = vc.deltaRof; verPar.allowSingleContribClusters = vc.allowSingleContribClusters; verPar.zCut = vc.zCut; verPar.phiCut = vc.phiCut; diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx index ceb5e4d0d44b7..3cc65fcf2c768 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx @@ -44,16 +44,18 @@ float smallestAngleDifference(float a, float b) return (diff < -constants::math::Pi) ? diff + constants::math::TwoPi : ((diff > constants::math::Pi) ? diff - constants::math::TwoPi : diff); } -template +template void trackleterKernelHost( const gsl::span& clustersNextLayer, // 0 2 const gsl::span& clustersCurrentLayer, // 1 1 + gsl::span usedClusters, int* indexTableNext, const float phiCut, std::vector& tracklets, gsl::span foundTracklets, const IndexTableUtils& utils, - const int rof, + const int pivotRof, + const int targetRof, const int rofFoundTrackletsOffset, const int maxTrackletsPerCluster = static_cast(2e3)) { @@ -78,13 +80,17 @@ void trackleterKernelHost( // loop on clusters next layer for (int iNextLayerClusterIndex{firstRowClusterIndex}; iNextLayerClusterIndex < maxRowClusterIndex && iNextLayerClusterIndex < static_cast(clustersNextLayer.size()); ++iNextLayerClusterIndex) { const Cluster& nextCluster{clustersNextLayer[iNextLayerClusterIndex]}; + if (usedClusters[nextCluster.clusterId]) { + continue; + } if (o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(currentCluster.phi, nextCluster.phi)) < phiCut) { if (storedTracklets < maxTrackletsPerCluster) { - if constexpr (!DryRun) { + if constexpr (!EvalRun) { + usedClusters[nextCluster.clusterId] = true; if constexpr (Mode == TrackletMode::Layer0Layer1) { - tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, rof, rof}; + tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof}; } else { - tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, rof, rof}; + tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof}; } } ++storedTracklets; @@ -93,7 +99,7 @@ void trackleterKernelHost( } } } - if constexpr (DryRun) { + if constexpr (EvalRun) { foundTracklets[iCurrentLayerClusterIndex] = storedTracklets; } else { cumulativeStoredTracklets += storedTracklets; @@ -173,31 +179,40 @@ void VertexerTraits::computeTracklets() #pragma omp parallel num_threads(mNThreads) { #pragma omp for schedule(dynamic) - for (int rofId = 0; rofId < mTimeFrame->getNrof(); ++rofId) { - trackleterKernelHost( - mTimeFrame->getClustersOnLayer(rofId, 0), - mTimeFrame->getClustersOnLayer(rofId, 1), - mTimeFrame->getIndexTable(rofId, 0).data(), - mVrtParams.phiCut, - mTimeFrame->getTracklets()[0], - mTimeFrame->getNTrackletsCluster(rofId, 0), - mIndexTableUtils, - rofId, - 0, - mVrtParams.maxTrackletsPerCluster); - trackleterKernelHost( - mTimeFrame->getClustersOnLayer(rofId, 2), - mTimeFrame->getClustersOnLayer(rofId, 1), - mTimeFrame->getIndexTable(rofId, 2).data(), - mVrtParams.phiCut, - mTimeFrame->getTracklets()[1], - mTimeFrame->getNTrackletsCluster(rofId, 1), - mIndexTableUtils, - rofId, - 0, - mVrtParams.maxTrackletsPerCluster); - mTimeFrame->getNTrackletsROf(rofId, 0) = std::accumulate(mTimeFrame->getNTrackletsCluster(rofId, 0).begin(), mTimeFrame->getNTrackletsCluster(rofId, 0).end(), 0); - mTimeFrame->getNTrackletsROf(rofId, 1) = std::accumulate(mTimeFrame->getNTrackletsCluster(rofId, 1).begin(), mTimeFrame->getNTrackletsCluster(rofId, 1).end(), 0); + for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { // Pivot rofId: the rof for which the tracklets are computed + int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; + int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; + std::cout << "pivoting on: " << pivotRofId << " start: " << startROF << " end: " << endROF << std::endl; + for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { + trackleterKernelHost( + mTimeFrame->getClustersOnLayer(targetRofId, 0), // Clusters to be matched with the next layer in target rof + mTimeFrame->getClustersOnLayer(pivotRofId, 1), // Clusters to be matched with the current layer in pivot rof + mTimeFrame->getUsedClusters(0), + mTimeFrame->getIndexTable(targetRofId, 0).data(), // Index table to access the data on the next layer in target rof + mVrtParams.phiCut, + mTimeFrame->getTracklets()[0], // Flat tracklet buffer + mTimeFrame->getNTrackletsCluster(pivotRofId, 0), // Span of the number of tracklets per each cluster in pivot rof + mIndexTableUtils, + pivotRofId, + targetRofId, + 0, // Offset in the tracklet buffer + mVrtParams.maxTrackletsPerCluster); + trackleterKernelHost( + mTimeFrame->getClustersOnLayer(targetRofId, 2), + mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClusters(2), + mTimeFrame->getIndexTable(targetRofId, 2).data(), + mVrtParams.phiCut, + mTimeFrame->getTracklets()[1], + mTimeFrame->getNTrackletsCluster(pivotRofId, 1), // Span of the number of tracklets per each cluster in pivot rof + mIndexTableUtils, + pivotRofId, + targetRofId, + 0, // Offset in the tracklet buffer + mVrtParams.maxTrackletsPerCluster); + mTimeFrame->getNTrackletsROf(pivotRofId, 0) += std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 0).end(), 0); + mTimeFrame->getNTrackletsROf(pivotRofId, 1) += std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 1).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 1).end(), 0); + } } #pragma omp single mTimeFrame->computeTrackletsScans(mNThreads); @@ -207,29 +222,37 @@ void VertexerTraits::computeTracklets() mTimeFrame->getTracklets()[1].resize(mTimeFrame->getTotalTrackletsTF(1)); #pragma omp for schedule(dynamic) - for (int rofId = 0; rofId < mTimeFrame->getNrof(); ++rofId) { - trackleterKernelHost( - mTimeFrame->getClustersOnLayer(rofId, 0), - mTimeFrame->getClustersOnLayer(rofId, 1), - mTimeFrame->getIndexTable(rofId, 0).data(), - mVrtParams.phiCut, - mTimeFrame->getTracklets()[0], - mTimeFrame->getNTrackletsCluster(rofId, 0), - mIndexTableUtils, - rofId, - mTimeFrame->getNTrackletsROf(rofId, 0), - mVrtParams.maxTrackletsPerCluster); - trackleterKernelHost( - mTimeFrame->getClustersOnLayer(rofId, 2), - mTimeFrame->getClustersOnLayer(rofId, 1), - mTimeFrame->getIndexTable(rofId, 2).data(), - mVrtParams.phiCut, - mTimeFrame->getTracklets()[1], - mTimeFrame->getNTrackletsCluster(rofId, 1), - mIndexTableUtils, - rofId, - mTimeFrame->getNTrackletsROf(rofId, 1), - mVrtParams.maxTrackletsPerCluster); + for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { + int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; + int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; + for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { + trackleterKernelHost( + mTimeFrame->getClustersOnLayer(targetRofId, 0), + mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClusters(0), + mTimeFrame->getIndexTable(pivotRofId, 0).data(), + mVrtParams.phiCut, + mTimeFrame->getTracklets()[0], + mTimeFrame->getNTrackletsCluster(pivotRofId, 0), + mIndexTableUtils, + pivotRofId, + targetRofId, + mTimeFrame->getNTrackletsROf(pivotRofId, 0), + mVrtParams.maxTrackletsPerCluster); + trackleterKernelHost( + mTimeFrame->getClustersOnLayer(targetRofId, 2), + mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClusters(2), + mTimeFrame->getIndexTable(targetRofId, 2).data(), + mVrtParams.phiCut, + mTimeFrame->getTracklets()[1], + mTimeFrame->getNTrackletsCluster(pivotRofId, 1), + mIndexTableUtils, + pivotRofId, + targetRofId, + mTimeFrame->getNTrackletsROf(pivotRofId, 1), + mVrtParams.maxTrackletsPerCluster); + } } } @@ -238,7 +261,7 @@ void VertexerTraits::computeTracklets() for (auto& trk : mTimeFrame->getTracklets()[0]) { MCCompLabel label; int sortedId0{mTimeFrame->getSortedIndex(trk.rof[0], 0, trk.firstClusterIndex)}; - int sortedId1{mTimeFrame->getSortedIndex(trk.rof[0], 1, trk.secondClusterIndex)}; + int sortedId1{mTimeFrame->getSortedIndex(trk.rof[1], 1, trk.secondClusterIndex)}; for (auto& lab0 : mTimeFrame->getClusterLabels(0, mTimeFrame->getClusters()[0][sortedId0].clusterId)) { for (auto& lab1 : mTimeFrame->getClusterLabels(1, mTimeFrame->getClusters()[1][sortedId1].clusterId)) { if (lab0 == lab1 && lab0.isValid()) { From f47789b3a3bf0b053068920661a261bed0b79435 Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Thu, 25 Apr 2024 16:33:37 +0200 Subject: [PATCH 2/6] TrackletFinder extended, slection is still bugged --- .../tracking/include/ITStracking/TimeFrame.h | 106 +++++++++--------- .../tracking/include/ITStracking/Tracklet.h | 7 +- .../tracking/include/ITStracking/Vertexer.h | 5 +- .../ITSMFT/ITS/tracking/src/TimeFrame.cxx | 44 ++++---- Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx | 8 +- .../ITSMFT/ITS/tracking/src/Vertexer.cxx | 47 +++++--- .../ITS/tracking/src/VertexerTraits.cxx | 96 ++++++++-------- .../ITS3/reconstruction/src/IOUtils.cxx | 2 +- 8 files changed, 173 insertions(+), 142 deletions(-) diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h index be07f3cd8d351..33aa68aaf2438 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h @@ -78,16 +78,16 @@ class TimeFrame friend class TimeFrameGPU; TimeFrame(int nLayers = 7); const Vertex& getPrimaryVertex(const int) const; - gsl::span getPrimaryVertices(int tf) const; + gsl::span getPrimaryVertices(int rof) const; gsl::span getPrimaryVertices(int romin, int romax) const; gsl::span> getPrimaryVerticesLabels(const int rof) const; - gsl::span> getPrimaryVerticesXAlpha(int tf) const; + gsl::span> getPrimaryVerticesXAlpha(int rof) const; void fillPrimaryVerticesXandAlpha(); int getPrimaryVerticesNum(int rofID = -1) const; void addPrimaryVertices(const std::vector& vertices); void addPrimaryVertices(const gsl::span& vertices); void addPrimaryVertices(const std::vector&); - void removePrimaryVerticesInROf(const int rofId); + void removePrimaryVerticesInROF(const int rofId); int loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span clusters, const dataformats::MCTruthContainer* mcLabels = nullptr); @@ -124,8 +124,8 @@ class TimeFrame gsl::span getClustersOnLayer(int rofId, int layerId) const; gsl::span getClustersPerROFrange(int rofMin, int range, int layerId) const; gsl::span getUnsortedClustersOnLayer(int rofId, int layerId) const; - gsl::span getROframesClustersPerROFrange(int rofMin, int range, int layerId) const; - gsl::span getROframeClusters(int layerId) const; + gsl::span getROFramesClustersPerROFrange(int rofMin, int range, int layerId) const; + gsl::span getROFrameClusters(int layerId) const; gsl::span getNClustersROFrange(int rofMin, int range, int layerId) const; gsl::span getIndexTablePerROFrange(int rofMin, int range, int layerId) const; gsl::span getIndexTable(int rofId, int layerId); @@ -146,7 +146,7 @@ class TimeFrame void resetRofPV() { mPrimaryVertices.clear(); - mROframesPV.resize(1, 0); + mROFramesPV.resize(1, 0); }; bool isClusterUsed(int layer, int clusterId) const; @@ -179,15 +179,19 @@ class TimeFrame bool checkMemory(unsigned long max) { return getArtefactsMemory() < max; } unsigned long getArtefactsMemory(); - int getROfCutClusterMult() const { return mCutClusterMult; }; - int getROfCutVertexMult() const { return mCutVertexMult; }; - int getROfCutAllMult() const { return mCutClusterMult + mCutVertexMult; } + int getROFCutClusterMult() const { return mCutClusterMult; }; + int getROFCutVertexMult() const { return mCutVertexMult; }; + int getROFCutAllMult() const { return mCutClusterMult + mCutVertexMult; } // Vertexer void computeTrackletsScans(const int nThreads = 1); - int& getNTrackletsROf(int tf, int combId); - std::vector& getLines(int tf); - std::vector& getTrackletClusters(int tf); + int& getNTrackletsROF(int rof, int combId); + std::vector& getLines(int rof); + int getNLinesTotal() const + { + return std::accumulate(mLines.begin(), mLines.end(), 0, [](int sum, const auto& l) { return sum + l.size(); }); + } + std::vector& getTrackletClusters(int rof); gsl::span getFoundTracklets(int rofId, int combId) const; gsl::span getFoundTracklets(int rofId, int combId); gsl::span getLabelsFoundTracklets(int rofId, int combId) const; @@ -250,15 +254,15 @@ class TimeFrame std::vector> mClusters; std::vector> mTrackingFrameInfo; std::vector> mClusterExternalIndices; - std::vector> mROframesClusters; + std::vector> mROFramesClusters; const dataformats::MCTruthContainer* mClusterLabels = nullptr; - std::array, 2> mNTrackletsPerCluster; // TODO: remove in favour of mNTrackletsPerROf + std::array, 2> mNTrackletsPerCluster; // TODO: remove in favour of mNTrackletsPerROF std::vector> mNClustersPerROF; std::vector> mIndexTables; std::vector> mTrackletsLookupTable; std::vector> mUsedClusters; int mNrof = 0; - std::vector mROframesPV = {0}; + std::vector mROFramesPV = {0}; std::vector mPrimaryVertices; // State if memory will be externally managed. @@ -306,10 +310,10 @@ class TimeFrame int mCutVertexMult; // Vertexer - std::vector> mNTrackletsPerROf; + std::vector> mNTrackletsPerROF; std::vector> mLines; std::vector> mTrackletClusters; - std::vector> mTrackletsIndexROf; + std::vector> mTrackletsIndexROF; std::vector> mLinesLabels; std::vector> mVerticesLabels; std::array mTotalTracklets = {0, 0}; @@ -320,43 +324,43 @@ inline const Vertex& TimeFrame::getPrimaryVertex(const int vertexIndex) const { inline gsl::span TimeFrame::getPrimaryVertices(int rof) const { - const int start = mROframesPV[rof]; + const int start = mROFramesPV[rof]; const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROframesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&mPrimaryVertices[start], static_cast::size_type>(delta)}; } inline gsl::span> TimeFrame::getPrimaryVerticesLabels(const int rof) const { - const int start = mROframesPV[rof]; + const int start = mROFramesPV[rof]; const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROframesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&mVerticesLabels[start], static_cast::size_type>(delta)}; } inline gsl::span TimeFrame::getPrimaryVertices(int romin, int romax) const { - return {&mPrimaryVertices[mROframesPV[romin]], static_cast::size_type>(mROframesPV[romax + 1] - mROframesPV[romin])}; + return {&mPrimaryVertices[mROFramesPV[romin]], static_cast::size_type>(mROFramesPV[romax + 1] - mROFramesPV[romin])}; } inline gsl::span> TimeFrame::getPrimaryVerticesXAlpha(int rof) const { - const int start = mROframesPV[rof]; + const int start = mROFramesPV[rof]; const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROframesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&(mPValphaX[start]), static_cast>::size_type>(delta)}; } inline int TimeFrame::getPrimaryVerticesNum(int rofID) const { - return rofID < 0 ? mPrimaryVertices.size() : mROframesPV[rofID + 1] - mROframesPV[rofID]; + return rofID < 0 ? mPrimaryVertices.size() : mROFramesPV[rofID + 1] - mROFramesPV[rofID]; } inline bool TimeFrame::empty() const { return getTotalClusters() == 0; } -inline int TimeFrame::getSortedIndex(int rof, int layer, int index) const { return mROframesClusters[layer][rof] + index; } +inline int TimeFrame::getSortedIndex(int rof, int layer, int index) const { return mROFramesClusters[layer][rof] + index; } -inline int TimeFrame::getSortedStartIndex(const int rof, const int layer) const { return mROframesClusters[layer][rof]; } +inline int TimeFrame::getSortedStartIndex(const int rof, const int layer) const { return mROFramesClusters[layer][rof]; } inline int TimeFrame::getNrof() const { return mNrof; } @@ -371,9 +375,9 @@ inline float TimeFrame::getBeamX() const { return mBeamPos[0]; } inline float TimeFrame::getBeamY() const { return mBeamPos[1]; } -inline gsl::span TimeFrame::getROframeClusters(int layerId) const +inline gsl::span TimeFrame::getROFrameClusters(int layerId) const { - return {&mROframesClusters[layerId][0], static_cast::size_type>(mROframesClusters[layerId].size())}; + return {&mROFramesClusters[layerId][0], static_cast::size_type>(mROFramesClusters[layerId].size())}; } inline gsl::span TimeFrame::getClustersOnLayer(int rofId, int layerId) @@ -381,8 +385,8 @@ inline gsl::span TimeFrame::getClustersOnLayer(int rofId, int layerId) if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - int startIdx{mROframesClusters[layerId][rofId]}; - return {&mClusters[layerId][startIdx], static_cast::size_type>(mROframesClusters[layerId][rofId + 1] - startIdx)}; + int startIdx{mROFramesClusters[layerId][rofId]}; + return {&mClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; } inline gsl::span TimeFrame::getClustersOnLayer(int rofId, int layerId) const @@ -390,8 +394,8 @@ inline gsl::span TimeFrame::getClustersOnLayer(int rofId, int lay if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - int startIdx{mROframesClusters[layerId][rofId]}; - return {&mClusters[layerId][startIdx], static_cast::size_type>(mROframesClusters[layerId][rofId + 1] - startIdx)}; + int startIdx{mROFramesClusters[layerId][rofId]}; + return {&mClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; } inline gsl::span TimeFrame::getClustersPerROFrange(int rofMin, int range, int layerId) const @@ -399,15 +403,15 @@ inline gsl::span TimeFrame::getClustersPerROFrange(int rofMin, in if (rofMin < 0 || rofMin >= mNrof) { return gsl::span(); } - int startIdx{mROframesClusters[layerId][rofMin]}; // First cluster of rofMin - int endIdx{mROframesClusters[layerId][std::min(rofMin + range, mNrof)]}; + int startIdx{mROFramesClusters[layerId][rofMin]}; // First cluster of rofMin + int endIdx{mROFramesClusters[layerId][std::min(rofMin + range, mNrof)]}; return {&mClusters[layerId][startIdx], static_cast::size_type>(endIdx - startIdx)}; } -inline gsl::span TimeFrame::getROframesClustersPerROFrange(int rofMin, int range, int layerId) const +inline gsl::span TimeFrame::getROFramesClustersPerROFrange(int rofMin, int range, int layerId) const { int chkdRange{std::min(range, mNrof - rofMin)}; - return {&mROframesClusters[layerId][rofMin], static_cast::size_type>(chkdRange)}; + return {&mROFramesClusters[layerId][rofMin], static_cast::size_type>(chkdRange)}; } inline gsl::span TimeFrame::getNClustersROFrange(int rofMin, int range, int layerId) const @@ -420,7 +424,7 @@ inline int TimeFrame::getTotalClustersPerROFrange(int rofMin, int range, int lay { int startIdx{rofMin}; // First cluster of rofMin int endIdx{std::min(rofMin + range, mNrof)}; - return mROframesClusters[layerId][endIdx] - mROframesClusters[layerId][startIdx]; + return mROFramesClusters[layerId][endIdx] - mROFramesClusters[layerId][startIdx]; } inline gsl::span TimeFrame::getIndexTablePerROFrange(int rofMin, int range, int layerId) const @@ -432,7 +436,7 @@ inline gsl::span TimeFrame::getIndexTablePerROFrange(int rofMin, int inline int TimeFrame::getClusterROF(int iLayer, int iCluster) { - return std::lower_bound(mROframesClusters[iLayer].begin(), mROframesClusters[iLayer].end(), iCluster + 1) - mROframesClusters[iLayer].begin() - 1; + return std::lower_bound(mROFramesClusters[iLayer].begin(), mROFramesClusters[iLayer].end(), iCluster + 1) - mROFramesClusters[iLayer].begin() - 1; } inline gsl::span TimeFrame::getUnsortedClustersOnLayer(int rofId, int layerId) const @@ -440,8 +444,8 @@ inline gsl::span TimeFrame::getUnsortedClustersOnLayer(int rofId, if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - int startIdx{mROframesClusters[layerId][rofId]}; - return {&mUnsortedClusters[layerId][startIdx], static_cast::size_type>(mROframesClusters[layerId][rofId + 1] - startIdx)}; + int startIdx{mROFramesClusters[layerId][rofId]}; + return {&mUnsortedClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; } inline const std::vector& TimeFrame::getTrackingFrameInfoOnLayer(int layerId) const @@ -559,13 +563,13 @@ inline gsl::span TimeFrame::getNTrackletsCluster(int rofId, int combId) if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - auto startIdx{mROframesClusters[1][rofId]}; - return {&mNTrackletsPerCluster[combId][startIdx], static_cast::size_type>(mROframesClusters[1][rofId + 1] - startIdx)}; + auto startIdx{mROFramesClusters[1][rofId]}; + return {&mNTrackletsPerCluster[combId][startIdx], static_cast::size_type>(mROFramesClusters[1][rofId + 1] - startIdx)}; } -inline int& TimeFrame::getNTrackletsROf(int rof, int combId) +inline int& TimeFrame::getNTrackletsROF(int rof, int combId) { - return mNTrackletsPerROf[combId][rof]; + return mNTrackletsPerROF[combId][rof]; } inline bool TimeFrame::isRoadFake(int i) const @@ -600,8 +604,8 @@ inline gsl::span TimeFrame::getFoundTracklets(int rofId, int combId) if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - auto startIdx{mNTrackletsPerROf[combId][rofId]}; - return {&mTracklets[combId][startIdx], static_cast::size_type>(mNTrackletsPerROf[combId][rofId + 1] - startIdx)}; + auto startIdx{mNTrackletsPerROF[combId][rofId]}; + return {&mTracklets[combId][startIdx], static_cast::size_type>(mNTrackletsPerROF[combId][rofId + 1] - startIdx)}; } inline gsl::span TimeFrame::getFoundTracklets(int rofId, int combId) const @@ -609,8 +613,8 @@ inline gsl::span TimeFrame::getFoundTracklets(int rofId, int com if (rofId < 0 || rofId >= mNrof) { return gsl::span(); } - auto startIdx{mNTrackletsPerROf[combId][rofId]}; - return {&mTracklets[combId][startIdx], static_cast::size_type>(mNTrackletsPerROf[combId][rofId + 1] - startIdx)}; + auto startIdx{mNTrackletsPerROF[combId][rofId]}; + return {&mTracklets[combId][startIdx], static_cast::size_type>(mNTrackletsPerROF[combId][rofId + 1] - startIdx)}; } inline gsl::span TimeFrame::getLabelsFoundTracklets(int rofId, int combId) const @@ -618,8 +622,8 @@ inline gsl::span TimeFrame::getLabelsFoundTracklets(int rofId if (rofId < 0 || rofId >= mNrof || !hasMCinformation()) { return gsl::span(); } - auto startIdx{mNTrackletsPerROf[combId][rofId]}; - return {&mTrackletLabels[combId][startIdx], static_cast::size_type>(mNTrackletsPerROf[combId][rofId + 1] - startIdx)}; + auto startIdx{mNTrackletsPerROF[combId][rofId]}; + return {&mTrackletLabels[combId][startIdx], static_cast::size_type>(mNTrackletsPerROF[combId][rofId + 1] - startIdx)}; } inline int TimeFrame::getNumberOfClusters() const diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h index a7f7adbb8abad..f37dfb67da6ef 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h @@ -36,6 +36,7 @@ struct Tracklet final { { return firstClusterIndex < 0 || secondClusterIndex < 0; } + GPUhdi() auto getDeltaRof() const { return rof[1] - rof[0]; } GPUhdi() void dump(); GPUhdi() void dump() const; GPUhdi() void dump(const int, const int); @@ -78,7 +79,7 @@ GPUhdi() Tracklet::Tracklet(const int idx0, const int idx1, float tanL, float ph // Nothing to do } -GPUhdi() bool Tracklet::operator==(const Tracklet& rhs) const +GPUhdi() bool Tracklet::operator==(const Tracklet & rhs) const { return this->firstClusterIndex == rhs.firstClusterIndex && this->secondClusterIndex == rhs.secondClusterIndex && @@ -88,7 +89,7 @@ GPUhdi() bool Tracklet::operator==(const Tracklet& rhs) const this->rof[1] == rhs.rof[1]; } -GPUhdi() bool Tracklet::operator!=(const Tracklet& rhs) const +GPUhdi() bool Tracklet::operator!=(const Tracklet & rhs) const { return this->firstClusterIndex != rhs.firstClusterIndex || this->secondClusterIndex != rhs.secondClusterIndex || @@ -96,7 +97,7 @@ GPUhdi() bool Tracklet::operator!=(const Tracklet& rhs) const this->phi != rhs.phi; } -GPUhdi() unsigned char Tracklet::operator<(const Tracklet& t) const +GPUhdi() unsigned char Tracklet::operator<(const Tracklet & t) const { if (isEmpty()) { return false; diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h index eb46c60305bda..544eb872042f7 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h @@ -88,7 +88,10 @@ class Vertexer void dumpTraits(); template float evaluateTask(void (Vertexer::*)(T...), const char*, std::function logger, T&&... args); - void printEpilog(std::function logger, const float total); + void printEpilog(std::function logger, + bool isHybrid, + const unsigned int trackletN01, const unsigned int trackletN12, const unsigned selectedN, + const float initT, const float trackletT, const float selecT, const float vertexT); private: std::uint32_t mTimeFrameCounter = 0; diff --git a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx index 6e7114d9ca54d..f314a3871b053 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx @@ -73,9 +73,9 @@ TimeFrame::TimeFrame(int nLayers) mTrackingFrameInfo.resize(nLayers); mClusterExternalIndices.resize(nLayers); mUsedClusters.resize(nLayers); - mROframesClusters.resize(nLayers, {0}); /// TBC: if resetting the timeframe is required, then this has to be done + mROFramesClusters.resize(nLayers, {0}); /// TBC: if resetting the timeframe is required, then this has to be done mNClustersPerROF.resize(nLayers); - mTrackletsIndexROf.resize(2, {0}); + mTrackletsIndexROF.resize(2, {0}); } void TimeFrame::addPrimaryVertices(const std::vector& vertices) @@ -89,7 +89,7 @@ void TimeFrame::addPrimaryVertices(const std::vector& vertices) mBeamPosWeight += w; } } - mROframesPV.push_back(mPrimaryVertices.size()); + mROFramesPV.push_back(mPrimaryVertices.size()); } void TimeFrame::addPrimaryVertices(const gsl::span& vertices) @@ -103,7 +103,7 @@ void TimeFrame::addPrimaryVertices(const gsl::span& vertices) mBeamPosWeight += w; } } - mROframesPV.push_back(mPrimaryVertices.size()); + mROFramesPV.push_back(mPrimaryVertices.size()); } void TimeFrame::addPrimaryVertices(const std::vector& lVertices) @@ -141,10 +141,10 @@ int TimeFrame::loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span rofs, deepVectorClear(mUnsortedClusters[iLayer]); deepVectorClear(mTrackingFrameInfo[iLayer]); deepVectorClear(mClusterExternalIndices[iLayer]); - mROframesClusters[iLayer].resize(1, 0); + mROFramesClusters[iLayer].resize(1, 0); if (iLayer < 2) { - deepVectorClear(mTrackletsIndexROf[iLayer]); + deepVectorClear(mTrackletsIndexROF[iLayer]); } } @@ -223,8 +223,8 @@ int TimeFrame::loadROFrameData(gsl::span rofs, addClusterExternalIndexToLayer(layer, clusterId); } for (unsigned int iL{0}; iL < mUnsortedClusters.size(); ++iL) { - // mNClustersPerROF[iL].push_back(mUnsortedClusters[iL].size() - mROframesClusters[iL].back()); - mROframesClusters[iL].push_back(mUnsortedClusters[iL].size()); + // mNClustersPerROF[iL].push_back(mUnsortedClusters[iL].size() - mROFramesClusters[iL].back()); + mROFramesClusters[iL].push_back(mUnsortedClusters[iL].size()); } mNrof++; } @@ -285,8 +285,8 @@ void TimeFrame::initialise(const int iteration, const TrackingParameters& trkPar mIndexTables.resize(mClusters.size(), std::vector(mNrof * (trkParam.ZBins * trkParam.PhiBins + 1), 0)); mLines.resize(mNrof); mTrackletClusters.resize(mNrof); - mNTrackletsPerROf.resize(2); - for (auto& v : mNTrackletsPerROf) { + mNTrackletsPerROF.resize(2); + for (auto& v : mNTrackletsPerROF) { v = std::vector(mNrof + 1, 0); } @@ -474,9 +474,9 @@ void TimeFrame::resizeVectors(int nLayers) mTrackingFrameInfo.resize(nLayers); mClusterExternalIndices.resize(nLayers); mUsedClusters.resize(nLayers); - mROframesClusters.resize(nLayers, {0}); + mROFramesClusters.resize(nLayers, {0}); mNClustersPerROF.resize(nLayers); - mTrackletsIndexROf.resize(2, {0}); + mTrackletsIndexROF.resize(2, {0}); } void TimeFrame::printTrackletLUTonLayer(int i) @@ -517,9 +517,9 @@ void TimeFrame::printCellLUTs() void TimeFrame::printVertices() { - std::cout << "Vertices in ROF (nROF = " << mNrof << ", lut size = " << mROframesPV.size() << ")" << std::endl; - for (unsigned int iR{0}; iR < mROframesPV.size(); ++iR) { - std::cout << mROframesPV[iR] << "\t"; + std::cout << "Vertices in ROF (nROF = " << mNrof << ", lut size = " << mROFramesPV.size() << ")" << std::endl; + for (unsigned int iR{0}; iR < mROFramesPV.size(); ++iR) { + std::cout << mROFramesPV[iR] << "\t"; } std::cout << "\n\n Vertices:" << std::endl; for (unsigned int iV{0}; iV < mPrimaryVertices.size(); ++iV) { @@ -531,9 +531,9 @@ void TimeFrame::printVertices() void TimeFrame::printROFoffsets() { std::cout << "--------" << std::endl; - for (unsigned int iLayer{0}; iLayer < mROframesClusters.size(); ++iLayer) { + for (unsigned int iLayer{0}; iLayer < mROFramesClusters.size(); ++iLayer) { std::cout << "Layer " << iLayer << std::endl; - for (auto value : mROframesClusters[iLayer]) { + for (auto value : mROFramesClusters[iLayer]) { std::cout << value << "\t"; } std::cout << std::endl; @@ -556,8 +556,8 @@ void TimeFrame::computeTrackletsScans(const int nThreads) { // #pragma omp parallel for num_threads(nThreads > 1 ? 2 : nThreads) for (ushort iLayer = 0; iLayer < 2; ++iLayer) { - mTotalTracklets[iLayer] = std::accumulate(mNTrackletsPerROf[iLayer].begin(), mNTrackletsPerROf[iLayer].end(), 0); - std::exclusive_scan(mNTrackletsPerROf[iLayer].begin(), mNTrackletsPerROf[iLayer].end(), mNTrackletsPerROf[iLayer].begin(), 0); + mTotalTracklets[iLayer] = std::accumulate(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), 0); + std::exclusive_scan(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), mNTrackletsPerROF[iLayer].begin(), 0); } } diff --git a/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx b/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx index 42df15b24f052..c9d5dc9451f22 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx @@ -102,10 +102,10 @@ void Tracker::clustersToTracks(std::function logger, std::f } iVertex++; } while (iVertex < maxNvertices); - logger(fmt::format("- Tracklet finding: {} tracklets found in {:.2f} ms", nTracklets, timeTracklets)); - logger(fmt::format("- Cell finding: {} cells found in {:.2f} ms", nCells, timeCells)); - logger(fmt::format("- Neighbours finding: {} neighbours found in {:.2f} ms", nNeighbours, timeNeighbours)); - logger(fmt::format("- Track finding: {} tracks found in {:.2f} ms", nTracks + mTimeFrame->getNumberOfTracks(), timeRoads)); + logger(fmt::format(" - Tracklet finding: {} tracklets found in {:.2f} ms", nTracklets, timeTracklets)); + logger(fmt::format(" - Cell finding: {} cells found in {:.2f} ms", nCells, timeCells)); + logger(fmt::format(" - Neighbours finding: {} neighbours found in {:.2f} ms", nNeighbours, timeNeighbours)); + logger(fmt::format(" - Track finding: {} tracks found in {:.2f} ms", nTracks + mTimeFrame->getNumberOfTracks(), timeRoads)); total += timeTracklets + timeCells + timeNeighbours + timeRoads; total += evaluateTask(&Tracker::extendTracks, "Extending tracks", logger, iteration); } diff --git a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx index 8a4a2cf7664a2..40e56f3b42f1a 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx @@ -38,30 +38,40 @@ Vertexer::Vertexer(VertexerTraits* traits) float Vertexer::clustersToVertices(std::function logger) { - float total{0.f}; + float timeTracklet, timeSelection, timeVertexing, timeInit; + unsigned int nTracklets01, nTracklets12, nSelectedTracklets; TrackingParameters trkPars; trkPars.PhiBins = mTraits->getVertexingParameters().PhiBins; trkPars.ZBins = mTraits->getVertexingParameters().ZBins; - total += evaluateTask(&Vertexer::initialiseVertexer, "Vertexer initialisation", logger, trkPars); - total += evaluateTask(&Vertexer::findTracklets, "Vertexer tracklet finding", logger); - total += evaluateTask(&Vertexer::validateTracklets, "Vertexer adjacent tracklets validation", logger); - total += evaluateTask(&Vertexer::findVertices, "Vertexer vertex finding", logger); - printEpilog(logger, total); - return total; + timeInit += evaluateTask( + &Vertexer::initialiseVertexer, "Vertexer initialisation", [](std::string) {}, trkPars); + timeTracklet = evaluateTask(&Vertexer::findTracklets, "Vertexer tracklet finding", [](std::string) {}); + nTracklets01 = mTimeFrame->getTotalTrackletsTF(0); + nTracklets12 = mTimeFrame->getTotalTrackletsTF(1); + timeSelection = evaluateTask(&Vertexer::validateTracklets, "Vertexer adjacent tracklets validation", [](std::string) {}); + nSelectedTracklets = mTimeFrame->getNLinesTotal(); + timeVertexing = evaluateTask(&Vertexer::findVertices, "Vertexer vertex finding", [](std::string) {}); + printEpilog(logger, false, nTracklets01, nTracklets12, nSelectedTracklets, timeInit, timeTracklet, timeSelection, timeVertexing); + return timeInit + timeTracklet + timeSelection + timeVertexing; } float Vertexer::clustersToVerticesHybrid(std::function logger) { - float total{0.f}; + float timeTracklet, timeSelection, timeVertexing, timeInit; + unsigned int nTracklets01, nTracklets12, nSelectedTracklets; TrackingParameters trkPars; trkPars.PhiBins = mTraits->getVertexingParameters().PhiBins; trkPars.ZBins = mTraits->getVertexingParameters().ZBins; - total += evaluateTask(&Vertexer::initialiseVertexerHybrid, "Hybrid Vertexer initialisation", logger, trkPars); - total += evaluateTask(&Vertexer::findTrackletsHybrid, "Hybrid Vertexer tracklet finding", logger); - total += evaluateTask(&Vertexer::validateTrackletsHybrid, "Hybrid Vertexer adjacent tracklets validation", logger); - total += evaluateTask(&Vertexer::findVerticesHybrid, "Hybrid Vertexer vertex finding", logger); - printEpilog(logger, total); - return total; + timeInit += evaluateTask( + &Vertexer::initialiseVertexerHybrid, "Hybrid Vertexer initialisation", [](std::string) {}, trkPars); + timeTracklet = evaluateTask(&Vertexer::findTrackletsHybrid, "Hybrid Vertexer tracklet finding", [](std::string) {}); + nTracklets01 = mTimeFrame->getTotalTrackletsTF(0); + nTracklets12 = mTimeFrame->getTotalTrackletsTF(1); + timeSelection = evaluateTask(&Vertexer::validateTrackletsHybrid, "Hybrid Vertexer adjacent tracklets validation", [](std::string) {}); + nSelectedTracklets = mTimeFrame->getNLinesTotal(); + timeVertexing = evaluateTask(&Vertexer::findVerticesHybrid, "Hybrid Vertexer vertex finding", [](std::string) {}); + printEpilog(logger, false, nTracklets01, nTracklets12, nSelectedTracklets, timeInit, timeTracklet, timeSelection, timeVertexing); + return timeInit + timeTracklet + timeSelection + timeVertexing; } void Vertexer::getGlobalConfiguration() @@ -103,8 +113,15 @@ void Vertexer::adoptTimeFrame(TimeFrame& tf) mTraits->adoptTimeFrame(&tf); } -void Vertexer::printEpilog(std::function logger, const float total) +void Vertexer::printEpilog(std::function logger, + bool isHybrid, + const unsigned int trackletN01, const unsigned int trackletN12, const unsigned selectedN, + const float initT, const float trackletT, const float selecT, const float vertexT) { + float total = initT + trackletT + selecT + vertexT; + logger(fmt::format(" - {}Vertexer: trackleting found {}|{} tracklets, completed in: {} ms", isHybrid ? "Hybrid" : "", trackletN01, trackletN12, trackletT)); + logger(fmt::format(" - {}Vertexer: selected {} tracklets, completed in: {} ms", isHybrid ? "Hybrid" : "", selectedN, selecT)); + logger(fmt::format(" - {}Vertexer: vertexing completed in: {} ms", isHybrid ? "Hybrid" : "", vertexT)); logger(fmt::format(" - Timeframe {} vertexing completed in: {} ms, using {} thread(s).", mTimeFrameCounter++, total, mTraits->getNThreads())); } diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx index 3cc65fcf2c768..6f26a4109fd65 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx @@ -19,7 +19,7 @@ #include "ITStracking/VertexerTraits.h" #include "ITStracking/ClusterLines.h" #include "ITStracking/Tracklet.h" - +#define VTX_DEBUG #ifdef VTX_DEBUG #include "TTree.h" #include "TFile.h" @@ -48,7 +48,6 @@ template void trackleterKernelHost( const gsl::span& clustersNextLayer, // 0 2 const gsl::span& clustersCurrentLayer, // 1 1 - gsl::span usedClusters, int* indexTableNext, const float phiCut, std::vector& tracklets, @@ -56,13 +55,13 @@ void trackleterKernelHost( const IndexTableUtils& utils, const int pivotRof, const int targetRof, - const int rofFoundTrackletsOffset, + int* rofFoundTrackletsOffset, // we want to change that, to keep track of the offset in deltaRof>0 const int maxTrackletsPerCluster = static_cast(2e3)) { const int PhiBins{utils.getNphiBins()}; const int ZBins{utils.getNzBins()}; // loop on layer1 clusters - int cumulativeStoredTracklets{0}; + // int cumulativeStoredTracklets{rofFoundTrackletsOffset}; for (int iCurrentLayerClusterIndex = 0; iCurrentLayerClusterIndex < clustersCurrentLayer.size(); ++iCurrentLayerClusterIndex) { int storedTracklets{0}; const Cluster& currentCluster{clustersCurrentLayer[iCurrentLayerClusterIndex]}; @@ -80,17 +79,13 @@ void trackleterKernelHost( // loop on clusters next layer for (int iNextLayerClusterIndex{firstRowClusterIndex}; iNextLayerClusterIndex < maxRowClusterIndex && iNextLayerClusterIndex < static_cast(clustersNextLayer.size()); ++iNextLayerClusterIndex) { const Cluster& nextCluster{clustersNextLayer[iNextLayerClusterIndex]}; - if (usedClusters[nextCluster.clusterId]) { - continue; - } if (o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(currentCluster.phi, nextCluster.phi)) < phiCut) { if (storedTracklets < maxTrackletsPerCluster) { if constexpr (!EvalRun) { - usedClusters[nextCluster.clusterId] = true; if constexpr (Mode == TrackletMode::Layer0Layer1) { - tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof}; + tracklets[*rofFoundTrackletsOffset + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof}; } else { - tracklets[rofFoundTrackletsOffset + cumulativeStoredTracklets + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof}; + tracklets[*rofFoundTrackletsOffset + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof}; } } ++storedTracklets; @@ -100,9 +95,9 @@ void trackleterKernelHost( } } if constexpr (EvalRun) { - foundTracklets[iCurrentLayerClusterIndex] = storedTracklets; + foundTracklets[iCurrentLayerClusterIndex] += storedTracklets; } else { - cumulativeStoredTracklets += storedTracklets; + *rofFoundTrackletsOffset += storedTracklets; } } } @@ -112,26 +107,32 @@ void trackletSelectionKernelHost( const gsl::span clusters1, // 1 const gsl::span& tracklets01, const gsl::span& tracklets12, + std::vector& usedTracklets, const gsl::span foundTracklets01, const gsl::span foundTracklets12, - std::vector& destTracklets, + std::vector& lines, const gsl::span& trackletLabels, std::vector& linesLabels, + const int rofDist = 0, // signed rof distance to consider for the tracklet to be promoted, so to be sure taht clusters0 and clusters1 are consistent with tracklet indices const float tanLambdaCut = 0.025f, const float phiCut = 0.005f, const int maxTracklets = static_cast(1e2)) { int offset01{0}, offset12{0}; - std::vector usedTracklets(tracklets01.size(), false); for (unsigned int iCurrentLayerClusterIndex{0}; iCurrentLayerClusterIndex < clusters1.size(); ++iCurrentLayerClusterIndex) { int validTracklets{0}; for (int iTracklet12{offset12}; iTracklet12 < offset12 + foundTracklets12[iCurrentLayerClusterIndex]; ++iTracklet12) { for (int iTracklet01{offset01}; iTracklet01 < offset01 + foundTracklets01[iCurrentLayerClusterIndex]; ++iTracklet01) { - const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(tracklets01[iTracklet01].tanLambda - tracklets12[iTracklet12].tanLambda)}; - const float deltaPhi{o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(tracklets01[iTracklet01].phi, tracklets12[iTracklet12].phi))}; + const auto& tracklet01{tracklets01[iTracklet01]}; + if (tracklet01.getDeltaRof() != rofDist) { + continue; + } + const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(tracklet01.tanLambda - tracklets12[iTracklet12].tanLambda)}; + const float deltaPhi{o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(tracklet01.phi, tracklets12[iTracklet12].phi))}; if (!usedTracklets[iTracklet01] && deltaTanLambda < tanLambdaCut && deltaPhi < phiCut && validTracklets != maxTracklets) { usedTracklets[iTracklet01] = true; - destTracklets.emplace_back(tracklets01[iTracklet01], clusters0.data(), clusters1.data()); + LOGP(info, "\t validated a tracklet over {}/{} with deltarof={}", tracklets01.size(), tracklets12.size(), tracklet01.getDeltaRof()); + lines.emplace_back(tracklet01, clusters0.data(), clusters1.data()); if (trackletLabels.size()) { linesLabels.emplace_back(trackletLabels[iTracklet01]); } @@ -182,12 +183,10 @@ void VertexerTraits::computeTracklets() for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { // Pivot rofId: the rof for which the tracklets are computed int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; - std::cout << "pivoting on: " << pivotRofId << " start: " << startROF << " end: " << endROF << std::endl; for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { trackleterKernelHost( - mTimeFrame->getClustersOnLayer(targetRofId, 0), // Clusters to be matched with the next layer in target rof - mTimeFrame->getClustersOnLayer(pivotRofId, 1), // Clusters to be matched with the current layer in pivot rof - mTimeFrame->getUsedClusters(0), + mTimeFrame->getClustersOnLayer(targetRofId, 0), // Clusters to be matched with the next layer in target rof + mTimeFrame->getClustersOnLayer(pivotRofId, 1), // Clusters to be matched with the current layer in pivot rof mTimeFrame->getIndexTable(targetRofId, 0).data(), // Index table to access the data on the next layer in target rof mVrtParams.phiCut, mTimeFrame->getTracklets()[0], // Flat tracklet buffer @@ -195,12 +194,11 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - 0, // Offset in the tracklet buffer + nullptr, // Offset in the tracklet buffer mVrtParams.maxTrackletsPerCluster); trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), mTimeFrame->getClustersOnLayer(pivotRofId, 1), - mTimeFrame->getUsedClusters(2), mTimeFrame->getIndexTable(targetRofId, 2).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[1], @@ -208,11 +206,11 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - 0, // Offset in the tracklet buffer + nullptr, // Offset in the tracklet buffer mVrtParams.maxTrackletsPerCluster); - mTimeFrame->getNTrackletsROf(pivotRofId, 0) += std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 0).end(), 0); - mTimeFrame->getNTrackletsROf(pivotRofId, 1) += std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 1).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 1).end(), 0); } + mTimeFrame->getNTrackletsROF(pivotRofId, 0) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 0).end(), 0); + mTimeFrame->getNTrackletsROF(pivotRofId, 1) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 1).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 1).end(), 0); } #pragma omp single mTimeFrame->computeTrackletsScans(mNThreads); @@ -225,24 +223,24 @@ void VertexerTraits::computeTracklets() for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; + auto mobileOffset0 = mTimeFrame->getNTrackletsROF(pivotRofId, 0); + auto mobileOffset1 = mTimeFrame->getNTrackletsROF(pivotRofId, 1); for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), mTimeFrame->getClustersOnLayer(pivotRofId, 1), - mTimeFrame->getUsedClusters(0), - mTimeFrame->getIndexTable(pivotRofId, 0).data(), + mTimeFrame->getIndexTable(targetRofId, 0).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[0], mTimeFrame->getNTrackletsCluster(pivotRofId, 0), mIndexTableUtils, pivotRofId, targetRofId, - mTimeFrame->getNTrackletsROf(pivotRofId, 0), + &mobileOffset0, mVrtParams.maxTrackletsPerCluster); trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), mTimeFrame->getClustersOnLayer(pivotRofId, 1), - mTimeFrame->getUsedClusters(2), mTimeFrame->getIndexTable(targetRofId, 2).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[1], @@ -250,7 +248,7 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - mTimeFrame->getNTrackletsROf(pivotRofId, 1), + &mobileOffset1, mVrtParams.maxTrackletsPerCluster); } } @@ -318,20 +316,28 @@ void VertexerTraits::computeTracklets() void VertexerTraits::computeTrackletMatching() { #pragma omp parallel for num_threads(mNThreads) schedule(dynamic) - for (int rofId = 0; rofId < mTimeFrame->getNrof(); ++rofId) { - mTimeFrame->getLines(rofId).reserve(mTimeFrame->getNTrackletsCluster(rofId, 0).size()); - trackletSelectionKernelHost( - mTimeFrame->getClustersOnLayer(rofId, 0), - mTimeFrame->getClustersOnLayer(rofId, 1), - mTimeFrame->getFoundTracklets(rofId, 0), - mTimeFrame->getFoundTracklets(rofId, 1), - mTimeFrame->getNTrackletsCluster(rofId, 0), - mTimeFrame->getNTrackletsCluster(rofId, 1), - mTimeFrame->getLines(rofId), - mTimeFrame->getLabelsFoundTracklets(rofId, 0), - mTimeFrame->getLinesLabel(rofId), - mVrtParams.tanLambdaCut, - mVrtParams.phiCut); + for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { + mTimeFrame->getLines(pivotRofId).reserve(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).size()); + std::vector usedTracklets(mTimeFrame->getFoundTracklets(pivotRofId, 0).size(), false); + int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; + int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; + for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { + LOGP(info, "Matching tracklets from rof {} to rof {}, delta={}", pivotRofId, targetRofId, pivotRofId - targetRofId); + trackletSelectionKernelHost( + mTimeFrame->getClustersOnLayer(targetRofId, 0), + mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getFoundTracklets(pivotRofId, 0), + mTimeFrame->getFoundTracklets(pivotRofId, 1), + usedTracklets, + mTimeFrame->getNTrackletsCluster(pivotRofId, 0), + mTimeFrame->getNTrackletsCluster(pivotRofId, 1), + mTimeFrame->getLines(pivotRofId), + mTimeFrame->getLabelsFoundTracklets(pivotRofId, 0), + mTimeFrame->getLinesLabel(pivotRofId), + pivotRofId - targetRofId, + mVrtParams.tanLambdaCut, + mVrtParams.phiCut); + } } #ifdef VTX_DEBUG diff --git a/Detectors/Upgrades/ITS3/reconstruction/src/IOUtils.cxx b/Detectors/Upgrades/ITS3/reconstruction/src/IOUtils.cxx index fd8a7d333037e..f930ae00f22c7 100644 --- a/Detectors/Upgrades/ITS3/reconstruction/src/IOUtils.cxx +++ b/Detectors/Upgrades/ITS3/reconstruction/src/IOUtils.cxx @@ -91,7 +91,7 @@ int loadROFrameDataITS3(its::TimeFrame* tf, tf->addClusterExternalIndexToLayer(layer, clusterId); } for (unsigned int iL{0}; iL < tf->getUnsortedClusters().size(); ++iL) { - tf->mROframesClusters[iL].push_back(tf->getUnsortedClusters()[iL].size()); + tf->mROFramesClusters[iL].push_back(tf->getUnsortedClusters()[iL].size()); } tf->mNrof++; } From 0be7198a970de6051b16e0122d8bdd007c3f8263 Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Wed, 8 May 2024 18:36:15 +0200 Subject: [PATCH 3/6] More vertices is better than less --- .../tracking/include/ITStracking/TimeFrame.h | 16 +++++++- .../tracking/include/ITStracking/Tracklet.h | 37 +++++++------------ .../ITSMFT/ITS/tracking/src/TimeFrame.cxx | 24 ++++++------ .../ITS/tracking/src/VertexerTraits.cxx | 32 ++++++++-------- 4 files changed, 58 insertions(+), 51 deletions(-) diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h index 33aa68aaf2438..a554417b09726 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h @@ -184,7 +184,8 @@ class TimeFrame int getROFCutAllMult() const { return mCutClusterMult + mCutVertexMult; } // Vertexer - void computeTrackletsScans(const int nThreads = 1); + void computeTrackletsPerROFScans(); + void computeTracletsPerClusterScans(); int& getNTrackletsROF(int rof, int combId); std::vector& getLines(int rof); int getNLinesTotal() const @@ -196,6 +197,7 @@ class TimeFrame gsl::span getFoundTracklets(int rofId, int combId); gsl::span getLabelsFoundTracklets(int rofId, int combId) const; gsl::span getNTrackletsCluster(int rofId, int combId); + gsl::span getExclusiveNTrackletsCluster(int rofId, int combId); uint32_t getTotalTrackletsTF(const int iLayer) { return mTotalTracklets[iLayer]; } int getTotalClustersPerROFrange(int rofMin, int range, int layerId) const; std::array& getBeamXY() { return mBeamPos; } @@ -225,6 +227,7 @@ class TimeFrame } virtual void setDevicePropagator(const o2::base::PropagatorImpl*){}; + const o2::base::PropagatorImpl* getDevicePropagator() const { return mPropagatorDevice; } template @@ -257,6 +260,7 @@ class TimeFrame std::vector> mROFramesClusters; const dataformats::MCTruthContainer* mClusterLabels = nullptr; std::array, 2> mNTrackletsPerCluster; // TODO: remove in favour of mNTrackletsPerROF + std::array, 2> mNTrackletsPerClusterSum; std::vector> mNClustersPerROF; std::vector> mIndexTables; std::vector> mTrackletsLookupTable; @@ -567,6 +571,16 @@ inline gsl::span TimeFrame::getNTrackletsCluster(int rofId, int combId) return {&mNTrackletsPerCluster[combId][startIdx], static_cast::size_type>(mROFramesClusters[1][rofId + 1] - startIdx)}; } +inline gsl::span TimeFrame::getExclusiveNTrackletsCluster(int rofId, int combId) +{ + if (rofId < 0 || rofId >= mNrof) { + return gsl::span(); + } + auto clusStartIdx{mROFramesClusters[1][rofId]}; + + return {&mNTrackletsPerClusterSum[combId][clusStartIdx], static_cast::size_type>(mROFramesClusters[1][rofId + 1] - clusStartIdx)}; +} + inline int& TimeFrame::getNTrackletsROF(int rof, int combId) { return mNTrackletsPerROF[combId][rof]; diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h index f37dfb67da6ef..51732f3063420 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h @@ -21,6 +21,10 @@ #include "GPUCommonMath.h" #include "GPUCommonDef.h" +#ifndef GPUCA_GPUCODE_DEVICE +#include +#endif + namespace o2 { namespace its @@ -42,6 +46,13 @@ struct Tracklet final { GPUhdi() void dump(const int, const int); GPUhdi() void dump(const int, const int) const; GPUhdi() unsigned char operator<(const Tracklet&) const; +#ifndef GPUCA_GPUCODE_DEVICE + std::string asString() const + { + return "fClIdx: " + std::to_string(firstClusterIndex) + " sClIdx: " + std::to_string(secondClusterIndex) + + " rof1: " + std::to_string(rof[0]) + " rof2: " + std::to_string(rof[1]); + } +#endif int firstClusterIndex; int secondClusterIndex; @@ -79,7 +90,7 @@ GPUhdi() Tracklet::Tracklet(const int idx0, const int idx1, float tanL, float ph // Nothing to do } -GPUhdi() bool Tracklet::operator==(const Tracklet & rhs) const +GPUhdi() bool Tracklet::operator==(const Tracklet& rhs) const { return this->firstClusterIndex == rhs.firstClusterIndex && this->secondClusterIndex == rhs.secondClusterIndex && @@ -89,7 +100,7 @@ GPUhdi() bool Tracklet::operator==(const Tracklet & rhs) const this->rof[1] == rhs.rof[1]; } -GPUhdi() bool Tracklet::operator!=(const Tracklet & rhs) const +GPUhdi() bool Tracklet::operator!=(const Tracklet& rhs) const { return this->firstClusterIndex != rhs.firstClusterIndex || this->secondClusterIndex != rhs.secondClusterIndex || @@ -97,7 +108,7 @@ GPUhdi() bool Tracklet::operator!=(const Tracklet & rhs) const this->phi != rhs.phi; } -GPUhdi() unsigned char Tracklet::operator<(const Tracklet & t) const +GPUhdi() unsigned char Tracklet::operator<(const Tracklet& t) const { if (isEmpty()) { return false; @@ -105,26 +116,6 @@ GPUhdi() unsigned char Tracklet::operator<(const Tracklet & t) const return true; } -// GPUhdi() void Tracklet::dump() -// { -// printf("fClIdx: %d sClIdx: %d rof1: %hu rof2: %hu phi: %f tl: %f \n", firstClusterIndex, secondClusterIndex, rof[0], rof[1], phi, tanLambda); -// } - -// GPUhdi() void Tracklet::dump() const -// { -// printf("fClIdx: %d sClIdx: %d rof1: %hu rof2: %hu phi: %f tl: %f \n", firstClusterIndex, secondClusterIndex, rof[0], rof[1], phi, tanLambda); -// } - -// GPUhdi() void Tracklet::dump(const int offsetFirst, const int offsetSecond) -// { -// printf("fClIdx: %d sClIdx: %d rof1: %hu rof2: %hu phi: %f tl: %f \n", firstClusterIndex + offsetFirst, secondClusterIndex + offsetSecond, rof[0], rof[1], phi, tanLambda); -// } - -// GPUhdi() void Tracklet::dump(const int offsetFirst, const int offsetSecond) const -// { -// printf("fClIdx: %d sClIdx: %d rof1: %hu rof2: %hu phi: %f tl: %f \n", firstClusterIndex + offsetFirst, secondClusterIndex + offsetSecond, rof[0], rof[1], phi, tanLambda); -// } - GPUhdi() void Tracklet::dump(const int offsetFirst, const int offsetSecond) { printf("fClIdx: %d sClIdx: %d rof1: %hu rof2: %hu\n", firstClusterIndex + offsetFirst, secondClusterIndex + offsetSecond, rof[0], rof[1]); diff --git a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx index f314a3871b053..e84dc21b05ded 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx @@ -229,8 +229,9 @@ int TimeFrame::loadROFrameData(gsl::span rofs, mNrof++; } - for (auto& v : mNTrackletsPerCluster) { - v.resize(mUnsortedClusters[1].size()); + for (auto i = 0; i < mNTrackletsPerCluster.size(); ++i) { + mNTrackletsPerCluster[i].resize(mUnsortedClusters[1].size()); + mNTrackletsPerClusterSum[i].resize(mUnsortedClusters[1].size() + 1); // Exc sum "prepends" a 0 } if (mcLabels) { @@ -433,6 +434,15 @@ void TimeFrame::fillPrimaryVerticesXandAlpha() } } +void TimeFrame::computeTrackletsPerROFScans() +{ + for (ushort iLayer = 0; iLayer < 2; ++iLayer) { + mTotalTracklets[iLayer] = std::accumulate(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), 0); + std::exclusive_scan(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), mNTrackletsPerROF[iLayer].begin(), 0); + std::exclusive_scan(mNTrackletsPerCluster[iLayer].begin(), mNTrackletsPerCluster[iLayer].end(), mNTrackletsPerClusterSum[iLayer].begin(), 0); + } +} + void TimeFrame::checkTrackletLUTs() { for (uint32_t iLayer{0}; iLayer < getTracklets().size(); ++iLayer) { @@ -551,15 +561,5 @@ void TimeFrame::printNClsPerROF() std::cout << std::endl; } } - -void TimeFrame::computeTrackletsScans(const int nThreads) -{ - // #pragma omp parallel for num_threads(nThreads > 1 ? 2 : nThreads) - for (ushort iLayer = 0; iLayer < 2; ++iLayer) { - mTotalTracklets[iLayer] = std::accumulate(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), 0); - std::exclusive_scan(mNTrackletsPerROF[iLayer].begin(), mNTrackletsPerROF[iLayer].end(), mNTrackletsPerROF[iLayer].begin(), 0); - } -} - } // namespace its } // namespace o2 diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx index 6f26a4109fd65..217352360c069 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx @@ -19,7 +19,7 @@ #include "ITStracking/VertexerTraits.h" #include "ITStracking/ClusterLines.h" #include "ITStracking/Tracklet.h" -#define VTX_DEBUG + #ifdef VTX_DEBUG #include "TTree.h" #include "TFile.h" @@ -55,13 +55,12 @@ void trackleterKernelHost( const IndexTableUtils& utils, const int pivotRof, const int targetRof, - int* rofFoundTrackletsOffset, // we want to change that, to keep track of the offset in deltaRof>0 + gsl::span rofFoundTrackletsOffsets, // we want to change those, to keep track of the offset in deltaRof>0 const int maxTrackletsPerCluster = static_cast(2e3)) { const int PhiBins{utils.getNphiBins()}; const int ZBins{utils.getNzBins()}; // loop on layer1 clusters - // int cumulativeStoredTracklets{rofFoundTrackletsOffset}; for (int iCurrentLayerClusterIndex = 0; iCurrentLayerClusterIndex < clustersCurrentLayer.size(); ++iCurrentLayerClusterIndex) { int storedTracklets{0}; const Cluster& currentCluster{clustersCurrentLayer[iCurrentLayerClusterIndex]}; @@ -83,9 +82,9 @@ void trackleterKernelHost( if (storedTracklets < maxTrackletsPerCluster) { if constexpr (!EvalRun) { if constexpr (Mode == TrackletMode::Layer0Layer1) { - tracklets[*rofFoundTrackletsOffset + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof}; + tracklets[rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof}; } else { - tracklets[*rofFoundTrackletsOffset + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof}; + tracklets[rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof}; } } ++storedTracklets; @@ -97,7 +96,7 @@ void trackleterKernelHost( if constexpr (EvalRun) { foundTracklets[iCurrentLayerClusterIndex] += storedTracklets; } else { - *rofFoundTrackletsOffset += storedTracklets; + rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] += storedTracklets; } } } @@ -124,6 +123,7 @@ void trackletSelectionKernelHost( for (int iTracklet12{offset12}; iTracklet12 < offset12 + foundTracklets12[iCurrentLayerClusterIndex]; ++iTracklet12) { for (int iTracklet01{offset01}; iTracklet01 < offset01 + foundTracklets01[iCurrentLayerClusterIndex]; ++iTracklet01) { const auto& tracklet01{tracklets01[iTracklet01]}; + const auto& tracklet12{tracklets12[iTracklet12]}; if (tracklet01.getDeltaRof() != rofDist) { continue; } @@ -131,7 +131,6 @@ void trackletSelectionKernelHost( const float deltaPhi{o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(tracklet01.phi, tracklets12[iTracklet12].phi))}; if (!usedTracklets[iTracklet01] && deltaTanLambda < tanLambdaCut && deltaPhi < phiCut && validTracklets != maxTracklets) { usedTracklets[iTracklet01] = true; - LOGP(info, "\t validated a tracklet over {}/{} with deltarof={}", tracklets01.size(), tracklets12.size(), tracklet01.getDeltaRof()); lines.emplace_back(tracklet01, clusters0.data(), clusters1.data()); if (trackletLabels.size()) { linesLabels.emplace_back(trackletLabels[iTracklet01]); @@ -194,7 +193,7 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - nullptr, // Offset in the tracklet buffer + gsl::span(), // Offset in the tracklet buffer mVrtParams.maxTrackletsPerCluster); trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), @@ -206,14 +205,14 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - nullptr, // Offset in the tracklet buffer + gsl::span(), // Offset in the tracklet buffer mVrtParams.maxTrackletsPerCluster); } mTimeFrame->getNTrackletsROF(pivotRofId, 0) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 0).end(), 0); mTimeFrame->getNTrackletsROF(pivotRofId, 1) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 1).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 1).end(), 0); } #pragma omp single - mTimeFrame->computeTrackletsScans(mNThreads); + mTimeFrame->computeTrackletsPerROFScans(); #pragma omp single mTimeFrame->getTracklets()[0].resize(mTimeFrame->getTotalTrackletsTF(0)); #pragma omp single @@ -223,8 +222,6 @@ void VertexerTraits::computeTracklets() for (int pivotRofId = 0; pivotRofId < mTimeFrame->getNrof(); ++pivotRofId) { int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; - auto mobileOffset0 = mTimeFrame->getNTrackletsROF(pivotRofId, 0); - auto mobileOffset1 = mTimeFrame->getNTrackletsROF(pivotRofId, 1); for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), @@ -236,7 +233,7 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - &mobileOffset0, + mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 0), mVrtParams.maxTrackletsPerCluster); trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), @@ -248,7 +245,7 @@ void VertexerTraits::computeTracklets() mIndexTableUtils, pivotRofId, targetRofId, - &mobileOffset1, + mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 1), mVrtParams.maxTrackletsPerCluster); } } @@ -303,8 +300,14 @@ void VertexerTraits::computeTracklets() std::ofstream out01("NTC01_cpu.txt"), out12("NTC12_cpu.txt"); for (int iRof{0}; iRof < mTimeFrame->getNrof(); ++iRof) { + out01 << "ROF: " << iRof << std::endl; + out12 << "ROF: " << iRof << std::endl; std::copy(mTimeFrame->getNTrackletsCluster(iRof, 0).begin(), mTimeFrame->getNTrackletsCluster(iRof, 0).end(), std::ostream_iterator(out01, "\t")); + out01 << std::endl; + std::copy(mTimeFrame->getExclusiveNTrackletsCluster(iRof, 0).begin(), mTimeFrame->getExclusiveNTrackletsCluster(iRof, 0).end(), std::ostream_iterator(out01, "\t")); std::copy(mTimeFrame->getNTrackletsCluster(iRof, 1).begin(), mTimeFrame->getNTrackletsCluster(iRof, 1).end(), std::ostream_iterator(out12, "\t")); + out12 << std::endl; + std::copy(mTimeFrame->getExclusiveNTrackletsCluster(iRof, 1).begin(), mTimeFrame->getExclusiveNTrackletsCluster(iRof, 1).end(), std::ostream_iterator(out12, "\t")); out01 << std::endl; out12 << std::endl; } @@ -322,7 +325,6 @@ void VertexerTraits::computeTrackletMatching() int startROF{std::max(0, pivotRofId - mVrtParams.deltaRof)}; int endROF{std::min(mTimeFrame->getNrof(), pivotRofId + mVrtParams.deltaRof + 1)}; for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) { - LOGP(info, "Matching tracklets from rof {} to rof {}, delta={}", pivotRofId, targetRofId, pivotRofId - targetRofId); trackletSelectionKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), mTimeFrame->getClustersOnLayer(pivotRofId, 1), From 826b395ad1957da7f8ae3bc462e7af62b541165c Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Wed, 15 May 2024 16:09:13 +0200 Subject: [PATCH 4/6] Improve tracklet validation and manage late vertices --- .../include/ITStracking/ClusterLines.h | 43 +++++++---- .../tracking/include/ITStracking/TimeFrame.h | 76 +++++++++++-------- .../ITSMFT/ITS/tracking/src/ClusterLines.cxx | 23 +++++- .../ITSMFT/ITS/tracking/src/TimeFrame.cxx | 11 ++- .../ITS/tracking/src/VertexerTraits.cxx | 16 ++-- 5 files changed, 111 insertions(+), 58 deletions(-) diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h index 550b8405770dc..08102df695abd 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h @@ -9,8 +9,8 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -#ifndef O2_ITSMFT_TRACKING_LINE_H_ -#define O2_ITSMFT_TRACKING_LINE_H_ +#ifndef O2_ITS_CLUSTERLINES_H +#define O2_ITS_CLUSTERLINES_H #include #include @@ -19,11 +19,8 @@ #include "ITStracking/Tracklet.h" #include "GPUCommonMath.h" -namespace o2 +namespace o2::its { -namespace its -{ - struct Line final { GPUhd() Line(); GPUhd() Line(const Line&); @@ -31,7 +28,7 @@ struct Line final { GPUhd() Line(const float firstPoint[3], const float secondPoint[3]); GPUhd() Line(const Tracklet&, const Cluster*, const Cluster*); - inline static float getDistanceFromPoint(const Line& line, const std::array& point); + static float getDistanceFromPoint(const Line& line, const std::array& point); GPUhd() static float getDistanceFromPoint(const Line& line, const float point[3]); static std::array getDCAComponents(const Line& line, const std::array point); GPUhd() static void getDCAComponents(const Line& line, const float point[3], float destArray[6]); @@ -39,12 +36,14 @@ struct Line final { static bool areParallel(const Line&, const Line&, const float precision = 1e-14); GPUhd() unsigned char isEmpty() const { return (originPoint[0] == 0.f && originPoint[1] == 0.f && originPoint[2] == 0.f) && (cosinesDirector[0] == 0.f && cosinesDirector[1] == 0.f && cosinesDirector[2] == 0.f); } + GPUhdi() auto getDeltaROF() const { return rof[1] - rof[0]; } GPUhd() void print() const; bool operator==(const Line&) const; bool operator!=(const Line&) const; + const short getMinROF() const { return rof[0] < rof[1] ? rof[0] : rof[1]; } - float originPoint[3], cosinesDirector[3]; // std::array originPoint, cosinesDirector; - float weightMatrix[6] = {1., 0., 0., 1., 0., 1.}; // std::array weightMatrix; + float originPoint[3], cosinesDirector[3]; + float weightMatrix[6] = {1., 0., 0., 1., 0., 1.}; // weightMatrix is a symmetric matrix internally stored as // 0 --> row = 0, col = 0 // 1 --> 0,1 @@ -52,11 +51,13 @@ struct Line final { // 3 --> 1,1 // 4 --> 1,2 // 5 --> 2,2 - // Debug quantities + short rof[2]; }; GPUhdi() Line::Line() : weightMatrix{1., 0., 0., 1., 0., 1.} { + rof[0] = 0; + rof[1] = 0; } GPUhdi() Line::Line(const Line& other) @@ -68,6 +69,9 @@ GPUhdi() Line::Line(const Line& other) for (int i{0}; i < 6; ++i) { weightMatrix[i] = other.weightMatrix[i]; } + for (int i{0}; i < 2; ++i) { + rof[i] = other.rof[i]; + } } GPUhdi() Line::Line(const float firstPoint[3], const float secondPoint[3]) @@ -83,6 +87,9 @@ GPUhdi() Line::Line(const float firstPoint[3], const float secondPoint[3]) for (int index{0}; index < 3; ++index) { cosinesDirector[index] *= inverseNorm; } + + rof[0] = 0; + rof[1] = 0; } GPUhdi() Line::Line(const Tracklet& tracklet, const Cluster* innerClusters, const Cluster* outerClusters) @@ -101,6 +108,9 @@ GPUhdi() Line::Line(const Tracklet& tracklet, const Cluster* innerClusters, cons for (int index{0}; index < 3; ++index) { cosinesDirector[index] *= inverseNorm; } + + rof[0] = tracklet.rof[0]; + rof[1] = tracklet.rof[1]; } // static functions: @@ -198,8 +208,8 @@ inline bool Line::operator!=(const Line& rhs) const GPUhdi() void Line::print() const { - printf("Line: originPoint = (%f, %f, %f), cosinesDirector = (%f, %f, %f)\n", - originPoint[0], originPoint[1], originPoint[2], cosinesDirector[0], cosinesDirector[1], cosinesDirector[2]); + printf("Line: originPoint = (%f, %f, %f), cosinesDirector = (%f, %f, %f), rofs = (%u, %u)\n", + originPoint[0], originPoint[1], originPoint[2], cosinesDirector[0], cosinesDirector[1], cosinesDirector[2], rof[0], rof[1]); } class ClusterLines final @@ -211,11 +221,13 @@ class ClusterLines final ClusterLines(const Line& firstLine, const Line& secondLine); void add(const int& lineLabel, const Line& line, const bool& weight = false); void computeClusterCentroid(); + void updateROFPoll(const Line&); inline std::vector& getLabels() { return mLabels; } inline int getSize() const { return mLabels.size(); } + inline short getROF() const { return mROF; } inline std::array getVertex() const { return mVertex; } inline std::array getRMS2() const { return mRMS2; } inline float getAvgDistance2() const { return mAvgDistance2; } @@ -230,8 +242,9 @@ class ClusterLines final std::array mVertex = {0.f}; // cluster centroid position std::array mRMS2 = {0.f}; // symmetric matrix: diagonal is RMS2 float mAvgDistance2 = 0.f; // substitute for chi2 + int mROFWeight = 0; // rof weight for voting + short mROF = -1; // rof }; -} // namespace its -} // namespace o2 -#endif /* O2_ITSMFT_TRACKING_LINE_H_ */ +} // namespace o2::its +#endif /* O2_ITS_CLUSTERLINES_H */ diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h index a554417b09726..6bfba0ebedaf3 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h @@ -78,12 +78,12 @@ class TimeFrame friend class TimeFrameGPU; TimeFrame(int nLayers = 7); const Vertex& getPrimaryVertex(const int) const; - gsl::span getPrimaryVertices(int rof) const; + gsl::span getPrimaryVertices(int rofId) const; gsl::span getPrimaryVertices(int romin, int romax) const; - gsl::span> getPrimaryVerticesLabels(const int rof) const; - gsl::span> getPrimaryVerticesXAlpha(int rof) const; + gsl::span> getPrimaryVerticesLabels(const int rofId) const; + gsl::span> getPrimaryVerticesXAlpha(int rofId) const; void fillPrimaryVerticesXandAlpha(); - int getPrimaryVerticesNum(int rofID = -1) const; + int getPrimaryVerticesNum(int rofId = -1) const; void addPrimaryVertices(const std::vector& vertices); void addPrimaryVertices(const gsl::span& vertices); void addPrimaryVertices(const std::vector&); @@ -100,7 +100,7 @@ class TimeFrame int getTotalClusters() const; bool empty() const; bool isGPU() const { return mIsGPU; } - int getSortedIndex(int rof, int layer, int i) const; + int getSortedIndex(int rofId, int layer, int i) const; int getSortedStartIndex(const int, const int) const; int getNrof() const; @@ -165,9 +165,9 @@ class TimeFrame std::vector>& getCellsNeighbours(); std::vector>& getCellsNeighboursLUT(); std::vector>& getRoads(); - std::vector& getTracks(int rof) { return mTracks[rof]; } - std::vector& getTracksLabel(const int rof) { return mTracksLabel[rof]; } - std::vector& getLinesLabel(const int rof) { return mLinesLabels[rof]; } + std::vector& getTracks(int rofId) { return mTracks[rofId]; } + std::vector& getTracksLabel(const int rofId) { return mTracksLabel[rofId]; } + std::vector& getLinesLabel(const int rofId) { return mLinesLabels[rofId]; } std::vector>& getVerticesLabels() { return mVerticesLabels; } int getNumberOfClusters() const; @@ -186,13 +186,13 @@ class TimeFrame // Vertexer void computeTrackletsPerROFScans(); void computeTracletsPerClusterScans(); - int& getNTrackletsROF(int rof, int combId); - std::vector& getLines(int rof); + int& getNTrackletsROF(int rofId, int combId); + std::vector& getLines(int rofId); int getNLinesTotal() const { return std::accumulate(mLines.begin(), mLines.end(), 0, [](int sum, const auto& l) { return sum + l.size(); }); } - std::vector& getTrackletClusters(int rof); + std::vector& getTrackletClusters(int rofId); gsl::span getFoundTracklets(int rofId, int combId) const; gsl::span getFoundTracklets(int rofId, int combId); gsl::span getLabelsFoundTracklets(int rofId, int combId) const; @@ -201,6 +201,7 @@ class TimeFrame uint32_t getTotalTrackletsTF(const int iLayer) { return mTotalTracklets[iLayer]; } int getTotalClustersPerROFrange(int rofMin, int range, int layerId) const; std::array& getBeamXY() { return mBeamPos; } + void insertLateVertex(const Vertex& vertex); // \Vertexer void initialiseRoadLabels(); @@ -326,19 +327,19 @@ class TimeFrame inline const Vertex& TimeFrame::getPrimaryVertex(const int vertexIndex) const { return mPrimaryVertices[vertexIndex]; } -inline gsl::span TimeFrame::getPrimaryVertices(int rof) const +inline gsl::span TimeFrame::getPrimaryVertices(int rofId) const { - const int start = mROFramesPV[rof]; - const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + const int start = mROFramesPV[rofId]; + const int stop_idx = rofId >= mNrof - 1 ? mNrof : rofId + 1; + int delta = mMultiplicityCutMask[rofId] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&mPrimaryVertices[start], static_cast::size_type>(delta)}; } -inline gsl::span> TimeFrame::getPrimaryVerticesLabels(const int rof) const +inline gsl::span> TimeFrame::getPrimaryVerticesLabels(const int rofId) const { - const int start = mROFramesPV[rof]; - const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + const int start = mROFramesPV[rofId]; + const int stop_idx = rofId >= mNrof - 1 ? mNrof : rofId + 1; + int delta = mMultiplicityCutMask[rofId] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&mVerticesLabels[start], static_cast::size_type>(delta)}; } @@ -347,24 +348,24 @@ inline gsl::span TimeFrame::getPrimaryVertices(int romin, int roma return {&mPrimaryVertices[mROFramesPV[romin]], static_cast::size_type>(mROFramesPV[romax + 1] - mROFramesPV[romin])}; } -inline gsl::span> TimeFrame::getPrimaryVerticesXAlpha(int rof) const +inline gsl::span> TimeFrame::getPrimaryVerticesXAlpha(int rofId) const { - const int start = mROFramesPV[rof]; - const int stop_idx = rof >= mNrof - 1 ? mNrof : rof + 1; - int delta = mMultiplicityCutMask[rof] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded + const int start = mROFramesPV[rofId]; + const int stop_idx = rofId >= mNrof - 1 ? mNrof : rofId + 1; + int delta = mMultiplicityCutMask[rofId] ? mROFramesPV[stop_idx] - start : 0; // return empty span if Rof is excluded return {&(mPValphaX[start]), static_cast>::size_type>(delta)}; } -inline int TimeFrame::getPrimaryVerticesNum(int rofID) const +inline int TimeFrame::getPrimaryVerticesNum(int rofId) const { - return rofID < 0 ? mPrimaryVertices.size() : mROFramesPV[rofID + 1] - mROFramesPV[rofID]; + return rofId < 0 ? mPrimaryVertices.size() : mROFramesPV[rofId + 1] - mROFramesPV[rofId]; } inline bool TimeFrame::empty() const { return getTotalClusters() == 0; } -inline int TimeFrame::getSortedIndex(int rof, int layer, int index) const { return mROFramesClusters[layer][rof] + index; } +inline int TimeFrame::getSortedIndex(int rofId, int layer, int index) const { return mROFramesClusters[layer][rofId] + index; } -inline int TimeFrame::getSortedStartIndex(const int rof, const int layer) const { return mROFramesClusters[layer][rof]; } +inline int TimeFrame::getSortedStartIndex(const int rofId, const int layer) const { return mROFramesClusters[layer][rofId]; } inline int TimeFrame::getNrof() const { return mNrof; } @@ -491,14 +492,14 @@ inline gsl::span TimeFrame::getIndexTable(int rofId, int layer) static_cast::size_type>(mIndexTableUtils.getNphiBins() * mIndexTableUtils.getNzBins() + 1)}; } -inline std::vector& TimeFrame::getLines(int rof) +inline std::vector& TimeFrame::getLines(int rofId) { - return mLines[rof]; + return mLines[rofId]; } -inline std::vector& TimeFrame::getTrackletClusters(int rof) +inline std::vector& TimeFrame::getTrackletClusters(int rofId) { - return mTrackletClusters[rof]; + return mTrackletClusters[rofId]; } template @@ -581,9 +582,9 @@ inline gsl::span TimeFrame::getExclusiveNTrackletsCluster(int rofId, int co return {&mNTrackletsPerClusterSum[combId][clusStartIdx], static_cast::size_type>(mROFramesClusters[1][rofId + 1] - clusStartIdx)}; } -inline int& TimeFrame::getNTrackletsROF(int rof, int combId) +inline int& TimeFrame::getNTrackletsROF(int rofId, int combId) { - return mNTrackletsPerROF[combId][rof]; + return mNTrackletsPerROF[combId][rofId]; } inline bool TimeFrame::isRoadFake(int i) const @@ -694,6 +695,15 @@ inline size_t TimeFrame::getNumberOfUsedClusters() const return nClusters; } +inline void TimeFrame::insertLateVertex(const Vertex& vertex) +{ + int rofId = vertex.getTimeStamp().getTimeStamp(); + mPrimaryVertices.insert(mPrimaryVertices.begin() + mROFramesPV[rofId + 1], vertex); // inserted at the end of sequence of vertices belonging to rofId + for (int i = rofId + 1; i < mROFramesPV.size(); ++i) { + mROFramesPV[i] += 1; + } +} + } // namespace its } // namespace o2 diff --git a/Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx b/Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx index 8605bb0ecada2..32b78d4f2ef18 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx @@ -84,6 +84,9 @@ ClusterLines::ClusterLines(const int firstLabel, const Line& firstLine, const in const bool weight) { + updateROFPoll(firstLine); + updateROFPoll(secondLine); + mLabels.push_back(firstLabel); if (secondLabel > 0) { mLabels.push_back(secondLabel); // don't add info in case of beamline used @@ -188,7 +191,8 @@ ClusterLines::ClusterLines(const Line& firstLine, const Line& secondLine) std::array covarianceFirst{1., 1., 1.}; std::array covarianceSecond{1., 1., 1.}; - + updateROFPoll(firstLine); + updateROFPoll(secondLine); for (int i{0}; i < 6; ++i) { mWeightMatrix[i] = firstLine.weightMatrix[i] + secondLine.weightMatrix[i]; } @@ -274,6 +278,7 @@ ClusterLines::ClusterLines(const Line& firstLine, const Line& secondLine) void ClusterLines::add(const int& lineLabel, const Line& line, const bool& weight) { mLabels.push_back(lineLabel); + updateROFPoll(line); std::array covariance{1., 1., 1.}; for (int i{0}; i < 6; ++i) { @@ -361,5 +366,21 @@ bool ClusterLines::operator==(const ClusterLines& rhs) const } return retval && this->mAvgDistance2 == rhs.mAvgDistance2; } + +GPUhdi() void ClusterLines::updateROFPoll(const Line& line) +{ + // Boyer-Moore voting for rof label + if (mROFWeight == 0) { + mROF = line.getMinROF(); + mROFWeight = 1; + } else { + if (mROF == line.getMinROF()) { + mROFWeight++; + } else { + mROFWeight--; + } + } +} + } // namespace its } // namespace o2 diff --git a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx index e84dc21b05ded..1810991630ce2 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx @@ -80,8 +80,13 @@ TimeFrame::TimeFrame(int nLayers) void TimeFrame::addPrimaryVertices(const std::vector& vertices) { + auto rofId = mROFramesPV.size() - 1; for (const auto& vertex : vertices) { - mPrimaryVertices.emplace_back(vertex); + if (vertex.getTimeStamp() != rofId) { + insertLateVertex(vertex); + } else { + mPrimaryVertices.emplace_back(vertex); + } if (!isBeamPositionOverridden) { const int w{vertex.getNContributors()}; mBeamPos[0] = (mBeamPos[0] * mBeamPosWeight + vertex.getX() * w) / (mBeamPosWeight + w); @@ -94,7 +99,11 @@ void TimeFrame::addPrimaryVertices(const std::vector& vertices) void TimeFrame::addPrimaryVertices(const gsl::span& vertices) { + auto rofId = mROFramesPV.size(); for (const auto& vertex : vertices) { + if (vertex.getTimeStamp() != rofId) { + LOG(warning) << "Vertex timestamp mismatch: " << vertex.getTimeStamp() << " vs " << rofId; + } mPrimaryVertices.emplace_back(vertex); if (!isBeamPositionOverridden) { const int w{vertex.getNContributors()}; diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx index 217352360c069..d49044d079607 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx @@ -124,7 +124,7 @@ void trackletSelectionKernelHost( for (int iTracklet01{offset01}; iTracklet01 < offset01 + foundTracklets01[iCurrentLayerClusterIndex]; ++iTracklet01) { const auto& tracklet01{tracklets01[iTracklet01]}; const auto& tracklet12{tracklets12[iTracklet12]}; - if (tracklet01.getDeltaRof() != rofDist) { + if (tracklet01.getDeltaRof() != rofDist || o2::gpu::GPUCommonMath::Abs(tracklet01.rof[0] - tracklet12.rof[1]) > 1) { continue; } const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(tracklet01.tanLambda - tracklets12[iTracklet12].tanLambda)}; @@ -336,7 +336,7 @@ void VertexerTraits::computeTrackletMatching() mTimeFrame->getLines(pivotRofId), mTimeFrame->getLabelsFoundTracklets(pivotRofId, 0), mTimeFrame->getLinesLabel(pivotRofId), - pivotRofId - targetRofId, + pivotRofId - targetRofId, // delta-rof mVrtParams.tanLambdaCut, mVrtParams.phiCut); } @@ -439,7 +439,7 @@ void VertexerTraits::computeVertices() std::array vertex2{}; for (int iCluster2{iCluster1 + 1}; iCluster2 < noClustersVec[rofId]; ++iCluster2) { vertex2 = mTimeFrame->getTrackletClusters(rofId)[iCluster2].getVertex(); - if (std::abs(vertex1[2] - vertex2[2]) < mVrtParams.clusterCut) { + if (o2::gpu::GPUCommonMath::Abs(vertex1[2] - vertex2[2]) < mVrtParams.clusterCut) { float distance{(vertex1[0] - vertex2[0]) * (vertex1[0] - vertex2[0]) + (vertex1[1] - vertex2[1]) * (vertex1[1] - vertex2[1]) + (vertex1[2] - vertex2[2]) * (vertex1[2] - vertex2[2])}; @@ -479,7 +479,7 @@ void VertexerTraits::computeVertices() } } - if (beamDistance2 < nsigmaCut && std::abs(mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[2]) < mVrtParams.maxZPositionAllowed) { + if (beamDistance2 < nsigmaCut && o2::gpu::GPUCommonMath::Abs(mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[2]) < mVrtParams.maxZPositionAllowed) { atLeastOneFound = true; vertices.emplace_back(o2::math_utils::Point3D(mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[0], mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[1], @@ -489,7 +489,7 @@ void VertexerTraits::computeVertices() mTimeFrame->getTrackletClusters(rofId)[iCluster].getSize(), // Contributors mTimeFrame->getTrackletClusters(rofId)[iCluster].getAvgDistance2()); // In place of chi2 - vertices.back().setTimeStamp(rofId); + vertices.back().setTimeStamp(mTimeFrame->getTrackletClusters(rofId)[iCluster].getROF()); if (mTimeFrame->hasMCinformation()) { mTimeFrame->getVerticesLabels().emplace_back(); for (auto& index : mTimeFrame->getTrackletClusters(rofId)[iCluster].getLabels()) { @@ -600,7 +600,7 @@ void VertexerTraits::computeVerticesInRof(int rofId, std::array vertex2{}; for (int iCluster2{iCluster1 + 1}; iCluster2 < nClusters; ++iCluster2) { vertex2 = clusterLines[iCluster2].getVertex(); - if (std::abs(vertex1[2] - vertex2[2]) < mVrtParams.clusterCut) { + if (o2::gpu::GPUCommonMath::Abs(vertex1[2] - vertex2[2]) < mVrtParams.clusterCut) { float distance{(vertex1[0] - vertex2[0]) * (vertex1[0] - vertex2[0]) + (vertex1[1] - vertex2[1]) * (vertex1[1] - vertex2[1]) + (vertex1[2] - vertex2[2]) * (vertex1[2] - vertex2[2])}; @@ -633,7 +633,7 @@ void VertexerTraits::computeVerticesInRof(int rofId, continue; } } - if (beamDistance2 < nsigmaCut && std::abs(clusterLines[iCluster].getVertex()[2]) < mVrtParams.maxZPositionAllowed) { + if (beamDistance2 < nsigmaCut && o2::gpu::GPUCommonMath::Abs(clusterLines[iCluster].getVertex()[2]) < mVrtParams.maxZPositionAllowed) { atLeastOneFound = true; ++foundVertices; vertices.emplace_back(o2::math_utils::Point3D(clusterLines[iCluster].getVertex()[0], @@ -643,7 +643,7 @@ void VertexerTraits::computeVerticesInRof(int rofId, // off-diagonal: square mean of projections on planes. clusterLines[iCluster].getSize(), // Contributors clusterLines[iCluster].getAvgDistance2()); // In place of chi2 - vertices.back().setTimeStamp(rofId); + vertices.back().setTimeStamp(clusterLines[iCluster].getROF()); if (labels) { for (auto& index : clusterLines[iCluster].getLabels()) { labels->push_back(tf->getLinesLabel(rofId)[index]); // then we can use nContributors from vertices to get the labels From 21ed1032e2f2519366384ef5f7a709f8e4997354 Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Thu, 16 May 2024 11:17:19 +0200 Subject: [PATCH 5/6] Fix GPU compilation --- .../ReconstructionDataFormats/Vertex.h | 4 ++-- .../include/CommonDataFormat/TimeStamp.h | 2 +- .../GPU/ITStrackingGPU/TimeFrameGPU.h | 10 +++++----- .../ITS/tracking/GPU/cuda/TimeFrameGPU.cu | 8 ++++---- .../tracking/GPU/cuda/VertexerTraitsGPU.cu | 20 +++++++++---------- .../include/ITStracking/ClusterLines.h | 2 +- 6 files changed, 23 insertions(+), 23 deletions(-) diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h index bbc2c01359276..7f9332f9447d4 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h @@ -146,8 +146,8 @@ class Vertex : public VertexBase GPUd() void setChi2(float v) { mChi2 = v; } GPUd() float getChi2() const { return mChi2; } - GPUd() const Stamp& getTimeStamp() const { return mTimeStamp; } - GPUd() Stamp& getTimeStamp() { return mTimeStamp; } + GPUhd() const Stamp& getTimeStamp() const { return mTimeStamp; } + GPUhd() Stamp& getTimeStamp() { return mTimeStamp; } GPUd() void setTimeStamp(const Stamp& v) { mTimeStamp = v; } protected: diff --git a/DataFormats/common/include/CommonDataFormat/TimeStamp.h b/DataFormats/common/include/CommonDataFormat/TimeStamp.h index 90019de5122be..354c937c1a24b 100644 --- a/DataFormats/common/include/CommonDataFormat/TimeStamp.h +++ b/DataFormats/common/include/CommonDataFormat/TimeStamp.h @@ -28,7 +28,7 @@ class TimeStamp GPUhdDefault() TimeStamp() CON_DEFAULT; GPUhdDefault() ~TimeStamp() CON_DEFAULT; GPUdi() TimeStamp(T time) { mTimeStamp = time; } - GPUdi() T getTimeStamp() const { return mTimeStamp; } + GPUhdi() T getTimeStamp() const { return mTimeStamp; } GPUdi() void setTimeStamp(T t) { mTimeStamp = t; } GPUdi() bool operator==(const TimeStamp& t) const { return mTimeStamp == t.mTimeStamp; } diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h index 2a2fadb21caeb..a67ef03e0ef22 100644 --- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h +++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h @@ -210,7 +210,7 @@ class TimeFrameGPU : public TimeFrame /// interface int getNClustersInRofSpan(const int, const int, const int) const; IndexTableUtils* getDeviceIndexTableUtils() { return mIndexTableUtilsDevice; } - int* getDeviceROframesClusters(const int layer) { return mROframesClustersDevice[layer]; } + int* getDeviceROFramesClusters(const int layer) { return mROFramesClustersDevice[layer]; } std::vector>& getVerticesInChunks() { return mVerticesInChunks; } std::vector>& getNVerticesInChunks() { return mNVerticesInChunks; } std::vector& getTrackITSExt() { return mTrackITSExt; } @@ -218,7 +218,7 @@ class TimeFrameGPU : public TimeFrame int getNAllocatedROFs() const { return mNrof; } // Allocated means maximum nROF for each chunk while populated is the number of loaded ones. StaticTrackingParameters* getDeviceTrackingParameters() { return mTrackingParamsDevice; } Vertex* getDeviceVertices() { return mVerticesDevice; } - int* getDeviceROframesPV() { return mROframesPVDevice; } + int* getDeviceROFramesPV() { return mROFramesPVDevice; } unsigned char* getDeviceUsedClusters(const int); const o2::base::Propagator* getChainPropagator(); @@ -251,10 +251,10 @@ class TimeFrameGPU : public TimeFrame // Device pointers StaticTrackingParameters* mTrackingParamsDevice; IndexTableUtils* mIndexTableUtilsDevice; - std::array mROframesClustersDevice; + std::array mROFramesClustersDevice; std::array mUsedClustersDevice; Vertex* mVerticesDevice; - int* mROframesPVDevice; + int* mROFramesPVDevice; // Hybrid pref std::array mClustersDevice; @@ -314,7 +314,7 @@ size_t TimeFrameGPU::loadChunkData(const size_t chunk, const size_t off template inline int TimeFrameGPU::getNClustersInRofSpan(const int rofIdstart, const int rofSpanSize, const int layerId) const { - return static_cast(mROframesClusters[layerId][(rofIdstart + rofSpanSize) < mROframesClusters.size() ? rofIdstart + rofSpanSize : mROframesClusters.size() - 1] - mROframesClusters[layerId][rofIdstart]); + return static_cast(mROFramesClusters[layerId][(rofIdstart + rofSpanSize) < mROFramesClusters.size() ? rofIdstart + rofSpanSize : mROFramesClusters.size() - 1] - mROFramesClusters[layerId][rofIdstart]); } } // namespace gpu } // namespace its diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu index e758bf4990c4f..1d0fca2585730 100644 --- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu +++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu @@ -443,21 +443,21 @@ void TimeFrameGPU::initDevice(const int chunks, mMemChunks[iChunk].allocate(GpuTimeFrameChunk::computeRofPerChunk(mGpuParams, mAvailMemGB), mGpuStreams[iChunk]); } for (auto iLayer{0}; iLayer < nLayers; ++iLayer) { - checkGPUError(cudaMalloc(reinterpret_cast(&mROframesClustersDevice[iLayer]), mROframesClusters[iLayer].size() * sizeof(int))); + checkGPUError(cudaMalloc(reinterpret_cast(&mROFramesClustersDevice[iLayer]), mROFramesClusters[iLayer].size() * sizeof(int))); checkGPUError(cudaMalloc(reinterpret_cast(&(mUsedClustersDevice[iLayer])), sizeof(unsigned char) * mGpuParams.clustersPerROfCapacity * mNrof)); } checkGPUError(cudaMalloc(reinterpret_cast(&mVerticesDevice), sizeof(Vertex) * mGpuParams.maxVerticesCapacity)); - checkGPUError(cudaMalloc(reinterpret_cast(&mROframesPVDevice), sizeof(int) * (mNrof + 1))); + checkGPUError(cudaMalloc(reinterpret_cast(&mROFramesPVDevice), sizeof(int) * (mNrof + 1))); mFirstInit = false; } if (maxLayers < nLayers) { // Vertexer for (auto iLayer{0}; iLayer < nLayers; ++iLayer) { - checkGPUError(cudaMemcpy(mROframesClustersDevice[iLayer], mROframesClusters[iLayer].data(), mROframesClusters[iLayer].size() * sizeof(int), cudaMemcpyHostToDevice)); + checkGPUError(cudaMemcpy(mROFramesClustersDevice[iLayer], mROFramesClusters[iLayer].data(), mROFramesClusters[iLayer].size() * sizeof(int), cudaMemcpyHostToDevice)); } } else { // Tracker checkGPUError(cudaMemcpy(mVerticesDevice, mPrimaryVertices.data(), sizeof(Vertex) * mPrimaryVertices.size(), cudaMemcpyHostToDevice)); - checkGPUError(cudaMemcpy(mROframesPVDevice, mROframesPV.data(), sizeof(int) * mROframesPV.size(), cudaMemcpyHostToDevice)); + checkGPUError(cudaMemcpy(mROFramesPVDevice, mROFramesPV.data(), sizeof(int) * mROFramesPV.size(), cudaMemcpyHostToDevice)); if (!iteration) { for (auto iLayer{0}; iLayer < nLayers; ++iLayer) { checkGPUError(cudaMemset(mUsedClustersDevice[iLayer], 0, sizeof(unsigned char) * mGpuParams.clustersPerROfCapacity * mNrof)); diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu index 9a1ed507ae5a4..49cffe2a02042 100644 --- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu +++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu @@ -627,8 +627,8 @@ void VertexerTraitsGPU::computeTracklets() gpu::trackleterKernelMultipleRof<<getStream(chunkId).get()>>>( mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(0), // const Cluster* clustersNextLayer, // 0 2 mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(1), // const Cluster* clustersCurrentLayer, // 1 1 - mTimeFrameGPU->getDeviceROframesClusters(0), // const int* sizeNextLClusters, - mTimeFrameGPU->getDeviceROframesClusters(1), // const int* sizeCurrentLClusters, + mTimeFrameGPU->getDeviceROFramesClusters(0), // const int* sizeNextLClusters, + mTimeFrameGPU->getDeviceROFramesClusters(1), // const int* sizeCurrentLClusters, mTimeFrameGPU->getChunk(chunkId).getDeviceIndexTables(0), // const int* nextIndexTables, mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(0), // Tracklet* Tracklets, mTimeFrameGPU->getChunk(chunkId).getDeviceNTrackletCluster(0), // int* foundTracklets, @@ -641,8 +641,8 @@ void VertexerTraitsGPU::computeTracklets() gpu::trackleterKernelMultipleRof<<getStream(chunkId).get()>>>( mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(2), // const Cluster* clustersNextLayer, // 0 2 mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(1), // const Cluster* clustersCurrentLayer, // 1 1 - mTimeFrameGPU->getDeviceROframesClusters(2), // const int* sizeNextLClusters, - mTimeFrameGPU->getDeviceROframesClusters(1), // const int* sizeCurrentLClusters, + mTimeFrameGPU->getDeviceROFramesClusters(2), // const int* sizeNextLClusters, + mTimeFrameGPU->getDeviceROFramesClusters(1), // const int* sizeCurrentLClusters, mTimeFrameGPU->getChunk(chunkId).getDeviceIndexTables(2), // const int* nextIndexTables, mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(1), // Tracklet* Tracklets, mTimeFrameGPU->getChunk(chunkId).getDeviceNTrackletCluster(1), // int* foundTracklets, @@ -655,8 +655,8 @@ void VertexerTraitsGPU::computeTracklets() gpu::trackletSelectionKernelMultipleRof<<getStream(chunkId).get()>>>( mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(0), // const Cluster* clusters0, // Clusters on layer 0 mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(1), // const Cluster* clusters1, // Clusters on layer 1 - mTimeFrameGPU->getDeviceROframesClusters(0), // const int* sizeClustersL0, // Number of clusters on layer 0 per ROF - mTimeFrameGPU->getDeviceROframesClusters(1), // const int* sizeClustersL1, // Number of clusters on layer 1 per ROF + mTimeFrameGPU->getDeviceROFramesClusters(0), // const int* sizeClustersL0, // Number of clusters on layer 0 per ROF + mTimeFrameGPU->getDeviceROFramesClusters(1), // const int* sizeClustersL1, // Number of clusters on layer 1 per ROF mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(0), // Tracklet* tracklets01, // Tracklets on layer 0-1 mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(1), // Tracklet* tracklets12, // Tracklets on layer 1-2 mTimeFrameGPU->getChunk(chunkId).getDeviceNTrackletCluster(0), // const int* nFoundTracklets01, // Number of tracklets found on layers 0-1 @@ -688,8 +688,8 @@ void VertexerTraitsGPU::computeTracklets() gpu::trackletSelectionKernelMultipleRof<<getStream(chunkId).get()>>>( mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(0), // const Cluster* clusters0, // Clusters on layer 0 mTimeFrameGPU->getChunk(chunkId).getDeviceClusters(1), // const Cluster* clusters1, // Clusters on layer 1 - mTimeFrameGPU->getDeviceROframesClusters(0), // const int* sizeClustersL0, // Number of clusters on layer 0 per ROF - mTimeFrameGPU->getDeviceROframesClusters(1), // const int* sizeClustersL1, // Number of clusters on layer 1 per ROF + mTimeFrameGPU->getDeviceROFramesClusters(0), // const int* sizeClustersL0, // Number of clusters on layer 0 per ROF + mTimeFrameGPU->getDeviceROFramesClusters(1), // const int* sizeClustersL1, // Number of clusters on layer 1 per ROF mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(0), // Tracklet* tracklets01, // Tracklets on layer 0-1 mTimeFrameGPU->getChunk(chunkId).getDeviceTracklets(1), // Tracklet* tracklets12, // Tracklets on layer 1-2 mTimeFrameGPU->getChunk(chunkId).getDeviceNTrackletCluster(0), // const int* nFoundTracklets01, // Number of tracklets found on layers 0-1 @@ -723,8 +723,8 @@ void VertexerTraitsGPU::computeTracklets() std::vector usedLines; for (int rofId{0}; rofId < rofs; ++rofId) { auto rof = offset + rofId; - auto clustersL1offsetRof = mTimeFrameGPU->getROframeClusters(1)[rof] - mTimeFrameGPU->getROframeClusters(1)[offset]; // starting cluster offset for this ROF - auto nClustersL1Rof = mTimeFrameGPU->getROframeClusters(1)[rof + 1] - mTimeFrameGPU->getROframeClusters(1)[rof]; // number of clusters for this ROF + auto clustersL1offsetRof = mTimeFrameGPU->getROFrameClusters(1)[rof] - mTimeFrameGPU->getROFrameClusters(1)[offset]; // starting cluster offset for this ROF + auto nClustersL1Rof = mTimeFrameGPU->getROFrameClusters(1)[rof + 1] - mTimeFrameGPU->getROFrameClusters(1)[rof]; // number of clusters for this ROF auto linesOffsetRof = exclusiveFoundLinesHost[clustersL1offsetRof]; // starting line offset for this ROF auto nLinesRof = exclusiveFoundLinesHost[clustersL1offsetRof + nClustersL1Rof] - linesOffsetRof; gsl::span linesInRof(lines.data() + linesOffsetRof, static_cast::size_type>(nLinesRof)); diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h index 08102df695abd..8a4a2193ca2ef 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h @@ -40,7 +40,7 @@ struct Line final { GPUhd() void print() const; bool operator==(const Line&) const; bool operator!=(const Line&) const; - const short getMinROF() const { return rof[0] < rof[1] ? rof[0] : rof[1]; } + short getMinROF() const { return rof[0] < rof[1] ? rof[0] : rof[1]; } float originPoint[3], cosinesDirector[3]; float weightMatrix[6] = {1., 0., 0., 1., 0., 1.}; From b38166f9082dded4db30a87a727870b6257a6d84 Mon Sep 17 00:00:00 2001 From: Matteo Concas Date: Fri, 7 Jun 2024 22:51:37 +0200 Subject: [PATCH 6/6] Fix leak and bugs --- .../tracking/GPU/cuda/VertexerTraitsGPU.cu | 2 +- .../tracking/include/ITStracking/TimeFrame.h | 33 +++++++++++++++---- .../tracking/include/ITStracking/Tracklet.h | 2 +- .../ITSMFT/ITS/tracking/src/TimeFrame.cxx | 32 +++++++++--------- .../ITS/tracking/src/TrackingInterface.cxx | 2 +- .../ITSMFT/ITS/tracking/src/Vertexer.cxx | 3 +- .../ITS/tracking/src/VertexerTraits.cxx | 30 +++++++++++++---- .../ITS3/workflow/src/TrackerSpec.cxx | 2 +- 8 files changed, 70 insertions(+), 36 deletions(-) diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu index 49cffe2a02042..f121609f9b019 100644 --- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu +++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu @@ -756,7 +756,7 @@ void VertexerTraitsGPU::computeTracklets() int start{0}; for (int rofId{0}; rofId < mTimeFrameGPU->getNVerticesInChunks()[chunkId].size(); ++rofId) { gsl::span rofVerts{mTimeFrameGPU->getVerticesInChunks()[chunkId].data() + start, static_cast::size_type>(mTimeFrameGPU->getNVerticesInChunks()[chunkId][rofId])}; - mTimeFrameGPU->addPrimaryVertices(rofVerts); + mTimeFrameGPU->addPrimaryVertices(rofVerts, rofId); if (mTimeFrameGPU->hasMCinformation()) { mTimeFrameGPU->getVerticesLabels().emplace_back(); // TODO: add MC labels diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h index 6bfba0ebedaf3..f3acb125e28ee 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h @@ -84,9 +84,9 @@ class TimeFrame gsl::span> getPrimaryVerticesXAlpha(int rofId) const; void fillPrimaryVerticesXandAlpha(); int getPrimaryVerticesNum(int rofId = -1) const; - void addPrimaryVertices(const std::vector& vertices); - void addPrimaryVertices(const gsl::span& vertices); - void addPrimaryVertices(const std::vector&); + void addPrimaryVertices(const std::vector& vertices, const int rofId); + void addPrimaryVertices(const gsl::span& vertices, const int rofId); + // void addPrimaryVertices(const std::vector&); void removePrimaryVerticesInROF(const int rofId); int loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span clusters, const dataformats::MCTruthContainer* mcLabels = nullptr); @@ -124,6 +124,8 @@ class TimeFrame gsl::span getClustersOnLayer(int rofId, int layerId) const; gsl::span getClustersPerROFrange(int rofMin, int range, int layerId) const; gsl::span getUnsortedClustersOnLayer(int rofId, int layerId) const; + gsl::span getUsedClustersROF(int rofId, int layerId); + gsl::span getUsedClustersROF(int rofId, int layerId) const; gsl::span getROFramesClustersPerROFrange(int rofMin, int range, int layerId) const; gsl::span getROFrameClusters(int layerId) const; gsl::span getNClustersROFrange(int rofMin, int range, int layerId) const; @@ -145,7 +147,8 @@ class TimeFrame void initialise(const int iteration, const TrackingParameters& trkParam, const int maxLayers = 7); void resetRofPV() { - mPrimaryVertices.clear(); + // mPrimaryVertices.clear(); + deepVectorClear(mPrimaryVertices); mROFramesPV.resize(1, 0); }; @@ -227,7 +230,7 @@ class TimeFrame } } - virtual void setDevicePropagator(const o2::base::PropagatorImpl*){}; + virtual void setDevicePropagator(const o2::base::PropagatorImpl*) {}; const o2::base::PropagatorImpl* getDevicePropagator() const { return mPropagatorDevice; } @@ -260,7 +263,7 @@ class TimeFrame std::vector> mClusterExternalIndices; std::vector> mROFramesClusters; const dataformats::MCTruthContainer* mClusterLabels = nullptr; - std::array, 2> mNTrackletsPerCluster; // TODO: remove in favour of mNTrackletsPerROF + std::array, 2> mNTrackletsPerCluster; std::array, 2> mNTrackletsPerClusterSum; std::vector> mNClustersPerROF; std::vector> mIndexTables; @@ -403,6 +406,24 @@ inline gsl::span TimeFrame::getClustersOnLayer(int rofId, int lay return {&mClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; } +inline gsl::span TimeFrame::getUsedClustersROF(int rofId, int layerId) +{ + if (rofId < 0 || rofId >= mNrof) { + return gsl::span(); + } + int startIdx{mROFramesClusters[layerId][rofId]}; + return {&mUsedClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; +} + +inline gsl::span TimeFrame::getUsedClustersROF(int rofId, int layerId) const +{ + if (rofId < 0 || rofId >= mNrof) { + return gsl::span(); + } + int startIdx{mROFramesClusters[layerId][rofId]}; + return {&mUsedClusters[layerId][startIdx], static_cast::size_type>(mROFramesClusters[layerId][rofId + 1] - startIdx)}; +} + inline gsl::span TimeFrame::getClustersPerROFrange(int rofMin, int range, int layerId) const { if (rofMin < 0 || rofMin >= mNrof) { diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h index 51732f3063420..37c02b90b7510 100644 --- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h +++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h @@ -50,7 +50,7 @@ struct Tracklet final { std::string asString() const { return "fClIdx: " + std::to_string(firstClusterIndex) + " sClIdx: " + std::to_string(secondClusterIndex) + - " rof1: " + std::to_string(rof[0]) + " rof2: " + std::to_string(rof[1]); + " rof1: " + std::to_string(rof[0]) + " rof2: " + std::to_string(rof[1]) + " delta: " + std::to_string(getDeltaRof()); } #endif diff --git a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx index 1810991630ce2..b75a475df3577 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx @@ -78,11 +78,10 @@ TimeFrame::TimeFrame(int nLayers) mTrackletsIndexROF.resize(2, {0}); } -void TimeFrame::addPrimaryVertices(const std::vector& vertices) +void TimeFrame::addPrimaryVertices(const std::vector& vertices, const int rofId) { - auto rofId = mROFramesPV.size() - 1; for (const auto& vertex : vertices) { - if (vertex.getTimeStamp() != rofId) { + if (vertex.getTimeStamp().getTimeStamp() != rofId) { insertLateVertex(vertex); } else { mPrimaryVertices.emplace_back(vertex); @@ -97,13 +96,10 @@ void TimeFrame::addPrimaryVertices(const std::vector& vertices) mROFramesPV.push_back(mPrimaryVertices.size()); } -void TimeFrame::addPrimaryVertices(const gsl::span& vertices) +void TimeFrame::addPrimaryVertices(const gsl::span& vertices, const int rofId) { - auto rofId = mROFramesPV.size(); + // auto rofId = mROFramesPV.size(); for (const auto& vertex : vertices) { - if (vertex.getTimeStamp() != rofId) { - LOG(warning) << "Vertex timestamp mismatch: " << vertex.getTimeStamp() << " vs " << rofId; - } mPrimaryVertices.emplace_back(vertex); if (!isBeamPositionOverridden) { const int w{vertex.getNContributors()}; @@ -115,15 +111,15 @@ void TimeFrame::addPrimaryVertices(const gsl::span& vertices) mROFramesPV.push_back(mPrimaryVertices.size()); } -void TimeFrame::addPrimaryVertices(const std::vector& lVertices) -{ - std::vector vertices; - for (auto& vertex : lVertices) { - vertices.emplace_back(o2::math_utils::Point3D(vertex.mX, vertex.mY, vertex.mZ), vertex.mRMS2, vertex.mContributors, vertex.mAvgDistance2); - vertices.back().setTimeStamp(vertex.mTimeStamp); - } - addPrimaryVertices(vertices); -} +// void TimeFrame::addPrimaryVertices(const std::vector& lVertices) +// { +// std::vector vertices; +// for (auto& vertex : lVertices) { +// vertices.emplace_back(o2::math_utils::Point3D(vertex.mX, vertex.mY, vertex.mZ), vertex.mRMS2, vertex.mContributors, vertex.mAvgDistance2); +// vertices.back().setTimeStamp(vertex.mTimeStamp); +// } +// addPrimaryVertices(vertices); +// } int TimeFrame::loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span clusters, const dataformats::MCTruthContainer* mcLabels) @@ -177,6 +173,8 @@ int TimeFrame::loadROFrameData(gsl::span rofs, if (iLayer < 2) { deepVectorClear(mTrackletsIndexROF[iLayer]); + deepVectorClear(mNTrackletsPerCluster[iLayer]); + deepVectorClear(mNTrackletsPerClusterSum[iLayer]); } } diff --git a/Detectors/ITSMFT/ITS/tracking/src/TrackingInterface.cxx b/Detectors/ITSMFT/ITS/tracking/src/TrackingInterface.cxx index 492976d164227..73056595c644c 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/TrackingInterface.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/TrackingInterface.cxx @@ -208,7 +208,7 @@ void ITSTrackingInterface::run(framework::ProcessingContext& pc) for (auto& v : vtxVecLoc) { vertices.push_back(v); } - mTimeFrame->addPrimaryVertices(vtxVecLoc); + mTimeFrame->addPrimaryVertices(vtxVecLoc, iRof); } } LOG(info) << fmt::format(" - rejected {}/{} ROFs: random/mult.sel:{} (seed {}), vtx.sel:{}", cutRandomMult + cutVertexMult, rofspan.size(), cutRandomMult, multEst.lastRandomSeed, cutVertexMult); diff --git a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx index 40e56f3b42f1a..16777c66c0425 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx @@ -43,8 +43,7 @@ float Vertexer::clustersToVertices(std::function logger) TrackingParameters trkPars; trkPars.PhiBins = mTraits->getVertexingParameters().PhiBins; trkPars.ZBins = mTraits->getVertexingParameters().ZBins; - timeInit += evaluateTask( - &Vertexer::initialiseVertexer, "Vertexer initialisation", [](std::string) {}, trkPars); + timeInit += evaluateTask(&Vertexer::initialiseVertexer, "Vertexer initialisation", [](std::string) {}, trkPars); timeTracklet = evaluateTask(&Vertexer::findTracklets, "Vertexer tracklet finding", [](std::string) {}); nTracklets01 = mTimeFrame->getTotalTrackletsTF(0); nTracklets12 = mTimeFrame->getTotalTrackletsTF(1); diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx index d49044d079607..1e8cea8f722b1 100644 --- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx +++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx @@ -46,8 +46,9 @@ float smallestAngleDifference(float a, float b) template void trackleterKernelHost( - const gsl::span& clustersNextLayer, // 0 2 - const gsl::span& clustersCurrentLayer, // 1 1 + const gsl::span& clustersNextLayer, // 0 2 + const gsl::span& clustersCurrentLayer, // 1 1 + const gsl::span& usedClustersNextLayer, // 0 2 int* indexTableNext, const float phiCut, std::vector& tracklets, @@ -77,6 +78,10 @@ void trackleterKernelHost( const int maxRowClusterIndex{indexTableNext[firstBinIndex + ZBins]}; // loop on clusters next layer for (int iNextLayerClusterIndex{firstRowClusterIndex}; iNextLayerClusterIndex < maxRowClusterIndex && iNextLayerClusterIndex < static_cast(clustersNextLayer.size()); ++iNextLayerClusterIndex) { + if (usedClustersNextLayer[iNextLayerClusterIndex]) { + LOGP(warning, "skipping used cluster"); + continue; + } const Cluster& nextCluster{clustersNextLayer[iNextLayerClusterIndex]}; if (o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(currentCluster.phi, nextCluster.phi)) < phiCut) { if (storedTracklets < maxTrackletsPerCluster) { @@ -104,6 +109,8 @@ void trackleterKernelHost( void trackletSelectionKernelHost( const gsl::span clusters0, // 0 const gsl::span clusters1, // 1 + gsl::span usedClusters0, // Layer 0 + gsl::span usedClusters2, // Layer 2 const gsl::span& tracklets01, const gsl::span& tracklets12, std::vector& usedTracklets, @@ -112,7 +119,8 @@ void trackletSelectionKernelHost( std::vector& lines, const gsl::span& trackletLabels, std::vector& linesLabels, - const int rofDist = 0, // signed rof distance to consider for the tracklet to be promoted, so to be sure taht clusters0 and clusters1 are consistent with tracklet indices + const int pivotRofId, + const int targetRofId, const float tanLambdaCut = 0.025f, const float phiCut = 0.005f, const int maxTracklets = static_cast(1e2)) @@ -124,12 +132,14 @@ void trackletSelectionKernelHost( for (int iTracklet01{offset01}; iTracklet01 < offset01 + foundTracklets01[iCurrentLayerClusterIndex]; ++iTracklet01) { const auto& tracklet01{tracklets01[iTracklet01]}; const auto& tracklet12{tracklets12[iTracklet12]}; - if (tracklet01.getDeltaRof() != rofDist || o2::gpu::GPUCommonMath::Abs(tracklet01.rof[0] - tracklet12.rof[1]) > 1) { + if (tracklet01.rof[0] != targetRofId || tracklet12.rof[1] != targetRofId) { continue; } const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(tracklet01.tanLambda - tracklets12[iTracklet12].tanLambda)}; const float deltaPhi{o2::gpu::GPUCommonMath::Abs(smallestAngleDifference(tracklet01.phi, tracklets12[iTracklet12].phi))}; if (!usedTracklets[iTracklet01] && deltaTanLambda < tanLambdaCut && deltaPhi < phiCut && validTracklets != maxTracklets) { + usedClusters0[tracklet01.firstClusterIndex] = true; + usedClusters2[tracklet12.secondClusterIndex] = true; usedTracklets[iTracklet01] = true; lines.emplace_back(tracklet01, clusters0.data(), clusters1.data()); if (trackletLabels.size()) { @@ -186,6 +196,7 @@ void VertexerTraits::computeTracklets() trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), // Clusters to be matched with the next layer in target rof mTimeFrame->getClustersOnLayer(pivotRofId, 1), // Clusters to be matched with the current layer in pivot rof + mTimeFrame->getUsedClustersROF(targetRofId, 0), // Span of the used clusters in the target rof mTimeFrame->getIndexTable(targetRofId, 0).data(), // Index table to access the data on the next layer in target rof mVrtParams.phiCut, mTimeFrame->getTracklets()[0], // Flat tracklet buffer @@ -198,6 +209,7 @@ void VertexerTraits::computeTracklets() trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClustersROF(targetRofId, 2), mTimeFrame->getIndexTable(targetRofId, 2).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[1], @@ -226,6 +238,7 @@ void VertexerTraits::computeTracklets() trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClustersROF(targetRofId, 0), mTimeFrame->getIndexTable(targetRofId, 0).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[0], @@ -238,6 +251,7 @@ void VertexerTraits::computeTracklets() trackleterKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 2), mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClustersROF(targetRofId, 2), mTimeFrame->getIndexTable(targetRofId, 2).data(), mVrtParams.phiCut, mTimeFrame->getTracklets()[1], @@ -328,6 +342,8 @@ void VertexerTraits::computeTrackletMatching() trackletSelectionKernelHost( mTimeFrame->getClustersOnLayer(targetRofId, 0), mTimeFrame->getClustersOnLayer(pivotRofId, 1), + mTimeFrame->getUsedClustersROF(targetRofId, 0), + mTimeFrame->getUsedClustersROF(targetRofId, 2), mTimeFrame->getFoundTracklets(pivotRofId, 0), mTimeFrame->getFoundTracklets(pivotRofId, 1), usedTracklets, @@ -336,7 +352,8 @@ void VertexerTraits::computeTrackletMatching() mTimeFrame->getLines(pivotRofId), mTimeFrame->getLabelsFoundTracklets(pivotRofId, 0), mTimeFrame->getLinesLabel(pivotRofId), - pivotRofId - targetRofId, // delta-rof + pivotRofId, + targetRofId, mVrtParams.tanLambdaCut, mVrtParams.phiCut); } @@ -488,7 +505,6 @@ void VertexerTraits::computeVertices() // off-diagonal: square mean of projections on planes. mTimeFrame->getTrackletClusters(rofId)[iCluster].getSize(), // Contributors mTimeFrame->getTrackletClusters(rofId)[iCluster].getAvgDistance2()); // In place of chi2 - vertices.back().setTimeStamp(mTimeFrame->getTrackletClusters(rofId)[iCluster].getROF()); if (mTimeFrame->hasMCinformation()) { mTimeFrame->getVerticesLabels().emplace_back(); @@ -498,7 +514,7 @@ void VertexerTraits::computeVertices() } } } - mTimeFrame->addPrimaryVertices(vertices); + mTimeFrame->addPrimaryVertices(vertices, rofId); } #ifdef VTX_DEBUG TFile* dbg_file = TFile::Open("artefacts_tf.root", "update"); diff --git a/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx index b3f5a26f9d8be..d94f985bc32ab 100644 --- a/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx +++ b/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx @@ -246,7 +246,7 @@ void TrackerDPL::run(ProcessingContext& pc) for (auto& v : vtxVecLoc) { vertices.push_back(v); } - timeFrame->addPrimaryVertices(vtxVecLoc); + timeFrame->addPrimaryVertices(vtxVecLoc, iRof); } } // LOG(info) << fmt::format(" - rejected {}/{} ROFs: random/mult.sel:{} (seed {}), vtx.sel:{}", cutRandomMult + cutVertexMult, rofspan.size(), cutRandomMult, multEst.lastRandomSeed, cutVertexMult);