Skip to content

Commit

Permalink
More vertices is better than less
Browse files Browse the repository at this point in the history
  • Loading branch information
@mconcas
mconcas committed May 10, 2024
1 parent f926c03 commit 6d1d3b5
Show file tree
Hide file tree
Showing 4 changed files with 58 additions and 51 deletions.
16 changes: 15 additions & 1 deletion Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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<Line>& getLines(int rof);
int getNLinesTotal() const
Expand All @@ -196,6 +197,7 @@ class TimeFrame
gsl::span<Tracklet> getFoundTracklets(int rofId, int combId);
gsl::span<const MCCompLabel> getLabelsFoundTracklets(int rofId, int combId) const;
gsl::span<int> getNTrackletsCluster(int rofId, int combId);
gsl::span<int> 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<float, 2>& getBeamXY() { return mBeamPos; }
Expand Down Expand Up @@ -225,6 +227,7 @@ class TimeFrame
}

virtual void setDevicePropagator(const o2::base::PropagatorImpl<float>*){};

const o2::base::PropagatorImpl<float>* getDevicePropagator() const { return mPropagatorDevice; }

template <typename... T>
Expand Down Expand Up @@ -257,6 +260,7 @@ class TimeFrame
std::vector<std::vector<int>> mROFramesClusters;
const dataformats::MCTruthContainer<MCCompLabel>* mClusterLabels = nullptr;
std::array<std::vector<int>, 2> mNTrackletsPerCluster; // TODO: remove in favour of mNTrackletsPerROF
std::array<std::vector<int>, 2> mNTrackletsPerClusterSum;
std::vector<std::vector<int>> mNClustersPerROF;
std::vector<std::vector<int>> mIndexTables;
std::vector<std::vector<int>> mTrackletsLookupTable;
Expand Down Expand Up @@ -567,6 +571,16 @@ inline gsl::span<int> TimeFrame::getNTrackletsCluster(int rofId, int combId)
return {&mNTrackletsPerCluster[combId][startIdx], static_cast<gsl::span<int>::size_type>(mROFramesClusters[1][rofId + 1] - startIdx)};
}

inline gsl::span<int> TimeFrame::getExclusiveNTrackletsCluster(int rofId, int combId)
{
if (rofId < 0 || rofId >= mNrof) {
return gsl::span<int>();
}
auto clusStartIdx{mROFramesClusters[1][rofId]};

return {&mNTrackletsPerClusterSum[combId][clusStartIdx], static_cast<gsl::span<int>::size_type>(mROFramesClusters[1][rofId + 1] - clusStartIdx)};
}

inline int& TimeFrame::getNTrackletsROF(int rof, int combId)
{
return mNTrackletsPerROF[combId][rof];
Expand Down
37 changes: 14 additions & 23 deletions Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,10 @@
#include "GPUCommonMath.h"
#include "GPUCommonDef.h"

#ifndef GPUCA_GPUCODE_DEVICE
#include <string>
#endif

namespace o2
{
namespace its
Expand All @@ -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;
Expand Down Expand Up @@ -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 &&
Expand All @@ -89,42 +100,22 @@ 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 ||
this->tanLambda != rhs.tanLambda ||
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;
}
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]);
Expand Down
24 changes: 12 additions & 12 deletions Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -229,8 +229,9 @@ int TimeFrame::loadROFrameData(gsl::span<o2::itsmft::ROFRecord> 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) {
Expand Down Expand Up @@ -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) {
Expand Down Expand Up @@ -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
32 changes: 17 additions & 15 deletions Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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"
Expand Down Expand Up @@ -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<int> rofFoundTrackletsOffsets, // we want to change those, to keep track of the offset in deltaRof>0
const int maxTrackletsPerCluster = static_cast<int>(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]};
Expand All @@ -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;
Expand All @@ -97,7 +96,7 @@ void trackleterKernelHost(
if constexpr (EvalRun) {
foundTracklets[iCurrentLayerClusterIndex] += storedTracklets;
} else {
*rofFoundTrackletsOffset += storedTracklets;
rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] += storedTracklets;
}
}
}
Expand All @@ -124,14 +123,14 @@ 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;
}
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;
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]);
Expand Down Expand Up @@ -194,7 +193,7 @@ void VertexerTraits::computeTracklets()
mIndexTableUtils,
pivotRofId,
targetRofId,
nullptr, // Offset in the tracklet buffer
gsl::span<int>(), // Offset in the tracklet buffer
mVrtParams.maxTrackletsPerCluster);
trackleterKernelHost<TrackletMode::Layer1Layer2, true>(
mTimeFrame->getClustersOnLayer(targetRofId, 2),
Expand All @@ -206,14 +205,14 @@ void VertexerTraits::computeTracklets()
mIndexTableUtils,
pivotRofId,
targetRofId,
nullptr, // Offset in the tracklet buffer
gsl::span<int>(), // 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
Expand All @@ -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<TrackletMode::Layer0Layer1, false>(
mTimeFrame->getClustersOnLayer(targetRofId, 0),
Expand All @@ -236,7 +233,7 @@ void VertexerTraits::computeTracklets()
mIndexTableUtils,
pivotRofId,
targetRofId,
&mobileOffset0,
mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 0),
mVrtParams.maxTrackletsPerCluster);
trackleterKernelHost<TrackletMode::Layer1Layer2, false>(
mTimeFrame->getClustersOnLayer(targetRofId, 2),
Expand All @@ -248,7 +245,7 @@ void VertexerTraits::computeTracklets()
mIndexTableUtils,
pivotRofId,
targetRofId,
&mobileOffset1,
mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 1),
mVrtParams.maxTrackletsPerCluster);
}
}
Expand Down Expand Up @@ -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<double>(out01, "\t"));
out01 << std::endl;
std::copy(mTimeFrame->getExclusiveNTrackletsCluster(iRof, 0).begin(), mTimeFrame->getExclusiveNTrackletsCluster(iRof, 0).end(), std::ostream_iterator<double>(out01, "\t"));
std::copy(mTimeFrame->getNTrackletsCluster(iRof, 1).begin(), mTimeFrame->getNTrackletsCluster(iRof, 1).end(), std::ostream_iterator<double>(out12, "\t"));
out12 << std::endl;
std::copy(mTimeFrame->getExclusiveNTrackletsCluster(iRof, 1).begin(), mTimeFrame->getExclusiveNTrackletsCluster(iRof, 1).end(), std::ostream_iterator<double>(out12, "\t"));
out01 << std::endl;
out12 << std::endl;
}
Expand All @@ -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),
Expand Down

0 comments on commit 6d1d3b5

Please sign in to comment.