Skip to content

Commit

Permalink
Improve tracklet validation and manage late vertices
Browse files Browse the repository at this point in the history
  • Loading branch information
mconcas committed May 15, 2024
1 parent 6d1d3b5 commit 1b81ea3
Show file tree
Hide file tree
Showing 5 changed files with 111 additions and 58 deletions.
43 changes: 28 additions & 15 deletions Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h
Original file line number Diff line number Diff line change
Expand Up @@ -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 <array>
#include <vector>
Expand All @@ -19,44 +19,45 @@
#include "ITStracking/Tracklet.h"
#include "GPUCommonMath.h"

namespace o2
namespace o2::its
{
namespace its
{

struct Line final {
GPUhd() Line();
GPUhd() Line(const Line&);
Line(std::array<float, 3> firstPoint, std::array<float, 3> secondPoint);
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<float, 3>& point);
static float getDistanceFromPoint(const Line& line, const std::array<float, 3>& point);
GPUhd() static float getDistanceFromPoint(const Line& line, const float point[3]);
static std::array<float, 6> getDCAComponents(const Line& line, const std::array<float, 3> point);
GPUhd() static void getDCAComponents(const Line& line, const float point[3], float destArray[6]);
GPUhd() static float getDCA(const Line&, const Line&, const float precision = 1e-14);
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<float, 3> originPoint, cosinesDirector;
float weightMatrix[6] = {1., 0., 0., 1., 0., 1.}; // std::array<float, 6> 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
// 2 --> 0,2
// 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)
Expand All @@ -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])
Expand All @@ -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)
Expand All @@ -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:
Expand Down Expand Up @@ -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
Expand All @@ -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<int>& getLabels()
{
return mLabels;
}
inline int getSize() const { return mLabels.size(); }
inline short getROF() const { return mROF; }
inline std::array<float, 3> getVertex() const { return mVertex; }
inline std::array<float, 6> getRMS2() const { return mRMS2; }
inline float getAvgDistance2() const { return mAvgDistance2; }
Expand All @@ -230,8 +242,9 @@ class ClusterLines final
std::array<float, 3> mVertex = {0.f}; // cluster centroid position
std::array<float, 6> 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 */
76 changes: 43 additions & 33 deletions Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h
Original file line number Diff line number Diff line change
Expand Up @@ -78,12 +78,12 @@ class TimeFrame
friend class TimeFrameGPU;
TimeFrame(int nLayers = 7);
const Vertex& getPrimaryVertex(const int) const;
gsl::span<const Vertex> getPrimaryVertices(int rof) const;
gsl::span<const Vertex> getPrimaryVertices(int rofId) const;
gsl::span<const Vertex> getPrimaryVertices(int romin, int romax) const;
gsl::span<const std::vector<MCCompLabel>> getPrimaryVerticesLabels(const int rof) const;
gsl::span<const std::array<float, 2>> getPrimaryVerticesXAlpha(int rof) const;
gsl::span<const std::vector<MCCompLabel>> getPrimaryVerticesLabels(const int rofId) const;
gsl::span<const std::array<float, 2>> getPrimaryVerticesXAlpha(int rofId) const;
void fillPrimaryVerticesXandAlpha();
int getPrimaryVerticesNum(int rofID = -1) const;
int getPrimaryVerticesNum(int rofId = -1) const;
void addPrimaryVertices(const std::vector<Vertex>& vertices);
void addPrimaryVertices(const gsl::span<const Vertex>& vertices);
void addPrimaryVertices(const std::vector<lightVertex>&);
Expand All @@ -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;

Expand Down Expand Up @@ -165,9 +165,9 @@ class TimeFrame
std::vector<std::vector<int>>& getCellsNeighbours();
std::vector<std::vector<int>>& getCellsNeighboursLUT();
std::vector<Road<5>>& getRoads();
std::vector<TrackITSExt>& getTracks(int rof) { return mTracks[rof]; }
std::vector<MCCompLabel>& getTracksLabel(const int rof) { return mTracksLabel[rof]; }
std::vector<MCCompLabel>& getLinesLabel(const int rof) { return mLinesLabels[rof]; }
std::vector<TrackITSExt>& getTracks(int rofId) { return mTracks[rofId]; }
std::vector<MCCompLabel>& getTracksLabel(const int rofId) { return mTracksLabel[rofId]; }
std::vector<MCCompLabel>& getLinesLabel(const int rofId) { return mLinesLabels[rofId]; }
std::vector<std::vector<MCCompLabel>>& getVerticesLabels() { return mVerticesLabels; }

int getNumberOfClusters() const;
Expand All @@ -186,13 +186,13 @@ class TimeFrame
// Vertexer
void computeTrackletsPerROFScans();
void computeTracletsPerClusterScans();
int& getNTrackletsROF(int rof, int combId);
std::vector<Line>& getLines(int rof);
int& getNTrackletsROF(int rofId, int combId);
std::vector<Line>& 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<ClusterLines>& getTrackletClusters(int rof);
std::vector<ClusterLines>& getTrackletClusters(int rofId);
gsl::span<const Tracklet> getFoundTracklets(int rofId, int combId) const;
gsl::span<Tracklet> getFoundTracklets(int rofId, int combId);
gsl::span<const MCCompLabel> getLabelsFoundTracklets(int rofId, int combId) const;
Expand All @@ -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<float, 2>& getBeamXY() { return mBeamPos; }
void insertLateVertex(const Vertex& vertex);
// \Vertexer

void initialiseRoadLabels();
Expand Down Expand Up @@ -326,19 +327,19 @@ class TimeFrame

inline const Vertex& TimeFrame::getPrimaryVertex(const int vertexIndex) const { return mPrimaryVertices[vertexIndex]; }

inline gsl::span<const Vertex> TimeFrame::getPrimaryVertices(int rof) const
inline gsl::span<const Vertex> 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<gsl::span<const Vertex>::size_type>(delta)};
}

inline gsl::span<const std::vector<MCCompLabel>> TimeFrame::getPrimaryVerticesLabels(const int rof) const
inline gsl::span<const std::vector<MCCompLabel>> 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<gsl::span<const Vertex>::size_type>(delta)};
}

Expand All @@ -347,24 +348,24 @@ inline gsl::span<const Vertex> TimeFrame::getPrimaryVertices(int romin, int roma
return {&mPrimaryVertices[mROFramesPV[romin]], static_cast<gsl::span<const Vertex>::size_type>(mROFramesPV[romax + 1] - mROFramesPV[romin])};
}

inline gsl::span<const std::array<float, 2>> TimeFrame::getPrimaryVerticesXAlpha(int rof) const
inline gsl::span<const std::array<float, 2>> 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<gsl::span<const std::array<float, 2>>::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; }

Expand Down Expand Up @@ -491,14 +492,14 @@ inline gsl::span<int> TimeFrame::getIndexTable(int rofId, int layer)
static_cast<gsl::span<int>::size_type>(mIndexTableUtils.getNphiBins() * mIndexTableUtils.getNzBins() + 1)};
}

inline std::vector<Line>& TimeFrame::getLines(int rof)
inline std::vector<Line>& TimeFrame::getLines(int rofId)
{
return mLines[rof];
return mLines[rofId];
}

inline std::vector<ClusterLines>& TimeFrame::getTrackletClusters(int rof)
inline std::vector<ClusterLines>& TimeFrame::getTrackletClusters(int rofId)
{
return mTrackletClusters[rof];
return mTrackletClusters[rofId];
}

template <typename... T>
Expand Down Expand Up @@ -581,9 +582,9 @@ inline gsl::span<int> TimeFrame::getExclusiveNTrackletsCluster(int rofId, int co
return {&mNTrackletsPerClusterSum[combId][clusStartIdx], static_cast<gsl::span<int>::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
Expand Down Expand Up @@ -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

Expand Down
23 changes: 22 additions & 1 deletion Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -188,7 +191,8 @@ ClusterLines::ClusterLines(const Line& firstLine, const Line& secondLine)

std::array<float, 3> covarianceFirst{1., 1., 1.};
std::array<float, 3> covarianceSecond{1., 1., 1.};

updateROFPoll(firstLine);
updateROFPoll(secondLine);
for (int i{0}; i < 6; ++i) {
mWeightMatrix[i] = firstLine.weightMatrix[i] + secondLine.weightMatrix[i];
}
Expand Down Expand Up @@ -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<float, 3> covariance{1., 1., 1.};

for (int i{0}; i < 6; ++i) {
Expand Down Expand Up @@ -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
11 changes: 10 additions & 1 deletion Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -80,8 +80,13 @@ TimeFrame::TimeFrame(int nLayers)

void TimeFrame::addPrimaryVertices(const std::vector<Vertex>& 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);
Expand All @@ -94,7 +99,11 @@ void TimeFrame::addPrimaryVertices(const std::vector<Vertex>& vertices)

void TimeFrame::addPrimaryVertices(const gsl::span<const Vertex>& 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()};
Expand Down
Loading

0 comments on commit 1b81ea3

Please sign in to comment.