Skip to content

Commit

Permalink
alksjdlkjalksdjlkas
Browse files Browse the repository at this point in the history
  • Loading branch information
f3sch committed Oct 23, 2023
1 parent 967a319 commit 6a27aa4
Show file tree
Hide file tree
Showing 4 changed files with 72 additions and 53 deletions.
49 changes: 25 additions & 24 deletions Common/DCAFitter/include/DCAFitter/DCAFitterN.h
Original file line number Diff line number Diff line change
Expand Up @@ -220,16 +220,17 @@ class DCAFitterN
float getMasStep() const { return mMaxStep; }
float getMinXSeed() const { return mMinXSeed; }

void setDebug(int curId)
void setDebug(int curId, bool shutup)
{
mDebug = true;
LOG_IF(info, mDebug) << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Current failed Track = " << mCurID;
mShutup = shutup;
LOG_IF(info, mDebug && !mShutup) << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Current failed Track = " << mCurID;
mCurID = curId;
mCrossings.setDebug();
mCrossings.setDebug(shutup);
}
void unsetDebug()
{
LOG_IF(info, mDebug) << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~";
LOG_IF(info, mDebug && !mShutup) << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~";
mDebug = false;
mCrossings.unsetDebug();
}
Expand Down Expand Up @@ -365,7 +366,7 @@ class DCAFitterN
float mMaxStep = 2.0; // Max step for propagation with Propagator
int mFitterID = 0; // locat fitter ID (mostly for debugging)
size_t mCallID = 0;
bool mDebug{false};
bool mDebug{false}, mShutup{false};
int mCurID{0};
ULong64_t nCrossings{0}, nPropFail{0}, nZCut{0}, nInvFail{0}, nCorrectFails{0}, nCloserAlt{0}, nChi2{0},
nConcentric{0}, nDistXY{0}, nNotTouching{0}, nChi2Called{0}, nTotCrossings{0};
Expand All @@ -382,20 +383,20 @@ int DCAFitterN<N, Args...>::process(const Tr&... args)
static_assert(sizeof...(args) == N, "incorrect number of input tracks");
assign(0, args...);
clear();
LOG_IF(info, mDebug) << "------------ TracksStart";
LOG_IF(info, mDebug && !mShutup) << "------------ TracksStart";
for (int i = 0; i < N; i++) {
mTrAux[i].set(*mOrigTrPtr[i], mBz);
if (mDebug) {
LOG_IF(info, mDebug) << "+++++ i=" << i;
if (mDebug && !mShutup) {
LOG(info) << "+++++ i=" << i;
mOrigTrPtr[i]->print();
mTrAux[i].print();
LOG_IF(info, mDebug) << "-----";
LOG(info) << "-----";
}
}
LOG_IF(info, mDebug) << "------------ TracksEnd";
LOG_IF(info, mDebug && !mShutup) << "------------ TracksEnd";

if (!mCrossings.set(mTrAux[0], *mOrigTrPtr[0], mTrAux[1], *mOrigTrPtr[1], mMaxDXYIni)) { // even for N>2 it should be enough to test just 1 loop
LOG_IF(info, mDebug) << "no crossings";
LOG_IF(info, mDebug && !mShutup) << "no crossings";
if (mDebug) {
if (mCrossings.bNotTouching) {
++nNotTouching;
Expand All @@ -420,7 +421,7 @@ int DCAFitterN<N, Args...>::process(const Tr&... args)
calcRMatrices(); // needed for fast residuals derivatives calculation in case of abs. distance minimization
}
if (mCrossings.nDCA == MAXHYP) { // if there are 2 candidates and they are too close, chose their mean as a starting point
LOG_IF(info, mDebug) << " `-> 2 candidates";
LOG_IF(info, mDebug && !mShutup) << " `-> 2 candidates";
auto dst2 = (mCrossings.xDCA[0] - mCrossings.xDCA[1]) * (mCrossings.xDCA[0] - mCrossings.xDCA[1]) +
(mCrossings.yDCA[0] - mCrossings.yDCA[1]) * (mCrossings.yDCA[0] - mCrossings.yDCA[1]);
if (dst2 < mMaxDist2ToMergeSeeds) {
Expand All @@ -434,10 +435,10 @@ int DCAFitterN<N, Args...>::process(const Tr&... args)
nTotCrossings += mCrossings.nDCA;
// check all crossings
for (int ic = 0; ic < mCrossings.nDCA; ic++) {
LOG_IF(info, mDebug) << "iCandidate=" << ic << "/" << mCrossings.nDCA << ": xDCA=" << mCrossings.xDCA[ic] << " yDCA=" << mCrossings.yDCA[ic];
LOG_IF(info, mDebug && !mShutup) << "iCandidate=" << ic << "/" << mCrossings.nDCA << ": xDCA=" << mCrossings.xDCA[ic] << " yDCA=" << mCrossings.yDCA[ic];
// check if radius is acceptable
if (mCrossings.xDCA[ic] * mCrossings.xDCA[ic] + mCrossings.yDCA[ic] * mCrossings.yDCA[ic] > mMaxR2) {
LOG_IF(info, mDebug) << " `-> radius unaccepatble";
LOG_IF(info, mDebug && !mShutup) << " `-> radius unaccepatble";
cont = true;
++nCont;
continue;
Expand All @@ -451,18 +452,18 @@ int DCAFitterN<N, Args...>::process(const Tr&... args)
mPCA[mCurHyp][1] = mCrossings.yDCA[ic];
++nChi2Called;
if (mUseAbsDCA ? minimizeChi2NoErr() : minimizeChi2()) {
LOG_IF(info, mDebug) << "minimizeChi2 true";
LOG_IF(info, mDebug && !mShutup) << "minimizeChi2 true";
mOrder[mCurHyp] = mCurHyp;
if (mPropagateToPCA && !propagateTracksToVertex(mCurHyp)) {
LOG_IF(info, mDebug) << "propagateTracksToVertex failed";
LOG_IF(info, mDebug && !mShutup) << "propagateTracksToVertex failed";
cont = true;
++nCont;
continue; // discard candidate if failed to propagate to it
}
mCurHyp++;
}
}
LOG_IF(info, mDebug && nCont) << "Hit Continued: " << nCont;
LOG_IF(info, mDebug && nCont && !mShutup) << "Hit Continued: " << nCont;

for (int i = mCurHyp; i--;) { // order in quality
for (int j = i; j--;) {
Expand All @@ -476,7 +477,7 @@ int DCAFitterN<N, Args...>::process(const Tr&... args)
recalculatePCAWithErrors(i);
}
}
LOG_IF(info, mDebug) << "number of candidates=" << mCurHyp;
LOG_IF(info, mDebug && !mShutup) << "number of candidates=" << mCurHyp;
return mCurHyp;
}

Expand Down Expand Up @@ -996,7 +997,7 @@ bool DCAFitterN<N, Args...>::minimizeChi2NoErr()
mCandTr[mCurHyp][i] = *mOrigTrPtr[i];
auto x = mTrAux[i].c * mPCA[mCurHyp][0] + mTrAux[i].s * mPCA[mCurHyp][1]; // X of PCA in the track frame
if (x < mMinXSeed || !propagateParamToX(mCandTr[mCurHyp][i], x)) {
LOG_IF(info, mDebug) << " " << i << ": ---> prop failed x of PCA";
LOG_IF(info, mDebug && !mShutup) << " " << i << ": ---> prop failed x of PCA";
if (mDebug)
++nPropFail;
return false;
Expand Down Expand Up @@ -1026,14 +1027,14 @@ bool DCAFitterN<N, Args...>::minimizeChi2NoErr()
}
VecND dx = mD2Chi2Dx2 * mDChi2Dx;
if (!correctTracks(dx)) {
LOG_IF(info, mDebug) << " -> correctTracks failed";
LOG_IF(info, mDebug && !mShutup) << " -> correctTracks failed";
if (mDebug)
++nCorrectFails;
return false;
}
calcPCANoErr(); // updated PCA
if (mCrossIDAlt >= 0 && closerToAlternative()) {
LOG_IF(info, mDebug) << " -> closerToAlternative failed";
LOG_IF(info, mDebug && !mShutup) << " -> closerToAlternative failed";
if (mDebug)
++nCloserAlt;
mAllowAltPreference = false;
Expand All @@ -1043,15 +1044,15 @@ bool DCAFitterN<N, Args...>::minimizeChi2NoErr()
chi2Upd = calcChi2NoErr(); // updated chi2
if (getAbsMax(dx) < mMinParamChange || chi2Upd > chi2 * mMinRelChi2Change) {
chi2 = chi2Upd;
LOG_IF(info, mDebug) << " ~-> converged: chi2 = " << chi2;
LOG_IF(info, mDebug && !mShutup) << " ~-> converged: chi2 = " << chi2;
break; // converged
}
chi2 = chi2Upd;
} while (++mNIters[mCurHyp] < mMaxIter);
//
mChi2[mCurHyp] = chi2 * NInv;
bool ret = mChi2[mCurHyp] < mMaxChi2;
LOG_IF(info, !ret && mDebug) << " ~-~-~> mChi2(=" << mChi2[mCurHyp] << ")<mMaxChi2(=" << mMaxChi2 << ") ==> " << ret;
LOG_IF(info, !ret && mDebug && !mShutup) << " ~-~-~> mChi2(=" << mChi2[mCurHyp] << ")<mMaxChi2(=" << mMaxChi2 << ") ==> " << ret;
if (!ret && mDebug) {
++nChi2;
vChi2.push_back(mChi2[mCurHyp]);
Expand All @@ -1072,7 +1073,7 @@ bool DCAFitterN<N, Args...>::roughDZCut()
if (mDebug) {
vDistZ.push_back(zdist);
}
LOG_IF(info, mDebug) << " i=" << i << " j=" << j << " ~~~~> zCut failed: " << zdist;
LOG_IF(info, mDebug && !mShutup) << " i=" << i << " j=" << j << " ~~~~> zCut failed: " << zdist;
accept = false;
break;
}
Expand Down
30 changes: 16 additions & 14 deletions Common/DCAFitter/include/DCAFitter/HelixHelper.h
Original file line number Diff line number Diff line change
Expand Up @@ -64,9 +64,11 @@ struct CrossInfo {
float yDCA[2] = {};
int nDCA = 0;
bool mDebug = false;
void setDebug()
bool mShutup = false;
void setDebug(bool shutup)
{
mDebug = true;
mShutup = shutup;
bConcentric = false;
bDistXY = false;
bNotTouching = false;
Expand All @@ -83,17 +85,17 @@ struct CrossInfo {
const auto& trcB = trax0.rC > trax1.rC ? trax1 : trax0;
float xDist = trcB.xC - trcA.xC, yDist = trcB.yC - trcA.yC;
float dist2 = xDist * xDist + yDist * yDist, dist = std::sqrt(dist2), rsum = trcA.rC + trcB.rC;
LOG_IF(info, mDebug) << "xDist=" << xDist << ", yDist=" << yDist << ", dist2=" << dist2 << ", dist=" << dist << ", rsum=" << rsum;
LOG_IF(info, mDebug && !mShutup) << "xDist=" << xDist << ", yDist=" << yDist << ", dist2=" << dist2 << ", dist=" << dist << ", rsum=" << rsum;
if (std::abs(dist) < 1e-12) {
LOG_IF(info, mDebug) << " -> concentric=" << nDCA;
LOG_IF(info, mDebug && !mShutup) << " -> concentric=" << nDCA;
bConcentric = true;
return nDCA; // circles are concentric?
}
if (dist > rsum) { // circles don't touch, chose a point in between
// the parametric equation of lines connecting the centers is
// x = x0 + t/dist * (x1-x0), y = y0 + t/dist * (y1-y0)
if (dist - rsum > maxDistXY) { // too large distance
LOG_IF(info, mDebug) << " -> too large distance=" << nDCA;
LOG_IF(info, mDebug && !mShutup) << " -> too large distance=" << nDCA;
vDistXY = dist - rsum;
vDist = dist;
vRsum = rsum;
Expand All @@ -105,14 +107,14 @@ struct CrossInfo {
// select the point of closest approach of 2 circles
notTouchingXY(dist, xDist, yDist, trcA, -trcB.rC);
} else { // 2 intersection points
LOG_IF(info, mDebug) << "2 intersection points";
LOG_IF(info, mDebug && !mShutup) << "2 intersection points";
// to simplify calculations, we move to new frame x->x+Xc0, y->y+Yc0, so that
// the 1st one is centered in origin
if (std::abs(xDist) < std::abs(yDist)) {
float a = (trcA.rC * trcA.rC - trcB.rC * trcB.rC + dist2) / (2. * yDist), b = -xDist / yDist, ab = a * b, bb = b * b;
float det = ab * ab - (1. + bb) * (a * a - trcA.rC * trcA.rC);
LOG_IF(info, mDebug) << "Abs(xDist)<Abs(yDist)"
<< "a=" << a << ", b=" << b << ", det=" << det;
LOG_IF(info, mDebug && !mShutup) << "Abs(xDist)<Abs(yDist)"
<< "a=" << a << ", b=" << b << ", det=" << det;
if (det > 0.) {
det = std::sqrt(det);
xDCA[0] = (-ab + det) / (1. + b * b);
Expand All @@ -122,16 +124,16 @@ struct CrossInfo {
yDCA[1] = a + b * xDCA[1] + trcA.yC;
xDCA[1] += trcA.xC;
nDCA = 2;
LOG_IF(info, mDebug) << "xDCA[0]=" << xDCA[0] << ", yDCA[0]=" << yDCA[0] << ", xDCA[1]=" << xDCA[1] << ", yDCA[1]=" << yDCA[1];
LOG_IF(info, mDebug && !mShutup) << "xDCA[0]=" << xDCA[0] << ", yDCA[0]=" << yDCA[0] << ", xDCA[1]=" << xDCA[1] << ", yDCA[1]=" << yDCA[1];
} else { // due to the finite precision the det<=0, i.e. the circles are barely touching, fall back to this special case
LOG_IF(info, mDebug) << " --> barely touching";
LOG_IF(info, mDebug && !mShutup) << " --> barely touching";
notTouchingXY(dist, xDist, yDist, trcA, trcB.rC);
}
} else {
float a = (trcA.rC * trcA.rC - trcB.rC * trcB.rC + dist2) / (2. * xDist), b = -yDist / xDist, ab = a * b, bb = b * b;
float det = ab * ab - (1. + bb) * (a * a - trcA.rC * trcA.rC);
LOG_IF(info, mDebug) << "Abs(xDist)>Abs(yDist)"
<< "a=" << a << ", b=" << b << ", det=" << det;
LOG_IF(info, mDebug && !mShutup) << "Abs(xDist)>Abs(yDist)"
<< "a=" << a << ", b=" << b << ", det=" << det;
if (det > 0.) {
det = std::sqrt(det);
yDCA[0] = (-ab + det) / (1. + bb);
Expand All @@ -141,9 +143,9 @@ struct CrossInfo {
xDCA[1] = a + b * yDCA[1] + trcA.xC;
yDCA[1] += trcA.yC;
nDCA = 2;
LOG_IF(info, mDebug) << "xDCA[0]=" << xDCA[0] << ", yDCA[0]=" << yDCA[0] << ", xDCA[1]=" << xDCA[1] << ", yDCA[1]=" << yDCA[1];
LOG_IF(info, mDebug && !mShutup) << "xDCA[0]=" << xDCA[0] << ", yDCA[0]=" << yDCA[0] << ", xDCA[1]=" << xDCA[1] << ", yDCA[1]=" << yDCA[1];
} else { // due to the finite precision the det<=0, i.e. the circles are barely touching, fall back to this special case
LOG_IF(info, mDebug) << " --> barely touching";
LOG_IF(info, mDebug && !mShutup) << " --> barely touching";
notTouchingXY(dist, xDist, yDist, trcA, trcB.rC);
}
}
Expand All @@ -163,7 +165,7 @@ struct CrossInfo {
auto t2d = (dist + trcA.rC - rBSign) / dist;
xDCA[0] = trcA.xC + 0.5 * (xDist * t2d);
yDCA[0] = trcA.yC + 0.5 * (yDist * t2d);
LOG_IF(info, mDebug) << " --> notTouchingXY: xDCA=" << xDCA[0] << " yDCA=" << yDCA[0];
LOG_IF(info, mDebug && !mShutup) << " --> notTouchingXY: xDCA=" << xDCA[0] << " yDCA=" << yDCA[0];
bNotTouching = true;
}

Expand Down
14 changes: 7 additions & 7 deletions Detectors/Vertexing/include/DetectorsVertexing/SVertexerParams.h
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@ struct SVertexerParams : public o2::conf::ConfigurableParamHelper<SVertexerParam
float minParamChange = 1e-3; ///< stop when tracks X-params being minimized change by less that this value
float minRelChi2Change = 0.9; ///< stop when chi2 changes by less than this value
float maxDZIni = 5.; ///< don't consider as a seed (circles intersection) if Z distance exceeds this
float maxDXYIni = 4.; ///< don't consider as a seed (circles intersection) if XY distance exceeds this
float maxDXYIni = 4.; ///< don't consider as a seed (circles intersection) if XY distance exceeds this
float maxRIni = 150; ///< don't consider as a seed (circles intersection) if its R exceeds this
/// bb
//
Expand All @@ -66,8 +66,8 @@ struct SVertexerParams : public o2::conf::ConfigurableParamHelper<SVertexerParam
float causalityRTolerance = 1.; ///< V0 radius cannot exceed its contributors minR by more than this value
float maxV0ToProngsRDiff = 50.; ///< V0 radius cannot be lower than this ammount wrt minR of contributors

float minCosPAXYMeanVertex = 0.95; ///< min cos of PA to beam line (mean vertex) in tr. plane for prompt V0 candidates
float minCosPAXYMeanVertexCascV0 = 0.8; ///< min cos of PA to beam line (mean vertex) in tr. plane for V0 of cascade cand.
float minCosPAXYMeanVertex = 0.95; ///< min cos of PA to beam line (mean vertex) in tr. plane for prompt V0 candidates
float minCosPAXYMeanVertexCascV0 = 0.8; ///< min cos of PA to beam line (mean vertex) in tr. plane for V0 of cascade cand.
float minCosPAXYMeanVertex3bodyV0 = 0.9; ///< min cos of PA to beam line (mean vertex) in tr. plane for 3body V0 cand.

float maxRToMeanVertexCascV0 = 80; // don't consider as a cascade V0 seed if above this R
Expand All @@ -85,10 +85,10 @@ struct SVertexerParams : public o2::conf::ConfigurableParamHelper<SVertexerParam
float maxDCAZ3Body = 1.; // max DCA of 3 body decay to PV in Z
float minCosPACasc = 0.7; // min cos of PA to PV for cascade candidates
float minCosPA3body = 0.8; // min cos of PA to PV for 3body decay candidates
float minPtCasc = 0.01; // cascade minimum pT
float maxTglCasc = 2.; // maximum tgLambda of cascade
float minPt3Body = 0.01; // minimum pT of 3body V0
float maxTgl3Body = 2.; // maximum tgLambda of 3body V0
float minPtCasc = 0.01; // cascade minimum pT
float maxTglCasc = 2.; // maximum tgLambda of cascade
float minPt3Body = 0.01; // minimum pT of 3body V0
float maxTgl3Body = 2.; // maximum tgLambda of 3body V0

float maxRIni3body = 90.; // don't consider as a 3body seed (circles/line intersection) if its R exceeds this

Expand Down
32 changes: 24 additions & 8 deletions Detectors/Vertexing/src/SVertexer.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -334,6 +334,21 @@ void SVertexer::init()
<< "header=" << header
<< "R=" << R
<< "\n";
// if (comb == 54) {
// mcparticle.Print();
// int i = 0;
// auto mPtr = o2::mcutils::MCTrackNavigator::getMother(mcparticle, pcontainer);
// while (mPtr != nullptr) {
// auto const& mm = *mPtr;
// LOGP(info, "Level: {}", i);
// mm.Print();
// mPtr == o2::mcutils::MCTrackNavigator::getMother(mm, pcontainer);
// --i;
// }
// d0->Print();
// d1->Print();
// LOGP(info, "~~~");
// }
}
}
}
Expand All @@ -355,6 +370,7 @@ void SVertexer::init()
LOGP(info, "~~~~~~~~~~~~~~~~~~~V0 pairs generated~~~~~~~~~~");
mCounterMC.printMC2();
LOGP(info, "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~");
// LOGP(fatal, "just stop");
}
}
}
Expand Down Expand Up @@ -674,13 +690,13 @@ bool SVertexer::checkV0(const TrackCand& seedP, const TrackCand& seedN, int iP,
int nCand = fitterV0.process(seedP, seedN);
if (nCand == 0) { // discard this pair
if (check) {
fitterV0.setDebug(nProcess++);
LOGP(info, "Label: 0");
seedP.gid.print();
lbl0.print();
LOGP(info, "Label: 1");
lbl1.print();
seedN.gid.print();
fitterV0.setDebug(nProcess++, true);
// LOGP(info, "Label: 0");
// seedP.gid.print();
// lbl0.print();
// LOGP(info, "Label: 1");
// lbl1.print();
// seedN.gid.print();
fitterV0.process(seedP, seedN);
fitterV0.unsetDebug();
}
Expand Down Expand Up @@ -829,7 +845,7 @@ bool SVertexer::checkV0(const TrackCand& seedP, const TrackCand& seedN, int iP,
}
}
}
LOG_IF(info, check) << "[------] V0: dca2(=" << dca2 << ") > mMaxDCAXY2ToMeanVertex(=" << mMaxDCAXY2ToMeanVertex << "); cosPAXY(=" << cosPAXY << ") < minCosPAXYMeanVertex(=" << mSVParams->minCosPAXYMeanVertex << ")";
// LOG_IF(info, check) << "[------] V0: dca2(=" << dca2 << ") > mMaxDCAXY2ToMeanVertex(=" << mMaxDCAXY2ToMeanVertex << "); cosPAXY(=" << cosPAXY << ") < minCosPAXYMeanVertex(=" << mSVParams->minCosPAXYMeanVertex << ")";
mCounterV0.inc(CHECKV0::ACOSPAXY, {}, pVtx, pVtxLbl, fitterV0.getPCACandidate(cand), seedP, seedN, lbl0, lbl1, ok, mD0V0Map, mD1V0Map, mcReader, mDebugStream, true, check, fitterV0.isPropagationFailure(), -1, cosPAXY, dca2);

auto vlist = seedP.vBracket.getOverlap(seedN.vBracket); // indices of vertices shared by both seeds
Expand Down

0 comments on commit 6a27aa4

Please sign in to comment.