Skip to content

Commit 876ce90

Browse files
creetz16alibuild
creetz16
and
alibuild
authored
[PWGLF] Updates to trackedHypertritonRecoTask.cxx (#17134)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent 1baffb7 commit 876ce90

4 files changed

Lines changed: 260 additions & 117 deletions

File tree

‎PWGLF/DataModel/Vtx3BodyTables.h‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -140,6 +140,10 @@ DECLARE_SOA_COLUMN(IsSurvEvSel, isSurvEvSel, int); //! flag if reco collisio
140140
DECLARE_SOA_COLUMN(IsReco, isreco, int); //! flag if candidate was reconstructed
141141
DECLARE_SOA_COLUMN(MotherLabel, motherLabel, int); //! label of the mother particle (signal: MC index, bkg: negative value)
142142

143+
// strangeness tracking information
144+
DECLARE_SOA_COLUMN(ITSTrackDCAXYToSVXY, itsTrackDcaXYToSv, float); //! DCAxy of ITS track to SV
145+
DECLARE_SOA_COLUMN(ITSTrackDCAZToSVXY, itsTrackDcaZToSv, float); //! DCAz of ITS track to SV
146+
143147
// Derived expressions
144148
// Momenta
145149
DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, //! 3 body pT in GeV/c
@@ -258,6 +262,11 @@ DECLARE_SOA_TABLE(Vtx3BodyCovs, "AOD", "VTX3BODYCOV", //!
258262
vtx3body::CovProton, vtx3body::CovPion, vtx3body::CovDeuteron,
259263
vtx3body::VtxCovMat);
260264

265+
// strangeness tracked extra table
266+
DECLARE_SOA_TABLE(Vtx3BodyTrackedInfo, "AOD", "VTX3BODYTR", //!
267+
vtx3body::ITSTrackDCAXYToSVXY,
268+
vtx3body::ITSTrackDCAZToSVXY);
269+
261270
// MC candidate table for analysis
262271
DECLARE_SOA_TABLE(McVtx3BodyDatas, "AOD", "MC3BODYDATA", //!
263272
o2::soa::Index<>,

‎PWGLF/TableProducer/Nuspex/decay3bodybuilder.cxx‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1115,10 +1115,10 @@ struct decay3bodyBuilder {
11151115
registry.fill(HIST("Counters/hTableBuildingStatistics"), kVtx3BodyDatas);
11161116
}
11171117
if (mEnabledTables[kVtx3BodyCovs]) {
1118-
products.vtx3bodycovs(helper.decay3body.covProton,
1119-
helper.decay3body.covPion,
1120-
helper.decay3body.covDeuteron,
1121-
helper.decay3body.covariance);
1118+
products.vtx3bodycovs(helper.decay3body.covProton.data(),
1119+
helper.decay3body.covPion.data(),
1120+
helper.decay3body.covDeuteron.data(),
1121+
helper.decay3body.covariance.data());
11221122
registry.fill(HIST("Counters/hTableBuildingStatistics"), kVtx3BodyCovs);
11231123
}
11241124
if (mEnabledTables[kMcVtx3BodyDatas]) {

‎PWGLF/TableProducer/Nuspex/trackedHypertritonRecoTask.cxx‎

Lines changed: 185 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -70,6 +70,8 @@
7070
using namespace o2;
7171
using namespace o2::framework;
7272

73+
o2::common::core::MetadataHelper metadataInfo{};
74+
7375
namespace
7476
{
7577
using Collisions = soa::Join<aod::Collisions, aod::EvSels, aod::CentFT0As, aod::CentFT0Cs, aod::CentFT0Ms, aod::EvTimeTOFFT0>;
@@ -94,12 +96,11 @@ enum ZorroTrigger : std::size_t {
9496
} // namespace
9597

9698
struct TrackedHypertritonRecoTask {
97-
o2::common::core::MetadataHelper metadataInfo{};
98-
9999
Produces<aod::DataHypCands> dataHypCands;
100100
Produces<aod::MCHypCands> mcHypCands;
101101
Produces<aod::Vtx3BodyDatas> vtx3BodyDatas;
102102
Produces<aod::Vtx3BodyCovs> vtx3BodyCovs;
103+
Produces<aod::Vtx3BodyTrackedInfo> vtx3BodyTrackedInfo;
103104
Produces<aod::McVtx3BodyDatas> mcVtx3BodyDatas;
104105

105106
Service<o2::ccdb::BasicCCDBManager> ccdb{};
@@ -124,7 +125,9 @@ struct TrackedHypertritonRecoTask {
124125

125126
struct : ConfigurableGroup {
126127
std::string prefix = "twoBody";
128+
Configurable<bool> useKFParticle{"useKFParticle", false, "Use KFParticle to build the two-body candidate"};
127129
Configurable<bool> useSelections{"useSelections", false, "Apply the two-body candidate selections"};
130+
Configurable<bool> kfSetTopologicalConstraint{"kfSetTopologicalConstraint", false, "Set topological vertex constraint in case of KFParticle reconstruction"};
128131
Configurable<float> maxEtaDaughters{"maxEtaDaughters", 1.e10f, "Maximum absolute daughter eta"};
129132
Configurable<float> minTPCNClsHe{"minTPCNClsHe", -1.f, "Minimum He TPC clusters"};
130133
Configurable<float> minTPCNClsPi{"minTPCNClsPi", -1.f, "Minimum pion TPC clusters"};
@@ -133,7 +136,7 @@ struct TrackedHypertritonRecoTask {
133136
Configurable<float> minTPCInnerParamHe{"minTPCInnerParamHe", -1.f, "Minimum He TPC rigidity"};
134137
Configurable<float> minPt{"minPt", -1.f, "Minimum candidate transverse momentum"};
135138
Configurable<float> massWindow{"massWindow", 1.e10f, "Half-width of the hypertriton mass window"};
136-
Configurable<float> maxDcaDaughters{"maxDcaDaughters", 1.e10f, "Maximum DCA between daughters"};
139+
Configurable<float> maxChi2{"maxChi2", 1.e10f, "KFParticle: Maximum SV chi2, DCA fitter: Maximum DCA between daughters"};
137140
Configurable<float> minCosPA{"minCosPA", -2.f, "Minimum cosine of the pointing angle"};
138141
Configurable<float> minDcaHeToPV{"minDcaHeToPV", -1.f, "Minimum absolute He DCA to PV"};
139142
Configurable<float> minDcaPiToPV{"minDcaPiToPV", -1.f, "Minimum absolute pion DCA to PV"};
@@ -179,6 +182,8 @@ struct TrackedHypertritonRecoTask {
179182
Configurable<float> maxCtau{"maxCtau", 100.f, "Maximum candidate c tau"};
180183
Configurable<float> minCosPA{"minCosPA", 0.9f, "Minimum candidate cosine of the pointing angle"};
181184
Configurable<float> maxChi2{"maxChi2", 100.f, "Maximum candidate chi2"};
185+
Configurable<float> maxITSDCAxytrackToSV{"maxITSDCAxytrackToSV", 10.0, "Maximum distance of ITS matched track to SV in xy"};
186+
Configurable<float> maxITSDCAztrackToSV{"maxITSDCAztrackToSV", 10.0, "Maximum distance of ITS matched track to SV in z"};
182187
} threeBody;
183188

184189
o2::vertexing::DCAFitterN<2> fitter2Body;
@@ -196,6 +201,22 @@ struct TrackedHypertritonRecoTask {
196201
std::vector<int> recoCollisionForMC;
197202
std::vector<bool> survivedMCEventSelection;
198203

204+
struct v0Candidate {
205+
// daughter properties
206+
std::array<float, 3> momHelium{};
207+
std::array<float, 3> momPion{};
208+
std::array<float, 3> posHelium{};
209+
std::array<float, 3> posPion{};
210+
// vertex properties
211+
float mass{};
212+
float chi2{};
213+
float cosPA{};
214+
std::array<float, 3> decayVertex{};
215+
std::array<float, 3> momentum{};
216+
};
217+
218+
v0Candidate v0;
219+
199220
struct TwoBodyMCInfo {
200221
float genPt = -1.f;
201222
float genPhi = -1.f;
@@ -536,52 +557,166 @@ struct TrackedHypertritonRecoTask {
536557
return (track.tpcSignal() - expected) / (expected * bbParamsHe[5]);
537558
}
538559

539-
template <typename TTrack, typename TCollision, typename TFillCandidate>
540-
void buildTwoBody(TTrack const& heTrack, TTrack const& piTrack, TCollision const& collision, float trackedClSize, TFillCandidate const& fillCandidate)
560+
template <typename TCollision, typename TTrack, typename TTrackParCov>
561+
bool fit2BodyWithKF(TCollision const& collision,
562+
TTrack const& trackHelium,
563+
TTrack const& trackPion,
564+
TTrackParCov const& trackHeliumCov,
565+
TTrackParCov const& trackPionCov)
566+
{
567+
// initialise KF primary vertex
568+
KFParticle kfpv(createKFPVertexFromCollision(collision));
569+
570+
// create KFParticle objects
571+
KFParticle kfpHelium, kfpPion;
572+
// helium
573+
std::array<float, 3> xyz{};
574+
std::array<float, 3> pxpypz{};
575+
std::array<float, 6> xyzpxpypz{};
576+
trackHeliumCov.getPxPyPzGlo(pxpypz);
577+
trackHeliumCov.getXYZGlo(xyz);
578+
for (int i = 0; i < 3; ++i) {
579+
xyzpxpypz[i] = xyz[i];
580+
xyzpxpypz[i + 3] = pxpypz[i] * 2;
581+
}
582+
std::array<float, 21> cv{};
583+
trackHeliumCov.getCovXYZPxPyPzGlo(cv);
584+
KFParticle kfHelium;
585+
kfHelium.Create(xyzpxpypz.data(), cv.data(), trackHelium.sign() * 2, constants::physics::MassHelium3);
586+
// pion
587+
kfpPion = createKFParticleFromTrackParCov(trackPionCov, trackPion.sign(), constants::physics::MassPionCharged);
588+
589+
// construct V0 vertex
590+
KFParticle KFV0;
591+
int nDaughtersV0 = 2;
592+
const KFParticle* DaughtersV0[2] = {&kfpHelium, &kfpPion};
593+
KFV0.SetConstructMethod(2);
594+
try {
595+
KFV0.Construct(DaughtersV0, nDaughtersV0);
596+
} catch (std::runtime_error& e) {
597+
LOG(debug) << "Failed to create V0 vertex." << e.what();
598+
return false;
599+
}
600+
601+
// topological constraint
602+
if (twoBody.kfSetTopologicalConstraint) {
603+
KFV0.SetProductionVertex(kfpv);
604+
KFV0.TransportToDecayVertex();
605+
}
606+
607+
// get vertex position and momentum
608+
v0.decayVertex[0] = KFV0.GetX();
609+
v0.decayVertex[1] = KFV0.GetY();
610+
v0.decayVertex[2] = KFV0.GetZ();
611+
v0.momentum[0] = KFV0.GetPx();
612+
v0.momentum[1] = KFV0.GetPy();
613+
v0.momentum[2] = KFV0.GetPz();
614+
615+
// transport all daughter tracks to hypertriton vertex
616+
// float position[3];
617+
// for (int i; i < 3; i++) {
618+
// position[i] = v0.decayVertex[i];
619+
// }
620+
kfpHelium.TransportToPoint(v0.decayVertex.data());
621+
kfpPion.TransportToPoint(v0.decayVertex.data());
622+
623+
// daughter positions
624+
v0.posHelium[0] = kfpHelium.GetX();
625+
v0.posHelium[1] = kfpHelium.GetY();
626+
v0.posHelium[2] = kfpHelium.GetZ();
627+
v0.posPion[0] = kfpPion.GetX();
628+
v0.posPion[1] = kfpPion.GetY();
629+
v0.posPion[2] = kfpPion.GetZ();
630+
631+
// daughter momenta
632+
v0.momHelium[0] = kfpHelium.GetPx();
633+
v0.momHelium[1] = kfpHelium.GetPy();
634+
v0.momHelium[2] = kfpHelium.GetPz();
635+
v0.momPion[0] = kfpPion.GetPx();
636+
v0.momPion[1] = kfpPion.GetPy();
637+
v0.momPion[2] = kfpPion.GetPz();
638+
639+
// candidate mass
640+
float mass{};
641+
float massErr{};
642+
KFV0.GetMass(mass, massErr);
643+
v0.mass = mass;
644+
645+
// vertex chi2
646+
v0.chi2 = KFV0.GetChi2() / KFV0.GetNDF();
647+
648+
return true;
649+
}
650+
651+
template <typename TTrackParCov>
652+
bool fit2bodyWithDCAFitter(TTrackParCov const& trackHeliumCov,
653+
TTrackParCov const& trackPionCov)
541654
{
542-
auto heTrackCov = getTrackParCov(heTrack);
543-
auto piTrackCov = getTrackParCov(piTrack);
544655
int nCandidates = 0;
545656
try {
546-
nCandidates = fitter2Body.process(heTrackCov, piTrackCov);
657+
nCandidates = fitter2Body.process(trackHeliumCov, trackPionCov);
547658
} catch (...) {
548659
LOG(error) << "Exception while fitting a tracked two-body candidate";
549-
return;
660+
return false;
550661
}
551662
if (nCandidates == 0) {
552-
return;
663+
return false;
553664
}
554665

555-
std::array<float, 3> heMomentum{};
556-
std::array<float, 3> piMomentum{};
557-
fitter2Body.getTrack(0).getPxPyPzGlo(heMomentum);
558-
fitter2Body.getTrack(1).getPxPyPzGlo(piMomentum);
559-
for (std::size_t i = 0; i < heMomentum.size(); ++i) {
560-
heMomentum[i] *= 2.f;
561-
}
562-
std::array<float, 3> momentum{heMomentum[0] + piMomentum[0], heMomentum[1] + piMomentum[1], heMomentum[2] + piMomentum[2]};
563-
if (twoBody.useSelections && std::hypot(momentum[0], momentum[1]) < twoBody.minPt) {
564-
return;
666+
// get daughter momenta
667+
fitter2Body.getTrack(0).getPxPyPzGlo(v0.momHelium);
668+
fitter2Body.getTrack(1).getPxPyPzGlo(v0.momPion);
669+
670+
for (std::size_t i = 0; i < v0.momHelium.size(); ++i) {
671+
v0.momHelium[i] *= 2.f;
565672
}
566673

567-
const float heMomentum2 = RecoDecay::sumOfSquares(heMomentum[0], heMomentum[1], heMomentum[2]);
568-
const float piMomentum2 = RecoDecay::sumOfSquares(piMomentum[0], piMomentum[1], piMomentum[2]);
569-
const float candidateMomentum2 = RecoDecay::sumOfSquares(momentum[0], momentum[1], momentum[2]);
674+
// compute candidate momentum
675+
v0.momentum = {v0.momHelium[0] + v0.momPion[0], v0.momHelium[1] + v0.momPion[1], v0.momHelium[2] + v0.momPion[2]};
676+
677+
// compute candidate mass
678+
const float heMomentum2 = RecoDecay::sumOfSquares(v0.momHelium[0], v0.momHelium[1], v0.momHelium[2]);
679+
const float piMomentum2 = RecoDecay::sumOfSquares(v0.momPion[0], v0.momPion[1], v0.momPion[2]);
680+
const float candidateMomentum2 = RecoDecay::sumOfSquares(v0.momentum[0], v0.momentum[1], v0.momentum[2]);
570681
const float heEnergy = std::sqrt(heMomentum2 + constants::physics::MassHelium3 * constants::physics::MassHelium3);
571682
const float piEnergy = std::sqrt(piMomentum2 + constants::physics::MassPionCharged * constants::physics::MassPionCharged);
572-
const float mass = std::sqrt((heEnergy + piEnergy) * (heEnergy + piEnergy) - candidateMomentum2);
573-
if (twoBody.useSelections && std::abs(mass - constants::physics::MassHyperTriton) > twoBody.massWindow) {
574-
return;
575-
}
683+
v0.mass = std::sqrt((heEnergy + piEnergy) * (heEnergy + piEnergy) - candidateMomentum2);
576684

685+
// get SV position
577686
const auto& secondaryVertex = fitter2Body.getPCACandidate();
687+
for (int i = 0; i < 3; i++) {
688+
v0.decayVertex[i] = secondaryVertex[i];
689+
}
690+
v0.chi2 = std::sqrt(fitter2Body.getChi2AtPCACandidate());
691+
692+
return true;
693+
}
694+
695+
template <typename TTrack, typename TCollision, typename TFillCandidate>
696+
void buildTwoBody(TTrack const& heTrack, TTrack const& piTrack, TCollision const& collision, float trackedClSize, TFillCandidate const& fillCandidate)
697+
{
578698
const std::array<float, 3> primaryVertex{collision.posX(), collision.posY(), collision.posZ()};
579-
const std::array<float, 3> decayVertex{static_cast<float>(secondaryVertex[0]), static_cast<float>(secondaryVertex[1]), static_cast<float>(secondaryVertex[2])};
580-
const float dcaDaughters = std::sqrt(fitter2Body.getChi2AtPCACandidate());
581-
if (twoBody.useSelections && (dcaDaughters > twoBody.maxDcaDaughters || RecoDecay::cpa(primaryVertex, decayVertex, momentum) < twoBody.minCosPA)) {
699+
700+
auto heTrackCov = getTrackParCov(heTrack);
701+
auto piTrackCov = getTrackParCov(piTrack);
702+
703+
const bool fitSucceeded = twoBody.useKFParticle
704+
? fit2BodyWithKF(collision, heTrack, piTrack, heTrackCov, piTrackCov)
705+
: fit2bodyWithDCAFitter(heTrackCov, piTrackCov);
706+
if (!fitSucceeded) {
707+
return;
708+
}
709+
710+
v0.cosPA = RecoDecay::cpa(primaryVertex, v0.decayVertex, v0.momentum);
711+
if (twoBody.useSelections &&
712+
(std::hypot(v0.momentum[0], v0.momentum[1]) < twoBody.minPt ||
713+
std::abs(v0.mass - constants::physics::MassHyperTriton) > twoBody.massWindow ||
714+
v0.chi2 > twoBody.maxChi2 ||
715+
v0.cosPA < twoBody.minCosPA)) {
582716
return;
583717
}
584718

719+
// Do propagation with Propagator including material interactions in all cases (KFParticle propagation would not include material)
585720
std::array<float, 2> dcaInfo{};
586721
o2::base::Propagator::Instance()->propagateToDCABxByBz({collision.posX(), collision.posY(), collision.posZ()}, heTrackCov, 2.f, fitter2Body.getMatCorrType(), &dcaInfo);
587722
const float dcaHe = dcaInfo[0];
@@ -605,10 +740,10 @@ struct TrackedHypertritonRecoTask {
605740
fillCandidate(collision.centFT0A(), collision.centFT0C(), collision.centFT0M(),
606741
collision.posX(), collision.posY(), collision.posZ(),
607742
runNumber, heTrack.sign() > 0,
608-
std::hypot(heMomentum[0], heMomentum[1]), std::atan2(heMomentum[1], heMomentum[0]), RecoDecay::eta(heMomentum),
609-
std::hypot(piMomentum[0], piMomentum[1]), std::atan2(piMomentum[1], piMomentum[0]), RecoDecay::eta(piMomentum),
610-
secondaryVertex[0], secondaryVertex[1], secondaryVertex[2],
611-
dcaDaughters, dcaHe, dcaPi,
743+
std::hypot(v0.momHelium[0], v0.momHelium[1]), std::atan2(v0.momHelium[1], v0.momHelium[0]), RecoDecay::eta(v0.momHelium),
744+
std::hypot(v0.momPion[0], v0.momPion[1]), std::atan2(v0.momPion[1], v0.momPion[0]), RecoDecay::eta(v0.momPion),
745+
v0.decayVertex[0], v0.decayVertex[1], v0.decayVertex[2],
746+
v0.chi2, dcaHe, dcaPi,
612747
nSigmaHe3(heTrack), heTrack.tpcNClsFound(), piTrack.tpcNClsFound(),
613748
static_cast<int16_t>(heTrack.tpcNClsFindable()) - heTrack.tpcNClsFindableMinusPID(),
614749
static_cast<int16_t>(piTrack.tpcNClsFindable()) - piTrack.tpcNClsFindableMinusPID(),
@@ -661,7 +796,8 @@ struct TrackedHypertritonRecoTask {
661796
candidate.averageITSClSize[0], candidate.averageITSClSize[1], candidate.averageITSClSize[2],
662797
static_cast<int>(candidate.tpcNCl[0]), static_cast<int>(candidate.tpcNCl[1]), static_cast<int>(candidate.tpcNCl[2]),
663798
static_cast<uint32_t>(candidate.pidForTrackingDeuteron));
664-
vtx3BodyCovs(candidate.covProton, candidate.covPion, candidate.covDeuteron, candidate.covariance);
799+
vtx3BodyCovs(candidate.covProton.data(), candidate.covPion.data(), candidate.covDeuteron.data(), candidate.covariance.data());
800+
vtx3BodyTrackedInfo(candidate.itsTrackDCAToSV[0], candidate.itsTrackDCAToSV[1]);
665801
}
666802

667803
void fillThreeBodyMCTable(ThreeBodyMCInfo const& info)
@@ -776,17 +912,30 @@ struct TrackedHypertritonRecoTask {
776912
if (decay3Body.collisionId() < 0 || !goodCollision[decay3Body.collisionId()] || (skimmedProcessing && !zorroDecision[decay3Body.collisionId()][kTracked3Body])) {
777913
continue;
778914
}
915+
779916
const auto collision = decay3Body.collision_as<Collisions>();
780917
const auto trackPositive = decay3Body.track0_as<Tracks>();
781918
const auto trackNegative = decay3Body.track1_as<Tracks>();
782919
const auto trackDeuteron = decay3Body.track2_as<Tracks>();
783920
const auto trackProton = trackDeuteron.sign() > 0 ? trackPositive : trackNegative;
784921
const auto trackPion = trackDeuteron.sign() > 0 ? trackNegative : trackPositive;
922+
785923
if (builder3Body.buildDecay3BodyCandidate(collision, trackProton, trackPion, trackDeuteron,
786924
decay3Body.globalIndex(), deuteronTOFNSigma(collision, trackDeuteron), tracked3Body.itsClsSize(),
787925
threeBody.useKFParticle, threeBody.setTopologicalConstraint,
788926
threeBody.useSelections, threeBody.useChi2Selection, threeBody.useTPCforPion,
789927
threeBody.acceptTPCOnly, threeBody.askOnlyITSMatch, threeBody.calculateCovariance)) {
928+
// get DCA of ITS track to SV
929+
const auto itsTrack = tracked3Body.itsTrack_as<Tracks>();
930+
auto itsTrackParCov = getTrackParCov(itsTrack);
931+
std::array<float, 2> dcaInfoItsTrack{};
932+
o2::base::Propagator::Instance()->propagateToDCABxByBz({builder3Body.decay3body.position[0], builder3Body.decay3body.position[1], builder3Body.decay3body.position[2]}, itsTrackParCov, 2.f, fitter2Body.getMatCorrType(), &dcaInfoItsTrack);
933+
builder3Body.decay3body.itsTrackDCAToSV[0] = dcaInfoItsTrack[0];
934+
builder3Body.decay3body.itsTrackDCAToSV[1] = dcaInfoItsTrack[1];
935+
if (threeBody.useSelections && (std::abs(builder3Body.decay3body.itsTrackDCAToSV[0]) > threeBody.maxITSDCAxytrackToSV || std::abs(builder3Body.decay3body.itsTrackDCAToSV[1]) > threeBody.maxITSDCAztrackToSV)) {
936+
continue;
937+
}
938+
790939
fillThreeBodyTables();
791940
}
792941
}
@@ -987,7 +1136,6 @@ struct TrackedHypertritonRecoTask {
9871136

9881137
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
9891138
{
990-
auto metadataInfo = o2::common::core::MetadataHelper{};
9911139
metadataInfo.initMetadata(cfgc);
992-
return WorkflowSpec{adaptAnalysisTask<TrackedHypertritonRecoTask>(cfgc, std::move(metadataInfo))};
1140+
return WorkflowSpec{adaptAnalysisTask<TrackedHypertritonRecoTask>(cfgc)};
9931141
}

0 commit comments

Comments
 (0)