Skip to content

Use DBSCAN for Time-Z clustering in PVertexing, debris reduction #5104

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
shahor02 merged 2 commits into from
Jan 4, 2021
Merged

Use DBSCAN for Time-Z clustering in PVertexing, debris reduction #5104

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
86c212a
Fix in MCEventLabel::getIDOnly
shahor02 Dec 25, 2020
bb959e7
Use DBSCAN for Time-Z clustering in PVertexing + debris reduction
shahor02 Dec 22, 2020
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
4 changes: 3 additions & 1 deletion DataFormats/simulation/include/SimulationDataFormat/MCEventLabel.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,8 @@ class MCEventLabel
static constexpr uint32_t MaskEvID = (0x1 << nbitsEvID) - 1;
// Mask to extract MC source ID
static constexpr uint32_t MaskSrcID = (0x1 << nbitsSrcID) - 1;
// Mask to extract MC source and event ID only
static constexpr uint32_t MaskSrcEvID = MaskSrcID | MaskEvID;
// Mask to extract MC correct contribution weight
static constexpr uint32_t MaskCorrW = (0x1 << nbitsCorrW) - 1;
static constexpr float WeightNorm = 1. / float(MaskCorrW);
Expand All @@ -53,7 +55,7 @@ class MCEventLabel
uint32_t getRawValue() const { return mLabel; }

// get only combined identifier, discarding weight info
uint32_t getIDOnly() const { return mLabel & (0x1 << (nbitsEvID + nbitsSrcID)); }
uint32_t getIDOnly() const { return mLabel & MaskSrcEvID; }

// compare
bool compare(const MCEventLabel& other, bool strict = false) const
Expand Down
12 changes: 7 additions & 5 deletions Detectors/Vertexing/include/DetectorsVertexing/PVertexer.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -80,7 +80,7 @@ class PVertexer
bool getValidateWithFT0() const { return mValidateWithFT0; }

auto& getTracksPool() const { return mTracksPool; }
auto& getTimeClusters() const { return mTimesClusters; }
auto& getTimeZClusters() const { return mTimeZClusters; }
auto& getSortedTrackIndices() const { return mSortedTrackID; }

auto& getMeanVertex() const { return mMeanVertex; }
Expand All @@ -92,7 +92,6 @@ class PVertexer

float estimateScale2()
{
float minrange = std::min(mPVParams->zHistoBinSize, mPVParams->minZSeedRange);
auto sc = mPVParams->zHistoBinSize * mPVParams->zHistoBinSize * mTukey2I / (mStatZErr.getMean() * mStatZErr.getMean());
return sc;
}
Expand All @@ -108,11 +107,12 @@ class PVertexer
void applyConstraint(VertexSeed& vtxSeed) const;
bool upscaleSigma(VertexSeed& vtxSeed) const;
void createTracksPool(gsl::span<const o2d::TrackTPCITS> tracksITSTPC);
void clusterizeTimeBruteForce(float margin = 0.1, float cut = 25);
void clusterizeTime(float binSize = 0.1, float maxTDist = 0.6);
int findVertices(const VertexingInput& input, std::vector<PVertex>& vertices, std::vector<int>& vertexTrackIDs, std::vector<V2TRef>& v2tRefs);
std::pair<int, int> getBestFT0Trigger(const PVertex& vtx, gsl::span<const o2::ft0::RecPoints> ft0Data, int& currEntry) const;

int dbscan_RangeQuery(int idxs, std::vector<int>& cand, const std::vector<int>& status);
void dbscan_clusterize();

o2::BunchFilling mBunchFilling;
std::array<int16_t, o2::constants::lhc::LHCMaxBunches> mClosestBunchAbove; // closest filled bunch from above
std::array<int16_t, o2::constants::lhc::LHCMaxBunches> mClosestBunchBelow; // closest filled bunch from below
Expand All @@ -123,7 +123,9 @@ class PVertexer
o2::math_utils::StatAccumulator mStatTErr;
std::vector<TrackVF> mTracksPool; ///< tracks in internal representation used for vertexing
std::vector<int> mSortedTrackID; ///< indices of tracks sorted in time
std::vector<TimeCluster> mTimesClusters; ///< set of time clusters
std::vector<TimeZCluster> mTimeZClusters; ///< set of time clusters
std::vector<int> mClusterTrackIDs; ///< IDs of tracks making the clusters

float mBz = 0.; ///< mag.field at beam line
bool mValidateWithFT0 = false; ///< require vertex validation with FT0 (if available)

Expand Down
22 changes: 19 additions & 3 deletions Detectors/Vertexing/include/DetectorsVertexing/PVertexerHelpers.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -116,6 +116,23 @@ struct TrackVF {
return z + tgL * (vx * cosAlp + vy * sinAlp - x);
}

// weighted distance^2 to other track (accounting for own errors only)
float getDist2(const TrackVF& o) const
{
auto dtnorm2 = (timeEst.getTimeStamp() - o.timeEst.getTimeStamp()) / timeEst.getTimeStampError();
auto dz = z - o.z;
return dtnorm2 + dz * dz * sig2ZI;
}

// weighted distance^2 to other track (accounting for both track errors errors only)
float getDist2Sym(const TrackVF& o) const
{
auto dt = timeEst.getTimeStamp() - o.timeEst.getTimeStamp();
auto dz = z - o.z;
float dte2 = o.timeEst.getTimeStampError() * o.timeEst.getTimeStampError() + timeEst.getTimeStampError() * timeEst.getTimeStampError();
return dt / dte2 + dz * dz / (1. / sig2ZI + 1. / o.sig2ZI);
}

float getResiduals(const PVertex& vtx, float& dy, float& dz) const
{
// get residuals (Y and Z DCA in track frame) and calculate chi2
Expand All @@ -139,7 +156,6 @@ struct TrackVF {
};

struct VertexingInput {
gsl::span<TrackVF> tracks;
gsl::span<int> idRange;
TimeEst timeEst{0, -1.}; // negative error means don't use time info
float scaleSigma2 = 10;
Expand Down Expand Up @@ -227,7 +243,7 @@ struct SeedHisto {
}
};

struct TimeCluster {
struct TimeZCluster {
TimeEst timeEst;
int first = -1;
int last = -1;
Expand Down Expand Up @@ -275,7 +291,7 @@ struct TimeCluster {
}
}

void merge(TimeCluster& c)
void merge(TimeZCluster& c)
{
if (c.first < last) {
first = c.first;
Expand Down
13 changes: 12 additions & 1 deletion Detectors/Vertexing/include/DetectorsVertexing/PVertexerParams.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,8 +25,19 @@ namespace vertexing
struct PVertexerParams : public o2::conf::ConfigurableParamHelper<PVertexerParams> {
static constexpr float kDefTukey = 5.0f; ///< def.value for tukey constant

// DBSCAN clustering settings
float dbscanMaxDist2 = 9.; ///< distance^2 cut (eps^2).
float dbscanDeltaT = 10.; ///< abs. time difference cut
float dbscanAdaptCoef = 0.1; ///< adapt dbscan minPts for each cluster as minPts=max(minPts, currentSize*dbscanAdaptCoef).

int maxVerticesPerCluster = 1; ///< max vertices per time-z cluster to look for
int maxTrialsPerCluster = 3; ///< max unsucessful trials for vertex search per vertex

// track selection
float dcaTolerance = 1.3; ///< consider tracks within this abs DCA to mean vertex
float pullIniCut = 9; ///< cut on pull (n^2 sigma) on dca to mean vertex

// parameters
float minZSeedRange = 0.5; ///< min proximity of Zseed which should be covered by 1st iteration
float zHistoRange = 20.; ///< +-range of the Zseed histo
float zHistoBinSize = 0.5; ///< size of the Zseed histo bin
float tukey = kDefTukey; ///< 1./[Tukey parameter]^2
Expand Down
Loading