Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -38,8 +38,8 @@ class GeometryTGeo;

/// TRK TimeFrame class that extends ITS TimeFrame functionality
/// This allows for customization of tracking algorithms specific to the TRK detector
template <int nLayers = 11>
class TimeFrame : public o2::its::TimeFrame<nLayers>
template <int NLayers = 11>
class TimeFrame : public o2::its::TimeFrame<NLayers>
{
public:
TimeFrame() = default;
Expand All @@ -50,8 +50,6 @@ class TimeFrame : public o2::its::TimeFrame<nLayers>

/// Process hits from TTree to initialize ROFs
/// \param hitsTree Tree containing TRK hits
/// \param mcHeaderTree Tree containing MC event headers
/// \param nEvents Number of events to process
/// \param gman TRK geometry manager instance
/// \param config Configuration parameters for hit reconstruction
int loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman, const nlohmann::json& config);
Expand All @@ -61,7 +59,8 @@ class TimeFrame : public o2::its::TimeFrame<nLayers>
/// \param nRofs Number of ROFs (Read-Out Frames)
/// \param nEvents Number of events to process
/// \param inROFpileup Number of events per ROF
void getPrimaryVerticesFromMC(TTree* mcHeaderTree, int nRofs, Long64_t nEvents, int inROFpileup);
/// \param rofLength ROF length in BCs (must match what was used in loadROFsFromHitTree)
void getPrimaryVerticesFromMC(TTree* mcHeaderTree, int nRofs, Long64_t nEvents, int inROFpileup, uint32_t rofLength = 198);
};

} // namespace trk
Expand Down
91 changes: 61 additions & 30 deletions Detectors/Upgrades/ALICE3/TRK/reconstruction/src/TimeFrame.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -23,11 +23,13 @@
#include <vector>
#include <array>

using o2::its::clearResizeBoundedVector;

namespace o2::trk
{

template <int nLayers>
int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman, const nlohmann::json& config)
template <int NLayers>
int TimeFrame<NLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman, const nlohmann::json& config)
{
constexpr std::array<int, 2> startLayer{0, 3};
const Long64_t nEvents = hitsTree->GetEntries();
Expand All @@ -39,23 +41,39 @@ int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman,

const int inROFpileup{config.contains("inROFpileup") ? config["inROFpileup"].get<int>() : 1};

// Calculate number of ROFs and initialize data structures
this->mNrof = (nEvents + inROFpileup - 1) / inROFpileup;
// Calculate number of ROFs
const int nRofs = (nEvents + inROFpileup - 1) / inROFpileup;

// Set up ROF timing for all layers (no staggering in TRK simulation, all layers read out together)
constexpr uint32_t rofLength = 198; // ROF length in BC
o2::its::ROFOverlapTable<NLayers> overlapTable;
for (int iLayer = 0; iLayer < NLayers; ++iLayer) {
overlapTable.defineLayer(iLayer, nRofs, rofLength, 0, 0, 0);
}
overlapTable.init();
this->setROFOverlapTable(overlapTable);

// Set up the vertex lookup table timing (pre-allocate, vertices will be filled later)
o2::its::ROFVertexLookupTable<NLayers> vtxLookupTable;
for (int iLayer = 0; iLayer < NLayers; ++iLayer) {
vtxLookupTable.defineLayer(iLayer, nRofs, rofLength, 0, 0, 0);
}
vtxLookupTable.init(); // pre-allocate without vertices
this->setROFVertexLookupTable(vtxLookupTable);

// Reset and prepare ROF data structures
for (int iLayer{0}; iLayer < nLayers; ++iLayer) {
for (int iLayer{0}; iLayer < NLayers; ++iLayer) {
this->mMinR[iLayer] = std::numeric_limits<float>::max();
this->mMaxR[iLayer] = std::numeric_limits<float>::lowest();
this->mROFramesClusters[iLayer].clear();
this->mROFramesClusters[iLayer].resize(this->mNrof + 1, 0);
this->mROFramesClusters[iLayer].resize(nRofs + 1, 0);
this->mUnsortedClusters[iLayer].clear();
this->mTrackingFrameInfo[iLayer].clear();
this->mClusterExternalIndices[iLayer].clear();
}

// Pre-count hits to reserve memory efficiently
int totalNHits{0};
std::array<int, nLayers> clusterCountPerLayer{};
std::array<int, NLayers> clusterCountPerLayer{};
for (Long64_t iEvent = 0; iEvent < nEvents; ++iEvent) {
hitsTree->GetEntry(iEvent);
for (const auto& hit : *trkHit) {
Expand All @@ -64,35 +82,35 @@ int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman,
}
int subDetID = gman->getSubDetID(hit.GetDetectorID());
const int layer = startLayer[subDetID] + gman->getLayer(hit.GetDetectorID());
if (layer >= nLayers) {
if (layer >= NLayers) {
continue;
}
++clusterCountPerLayer[layer];
totalNHits++;
}
}

// Reserve memory for all layers
for (int iLayer{0}; iLayer < nLayers; ++iLayer) {
// Reserve memory for all layers (mClusterSize is now per-layer)
for (int iLayer{0}; iLayer < NLayers; ++iLayer) {
this->mUnsortedClusters[iLayer].reserve(clusterCountPerLayer[iLayer]);
this->mTrackingFrameInfo[iLayer].reserve(clusterCountPerLayer[iLayer]);
this->mClusterExternalIndices[iLayer].reserve(clusterCountPerLayer[iLayer]);
clearResizeBoundedVector(this->mClusterSize[iLayer], clusterCountPerLayer[iLayer], this->mMemoryPool.get());
}
clearResizeBoundedVector(this->mClusterSize, totalNHits, this->mMemoryPool.get());

std::array<float, 11> resolution{0.001, 0.001, 0.001, 0.001, 0.004, 0.004, 0.004, 0.004, 0.004, 0.004, 0.004};
if (config["geometry"]["pitch"].size() == nLayers) {
for (int iLayer{0}; iLayer < config["geometry"]["pitch"].size(); ++iLayer) {
if (config["geometry"]["pitch"].size() == static_cast<size_t>(NLayers)) {
for (size_t iLayer{0}; iLayer < config["geometry"]["pitch"].size(); ++iLayer) {
LOGP(info, "Setting resolution for layer {} from config", iLayer);
LOGP(info, "Layer {} pitch {} cm", iLayer, config["geometry"]["pitch"][iLayer].get<float>());
resolution[iLayer] = config["geometry"]["pitch"][iLayer].get<float>() / std::sqrt(12.f);
}
}
LOGP(info, "Number of active parts in VD: {}", gman->getNumberOfActivePartsVD());

int hitCounter{0};
auto labels = new dataformats::MCTruthContainer<MCCompLabel>();
// One shared MC label container for all layers
auto* labels = new dataformats::MCTruthContainer<MCCompLabel>();

int hitCounter{0};
int iRof{0}; // Current ROF index
for (Long64_t iEvent = 0; iEvent < nEvents; ++iEvent) {
hitsTree->GetEntry(iEvent);
Expand All @@ -108,7 +126,7 @@ int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman,
o2::math_utils::Point3D<float> gloXYZ;
o2::math_utils::Point3D<float> trkXYZ;
float r{0.f};
if (layer >= nLayers) {
if (layer >= NLayers) {
continue;
}
if (layer >= 3) {
Expand Down Expand Up @@ -139,11 +157,12 @@ int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman,
std::array<float, 2>{trkXYZ.y(), trkXYZ.z()},
std::array<float, 3>{resolution[layer] * resolution[layer], 0., resolution[layer] * resolution[layer]});
/// Rotate to the global frame
this->addClusterToLayer(layer, gloXYZ.x(), gloXYZ.y(), gloXYZ.z(), this->mUnsortedClusters[layer].size());
const int clusterIdxInLayer = this->mUnsortedClusters[layer].size();
this->addClusterToLayer(layer, gloXYZ.x(), gloXYZ.y(), gloXYZ.z(), clusterIdxInLayer);
this->addClusterExternalIndexToLayer(layer, hitCounter);
MCCompLabel label{hit.GetTrackID(), static_cast<int>(iEvent), 0};
labels->addElement(hitCounter, label);
this->mClusterSize[hitCounter] = 1; // For compatibility with cluster-based tracking, set cluster size to 1 for hits
this->mClusterSize[layer][clusterIdxInLayer] = 1;
hitCounter++;
}
trkHit->clear();
Expand All @@ -154,21 +173,23 @@ int TimeFrame<nLayers>::loadROFsFromHitTree(TTree* hitsTree, GeometryTGeo* gman,
for (unsigned int iLayer{0}; iLayer < this->mUnsortedClusters.size(); ++iLayer) {
this->mROFramesClusters[iLayer][iRof] = this->mUnsortedClusters[iLayer].size(); // effectively calculating an exclusive sum
}
// Update primary vertices ROF structure
}
this->mClusterLabels = labels;
}
return this->mNrof;

// Set the shared labels container for all layers
for (int iLayer = 0; iLayer < NLayers; ++iLayer) {
this->mClusterLabels[iLayer] = labels;
}

return nRofs;
}

template <int nLayers>
void TimeFrame<nLayers>::getPrimaryVerticesFromMC(TTree* mcHeaderTree, int nRofs, Long64_t nEvents, int inROFpileup)
template <int NLayers>
void TimeFrame<NLayers>::getPrimaryVerticesFromMC(TTree* mcHeaderTree, int nRofs, Long64_t nEvents, int inROFpileup, uint32_t rofLength)
{
auto mcheader = new o2::dataformats::MCEventHeader;
mcHeaderTree->SetBranchAddress("MCEventHeader.", &mcheader);

this->mROFramesPV.clear();
this->mROFramesPV.resize(nRofs + 1, 0);
this->mPrimaryVertices.clear();

int iRof{0};
Expand All @@ -178,14 +199,24 @@ void TimeFrame<nLayers>::getPrimaryVerticesFromMC(TTree* mcHeaderTree, int nRofs
vertex.setXYZ(mcheader->GetX(), mcheader->GetY(), mcheader->GetZ());
vertex.setNContributors(30);
vertex.setChi2(0.f);
LOGP(debug, "ROF {}: Added primary vertex at ({}, {}, {})", iRof, mcheader->GetX(), mcheader->GetY(), mcheader->GetZ());
this->mPrimaryVertices.push_back(vertex);

// Set proper BC timestamp for vertex-ROF compatibility
// The vertex timestamp is set to the center of its ROF with half-ROF as error
const uint32_t rofCenter = static_cast<uint32_t>(rofLength * iRof + rofLength / 2);
const uint16_t rofHalf = static_cast<uint16_t>(rofLength / 2);
vertex.setTimeStamp({rofCenter, rofHalf});

LOGP(debug, "ROF {}: Added primary vertex at ({}, {}, {}) with BC timestamp [{}, +/-{}]",
iRof, mcheader->GetX(), mcheader->GetY(), mcheader->GetZ(), rofCenter, rofHalf);
this->addPrimaryVertex(vertex);
if ((iEvent + 1) % inROFpileup == 0 || iEvent == nEvents - 1) {
iRof++;
this->mROFramesPV[iRof] = this->mPrimaryVertices.size(); // effectively calculating an exclusive sum
}
}
this->mMultiplicityCutMask.resize(nRofs, true); /// all ROFs are valid with MC primary vertices.

// Update the vertex lookup table with the newly added vertices
this->updateROFVertexLookupTable();
}

// Explicit template instantiation for TRK with 11 layers
Expand Down
92 changes: 36 additions & 56 deletions Detectors/Upgrades/ALICE3/TRK/workflow/src/TrackerSpec.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,6 @@ namespace o2
using namespace framework;
namespace trk
{
using Vertex = o2::dataformats::Vertex<o2::dataformats::TimeStamp<int>>;

TrackerDPL::TrackerDPL(std::shared_ptr<o2::base::GRPGeomRequest> gr,
bool isMC,
Expand Down Expand Up @@ -84,18 +83,12 @@ std::vector<o2::its::TrackingParameters> TrackerDPL::createTrackingParamsFromCon
if (paramConfig.contains("NLayers")) {
params.NLayers = paramConfig["NLayers"].get<int>();
}
if (paramConfig.contains("DeltaROF")) {
params.DeltaROF = paramConfig["DeltaROF"].get<int>();
}
if (paramConfig.contains("ZBins")) {
params.ZBins = paramConfig["ZBins"].get<int>();
}
if (paramConfig.contains("PhiBins")) {
params.PhiBins = paramConfig["PhiBins"].get<int>();
}
if (paramConfig.contains("nROFsPerIterations")) {
params.nROFsPerIterations = paramConfig["nROFsPerIterations"].get<int>();
}
if (paramConfig.contains("ClusterSharing")) {
params.ClusterSharing = paramConfig["ClusterSharing"].get<int>();
}
Expand All @@ -119,27 +112,21 @@ std::vector<o2::its::TrackingParameters> TrackerDPL::createTrackingParamsFromCon
if (paramConfig.contains("TrackletMinPt")) {
params.TrackletMinPt = paramConfig["TrackletMinPt"].get<float>();
}
if (paramConfig.contains("TrackletsPerClusterLimit")) {
params.TrackletsPerClusterLimit = paramConfig["TrackletsPerClusterLimit"].get<float>();
}
if (paramConfig.contains("CellDeltaTanLambdaSigma")) {
params.CellDeltaTanLambdaSigma = paramConfig["CellDeltaTanLambdaSigma"].get<float>();
}
if (paramConfig.contains("CellsPerClusterLimit")) {
params.CellsPerClusterLimit = paramConfig["CellsPerClusterLimit"].get<float>();
}
if (paramConfig.contains("MaxChi2ClusterAttachment")) {
params.MaxChi2ClusterAttachment = paramConfig["MaxChi2ClusterAttachment"].get<float>();
}
if (paramConfig.contains("MaxChi2NDF")) {
params.MaxChi2NDF = paramConfig["MaxChi2NDF"].get<float>();
}
if (paramConfig.contains("TrackFollowerNSigmaCutZ")) {
params.TrackFollowerNSigmaCutZ = paramConfig["TrackFollowerNSigmaCutZ"].get<float>();
}
if (paramConfig.contains("TrackFollowerNSigmaCutPhi")) {
params.TrackFollowerNSigmaCutPhi = paramConfig["TrackFollowerNSigmaCutPhi"].get<float>();
}
// if (paramConfig.contains("TrackFollowerNSigmaCutZ")) {
// params.TrackFollowerNSigmaCutZ = paramConfig["TrackFollowerNSigmaCutZ"].get<float>();
// }
// if (paramConfig.contains("TrackFollowerNSigmaCutPhi")) {
// params.TrackFollowerNSigmaCutPhi = paramConfig["TrackFollowerNSigmaCutPhi"].get<float>();
// }

// Parse boolean parameters
if (paramConfig.contains("UseDiamond")) {
Expand All @@ -154,9 +141,9 @@ std::vector<o2::its::TrackingParameters> TrackerDPL::createTrackingParamsFromCon
if (paramConfig.contains("ShiftRefToCluster")) {
params.ShiftRefToCluster = paramConfig["ShiftRefToCluster"].get<bool>();
}
if (paramConfig.contains("FindShortTracks")) {
params.FindShortTracks = paramConfig["FindShortTracks"].get<bool>();
}
// if (paramConfig.contains("FindShortTracks")) {
// params.FindShortTracks = paramConfig["FindShortTracks"].get<bool>();
// }
if (paramConfig.contains("PerPrimaryVertexProcessing")) {
params.PerPrimaryVertexProcessing = paramConfig["PerPrimaryVertexProcessing"].get<bool>();
}
Expand All @@ -169,18 +156,18 @@ std::vector<o2::its::TrackingParameters> TrackerDPL::createTrackingParamsFromCon
if (paramConfig.contains("FataliseUponFailure")) {
params.FataliseUponFailure = paramConfig["FataliseUponFailure"].get<bool>();
}
if (paramConfig.contains("UseTrackFollower")) {
params.UseTrackFollower = paramConfig["UseTrackFollower"].get<bool>();
}
if (paramConfig.contains("UseTrackFollowerTop")) {
params.UseTrackFollowerTop = paramConfig["UseTrackFollowerTop"].get<bool>();
}
if (paramConfig.contains("UseTrackFollowerBot")) {
params.UseTrackFollowerBot = paramConfig["UseTrackFollowerBot"].get<bool>();
}
if (paramConfig.contains("UseTrackFollowerMix")) {
params.UseTrackFollowerMix = paramConfig["UseTrackFollowerMix"].get<bool>();
}
// if (paramConfig.contains("UseTrackFollower")) {
// params.UseTrackFollower = paramConfig["UseTrackFollower"].get<bool>();
// }
// if (paramConfig.contains("UseTrackFollowerTop")) {
// params.UseTrackFollowerTop = paramConfig["UseTrackFollowerTop"].get<bool>();
// }
// if (paramConfig.contains("UseTrackFollowerBot")) {
// params.UseTrackFollowerBot = paramConfig["UseTrackFollowerBot"].get<bool>();
// }
// if (paramConfig.contains("UseTrackFollowerMix")) {
// params.UseTrackFollowerMix = paramConfig["UseTrackFollowerMix"].get<bool>();
// }
if (paramConfig.contains("createArtefactLabels")) {
params.createArtefactLabels = paramConfig["createArtefactLabels"].get<bool>();
}
Expand Down Expand Up @@ -297,44 +284,37 @@ void TrackerDPL::run(ProcessingContext& pc)
for (size_t iter{0}; iter < trackingParams.size(); ++iter) {
LOGP(info, "{}", trackingParams[iter].asString());
timeFrame.initialise(iter, trackingParams[iter], 11, false);
itsTrackerTraits.computeLayerTracklets(iter, -1, -1);
itsTrackerTraits.computeLayerTracklets(iter, -1);
LOGP(info, "Number of tracklets in iteration {}: {}", iter, timeFrame.getNumberOfTracklets());
itsTrackerTraits.computeLayerCells(iter);
LOGP(info, "Number of cells in iteration {}: {}", iter, timeFrame.getNumberOfCells());
itsTrackerTraits.findCellsNeighbours(iter);
LOGP(info, "Number of cell neighbours in iteration {}: {}", iter, timeFrame.getNumberOfNeighbours());
itsTrackerTraits.findRoads(iter);
LOGP(info, "Number of roads in iteration {}: {}", iter, timeFrame.getNumberOfTracks());
itsTrackerTraits.extendTracks(iter);
LOGP(info, "Number of tracks in iteration {}: {}", iter, timeFrame.getNumberOfTracks());
}
const auto trackingLoopElapsedMs = std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::steady_clock::now() - trackingLoopStart).count();
LOGP(info, "Tracking iterations block took {} ms", trackingLoopElapsedMs);

itsTracker.computeTracksMClabels();

// Stream tracks and their MC labels to the output
// Collect all tracks and labels from all ROFs
std::vector<o2::its::TrackITS> allTracks;
std::vector<o2::MCCompLabel> allLabels;
// Collect tracks and labels (flat vectors in the new interface)
const auto& tracks = timeFrame.getTracks();
const auto& labels = timeFrame.getTracksLabel();

int totalTracks = 0;
// Copy to output vectors (TrackITSExt -> TrackITS slicing for output compatibility)
std::vector<o2::its::TrackITS> allTracks(tracks.begin(), tracks.end());
std::vector<o2::MCCompLabel> allLabels(labels.begin(), labels.end());

int totalTracks = allTracks.size();
int goodTracks = 0;
int fakeTracks = 0;

for (int iRof = 0; iRof < nRofs; ++iRof) {
const auto& rofTracks = timeFrame.getTracks(iRof);
const auto& rofLabels = timeFrame.getTracksLabel(iRof);

allTracks.insert(allTracks.end(), rofTracks.begin(), rofTracks.end());
allLabels.insert(allLabels.end(), rofLabels.begin(), rofLabels.end());

totalTracks += rofTracks.size();
for (const auto& label : rofLabels) {
if (label.isFake()) {
fakeTracks++;
} else {
goodTracks++;
}
for (const auto& label : allLabels) {
if (label.isFake()) {
fakeTracks++;
} else {
goodTracks++;
}
}

Expand Down
Loading