Project
Loading...
Searching...
No Matches
VertexerTraits.cxx
Go to the documentation of this file.
1// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3// All rights not expressly granted are reserved.
4//
5// This software is distributed under the terms of the GNU General Public
6// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7//
8// In applying this license CERN does not waive the privileges and immunities
9// granted to it by virtue of its status as an Intergovernmental Organization
10// or submit itself to any jurisdiction.
12
13#include <memory>
14#include <ranges>
15#include <map>
16#include <algorithm>
17
18#include <oneapi/tbb/blocked_range.h>
19#include <oneapi/tbb/parallel_for.h>
20
30
31using namespace o2::its;
32
33template <TrackletMode Mode, bool EvalRun>
34static void trackleterKernelHost(
35 const gsl::span<const Cluster>& clustersNextLayer, // 0 2
36 const gsl::span<const Cluster>& clustersCurrentLayer, // 1 1
37 const gsl::span<uint8_t>& usedClustersNextLayer, // 0 2
38 int* indexTableNext,
39 const float phiCut,
41 gsl::span<int> foundTracklets,
43 const short pivotRof,
44 const short targetRof,
45 gsl::span<int> rofFoundTrackletsOffsets, // we want to change those, to keep track of the offset in deltaRof>0
46 const int maxTrackletsPerCluster = static_cast<int>(2e3))
47{
48 const int PhiBins{utils.getNphiBins()};
49 const int ZBins{utils.getNzBins()};
50 // loop on layer1 clusters
51 for (int iCurrentLayerClusterIndex = 0; iCurrentLayerClusterIndex < clustersCurrentLayer.size(); ++iCurrentLayerClusterIndex) {
52 int storedTracklets{0};
53 const Cluster& currentCluster{clustersCurrentLayer[iCurrentLayerClusterIndex]};
54 const int4 selectedBinsRect{VertexerTraits::getBinsRect(currentCluster, (int)Mode, 0.f, 50.f, phiCut / 2, utils)};
55 if (selectedBinsRect.x != 0 || selectedBinsRect.y != 0 || selectedBinsRect.z != 0 || selectedBinsRect.w != 0) {
56 int phiBinsNum{selectedBinsRect.w - selectedBinsRect.y + 1};
57 if (phiBinsNum < 0) {
58 phiBinsNum += PhiBins;
59 }
60 // loop on phi bins next layer
61 for (int iPhiBin{selectedBinsRect.y}, iPhiCount{0}; iPhiCount < phiBinsNum; iPhiBin = ++iPhiBin == PhiBins ? 0 : iPhiBin, iPhiCount++) {
62 const int firstBinIndex{utils.getBinIndex(selectedBinsRect.x, iPhiBin)};
63 const int firstRowClusterIndex{indexTableNext[firstBinIndex]};
64 const int maxRowClusterIndex{indexTableNext[firstBinIndex + ZBins]};
65 // loop on clusters next layer
66 for (int iNextLayerClusterIndex{firstRowClusterIndex}; iNextLayerClusterIndex < maxRowClusterIndex && iNextLayerClusterIndex < static_cast<int>(clustersNextLayer.size()); ++iNextLayerClusterIndex) {
67 if (usedClustersNextLayer[iNextLayerClusterIndex]) {
68 continue;
69 }
70 const Cluster& nextCluster{clustersNextLayer[iNextLayerClusterIndex]};
71 if (o2::gpu::GPUCommonMath::Abs(math_utils::smallestAngleDifference(currentCluster.phi, nextCluster.phi)) < phiCut) {
72 if (storedTracklets < maxTrackletsPerCluster) {
73 if constexpr (!EvalRun) {
74 if constexpr (Mode == TrackletMode::Layer0Layer1) {
75 tracklets[rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] + storedTracklets] = Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, targetRof, pivotRof};
76 } else {
77 tracklets[rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] + storedTracklets] = Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, pivotRof, targetRof};
78 }
79 }
80 ++storedTracklets;
81 }
82 }
83 }
84 }
85 }
86 if constexpr (EvalRun) {
87 foundTracklets[iCurrentLayerClusterIndex] += storedTracklets;
88 } else {
89 rofFoundTrackletsOffsets[iCurrentLayerClusterIndex] += storedTracklets;
90 }
91 }
92}
93
94static void trackletSelectionKernelHost(
95 const gsl::span<const Cluster> clusters0, // 0
96 const gsl::span<const Cluster> clusters1, // 1
97 gsl::span<unsigned char> usedClusters0, // Layer 0
98 gsl::span<unsigned char> usedClusters2, // Layer 2
99 const gsl::span<const Tracklet>& tracklets01,
100 const gsl::span<const Tracklet>& tracklets12,
101 bounded_vector<bool>& usedTracklets,
102 const gsl::span<int> foundTracklets01,
103 const gsl::span<int> foundTracklets12,
105 const gsl::span<const o2::MCCompLabel>& trackletLabels,
107 const short targetRofId0,
108 const short targetRofId2,
109 bool safeWrites = false,
110 const float tanLambdaCut = 0.025f,
111 const float phiCut = 0.005f,
112 const int maxTracklets = static_cast<int>(1e2))
113{
114 int offset01{0}, offset12{0};
115 for (unsigned int iCurrentLayerClusterIndex{0}; iCurrentLayerClusterIndex < clusters1.size(); ++iCurrentLayerClusterIndex) {
116 int validTracklets{0};
117 for (int iTracklet12{offset12}; iTracklet12 < offset12 + foundTracklets12[iCurrentLayerClusterIndex]; ++iTracklet12) {
118 for (int iTracklet01{offset01}; iTracklet01 < offset01 + foundTracklets01[iCurrentLayerClusterIndex]; ++iTracklet01) {
119 if (usedTracklets[iTracklet01]) {
120 continue;
121 }
122
123 const auto& tracklet01{tracklets01[iTracklet01]};
124 const auto& tracklet12{tracklets12[iTracklet12]};
125
126 if (tracklet01.rof[0] != targetRofId0 || tracklet12.rof[1] != targetRofId2) {
127 continue;
128 }
129
130 const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(tracklet01.tanLambda - tracklet12.tanLambda)};
131 const float deltaPhi{o2::gpu::GPUCommonMath::Abs(math_utils::smallestAngleDifference(tracklet01.phi, tracklet12.phi))};
132 if (deltaTanLambda < tanLambdaCut && deltaPhi < phiCut && validTracklets != maxTracklets) {
133 if (safeWrites) {
134 __atomic_store_n(&usedClusters0[tracklet01.firstClusterIndex], 1, __ATOMIC_RELAXED);
135 __atomic_store_n(&usedClusters2[tracklet12.secondClusterIndex], 1, __ATOMIC_RELAXED);
136 } else {
137 usedClusters0[tracklet01.firstClusterIndex] = 1;
138 usedClusters2[tracklet12.secondClusterIndex] = 1;
139 }
140 usedTracklets[iTracklet01] = true;
141 lines.emplace_back(tracklet01, clusters0.data(), clusters1.data());
142 if (!trackletLabels.empty()) {
143 linesLabels.emplace_back(trackletLabels[iTracklet01]);
144 }
145 ++validTracklets;
146 }
147 }
148 }
149 offset01 += foundTracklets01[iCurrentLayerClusterIndex];
150 offset12 += foundTracklets12[iCurrentLayerClusterIndex];
151 }
152}
153
154void VertexerTraits::updateVertexingParameters(const std::vector<VertexingParameters>& vrtPar, const TimeFrameGPUParameters& tfPar)
155{
156 mVrtParams = vrtPar;
158 for (auto& par : mVrtParams) {
159 par.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * par.phiCut / o2::constants::math::TwoPI));
160 par.zSpan = static_cast<int>(std::ceil(par.zCut * mIndexTableUtils.getInverseZCoordinate(0)));
161 }
162}
163
164// Main functions
165void VertexerTraits::computeTracklets(const int iteration)
166{
167 mTaskArena->execute([&] {
168 tbb::parallel_for(
169 tbb::blocked_range<short>(0, (short)mTimeFrame->getNrof()),
170 [&](const tbb::blocked_range<short>& Rofs) {
171 for (short pivotRofId = Rofs.begin(); pivotRofId < Rofs.end(); ++pivotRofId) {
172 bool skipROF = iteration && (int)mTimeFrame->getPrimaryVertices(pivotRofId).size() > mVrtParams[iteration].vertPerRofThreshold;
173 short startROF{std::max((short)0, static_cast<short>(pivotRofId - mVrtParams[iteration].deltaRof))};
174 short endROF{std::min(static_cast<short>(mTimeFrame->getNrof()), static_cast<short>(pivotRofId + mVrtParams[iteration].deltaRof + 1))};
175 for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) {
176 trackleterKernelHost<TrackletMode::Layer0Layer1, true>(
177 !skipROF ? mTimeFrame->getClustersOnLayer(targetRofId, 0) : gsl::span<Cluster>(), // Clusters to be matched with the next layer in target rof
178 !skipROF ? mTimeFrame->getClustersOnLayer(pivotRofId, 1) : gsl::span<Cluster>(), // Clusters to be matched with the current layer in pivot rof
179 mTimeFrame->getUsedClustersROF(targetRofId, 0), // Span of the used clusters in the target rof
180 mTimeFrame->getIndexTable(targetRofId, 0).data(), // Index table to access the data on the next layer in target rof
181 mVrtParams[iteration].phiCut,
182 mTimeFrame->getTracklets()[0], // Flat tracklet buffer
183 mTimeFrame->getNTrackletsCluster(pivotRofId, 0), // Span of the number of tracklets per each cluster in pivot rof
184 mIndexTableUtils,
185 pivotRofId,
186 targetRofId,
187 gsl::span<int>(), // Offset in the tracklet buffer
188 mVrtParams[iteration].maxTrackletsPerCluster);
189 trackleterKernelHost<TrackletMode::Layer1Layer2, true>(
190 !skipROF ? mTimeFrame->getClustersOnLayer(targetRofId, 2) : gsl::span<Cluster>(),
191 !skipROF ? mTimeFrame->getClustersOnLayer(pivotRofId, 1) : gsl::span<Cluster>(),
192 mTimeFrame->getUsedClustersROF(targetRofId, 2),
193 mTimeFrame->getIndexTable(targetRofId, 2).data(),
194 mVrtParams[iteration].phiCut,
195 mTimeFrame->getTracklets()[1],
196 mTimeFrame->getNTrackletsCluster(pivotRofId, 1), // Span of the number of tracklets per each cluster in pivot rof
197 mIndexTableUtils,
198 pivotRofId,
199 targetRofId,
200 gsl::span<int>(), // Offset in the tracklet buffer
201 mVrtParams[iteration].maxTrackletsPerCluster);
202 }
203 mTimeFrame->getNTrackletsROF(pivotRofId, 0) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 0).end(), 0);
204 mTimeFrame->getNTrackletsROF(pivotRofId, 1) = std::accumulate(mTimeFrame->getNTrackletsCluster(pivotRofId, 1).begin(), mTimeFrame->getNTrackletsCluster(pivotRofId, 1).end(), 0);
205 }
206 });
207
209 if (auto tot0 = mTimeFrame->getTotalTrackletsTF(0), tot1 = mTimeFrame->getTotalTrackletsTF(1);
210 tot0 == 0 || tot1 == 0) {
211 return;
212 } else {
213 mTimeFrame->getTracklets()[0].resize(tot0);
214 mTimeFrame->getTracklets()[1].resize(tot1);
215 }
216
217 tbb::parallel_for(
218 tbb::blocked_range<short>(0, (short)mTimeFrame->getNrof()),
219 [&](const tbb::blocked_range<short>& Rofs) {
220 for (short pivotRofId = Rofs.begin(); pivotRofId < Rofs.end(); ++pivotRofId) {
221 bool skipROF = iteration && (int)mTimeFrame->getPrimaryVertices(pivotRofId).size() > mVrtParams[iteration].vertPerRofThreshold;
222 short startROF{std::max((short)0, static_cast<short>(pivotRofId - mVrtParams[iteration].deltaRof))};
223 short endROF{std::min(static_cast<short>(mTimeFrame->getNrof()), static_cast<short>(pivotRofId + mVrtParams[iteration].deltaRof + 1))};
224 auto mobileOffset0 = mTimeFrame->getNTrackletsROF(pivotRofId, 0);
225 auto mobileOffset1 = mTimeFrame->getNTrackletsROF(pivotRofId, 1);
226 for (auto targetRofId = startROF; targetRofId < endROF; ++targetRofId) {
227 trackleterKernelHost<TrackletMode::Layer0Layer1, false>(
228 !skipROF ? mTimeFrame->getClustersOnLayer(targetRofId, 0) : gsl::span<Cluster>(),
229 !skipROF ? mTimeFrame->getClustersOnLayer(pivotRofId, 1) : gsl::span<Cluster>(),
230 mTimeFrame->getUsedClustersROF(targetRofId, 0),
231 mTimeFrame->getIndexTable(targetRofId, 0).data(),
232 mVrtParams[iteration].phiCut,
233 mTimeFrame->getTracklets()[0],
234 mTimeFrame->getNTrackletsCluster(pivotRofId, 0),
235 mIndexTableUtils,
236 pivotRofId,
237 targetRofId,
238 mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 0),
239 mVrtParams[iteration].maxTrackletsPerCluster);
240 trackleterKernelHost<TrackletMode::Layer1Layer2, false>(
241 !skipROF ? mTimeFrame->getClustersOnLayer(targetRofId, 2) : gsl::span<Cluster>(),
242 !skipROF ? mTimeFrame->getClustersOnLayer(pivotRofId, 1) : gsl::span<Cluster>(),
243 mTimeFrame->getUsedClustersROF(targetRofId, 2),
244 mTimeFrame->getIndexTable(targetRofId, 2).data(),
245 mVrtParams[iteration].phiCut,
246 mTimeFrame->getTracklets()[1],
247 mTimeFrame->getNTrackletsCluster(pivotRofId, 1),
248 mIndexTableUtils,
249 pivotRofId,
250 targetRofId,
251 mTimeFrame->getExclusiveNTrackletsCluster(pivotRofId, 1),
252 mVrtParams[iteration].maxTrackletsPerCluster);
253 }
254 }
255 });
256 });
257
259 if (mTimeFrame->hasMCinformation()) {
260 for (const auto& trk : mTimeFrame->getTracklets()[0]) {
262 if (!trk.isEmpty()) {
263 int sortedId0{mTimeFrame->getSortedIndex(trk.rof[0], 0, trk.firstClusterIndex)};
264 int sortedId1{mTimeFrame->getSortedIndex(trk.rof[1], 1, trk.secondClusterIndex)};
265 for (const auto& lab0 : mTimeFrame->getClusterLabels(0, mTimeFrame->getClusters()[0][sortedId0].clusterId)) {
266 for (const auto& lab1 : mTimeFrame->getClusterLabels(1, mTimeFrame->getClusters()[1][sortedId1].clusterId)) {
267 if (lab0 == lab1 && lab0.isValid()) {
268 label = lab0;
269 break;
270 }
271 }
272 if (label.isValid()) {
273 break;
274 }
275 }
276 }
277 mTimeFrame->getTrackletsLabel(0).emplace_back(label);
278 }
279 }
280
281#ifdef VTX_DEBUG
282 debugComputeTracklets(iteration);
283#endif
284}
285
287{
288 mTaskArena->execute([&] {
289 tbb::parallel_for(
290 tbb::blocked_range<short>(0, (short)mTimeFrame->getNrof()),
291 [&](const tbb::blocked_range<short>& Rofs) {
292 for (short pivotRofId = Rofs.begin(); pivotRofId < Rofs.end(); ++pivotRofId) {
293 if (iteration && (int)mTimeFrame->getPrimaryVertices(pivotRofId).size() > mVrtParams[iteration].vertPerRofThreshold) {
294 continue;
295 }
296 if (mTimeFrame->getFoundTracklets(pivotRofId, 0).empty()) {
297 continue;
298 }
299 mTimeFrame->getLines(pivotRofId).reserve(mTimeFrame->getNTrackletsCluster(pivotRofId, 0).size());
300 bounded_vector<bool> usedTracklets(mTimeFrame->getFoundTracklets(pivotRofId, 0).size(), false, mMemoryPool.get());
301 short startROF{std::max((short)0, static_cast<short>(pivotRofId - mVrtParams[iteration].deltaRof))};
302 short endROF{std::min(static_cast<short>(mTimeFrame->getNrof()), static_cast<short>(pivotRofId + mVrtParams[iteration].deltaRof + 1))};
303
304 // needed only if multi-threaded using deltaRof and only at the overlap edges of the ranges
305 bool safeWrite = mTaskArena->max_concurrency() > 1 && mVrtParams[iteration].deltaRof != 0 && ((Rofs.begin() - startROF < 0) || (endROF - Rofs.end() > 0));
306
307 for (short targetRofId0 = startROF; targetRofId0 < endROF; ++targetRofId0) {
308 for (short targetRofId2 = startROF; targetRofId2 < endROF; ++targetRofId2) {
309 if (std::abs(targetRofId0 - targetRofId2) > mVrtParams[iteration].deltaRof) { // do not allow over 3 ROFs
310 continue;
311 }
312 trackletSelectionKernelHost(
313 mTimeFrame->getClustersOnLayer(targetRofId0, 0),
314 mTimeFrame->getClustersOnLayer(pivotRofId, 1),
315 mTimeFrame->getUsedClustersROF(targetRofId0, 0),
316 mTimeFrame->getUsedClustersROF(targetRofId2, 2),
317 mTimeFrame->getFoundTracklets(pivotRofId, 0),
318 mTimeFrame->getFoundTracklets(pivotRofId, 1),
319 usedTracklets,
320 mTimeFrame->getNTrackletsCluster(pivotRofId, 0),
321 mTimeFrame->getNTrackletsCluster(pivotRofId, 1),
322 mTimeFrame->getLines(pivotRofId),
323 mTimeFrame->getLabelsFoundTracklets(pivotRofId, 0),
324 mTimeFrame->getLinesLabel(pivotRofId),
325 targetRofId0,
326 targetRofId2,
327 safeWrite,
328 mVrtParams[iteration].tanLambdaCut,
329 mVrtParams[iteration].phiCut);
330 }
331 }
332 }
333 });
334 });
335
336#ifdef VTX_DEBUG
337 debugComputeTrackletMatching(iteration);
338#endif
339
340 // from here on we do not use tracklets from L1-2 anymore, so let's free them
341 deepVectorClear(mTimeFrame->getTracklets()[1]);
342}
343
344void VertexerTraits::computeVertices(const int iteration)
345{
346 auto nsigmaCut{std::min(mVrtParams[iteration].vertNsigmaCut * mVrtParams[iteration].vertNsigmaCut * (mVrtParams[iteration].vertRadiusSigma * mVrtParams[iteration].vertRadiusSigma + mVrtParams[iteration].trackletSigma * mVrtParams[iteration].trackletSigma), 1.98f)};
347 bounded_vector<Vertex> vertices(mMemoryPool.get());
348 bounded_vector<std::pair<o2::MCCompLabel, float>> polls(mMemoryPool.get());
349 bounded_vector<o2::MCCompLabel> contLabels(mMemoryPool.get());
350 bounded_vector<int> noClustersVec(mTimeFrame->getNrof(), 0, mMemoryPool.get());
351 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
352 if (iteration && (int)mTimeFrame->getPrimaryVertices(rofId).size() > mVrtParams[iteration].vertPerRofThreshold) {
353 continue;
354 }
355 const int numTracklets{static_cast<int>(mTimeFrame->getLines(rofId).size())};
356
357 bounded_vector<bool> usedTracklets(numTracklets, false, mMemoryPool.get());
358 for (int line1{0}; line1 < numTracklets; ++line1) {
359 if (usedTracklets[line1]) {
360 continue;
361 }
362 for (int line2{line1 + 1}; line2 < numTracklets; ++line2) {
363 if (usedTracklets[line2]) {
364 continue;
365 }
366 auto dca{Line::getDCA(mTimeFrame->getLines(rofId)[line1], mTimeFrame->getLines(rofId)[line2])};
367 if (dca < mVrtParams[iteration].pairCut) {
368 mTimeFrame->getTrackletClusters(rofId).emplace_back(line1, mTimeFrame->getLines(rofId)[line1], line2, mTimeFrame->getLines(rofId)[line2]);
369 std::array<float, 3> tmpVertex{mTimeFrame->getTrackletClusters(rofId).back().getVertex()};
370 if (tmpVertex[0] * tmpVertex[0] + tmpVertex[1] * tmpVertex[1] > 4.f) {
371 mTimeFrame->getTrackletClusters(rofId).pop_back();
372 break;
373 }
374 usedTracklets[line1] = true;
375 usedTracklets[line2] = true;
376 for (int tracklet3{0}; tracklet3 < numTracklets; ++tracklet3) {
377 if (usedTracklets[tracklet3]) {
378 continue;
379 }
380 if (Line::getDistanceFromPoint(mTimeFrame->getLines(rofId)[tracklet3], tmpVertex) < mVrtParams[iteration].pairCut) {
381 mTimeFrame->getTrackletClusters(rofId).back().add(tracklet3, mTimeFrame->getLines(rofId)[tracklet3]);
382 usedTracklets[tracklet3] = true;
383 tmpVertex = mTimeFrame->getTrackletClusters(rofId).back().getVertex();
384 }
385 }
386 break;
387 }
388 }
389 }
390 if (mVrtParams[iteration].allowSingleContribClusters) {
391 auto beamLine = Line{{mTimeFrame->getBeamX(), mTimeFrame->getBeamY(), -50.f}, {mTimeFrame->getBeamX(), mTimeFrame->getBeamY(), 50.f}}; // use beam position as contributor
392 for (size_t iLine{0}; iLine < numTracklets; ++iLine) {
393 if (!usedTracklets[iLine]) {
394 auto dca = Line::getDCA(mTimeFrame->getLines(rofId)[iLine], beamLine);
395 if (dca < mVrtParams[iteration].pairCut) {
396 mTimeFrame->getTrackletClusters(rofId).emplace_back(iLine, mTimeFrame->getLines(rofId)[iLine], -1, beamLine); // beamline must be passed as second line argument
397 }
398 }
399 }
400 }
401
402 // Cluster merging
403 std::sort(mTimeFrame->getTrackletClusters(rofId).begin(), mTimeFrame->getTrackletClusters(rofId).end(),
404 [](ClusterLines& cluster1, ClusterLines& cluster2) { return cluster1.getSize() > cluster2.getSize(); });
405 noClustersVec[rofId] = static_cast<int>(mTimeFrame->getTrackletClusters(rofId).size());
406 for (int iCluster1{0}; iCluster1 < noClustersVec[rofId]; ++iCluster1) {
407 std::array<float, 3> vertex1{mTimeFrame->getTrackletClusters(rofId)[iCluster1].getVertex()};
408 std::array<float, 3> vertex2{};
409 for (int iCluster2{iCluster1 + 1}; iCluster2 < noClustersVec[rofId]; ++iCluster2) {
410 vertex2 = mTimeFrame->getTrackletClusters(rofId)[iCluster2].getVertex();
411 if (o2::gpu::GPUCommonMath::Abs(vertex1[2] - vertex2[2]) < mVrtParams[iteration].clusterCut) {
412 float distance{(vertex1[0] - vertex2[0]) * (vertex1[0] - vertex2[0]) +
413 (vertex1[1] - vertex2[1]) * (vertex1[1] - vertex2[1]) +
414 (vertex1[2] - vertex2[2]) * (vertex1[2] - vertex2[2])};
415 if (distance < mVrtParams[iteration].pairCut * mVrtParams[iteration].pairCut) {
416 for (auto label : mTimeFrame->getTrackletClusters(rofId)[iCluster2].getLabels()) {
417 mTimeFrame->getTrackletClusters(rofId)[iCluster1].add(label, mTimeFrame->getLines(rofId)[label]);
418 vertex1 = mTimeFrame->getTrackletClusters(rofId)[iCluster1].getVertex();
419 }
420 mTimeFrame->getTrackletClusters(rofId).erase(mTimeFrame->getTrackletClusters(rofId).begin() + iCluster2);
421 --iCluster2;
422 --noClustersVec[rofId];
423 }
424 }
425 }
426 }
427 }
428 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
429 std::sort(mTimeFrame->getTrackletClusters(rofId).begin(), mTimeFrame->getTrackletClusters(rofId).end(),
430 [](const ClusterLines& cluster1, const ClusterLines& cluster2) { return cluster1.getSize() > cluster2.getSize(); }); // ensure clusters are ordered by contributors, so that we can cat after the first.
431 bool atLeastOneFound{false};
432 for (int iCluster{0}; iCluster < noClustersVec[rofId]; ++iCluster) {
433 bool lowMultCandidate{false};
434 double beamDistance2{(mTimeFrame->getBeamX() - mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[0]) * (mTimeFrame->getBeamX() - mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[0]) +
435 (mTimeFrame->getBeamY() - mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[1]) * (mTimeFrame->getBeamY() - mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[1])};
436 if (atLeastOneFound && (lowMultCandidate = mTimeFrame->getTrackletClusters(rofId)[iCluster].getSize() < mVrtParams[iteration].clusterContributorsCut)) { // We might have pile up with nContr > cut.
437 lowMultCandidate &= (beamDistance2 < mVrtParams[iteration].lowMultBeamDistCut * mVrtParams[iteration].lowMultBeamDistCut);
438 if (!lowMultCandidate) { // Not the first cluster and not a low multiplicity candidate, we can remove it
439 mTimeFrame->getTrackletClusters(rofId).erase(mTimeFrame->getTrackletClusters(rofId).begin() + iCluster);
440 noClustersVec[rofId]--;
441 continue;
442 }
443 }
444
445 if (beamDistance2 < nsigmaCut && o2::gpu::GPUCommonMath::Abs(mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[2]) < mVrtParams[iteration].maxZPositionAllowed) {
446 atLeastOneFound = true;
447 auto& vertex = vertices.emplace_back(o2::math_utils::Point3D<float>(mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[0],
448 mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[1],
449 mTimeFrame->getTrackletClusters(rofId)[iCluster].getVertex()[2]),
450 mTimeFrame->getTrackletClusters(rofId)[iCluster].getRMS2(), // Symm matrix. Diagonal: RMS2 components,
451 // off-diagonal: square mean of projections on planes.
452 mTimeFrame->getTrackletClusters(rofId)[iCluster].getSize(), // Contributors
453 mTimeFrame->getTrackletClusters(rofId)[iCluster].getAvgDistance2()); // In place of chi2
454
455 if (iteration) {
456 vertex.setFlags(Vertex::UPCMode);
457 }
458 vertex.setTimeStamp(mTimeFrame->getTrackletClusters(rofId)[iCluster].getROF());
460 bounded_vector<o2::MCCompLabel> labels(mMemoryPool.get());
461 for (auto& index : mTimeFrame->getTrackletClusters(rofId)[iCluster].getLabels()) {
462 labels.push_back(mTimeFrame->getLinesLabel(rofId)[index]); // then we can use nContributors from vertices to get the labels
463 }
464 polls.push_back(computeMain(labels));
465 if (mVrtParams[iteration].outputContLabels) {
466 contLabels.insert(contLabels.end(), labels.begin(), labels.end());
467 }
468 }
469 }
470 }
471 if (!iteration) {
472 mTimeFrame->addPrimaryVertices(vertices, iteration);
475 if (mVrtParams[iteration].outputContLabels) {
477 }
478 }
479 } else {
480 mTimeFrame->addPrimaryVerticesInROF(vertices, rofId, iteration);
483 if (mVrtParams[iteration].outputContLabels) {
485 }
486 }
487 }
488 if (vertices.empty() && !(iteration && (int)mTimeFrame->getPrimaryVertices(rofId).size() > mVrtParams[iteration].vertPerRofThreshold)) {
490 }
491 vertices.clear();
492 polls.clear();
493 }
494
495#ifdef VTX_DEBUG
496 debugComputeVertices(iteration);
497#endif
498}
499
501{
502 LOGP(info, "Using truth seeds as vertices; will skip computations");
504 const auto dc = o2::steer::DigitizationContext::loadFromFile("collisioncontext.root");
505 const auto irs = dc->getEventRecords();
506 int64_t roFrameBiasInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameBiasInBC;
507 int64_t roFrameLengthInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameLengthInBC;
509 struct VertInfo {
510 bounded_vector<Vertex> vertices;
512 bounded_vector<int> events;
513 };
514 std::map<int, VertInfo> vertices;
515 for (int iSrc{0}; iSrc < mcReader.getNSources(); ++iSrc) {
516 auto eveId2colId = dc->getCollisionIndicesForSource(iSrc);
517 for (int iEve{0}; iEve < mcReader.getNEvents(iSrc); ++iEve) {
518 const auto& ir = irs[eveId2colId[iEve]];
519 if (!ir.isDummy()) { // do we need this, is this for diffractive events?
520 const auto& eve = mcReader.getMCEventHeader(iSrc, iEve);
521 int rofId = ((ir - raw::HBFUtils::Instance().getFirstSampledTFIR()).toLong() - roFrameBiasInBC) / roFrameLengthInBC;
522 if (!vertices.contains(rofId)) {
523 vertices[rofId] = {
524 .vertices = bounded_vector<Vertex>(mMemoryPool.get()),
525 .srcs = bounded_vector<int>(mMemoryPool.get()),
526 .events = bounded_vector<int>(mMemoryPool.get()),
527 };
528 }
529 Vertex vert;
530 vert.setTimeStamp(rofId);
531 vert.setNContributors(std::ranges::count_if(mcReader.getTracks(iSrc, iEve), [](const auto& trk) {
532 return trk.isPrimary() && trk.GetPt() > 0.2 && std::abs(trk.GetEta()) < 1.3;
533 }));
534 vert.setXYZ((float)eve.GetX(), (float)eve.GetY(), (float)eve.GetZ());
535 vert.setChi2(1);
536 constexpr float cov = 50e-9;
537 vert.setCov(cov, cov, cov, cov, cov, cov);
538 vertices[rofId].vertices.push_back(vert);
539 vertices[rofId].srcs.push_back(iSrc);
540 vertices[rofId].events.push_back(iEve);
541 }
542 }
543 }
544 size_t nVerts{0};
545 for (int iROF{0}; iROF < mTimeFrame->getNrof(); ++iROF) {
546 bounded_vector<Vertex> verts(mMemoryPool.get());
547 bounded_vector<std::pair<o2::MCCompLabel, float>> polls(mMemoryPool.get());
548 if (vertices.contains(iROF)) {
549 const auto& vertInfo = vertices[iROF];
550 verts = vertInfo.vertices;
551 nVerts += verts.size();
552 for (size_t i{0}; i < verts.size(); ++i) {
553 o2::MCCompLabel lbl(o2::MCCompLabel::maxTrackID(), vertInfo.events[i], vertInfo.srcs[i], false);
554 polls.emplace_back(lbl, 1.f);
555 }
556 } else {
558 }
561 }
562 LOGP(info, "Found {}/{} ROFs with {} vertices -> <NV>={:.2f}", vertices.size(), mTimeFrame->getNrof(), nVerts, (float)nVerts / (float)vertices.size());
563}
564
565void VertexerTraits::setNThreads(int n, std::shared_ptr<tbb::task_arena>& arena)
566{
567#if defined(VTX_DEBUG)
568 LOGP(info, "Vertexer with debug output forcing single thread");
569 mTaskArena = std::make_shared<tbb::task_arena>(1);
570#else
571 if (arena == nullptr) {
572 mTaskArena = std::make_shared<tbb::task_arena>(std::abs(n));
573 LOGP(info, "Setting seeding vertexer with {} threads.", n);
574 } else {
575 mTaskArena = arena;
576 LOGP(info, "Attaching vertexer to calling thread's arena");
577 }
578#endif
579}
580
581void VertexerTraits::debugComputeTracklets(int iteration)
582{
583 auto stream = new utils::TreeStreamRedirector("artefacts_tf.root", "recreate");
584 LOGP(info, "writing debug output for computeTracklets");
585 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
586 const auto& strk0 = mTimeFrame->getFoundTracklets(rofId, 0);
587 std::vector<Tracklet> trk0(strk0.begin(), strk0.end());
588 const auto& strk1 = mTimeFrame->getFoundTracklets(rofId, 1);
589 std::vector<Tracklet> trk1(strk1.begin(), strk1.end());
590 (*stream) << "tracklets"
591 << "Tracklets0=" << trk0
592 << "Tracklets1=" << trk1
593 << "iteration=" << iteration
594 << "\n";
595 }
596 stream->Close();
597 delete stream;
598}
599
600void VertexerTraits::debugComputeTrackletMatching(int iteration)
601{
602 auto stream = new utils::TreeStreamRedirector("artefacts_tf.root", "update");
603 LOGP(info, "writing debug output for computeTrackletMatching");
604 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
605 (*stream) << "lines"
606 << "Lines=" << toSTDVector(mTimeFrame->getLines(rofId))
607 << "NTrackletCluster01=" << mTimeFrame->getNTrackletsCluster(rofId, 0)
608 << "NTrackletCluster12=" << mTimeFrame->getNTrackletsCluster(rofId, 1)
609 << "iteration=" << iteration
610 << "\n";
611 }
612
614 LOGP(info, "\tdumping also MC information");
615 const auto dc = o2::steer::DigitizationContext::loadFromFile("collisioncontext.root");
616 const auto irs = dc->getEventRecords();
617 int64_t roFrameBiasInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameBiasInBC;
618 int64_t roFrameLengthInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameLengthInBC;
620
621 std::map<int, int> eve2BcInROF, bcInRofNEve;
622 for (int iSrc{0}; iSrc < mcReader.getNSources(); ++iSrc) {
623 auto eveId2colId = dc->getCollisionIndicesForSource(iSrc);
624 for (int iEve{0}; iEve < mcReader.getNEvents(iSrc); ++iEve) {
625 const auto& ir = irs[eveId2colId[iEve]];
626 if (!ir.isDummy()) { // do we need this, is this for diffractive events?
627 const auto& eve = mcReader.getMCEventHeader(iSrc, iEve);
628 const int bcInROF = ((ir - raw::HBFUtils::Instance().getFirstSampledTFIR()).toLong() - roFrameBiasInBC) % roFrameLengthInBC;
629 eve2BcInROF[iEve] = bcInROF;
630 ++bcInRofNEve[bcInROF];
631 }
632 }
633 }
634
635 std::unordered_map<int, int> bcROFNTracklets01, bcROFNTracklets12;
636 std::vector<std::vector<int>> tracklet01BC, tracklet12BC;
637 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
638 { // 0-1
639 const auto& tracklet01 = mTimeFrame->getFoundTracklets(rofId, 0);
640 const auto& lbls01 = mTimeFrame->getLabelsFoundTracklets(rofId, 0);
641 auto& trkls01 = tracklet01BC.emplace_back();
642 for (int iTrklt{0}; iTrklt < (int)tracklet01.size(); ++iTrklt) {
643 const auto& tracklet = tracklet01[iTrklt];
644 const auto& lbl = lbls01[iTrklt];
645 if (lbl.isCorrect()) {
646 ++bcROFNTracklets01[eve2BcInROF[lbl.getEventID()]];
647 trkls01.push_back(eve2BcInROF[lbl.getEventID()]);
648 } else {
649 trkls01.push_back(-1);
650 }
651 }
652 }
653 { // 1-2 computed on the fly!
654 const auto& tracklet12 = mTimeFrame->getFoundTracklets(rofId, 1);
655 auto& trkls12 = tracklet12BC.emplace_back();
656 for (int iTrklt{0}; iTrklt < (int)tracklet12.size(); ++iTrklt) {
657 const auto& tracklet = tracklet12[iTrklt];
659
660 int sortedId1{mTimeFrame->getSortedIndex(tracklet.rof[0], 1, tracklet.firstClusterIndex)};
661 int sortedId2{mTimeFrame->getSortedIndex(tracklet.rof[1], 2, tracklet.secondClusterIndex)};
662 for (const auto& lab1 : mTimeFrame->getClusterLabels(1, mTimeFrame->getClusters()[1][sortedId1].clusterId)) {
663 for (const auto& lab2 : mTimeFrame->getClusterLabels(2, mTimeFrame->getClusters()[2][sortedId2].clusterId)) {
664 if (lab1 == lab2 && lab1.isValid()) {
665 label = lab1;
666 break;
667 }
668 }
669 if (label.isValid()) {
670 break;
671 }
672 }
673
674 if (label.isCorrect()) {
675 ++bcROFNTracklets12[eve2BcInROF[label.getEventID()]];
676 trkls12.push_back(eve2BcInROF[label.getEventID()]);
677 } else {
678 trkls12.push_back(-1);
679 }
680 }
681 }
682 }
683 LOGP(info, "\tdumping ntracklets/RofBC ({})", bcInRofNEve.size());
684 for (const auto& [bcInRof, neve] : bcInRofNEve) {
685 (*stream) << "ntracklets"
686 << "bcInROF=" << bcInRof
687 << "ntrkl01=" << bcROFNTracklets01[bcInRof]
688 << "ntrkl12=" << bcROFNTracklets12[bcInRof]
689 << "neve=" << neve
690 << "iteration=" << iteration
691 << "\n";
692 }
693
694 std::unordered_map<int, int> bcROFNLines;
695 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
696 const auto& lines = mTimeFrame->getLines(rofId);
697 const auto& lbls = mTimeFrame->getLinesLabel(rofId);
698 for (int iLine{0}; iLine < (int)lines.size(); ++iLine) {
699 const auto& line = lines[iLine];
700 const auto& lbl = lbls[iLine];
701 if (lbl.isCorrect()) {
702 ++bcROFNLines[eve2BcInROF[lbl.getEventID()]];
703 }
704 }
705 }
706
707 LOGP(info, "\tdumping nlines/RofBC");
708 for (const auto& [bcInRof, neve] : bcInRofNEve) {
709 (*stream) << "nlines"
710 << "bcInROF=" << bcInRof
711 << "nline=" << bcROFNLines[bcInRof]
712 << "neve=" << neve
713 << "iteration=" << iteration
714 << "\n";
715 }
716 }
717 stream->Close();
718 delete stream;
719}
720
721void VertexerTraits::debugComputeVertices(int iteration)
722{
723 auto stream = new utils::TreeStreamRedirector("artefacts_tf.root", "update");
724 LOGP(info, "writing debug output for computeVertices");
725 for (auto rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
726 (*stream) << "clusterlines"
727 << "clines_post=" << toSTDVector(mTimeFrame->getTrackletClusters(rofId))
728 << "iteration=" << iteration
729 << "\n";
730 }
731
733 LOGP(info, "\tdumping also MC information");
734 const auto dc = o2::steer::DigitizationContext::loadFromFile("collisioncontext.root");
735 const auto irs = dc->getEventRecords();
736 int64_t roFrameBiasInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameBiasInBC;
737 int64_t roFrameLengthInBC = o2::itsmft::DPLAlpideParam<o2::detectors::DetID::ITS>::Instance().roFrameLengthInBC;
739
740 std::map<int, int> eve2BcInROF, bcInRofNEve;
741 for (int iSrc{0}; iSrc < mcReader.getNSources(); ++iSrc) {
742 auto eveId2colId = dc->getCollisionIndicesForSource(iSrc);
743 for (int iEve{0}; iEve < mcReader.getNEvents(iSrc); ++iEve) {
744 const auto& ir = irs[eveId2colId[iEve]];
745 if (!ir.isDummy()) { // do we need this, is this for diffractive events?
746 const auto& eve = mcReader.getMCEventHeader(iSrc, iEve);
747 const int bcInROF = ((ir - raw::HBFUtils::Instance().getFirstSampledTFIR()).toLong() - roFrameBiasInBC) % roFrameLengthInBC;
748 eve2BcInROF[iEve] = bcInROF;
749 ++bcInRofNEve[bcInROF];
750 }
751 }
752 }
753
754 std::unordered_map<int, int> bcROFNVtx;
755 std::unordered_map<int, float> bcROFNPur;
756 std::unordered_map<o2::MCCompLabel, size_t> uniqueVertices;
757 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
758 const auto& pvs = mTimeFrame->getPrimaryVertices(rofId);
759 const auto& lblspv = mTimeFrame->getPrimaryVerticesMCRecInfo(rofId);
760 for (int i{0}; i < (int)pvs.size(); ++i) {
761 const auto& pv = pvs[i];
762 const auto& [lbl, pur] = lblspv[i];
763 if (lbl.isCorrect()) {
764 ++uniqueVertices[lbl];
765 ++bcROFNVtx[eve2BcInROF[lbl.getEventID()]];
766 bcROFNPur[eve2BcInROF[lbl.getEventID()]] += pur;
767 }
768 }
769 }
770
771 std::unordered_map<int, int> bcROFNUVtx, bcROFNCVtx;
772 for (const auto& [k, _] : eve2BcInROF) {
773 bcROFNUVtx[k] = bcROFNCVtx[k] = 0;
774 }
775
776 for (const auto& [lbl, c] : uniqueVertices) {
777 if (c <= 1) {
778 ++bcROFNUVtx[eve2BcInROF[lbl.getEventID()]];
779 } else {
780 ++bcROFNCVtx[eve2BcInROF[lbl.getEventID()]];
781 }
782 }
783
784 LOGP(info, "\tdumping nvtx/RofBC");
785 for (const auto& [bcInRof, neve] : bcInRofNEve) {
786 (*stream) << "nvtx"
787 << "bcInROF=" << bcInRof
788 << "nvtx=" << bcROFNVtx[bcInRof] // all vertices
789 << "nuvtx=" << bcROFNUVtx[bcInRof] // unique vertices
790 << "ncvtx=" << bcROFNCVtx[bcInRof] // cloned vertices
791 << "npur=" << bcROFNPur[bcInRof]
792 << "neve=" << neve
793 << "iteration=" << iteration
794 << "\n";
795 }
796
797 // check dist of clones
798 std::unordered_map<o2::MCCompLabel, std::vector<Vertex>> cVtx;
799 for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
800 const auto& pvs = mTimeFrame->getPrimaryVertices(rofId);
801 const auto& lblspv = mTimeFrame->getPrimaryVerticesMCRecInfo(rofId);
802 for (int i{0}; i < (int)pvs.size(); ++i) {
803 const auto& pv = pvs[i];
804 const auto& [lbl, pur] = lblspv[i];
805 if (lbl.isCorrect() && uniqueVertices.contains(lbl) && uniqueVertices[lbl] > 1) {
806 if (!cVtx.contains(lbl)) {
807 cVtx[lbl] = std::vector<Vertex>();
808 }
809 cVtx[lbl].push_back(pv);
810 }
811 }
812 }
813
814 for (auto& [_, vertices] : cVtx) {
815 std::sort(vertices.begin(), vertices.end(), [](const Vertex& a, const Vertex& b) { return a.getNContributors() > b.getNContributors(); });
816 for (int i{0}; i < (int)vertices.size(); ++i) {
817 const auto vtx = vertices[i];
818 (*stream) << "cvtx"
819 << "vertex=" << vtx
820 << "i=" << i
821 << "dx=" << vertices[0].getX() - vtx.getX()
822 << "dy=" << vertices[0].getY() - vtx.getY()
823 << "dz=" << vertices[0].getZ() - vtx.getZ()
824 << "drof=" << vertices[0].getTimeStamp().getTimeStamp() - vtx.getTimeStamp().getTimeStamp()
825 << "dnc=" << vertices[0].getNContributors() - vtx.getNContributors()
826 << "iteration=" << iteration
827 << "\n";
828 }
829 }
830 }
831 stream->Close();
832 delete stream;
833}
uint64_t vertex
Definition RawEventData.h:9
int32_t i
Mode
Definition Utils.h:89
uint32_t c
Definition RawData.h:2
Class to compute the primary vertex in ITS from tracklets.
static constexpr int maxTrackID()
void setTrackingParameters(const T &params)
float getInverseZCoordinate(const int layerIndex) const
virtual void computeTrackletMatching(const int iteration=0)
virtual void computeTracklets(const int iteration=0)
void setNThreads(int n, std::shared_ptr< tbb::task_arena > &arena)
virtual void computeVertices(const int iteration=0)
std::vector< VertexingParameters > mVrtParams
IndexTableUtils mIndexTableUtils
virtual void updateVertexingParameters(const std::vector< VertexingParameters > &vrtPar, const TimeFrameGPUParameters &gpuTfPar)
static std::pair< o2::MCCompLabel, float > computeMain(const bounded_vector< o2::MCCompLabel > &elements)
static DigitizationContext * loadFromFile(std::string_view filename="")
size_t getNEvents(int source) const
Get number of events.
o2::dataformats::MCEventHeader const & getMCEventHeader(int source, int event) const
retrieves the MCEventHeader for a given eventID and sourceID
size_t getNSources() const
Get number of sources.
std::vector< MCTrack > const & getTracks(int source, int event) const
variant returning all tracks for source and event at once
GLdouble n
Definition glcorearb.h:1982
GLuint index
Definition glcorearb.h:781
GLboolean GLboolean GLboolean b
Definition glcorearb.h:1233
GLsizei GLsizei GLfloat distance
Definition glcorearb.h:5506
GLuint GLsizei const GLchar * label
Definition glcorearb.h:2519
GLboolean GLboolean GLboolean GLboolean a
Definition glcorearb.h:1233
GLuint GLuint stream
Definition glcorearb.h:1806
constexpr float TwoPI
int32_t const char int32_t line
void deepVectorClear(std::vector< T > &vec)
std::pmr::vector< T > bounded_vector
std::vector< T > toSTDVector(const bounded_vector< T > &b)
std::vector< Cluster > getClusters(int event)
Common utility functions.
int32_t w
const Cluster const Cluster *static float getDistanceFromPoint(const Line &line, const std::array< float, 3 > &point)
gsl::span< const Tracklet > getFoundTracklets(int rofId, int combId) const
Definition TimeFrame.h:604
int getNrof() const
Definition TimeFrame.h:102
void addPrimaryVertices(const bounded_vector< Vertex > &vertices, const int iteration)
Definition TimeFrame.cxx:62
float getBeamX() const
Definition TimeFrame.h:111
bool hasMCinformation() const
Definition TimeFrame.h:148
unsigned int & getNoVertexROF()
Definition TimeFrame.h:217
gsl::span< const std::pair< MCCompLabel, float > > getPrimaryVerticesMCRecInfo(const int rofId) const
Definition TimeFrame.h:368
void addPrimaryVerticesInROF(const bounded_vector< Vertex > &vertices, const int rofId, const int iteration)
Definition TimeFrame.cxx:90
auto & getTracklets()
Definition TimeFrame.h:161
auto & getTrackletClusters(int rofId)
Definition TimeFrame.h:208
void computeTrackletsPerROFScans()
void addPrimaryVerticesLabelsInROF(const bounded_vector< std::pair< MCCompLabel, float > > &labels, const int rofId)
void addPrimaryVerticesContributorLabelsInROF(const bounded_vector< MCCompLabel > &labels, const int rofId)
void addPrimaryVerticesLabels(bounded_vector< std::pair< MCCompLabel, float > > &labels)
Definition TimeFrame.cxx:78
gsl::span< const MCCompLabel > getLabelsFoundTracklets(int rofId, int combId) const
Definition TimeFrame.h:614
gsl::span< const Vertex > getPrimaryVertices(int rofId) const
Definition TimeFrame.h:356
auto & getLinesLabel(const int rofId)
Definition TimeFrame.h:175
float getBeamY() const
Definition TimeFrame.h:112
gsl::span< int > getNTrackletsCluster(int rofId, int combId)
Definition TimeFrame.h:573
uint32_t getTotalTrackletsTF(const int iLayer)
Definition TimeFrame.h:214
void addPrimaryVerticesContributorLabels(bounded_vector< MCCompLabel > &labels)
Definition TimeFrame.cxx:84
auto & getLines(int rofId)
Definition TimeFrame.h:203
int getSortedIndex(int rofId, int layer, int idx) const
Definition TimeFrame.h:100
IR getFirstSampledTFIR() const
get TF and HB (abs) for this IR
Definition HBFUtils.h:74
o2::InteractionRecord ir(0, 0)
std::vector< Tracklet64 > tracklets