diff --git a/ALICE3/Core/Decayer.h b/ALICE3/Core/Decayer.h index 22cc46139cd..33696b19683 100644 --- a/ALICE3/Core/Decayer.h +++ b/ALICE3/Core/Decayer.h @@ -69,7 +69,7 @@ class Decayer vy = pos[1] + rxyz * (mom[1] / track.getP()); vz = pos[2] + rxyz * (mom[2] / track.getP()); px = mom[0]; - py = mom[2]; + py = mom[1]; } else { float sna, csa; o2::math_utils::CircleXYf_t circle; diff --git a/ALICE3/Core/TrackUtilities.h b/ALICE3/Core/TrackUtilities.h index dcfb3bfb245..5ddce03fa15 100644 --- a/ALICE3/Core/TrackUtilities.h +++ b/ALICE3/Core/TrackUtilities.h @@ -33,8 +33,10 @@ struct OTFParticle { float mE; float mVx, mVy, mVz; float mPx, mPy, mPz; + bool mIsAlive; // Setters + void setIsAlive(bool isAlive) { mIsAlive = isAlive; } void setPDG(int pdg) { mPdgCode = pdg; } void setVxVyVz(float vx, float vy, float vz) { @@ -52,6 +54,7 @@ struct OTFParticle { // Getters int pdgCode() const { return mPdgCode; } + int isAlive() const { return mIsAlive; } float vx() const { return mVx; } float vy() const { return mVy; } float vz() const { return mVz; } @@ -97,7 +100,9 @@ void convertTLorentzVectorToO2Track(int pdgCode, /// \param o2track the address of the resulting TrackParCov /// \param pdg the pdg service template -void convertOTFParticleToO2Track(const OTFParticle& particle, o2::track::TrackParCov& o2track, const PdgService& pdg) +void convertOTFParticleToO2Track(const OTFParticle& particle, + o2::track::TrackParCov& o2track, + const PdgService& pdg) { static TLorentzVector tlv; tlv.SetPxPyPzE(particle.px(), particle.py(), particle.pz(), particle.e()); diff --git a/ALICE3/DataModel/OTFMCParticle.h b/ALICE3/DataModel/OTFMCParticle.h index 92f5af37c3e..515693edfbe 100644 --- a/ALICE3/DataModel/OTFMCParticle.h +++ b/ALICE3/DataModel/OTFMCParticle.h @@ -32,6 +32,9 @@ DECLARE_SOA_COLUMN(Eta, eta, float); DECLARE_SOA_COLUMN(Pt, pt, float); DECLARE_SOA_COLUMN(P, p, float); DECLARE_SOA_COLUMN(Y, y, float); +DECLARE_SOA_COLUMN(IsAlive, isAlive, bool); +DECLARE_SOA_COLUMN(IsPrimary, isPrimary, bool); + } // namespace otfmcparticle DECLARE_SOA_TABLE_FULL(McParticlesWithDau, "McParticlesWithDau", "AOD", "MCPARTICLEWITHDAU", @@ -56,6 +59,8 @@ DECLARE_SOA_TABLE_FULL(McParticlesWithDau, "McParticlesWithDau", "AOD", "MCPARTI otfmcparticle::Pt, otfmcparticle::P, otfmcparticle::Y, + otfmcparticle::IsAlive, + otfmcparticle::IsPrimary, mcparticle::PVector, mcparticle::ProducedByGenerator, mcparticle::FromBackgroundEvent, diff --git a/ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx b/ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx index f61bab5605f..f8bcb033654 100644 --- a/ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx +++ b/ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx @@ -51,7 +51,7 @@ using namespace o2::framework; static constexpr int NumDecays = 7; static constexpr int NumParameters = 1; static constexpr int DefaultParameters[NumDecays][NumParameters]{{1}, {1}, {1}, {1}, {1}, {1}, {1}}; -static constexpr float VerySmall = 1e-7f; +static constexpr float Tolerance = 1e-7f; static const std::vector parameterNames{"enable"}; static const std::vector particleNames{"K0s", "Lambda", @@ -102,6 +102,12 @@ struct OnTheFlyDecayer { histos.add("AntiXi/hGenAntiXi", "hGenAntiXi", kTH1D, {axisPt}); histos.add("Omega/hGenOmega", "hGenOmega", kTH1D, {axisPt}); histos.add("AntiOmega/hGenAntiOmega", "hGenAntiOmega", kTH1D, {axisPt}); + + histos.add("GeneratedElectron/hGenEl", "hGenEl", kTH1D, {axisPt}); + histos.add("GeneratedMuon/hGenMu", "hGenMu", kTH1D, {axisPt}); + histos.add("GeneratedPion/hGenPi", "hGenPi", kTH1D, {axisPt}); + histos.add("GeneratedKaon/hGenKa", "hGenKa", kTH1D, {axisPt}); + histos.add("GeneratedProton/hGenPr", "hGenPr", kTH1D, {axisPt}); } bool canDecay(const int pdgCode) @@ -123,10 +129,11 @@ struct OnTheFlyDecayer { o2::upgrade::convertMCParticleToO2Track(particle, o2track, pdgDB); decayDaughters = decayer.decayParticle(pdgDB, o2track, particle.pdgCode()); for (size_t idau{0}; idau < decayDaughters.size(); ++idau) { - const auto& dau = decayDaughters[idau]; + o2::upgrade::OTFParticle& dau = decayDaughters[idau]; o2::track::TrackParCov dauTrack; o2::upgrade::convertOTFParticleToO2Track(dau, dauTrack, pdgDB); if (canDecay(dau.pdgCode())) { + dau.setIsAlive(false); std::vector cascadingDaughers = decayer.decayParticle(pdgDB, dauTrack, dau.pdgCode()); for (const auto& daudau : cascadingDaughers) { decayDaughters.push_back(daudau); @@ -134,6 +141,8 @@ struct OnTheFlyDecayer { LOG(error) << "Seemingly stuck trying to perpetually decay products from pdg: " << particle.pdgCode(); } } + } else { + dau.setIsAlive(true); } } @@ -141,22 +150,39 @@ struct OnTheFlyDecayer { LOG(error) << "Attempted to decay " << particle.pdgCode() << " but resulting vector of daugthers were empty"; continue; } - } - if (particle.pdgCode() == kK0Short) { - histos.fill(HIST("K0S/hGenK0S"), particle.pt()); - } else if (particle.pdgCode() == kLambda0) { - histos.fill(HIST("Lambda/hGenLambda"), particle.pt()); - } else if (particle.pdgCode() == kLambda0Bar) { - histos.fill(HIST("AntiLambda/hGenAntiLambda"), particle.pt()); - } else if (particle.pdgCode() == kXiMinus) { - histos.fill(HIST("Xi/hGenXi"), particle.pt()); - } else if (particle.pdgCode() == kXiPlusBar) { - histos.fill(HIST("AntiXi/hGenAntiXi"), particle.pt()); - } else if (particle.pdgCode() == kOmegaMinus) { - histos.fill(HIST("Omega/hGenOmega"), particle.pt()); - } else if (particle.pdgCode() == kOmegaPlusBar) { - histos.fill(HIST("AntiOmega/hGenAntiOmega"), particle.pt()); + switch (particle.pdgCode()) { + case kK0Short: + histos.fill(HIST("K0S/hGenK0S"), particle.pt()); + break; + + case kLambda0: + histos.fill(HIST("Lambda/hGenLambda"), particle.pt()); + break; + + case kLambda0Bar: + histos.fill(HIST("AntiLambda/hGenAntiLambda"), particle.pt()); + break; + + case kXiMinus: + histos.fill(HIST("Xi/hGenXi"), particle.pt()); + break; + + case kXiPlusBar: + histos.fill(HIST("AntiXi/hGenAntiXi"), particle.pt()); + break; + + case kOmegaMinus: + histos.fill(HIST("Omega/hGenOmega"), particle.pt()); + break; + + case kOmegaPlusBar: + histos.fill(HIST("AntiOmega/hGenAntiOmega"), particle.pt()); + break; + + default: + break; + } } int daughtersIdSlice[2]; @@ -178,13 +204,13 @@ struct OnTheFlyDecayer { float p = std::sqrt(particle.px() * particle.px() + particle.py() * particle.py() + particle.pz() * particle.pz()); float y; // As https://github.com/AliceO2Group/AliceO2/blob/dev/Framework/Core/include/Framework/AnalysisDataModel.h#L1943 - if ((p - particle.pz()) < VerySmall) { + if ((p - particle.pz()) < Tolerance) { eta = (particle.pz() < 0.0f) ? -100.0f : 100.0f; } else { eta = 0.5f * std::log((p + particle.pz()) / (p - particle.pz())); } - if ((particle.e() - particle.pz()) < VerySmall) { + if ((particle.e() - particle.pz()) < Tolerance) { y = (particle.pz() < 0.0f) ? -100.0f : 100.0f; } else { y = 0.5f * std::log((particle.e() + particle.pz()) / (particle.e() - particle.pz())); @@ -196,7 +222,7 @@ struct OnTheFlyDecayer { particle.flags(), motherSpan, daughtersIdSlice, particle.weight(), particle.px(), particle.py(), particle.pz(), particle.e(), particle.vx(), particle.vy(), particle.vz(), particle.vt(), - phi, eta, pt, p, y); + phi, eta, pt, p, y, !canDecay(particle.pdgCode()), true); } int daughtersIdSlice[2] = {-1, -1}; @@ -217,18 +243,43 @@ struct OnTheFlyDecayer { float p = std::sqrt(dau.px() * dau.px() + dau.py() * dau.py() + dau.pz() * dau.pz()); float y; // Conditional as https://github.com/AliceO2Group/AliceO2/blob/dev/Framework/Core/include/Framework/AnalysisDataModel.h#L1943 - if ((p - dau.pz()) < VerySmall) { + if ((p - dau.pz()) < Tolerance) { eta = (dau.pz() < 0.0f) ? -100.0f : 100.0f; } else { eta = 0.5f * std::log((p + dau.pz()) / (p - dau.pz())); } - if ((dau.e() - dau.pz()) < VerySmall) { + if ((dau.e() - dau.pz()) < Tolerance) { y = (dau.pz() < 0.0f) ? -100.0f : 100.0f; } else { y = 0.5f * std::log((dau.e() + dau.pz()) / (dau.e() - dau.pz())); } + switch (dau.pdgCode()) { + case kElectron: + histos.fill(HIST("GeneratedElectron/hGenEl"), pt); + break; + + case kMuonMinus: + histos.fill(HIST("GeneratedMuon/hGenMu"), pt); + break; + + case kPiPlus: + histos.fill(HIST("GeneratedPion/hGenPi"), pt); + break; + + case kKPlus: + histos.fill(HIST("GeneratedKaon/hGenKa"), pt); + break; + + case kProton: + histos.fill(HIST("GeneratedProton/hGenPr"), pt); + break; + + default: + break; + } + // TODO: Particle status code // TODO: Expression columns // TODO: vt @@ -236,7 +287,7 @@ struct OnTheFlyDecayer { mother.flags(), motherSpan, daughtersIdSlice, mother.weight(), dau.px(), dau.py(), dau.pz(), dau.e(), dau.vx(), dau.vy(), dau.vz(), mother.vt(), - phi, eta, pt, p, y); + phi, eta, pt, p, y, dau.isAlive(), false); } } } diff --git a/ALICE3/TableProducer/OTF/onTheFlyTracker.cxx b/ALICE3/TableProducer/OTF/onTheFlyTracker.cxx index 252a9c6cab9..fcaa1ec994c 100644 --- a/ALICE3/TableProducer/OTF/onTheFlyTracker.cxx +++ b/ALICE3/TableProducer/OTF/onTheFlyTracker.cxx @@ -27,6 +27,7 @@ #include "ALICE3/Core/DetLayer.h" #include "ALICE3/Core/FastTracker.h" #include "ALICE3/Core/TrackUtilities.h" +#include "ALICE3/DataModel/OTFMCParticle.h" #include "ALICE3/DataModel/OTFStrangeness.h" #include "ALICE3/DataModel/OTFTracks.h" #include "ALICE3/DataModel/collisionAlice3.h" @@ -266,6 +267,17 @@ struct OnTheFlyTracker { kLambda0, kLambda0Bar}; + static constexpr std::array longLivedHandledPDGs = {kElectron, + kMuonMinus, + kPiPlus, + kKPlus, + kProton}; + + static constexpr std::array nucleiPDGs = {o2::constants::physics::kDeuteron, + o2::constants::physics::kTriton, + o2::constants::physics::kHelium3, + o2::constants::physics::kAlpha}; + // necessary for particle charges Service pdgDB; @@ -690,17 +702,6 @@ struct OnTheFlyTracker { auto ir = irSampler.generateCollisionTime(); const float eventCollisionTimeNS = ir.timeInBCNS; - constexpr std::array longLivedHandledPDGs = {kElectron, - kMuonMinus, - kPiPlus, - kKPlus, - kProton}; - - constexpr std::array nucleiPDGs = {o2::constants::physics::kDeuteron, - o2::constants::physics::kTriton, - o2::constants::physics::kHelium3, - o2::constants::physics::kAlpha}; - // First we compute the number of charged particles in the event dNdEta = 0.f; LOG(debug) << "Processing " << mcParticles.size() << " MC particles to compute dNch/deta"; @@ -1588,13 +1589,278 @@ struct OnTheFlyTracker { LOG(debug) << " <- Finished processing OTF tracking with LUT configuration ID " << icfg; } // end process - void process(aod::McCollision const& mcCollision, aod::McParticles const& mcParticles) + void processOnTheFly(aod::McCollision const& mcCollision, aod::McParticles const& mcParticles) { for (size_t icfg = 0; icfg < mSmearer.size(); ++icfg) { LOG(debug) << " -> Processing OTF tracking with LUT configuration ID " << icfg; processWithLUTs(mcCollision, mcParticles, static_cast(icfg)); } } + + void processConfigurationDev(aod::McCollision const& mcCollision, aod::McParticlesWithDau const& mcParticles, int const& icfg) + { + // const int lastTrackIndex = tableStoredTracksCov.lastIndex() + 1; // bookkeep the last added track + const std::string histPath = "Configuration_" + std::to_string(icfg) + "/"; + tracksAlice3.clear(); + ghostTracksAlice3.clear(); + bcData.clear(); + cascadesAlice3.clear(); + v0sAlice3.clear(); + + o2::dataformats::DCA dcaInfo; + o2::dataformats::VertexBase vtx; + + // generate collision time + auto ir = irSampler.generateCollisionTime(); + const float eventCollisionTimeNS = ir.timeInBCNS; + + // First we compute the number of charged particles in the event + dNdEta = 0.f; + for (const auto& mcParticle : mcParticles) { + if (std::abs(mcParticle.eta()) > multEtaRange) { + continue; + } + + if (!mcParticle.isPhysicalPrimary()) { + continue; + } + + const auto pdg = std::abs(mcParticle.pdgCode()); + const bool longLivedToBeHandled = std::find(longLivedHandledPDGs.begin(), longLivedHandledPDGs.end(), pdg) != longLivedHandledPDGs.end(); + const bool nucleiToBeHandled = std::find(nucleiPDGs.begin(), nucleiPDGs.end(), pdg) != nucleiPDGs.end(); + const bool pdgsToBeHandled = longLivedToBeHandled || (enableNucleiSmearing && nucleiToBeHandled); + if (!pdgsToBeHandled) { + continue; + } + + const auto& pdgInfo = pdgDB->GetParticle(mcParticle.pdgCode()); + if (!pdgInfo) { + LOG(warning) << "PDG code " << mcParticle.pdgCode() << " not found in the database"; + continue; + } + if (pdgInfo->Charge() == 0) { + continue; + } + dNdEta += 1.f; + } + + dNdEta /= (multEtaRange * 2.0f); + uint32_t multiplicityCounter = 0; + getHist(TH1, histPath + "hLUTMultiplicity")->Fill(dNdEta); + + // Now that the multiplicity is known, we can process the particles to smear them + for (const auto& mcParticle : mcParticles) { + const bool longLivedToBeHandled = std::find(longLivedHandledPDGs.begin(), longLivedHandledPDGs.end(), std::abs(mcParticle.pdgCode())) != longLivedHandledPDGs.end(); + const bool nucleiToBeHandled = std::find(nucleiPDGs.begin(), nucleiPDGs.end(), std::abs(mcParticle.pdgCode())) != nucleiPDGs.end(); + const bool pdgsToBeHandled = longLivedToBeHandled || (enableNucleiSmearing && nucleiToBeHandled); + if (!pdgsToBeHandled) { + continue; + } + + if (std::fabs(mcParticle.eta()) > maxEta) { + continue; + } + + if (mcParticle.pt() < minPt) { + continue; + } + + if (enablePrimarySmearing) { + getHist(TH1, histPath + "hPtGenerated")->Fill(mcParticle.pt()); + getHist(TH1, histPath + "hPhiGenerated")->Fill(mcParticle.phi()); + switch (std::abs(mcParticle.pdgCode())) { + case kElectron: + getHist(TH1, histPath + "hPtGeneratedEl")->Fill(mcParticle.pt()); + break; + case kPiPlus: + getHist(TH1, histPath + "hPtGeneratedPi")->Fill(mcParticle.pt()); + break; + case kKPlus: + getHist(TH1, histPath + "hPtGeneratedKa")->Fill(mcParticle.pt()); + break; + case kProton: + getHist(TH1, histPath + "hPtGeneratedPr")->Fill(mcParticle.pt()); + break; + } + } + + multiplicityCounter++; + o2::track::TrackParCov trackParCov; + const bool isDecayDaughter = (mcParticle.getProcess() == TMCProcess::kPDecay); + const float timeResolutionNs = 100.f; // ns + const float nsToMus = 1e-3f; + const float timeResolutionUs = timeResolutionNs * nsToMus; // us + const float time = (eventCollisionTimeNS + gRandom->Gaus(0., timeResolutionNs)) * nsToMus; + + bool reconstructed = false; + if (enablePrimarySmearing && mcParticle.isPrimary()) { + o2::upgrade::convertMCParticleToO2Track(mcParticle, trackParCov, pdgDB); + reconstructed = mSmearer[icfg]->smearTrack(trackParCov, mcParticle.pdgCode(), dNdEta); + } else if (enableSecondarySmearing) { + o2::track::TrackParCov perfectTrackParCov; + o2::upgrade::convertMCParticleToO2Track(mcParticle, perfectTrackParCov, pdgDB); + const int nHits = fastTracker[icfg]->FastTrack(perfectTrackParCov, trackParCov, dNdEta); + if (nHits < fastTrackerSettings.minSiliconHits) { + reconstructed = false; + } else { + reconstructed = true; + } + } + + if (!reconstructed && processUnreconstructedTracks) { + continue; // failed to reconstruct track + } + + if (std::isnan(trackParCov.getZ())) { + histos.fill(HIST("hNaNBookkeeping"), 0.0f, 0.0f); + continue; // capture smearing mistakes + } + + histos.fill(HIST("hNaNBookkeeping"), 0.0f, 1.0f); + if (enablePrimarySmearing) { + getHist(TH1, histPath + "hPtReconstructed")->Fill(trackParCov.getPt()); + if (std::abs(mcParticle.pdgCode()) == kElectron) + getHist(TH1, histPath + "hPtReconstructedEl")->Fill(trackParCov.getPt()); + if (std::abs(mcParticle.pdgCode()) == kPiPlus) + getHist(TH1, histPath + "hPtReconstructedPi")->Fill(trackParCov.getPt()); + if (std::abs(mcParticle.pdgCode()) == kKPlus) + getHist(TH1, histPath + "hPtReconstructedKa")->Fill(trackParCov.getPt()); + if (std::abs(mcParticle.pdgCode()) == kProton) + getHist(TH1, histPath + "hPtReconstructedPr")->Fill(trackParCov.getPt()); + } + + if (reconstructed) { + tracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isDecayDaughter}); + } else { + ghostTracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isDecayDaughter}); + } + } + + o2::vertexing::PVertex primaryVertex; + if (enablePrimaryVertexing) { // disabled for now + primaryVertex.setXYZ(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ()); + } else { + primaryVertex.setXYZ(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ()); + } + + getHist(TH1, histPath + "hSimMultiplicity")->Fill(multiplicityCounter); + getHist(TH1, histPath + "hRecoMultiplicity")->Fill(tracksAlice3.size()); + getHist(TH1, histPath + "hPVz")->Fill(primaryVertex.getZ()); + + // *+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+* + // populate collisions + tableCollisions(-1, // BC is irrelevant in synthetic MC tests for now, could be adjusted in future + primaryVertex.getX(), primaryVertex.getY(), primaryVertex.getZ(), + primaryVertex.getSigmaX2(), primaryVertex.getSigmaXY(), primaryVertex.getSigmaY2(), + primaryVertex.getSigmaXZ(), primaryVertex.getSigmaYZ(), primaryVertex.getSigmaZ2(), + 0, primaryVertex.getChi2(), primaryVertex.getNContributors(), + eventCollisionTimeNS, 0.f); // For the moment the event collision time is taken as the "GEANT" time, the computation of the event time is done a posteriori from the tracks in the OTF TOF PID task + tableMcCollisionLabels(mcCollision.globalIndex(), 0); + tableCollisionsAlice3(dNdEta); + tableOTFLUTConfigId(icfg); // Populate OTF LUT configuration ID + + // *+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+* + // populate tracks + for (const auto& trackParCov : tracksAlice3) { + aod::track::TrackTypeEnum trackType = aod::track::Track; + + if (populateTracksDCA) { + float dcaXY = 1e+10, dcaZ = 1e+10; + o2::track::TrackParCov trackParametrization(trackParCov); + if (trackParametrization.propagateToDCA(primaryVertex, mMagneticField, &dcaInfo)) { + dcaXY = dcaInfo.getY(); + dcaZ = dcaInfo.getZ(); + } + + tableTracksDCA(dcaXY, dcaZ); + if (populateTracksDCACov) { + tableTracksDCACov(dcaInfo.getSigmaY2(), dcaInfo.getSigmaZ2()); + } + } + + tableStoredTracks(tableCollisions.lastIndex(), trackType, trackParCov.getX(), trackParCov.getAlpha(), trackParCov.getY(), trackParCov.getZ(), trackParCov.getSnp(), trackParCov.getTgl(), trackParCov.getQ2Pt()); + tableTracksExtension(trackParCov.getPt(), trackParCov.getP(), trackParCov.getEta(), trackParCov.getPhi()); + + // TODO do we keep the rho as 0? Also the sigma's are duplicated information + tableStoredTracksCov(std::sqrt(trackParCov.getSigmaY2()), std::sqrt(trackParCov.getSigmaZ2()), std::sqrt(trackParCov.getSigmaSnp2()), + std::sqrt(trackParCov.getSigmaTgl2()), std::sqrt(trackParCov.getSigma1Pt2()), 0, 0, 0, 0, 0, 0, 0, 0, 0, 0); + tableTracksCovExtension(trackParCov.getSigmaY2(), trackParCov.getSigmaZY(), trackParCov.getSigmaZ2(), trackParCov.getSigmaSnpY(), + trackParCov.getSigmaSnpZ(), trackParCov.getSigmaSnp2(), trackParCov.getSigmaTglY(), trackParCov.getSigmaTglZ(), trackParCov.getSigmaTglSnp(), + trackParCov.getSigmaTgl2(), trackParCov.getSigma1PtY(), trackParCov.getSigma1PtZ(), trackParCov.getSigma1PtSnp(), trackParCov.getSigma1PtTgl(), + trackParCov.getSigma1Pt2()); + tableMcTrackLabels(trackParCov.mcLabel, 0); + tableTracksExtraA3(trackParCov.nSiliconHits, trackParCov.nTPCHits); + + // populate extra tables if required to do so + if (populateTracksExtra) { + tableStoredTracksExtra(0.0f, static_cast(0), static_cast(0), static_cast(0), static_cast(0), + static_cast(0), static_cast(0), static_cast(0), static_cast(0), + 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, + 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f); + } + if (populateTrackSelection) { + tableTrackSelection(static_cast(0), false, false, false, false, false, false); + tableTrackSelectionExtension(false, false, false, false, false, false, false, false, false, false, false, false, false, false, false, false, false); + } + tableTracksAlice3(true); + } + + // *+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+*+~+* + // populate ghost tracks + for (const auto& trackParCov : ghostTracksAlice3) { + aod::track::TrackTypeEnum trackType = aod::track::Track; + + if (populateTracksDCA) { + float dcaXY = 1e+10, dcaZ = 1e+10; + o2::track::TrackParCov trackParametrization(trackParCov); + if (trackParametrization.propagateToDCA(primaryVertex, mMagneticField, &dcaInfo)) { + dcaXY = dcaInfo.getY(); + dcaZ = dcaInfo.getZ(); + } + + tableTracksDCA(dcaXY, dcaZ); + if (populateTracksDCACov) { + tableTracksDCACov(dcaInfo.getSigmaY2(), dcaInfo.getSigmaZ2()); + } + } + + tableStoredTracks(tableCollisions.lastIndex(), trackType, trackParCov.getX(), trackParCov.getAlpha(), trackParCov.getY(), trackParCov.getZ(), trackParCov.getSnp(), trackParCov.getTgl(), trackParCov.getQ2Pt()); + tableTracksExtension(trackParCov.getPt(), trackParCov.getP(), trackParCov.getEta(), trackParCov.getPhi()); + + // TODO do we keep the rho as 0? Also the sigma's are duplicated information + tableStoredTracksCov(std::sqrt(trackParCov.getSigmaY2()), std::sqrt(trackParCov.getSigmaZ2()), std::sqrt(trackParCov.getSigmaSnp2()), + std::sqrt(trackParCov.getSigmaTgl2()), std::sqrt(trackParCov.getSigma1Pt2()), 0, 0, 0, 0, 0, 0, 0, 0, 0, 0); + tableTracksCovExtension(trackParCov.getSigmaY2(), trackParCov.getSigmaZY(), trackParCov.getSigmaZ2(), trackParCov.getSigmaSnpY(), + trackParCov.getSigmaSnpZ(), trackParCov.getSigmaSnp2(), trackParCov.getSigmaTglY(), trackParCov.getSigmaTglZ(), trackParCov.getSigmaTglSnp(), + trackParCov.getSigmaTgl2(), trackParCov.getSigma1PtY(), trackParCov.getSigma1PtZ(), trackParCov.getSigma1PtSnp(), trackParCov.getSigma1PtTgl(), + trackParCov.getSigma1Pt2()); + tableMcTrackLabels(trackParCov.mcLabel, 0); + tableTracksExtraA3(trackParCov.nSiliconHits, trackParCov.nTPCHits); + + // populate extra tables if required to do so + if (populateTracksExtra) { + tableStoredTracksExtra(0.0f, static_cast(0), static_cast(0), static_cast(0), static_cast(0), + static_cast(0), static_cast(0), static_cast(0), static_cast(0), + 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, + 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f); + } + if (populateTrackSelection) { + tableTrackSelection(static_cast(0), false, false, false, false, false, false); + tableTrackSelectionExtension(false, false, false, false, false, false, false, false, false, false, false, false, false, false, false, false, false); + } + tableTracksAlice3(false); + } + } + + void processDecayer(aod::McCollision const& mcCollision, aod::McParticlesWithDau const& mcParticles) + { + for (size_t icfg = 0; icfg < mSmearer.size(); ++icfg) { + processConfigurationDev(mcCollision, mcParticles, static_cast(icfg)); + } + } + + PROCESS_SWITCH(OnTheFlyTracker, processOnTheFly, "Enable default workflow", true); + PROCESS_SWITCH(OnTheFlyTracker, processDecayer, "Enable experimental decayer workflow", false); }; /// Extends TracksExtra if necessary