Skip to content

Commit b3ee6aa

Browse files
f3schf3sch
committed
ITS: debug dump for cell seeds
1 parent 78bacab commit b3ee6aa

File tree

6 files changed

+151
-25
lines changed

6 files changed

+151
-25
lines changed

‎Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h‎

Lines changed: 11 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -79,16 +79,17 @@ struct TrackingParameters {
7979
float TrackFollowerNSigmaCutPhi = 1.f;
8080

8181
/// Seeding parameters
82-
float seedingDCATolerance{1.8f}; // maximum allowed DCA to meanvertex for track to enter pool
83-
float seedingDCAMaxPull{9.f}; // maximum allowed initial pull on DCA to meanvertex
84-
float seedingMaxChi2Iter{3.f}; // maximum chi2 change required to end iteration
85-
float seedingTukeyStartIter{5.f}; // start value for tukey scaling for iteration
86-
int seedingMaxFitIter{3}; // maximum iterations for fit
87-
int seedingMinTracksIter{3}; // minimum tracks needed for one iteration
88-
int seedingDBScanMinPt{3}; // DBSCAN: minimum number of cluster points
89-
float seedingDBScanEpsZ{0.1f}; // DBSCAN: maximum epsilon z
90-
float seedingDBScanEpsT{10.f}; // DBSCAN: maximum epsilon t (BC)
91-
float seedingVertexExtraErr2[6] = {0.f}; // impose additional errors on seeding vertices
82+
float SeedingDCATolerance{1.8f}; // maximum allowed DCA to meanvertex for track to enter pool
83+
float SeedingDCAMaxPull{9.f}; // maximum allowed initial pull on DCA to meanvertex
84+
float SeedingMaxChi2Iter{3.f}; // maximum chi2 change required to end iteration
85+
float SeedingTukeyStartIter{5.f}; // start value for tukey scaling for iteration
86+
float SeedingMinWghTrk{0.05f}; // minimum weight a track has to have to be accounted
87+
int SeedingMaxFitIter{3}; // maximum iterations for fit
88+
int SeedingMinTracksIter{50}; // minimum tracks needed for one iteration
89+
int SeedingDBScanMinPt{100}; // DBSCAN: minimum number of cluster points
90+
float SeedingDBScanEpsZ{0.1f}; // DBSCAN: maximum epsilon z
91+
float SeedingDBScanEpsT{10.f}; // DBSCAN: maximum epsilon t (BC)
92+
float SeedingVertexExtraErr2[6] = {0.f}; // impose additional errors on seeding vertices
9293

9394
bool createArtefactLabels{false};
9495

‎Detectors/ITSMFT/ITS/tracking/include/ITStracking/Seeding.h‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -95,11 +95,13 @@ struct VertexSeed : public Vertex {
9595
float wghChi2 = 0.f; // sum of tracks weighted chi2's
9696
float scaleSig2ITuk2I = 0.f; // inverted scaled tukey weight
9797
float tukeyC = 0.f; // tukey constant
98+
float wghMin = 0.f; // minimum weight to account track
9899
ROOT::Math::SMatrix<float, 3, 3, ROOT::Math::MatRepSym<float, 3>> C; // C matrix
99100
ROOT::Math::SVector<float, 3> b; // b vector
100101

101102
void updateTukeyScale(const LinearizedTrack*, const int*, size_t, const LTRef&);
102103
void resetForNewIteration();
104+
bool acceptTrack(const LinearizedTrack&);
103105
void accountTrack(const LinearizedTrack&);
104106
void solveVertex();
105107
};

‎Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h‎

Lines changed: 12 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -55,7 +55,18 @@ struct TrackerParamConfig : public o2::conf::ConfigurableParamHelper<TrackerPara
5555
bool doUPCIteration = false; // Perform an additional iteration for UPC events on tagged vertices. You want to combine this config with VertexerParamConfig.nIterations=2
5656

5757
/// seeding
58-
float seedingMeanVertexExtraErr2{0.f}; // additional error imposed on mean vertex cov.
58+
float seedingMeanVertexExtraErr2{0.f}; // additional error imposed on mean vertex cov.
59+
float seedingDCATolerance{1.8f}; // maximum allowed DCA to meanvertex for track to enter pool
60+
float seedingDCAMaxPull{9.f}; // maximum allowed initial pull on DCA to meanvertex
61+
float seedingMaxChi2Iter{3.f}; // maximum chi2 change required to end iteration
62+
float seedingTukeyStartIter{5.f}; // start value for tukey scaling for iteration
63+
float seedingMinWghTrk{0.05f}; // minimum weight a track has to have to be accounted
64+
int seedingMaxFitIter{3}; // maximum iterations for fit
65+
int seedingMinTracksIter{50}; // minimum tracks needed for one iteration
66+
int seedingDBScanMinPt{100}; // DBSCAN: minimum number of cluster points
67+
float seedingDBScanEpsZ{0.1f}; // DBSCAN: maximum epsilon z
68+
float seedingDBScanEpsT{10.f}; // DBSCAN: maximum epsilon t (BC)
69+
float seedingVertexExtraErr2[6] = {0.f}; // impose additional errors on seeding vertices
5970

6071
bool createArtefactLabels{false}; // create on-the-fly labels for the artefacts
6172

‎Detectors/ITSMFT/ITS/tracking/src/Configuration.cxx‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -199,6 +199,16 @@ std::vector<RecoIteration> TrackingMode::getRecoIterations(TrackingMode::Type mo
199199
recoIterations[0].params.NLayers = 3; // only do cell finding up until the third layer
200200
recoIterations[0].params.UseDiamond = true; // use the blown up diamond constrain to (e.g. luminous region)
201201
recoIterations[0].params.CorrType = o2::base::PropagatorF::MatCorrType::USEMatCorrNONE; // do not use material
202+
recoIterations[0].params.SeedingDCATolerance = tc.seedingDCATolerance;
203+
recoIterations[0].params.SeedingDCAMaxPull = tc.seedingDCAMaxPull;
204+
recoIterations[0].params.SeedingMaxChi2Iter = tc.seedingMaxChi2Iter;
205+
recoIterations[0].params.SeedingTukeyStartIter = tc.seedingTukeyStartIter;
206+
recoIterations[0].params.SeedingMinWghTrk = tc.seedingMinWghTrk;
207+
recoIterations[0].params.SeedingMaxFitIter = tc.seedingMaxFitIter;
208+
recoIterations[0].params.SeedingMinTracksIter = tc.seedingMinTracksIter;
209+
recoIterations[0].params.SeedingDBScanMinPt = tc.seedingDBScanMinPt;
210+
recoIterations[0].params.SeedingDBScanEpsZ = tc.seedingDBScanEpsZ;
211+
recoIterations[0].params.SeedingDBScanEpsT = tc.seedingDBScanEpsT;
202212
}
203213

204214
if (mode == TrackingMode::Async) {

‎Detectors/ITSMFT/ITS/tracking/src/Seeding.cxx‎

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -78,11 +78,18 @@ void VertexSeed::resetForNewIteration()
7878
std::fill(b.begin(), b.end(), 0.f);
7979
}
8080

81+
bool VertexSeed::acceptTrack(const LinearizedTrack& lt)
82+
{
83+
float chi2Red = lt.getChi2((const dataformats::VertexBase&)*this);
84+
float wghT = (1.f - chi2Red * scaleSig2ITuk2I); // weighted distance to vertex
85+
return wghT >= wghMin;
86+
}
87+
8188
void VertexSeed::accountTrack(const LinearizedTrack& lt)
8289
{
8390
float chi2Red = lt.getChi2((const dataformats::VertexBase&)*this);
8491
float wghT = (1.f - chi2Red * scaleSig2ITuk2I); // weighted distance to vertex
85-
if (wghT < constants::Tolerance) {
92+
if (wghT < wghMin) {
8693
return;
8794
}
8895
wghT *= wghT;

‎Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx‎

Lines changed: 108 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,6 @@
2121
#include <limits>
2222
#include <ranges>
2323
#include <type_traits>
24-
#include "ReconstructionDataFormats/Vertex.h"
2524

2625
#ifdef OPTIMISATION_OUTPUT
2726
#include <format>
@@ -42,6 +41,7 @@
4241
#include "ITStracking/IndexTableUtils.h"
4342
#include "ITStracking/Tracklet.h"
4443
#include "ReconstructionDataFormats/Track.h"
44+
#include "CommonUtils/TreeStreamRedirector.h"
4545

4646
using o2::base::PropagatorF;
4747

@@ -456,9 +456,56 @@ template <int NLayers>
456456
void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
457457
{
458458
const auto& propagator = o2::base::Propagator::Instance();
459+
auto& cells = mTimeFrame->getCells()[0];
460+
461+
/// DEBUG OUTPUT
462+
// dump cells and pot. their labels
463+
auto ts = new utils::TreeStreamRedirector("cells.root");
464+
auto& s = (*ts) << "seeding";
465+
for (int iCell{0}; iCell < cells.size(); ++iCell) {
466+
auto cell = cells[iCell];
467+
dataformats::VertexBase vtx;
468+
dataformats::DCA dca;
469+
float z = cell.getZAt(0., getBz());
470+
if (z < -999.f) {
471+
if (mTimeFrame->hasMeanVertex()) {
472+
z = mTimeFrame->getMeanVertex()->getZ();
473+
} else {
474+
continue;
475+
}
476+
}
477+
if (mTimeFrame->hasMeanVertex()) {
478+
mTimeFrame->getMeanVertex()->setMeanXYVertexAtZ(vtx, z);
479+
}
480+
const auto& cls = cell.getClusters();
481+
const int sta = cell.getUserField();
482+
int startBC = std::numeric_limits<int>::max(), endBC = std::numeric_limits<int>::min();
483+
for (int i{sta}; i < 3; ++i) {
484+
int rof = mTimeFrame->getClusterROF(i, cls[i]);
485+
int rofStartBC = mTimeFrame->getROFOverlapTableView().getLayer(i).getROFStartInBC(rof);
486+
int rofEndBC = mTimeFrame->getROFOverlapTableView().getLayer(i).getROFEndInBC(rof);
487+
startBC = o2::gpu::CAMath::Min(startBC, rofStartBC);
488+
endBC = o2::gpu::CAMath::Max(endBC, rofEndBC);
489+
}
490+
auto cellC = cell;
491+
if (propagator->propagateToDCA(vtx, cellC, getBz(), 2.0f, o2::base::Propagator::MatCorrType::USEMatCorrLUT, &dca, nullptr, 0, mRecoParams[iteration].params.SeedingDCATolerance)) {
492+
s << "cell=" << (const o2::track::TrackParCov&)cell
493+
<< "cellPV=" << (o2::track::TrackParCov)cellC
494+
<< "vtx=" << vtx
495+
<< "staBC=" << startBC
496+
<< "endBC=" << endBC;
497+
if (mTimeFrame->hasMCinformation() && mRecoParams[iteration].params.createArtefactLabels) {
498+
s << "cellLbl=" << mTimeFrame->getCellsLabel(0)[iCell];
499+
}
500+
s << "\n";
501+
}
502+
}
503+
ts->Close();
504+
delete ts;
505+
ts = nullptr;
506+
459507
mTaskArena->execute([&] {
460508
// we get cells from the first three layers
461-
auto& cells = mTimeFrame->getCells()[0];
462509
bounded_vector<LinearizedTrack> ltracks;
463510
ltracks.resize(cells.size());
464511
tbb::parallel_for(size_t(0), cells.size(), [&](size_t iCell) {
@@ -478,11 +525,11 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
478525
if (mTimeFrame->hasMeanVertex()) {
479526
mTimeFrame->getMeanVertex()->setMeanXYVertexAtZ(vtx, z);
480527
}
481-
if (!propagator->propagateToDCA(vtx, cell, getBz(), 2.0f, o2::base::Propagator::MatCorrType::USEMatCorrLUT, &dca, nullptr, 0, mRecoParams[iteration].params.seedingDCATolerance)) {
528+
if (!propagator->propagateToDCA(vtx, cell, getBz(), 2.0f, o2::base::Propagator::MatCorrType::USEMatCorrLUT, &dca, nullptr, 0, mRecoParams[iteration].params.SeedingDCATolerance)) {
482529
ltracks[iCell].markDead();
483530
return;
484531
}
485-
if (dca.getY() * dca.getY() / (dca.getSigmaY2()) >= mRecoParams[iteration].params.seedingDCAMaxPull) {
532+
if (dca.getY() * dca.getY() / (dca.getSigmaY2()) >= mRecoParams[iteration].params.SeedingDCAMaxPull) {
486533
ltracks[iCell].markDead();
487534
return;
488535
}
@@ -531,6 +578,11 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
531578
// scan gives association in time and z
532579
const auto dbRes = mDBScan.cluster(ltracks.data(), ltracks.size());
533580
bounded_vector<VertexSeed> vtxSeeds(dbRes.nClusters, mMemoryPool.get());
581+
bounded_vector<std::pair<o2::MCCompLabel, float>> vtxSeedsLbl(dbRes.nClusters, mMemoryPool.get());
582+
// we only care about the source&event of the tracks, not the trackId
583+
auto composeVtxLabel = [](const o2::MCCompLabel& lbl) -> o2::MCCompLabel {
584+
return {o2::MCCompLabel::maxTrackID(), lbl.getEventID(), lbl.getSourceID(), lbl.isFake()};
585+
};
534586
tbb::parallel_for(0, dbRes.nClusters, [&](const int32_t iCls) {
535587
// create vertex seeds based on meanvertex and scanned z which one can take as start
536588
auto& seed = vtxSeeds[iCls];
@@ -539,11 +591,11 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
539591
mTimeFrame->getMeanVertex()->setMeanXYVertexAtZ(seed, seed.getZ());
540592
}
541593
seed.idx = iCls;
542-
seed.tukeyC = mRecoParams[iteration].params.seedingTukeyStartIter;
594+
seed.tukeyC = mRecoParams[iteration].params.SeedingTukeyStartIter;
543595
seed.scaleSig2ITuk2I = 1.f / (seed.tukeyC * seed.tukeyC);
544596
const auto& range = dbRes.ranges[iCls];
545597
// fit iteratively
546-
for (int iter{0}; iter < mRecoParams[iteration].params.seedingMaxFitIter; ++iter) {
598+
for (int iter{0}; iter < mRecoParams[iteration].params.SeedingMaxFitIter; ++iter) {
547599
seed.iteration = iter;
548600
// 1. update tukey scaling
549601
seed.updateTukeyScale(ltracks.data(), dbRes.labels.data(), ltracks.size(), range);
@@ -558,7 +610,7 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
558610
seed.accountTrack(lt);
559611
}
560612
// 4. if needed relax scaling
561-
if (seed.getNContributors() < mRecoParams[iteration].params.seedingMinTracksIter) {
613+
if (seed.getNContributors() < mRecoParams[iteration].params.SeedingMinTracksIter) {
562614
seed.status = VertexSeed::kIterateRelaxScale;
563615
continue;
564616
}
@@ -577,7 +629,7 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
577629
seed.solveVertex();
578630
// 7. check for convergence
579631
float avgChi2 = seed.wghChi2 / o2::gpu::CAMath::Max(1.f, seed.wghSum);
580-
if (avgChi2 < mRecoParams[iteration].params.seedingMaxChi2Iter) {
632+
if (avgChi2 < mRecoParams[iteration].params.SeedingMaxChi2Iter) {
581633
seed.status = VertexSeed::kConverged;
582634
break;
583635
}
@@ -591,7 +643,7 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
591643
continue;
592644
}
593645
float chi2Red = lt.getChi2(seed);
594-
if (chi2Red > mRecoParams[iteration].params.seedingMaxChi2Iter) {
646+
if (chi2Red > mRecoParams[iteration].params.SeedingMaxChi2Iter) {
595647
continue;
596648
}
597649
startBC = o2::gpu::CAMath::Min(startBC, lt.time.getTimeStamp());
@@ -601,7 +653,47 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
601653
seed.getTimeStamp().setTimeStampError(endBC - startBC);
602654
// add additional errors
603655
for (int i{0}; i < dataformats::VertexBase::kNCov; ++i) {
604-
seed.setCov(seed.getCov(i) + mRecoParams[iteration].params.seedingVertexExtraErr2[i], i);
656+
seed.setCov(seed.getCov(i) + mRecoParams[iteration].params.SeedingVertexExtraErr2[i], i);
657+
}
658+
// if mc is present calculate labels and purity
659+
if (mTimeFrame->hasMCinformation()) {
660+
// use boyer-moore voting
661+
int accepted{0}, weight{0};
662+
o2::MCCompLabel lbl;
663+
for (uint32_t entry{range.getFirstEntry()}; entry < range.getEntriesBound(); ++entry) {
664+
const auto& lt = ltracks[entry];
665+
if (lt.isDead() || dbRes.labels[entry] != iCls) {
666+
continue;
667+
}
668+
if (seed.acceptTrack(lt)) {
669+
++accepted;
670+
// create cell label
671+
const auto& cell = mTimeFrame->getCells()[0][lt.cellIdx];
672+
o2::MCCompLabel cl;
673+
for (const auto& lab1 : mTimeFrame->getClusterLabels(0, cell.getFirstClusterIndex())) {
674+
for (const auto& lab2 : mTimeFrame->getClusterLabels(1, cell.getSecondClusterIndex())) {
675+
if (lab1 == lab2 && lab1.isValid()) {
676+
cl = lab1;
677+
break;
678+
}
679+
}
680+
if (cl.isValid()) {
681+
break;
682+
}
683+
}
684+
if (weight == 0) {
685+
lbl = cl;
686+
weight = 1;
687+
} else {
688+
(cl == lbl) ? ++weight : --weight;
689+
}
690+
}
691+
}
692+
if (accepted > 0) {
693+
vtxSeedsLbl[iCls] = {lbl, static_cast<float>(weight) / static_cast<float>(accepted)};
694+
} else {
695+
vtxSeedsLbl[iCls] = {lbl, 0};
696+
}
605697
}
606698
});
607699
// TODO reduce debris vertices (eliminate small mult vertices within vicinity of high multiplicity ones)
@@ -611,6 +703,9 @@ void TrackerTraits<NLayers>::findCellSeeds(const int iteration)
611703
continue;
612704
}
613705
mTimeFrame->addPrimaryVertex(vtxSeeds[i]);
706+
if (mTimeFrame->hasMCinformation()) {
707+
mTimeFrame->addPrimaryVertexLabel(vtxSeedsLbl[i]);
708+
}
614709
}
615710
});
616711
}
@@ -1413,9 +1508,9 @@ void TrackerTraits<NLayers>::updateTrackingParameters(const std::vector<RecoIter
14131508
if (static bool onceDone{false}; !onceDone) {
14141509
onceDone = true;
14151510
dbscan::DBSCANParams p;
1416-
p.minPts = mRecoParams[0].params.seedingDBScanMinPt;
1417-
p.eps[0] = mRecoParams[0].params.seedingDBScanEpsZ;
1418-
p.eps[1] = mRecoParams[0].params.seedingDBScanEpsT;
1511+
p.minPts = mRecoParams[0].params.SeedingDBScanMinPt;
1512+
p.eps[0] = mRecoParams[0].params.SeedingDBScanEpsZ;
1513+
p.eps[1] = mRecoParams[0].params.SeedingDBScanEpsT;
14191514
mDBScan = dbscan::DBSCAN(p);
14201515
}
14211516
}

0 commit comments

Comments
 (0)