Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
17 changes: 7 additions & 10 deletions ALICE3/Core/OTFParticle.h
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,7 @@ class OTFParticle
mVt = particle.vt();
mFlag = particle.flags();
mStatusCode = particle.statusCode();
mIsFromMcParticles = true;
setBitOff(DecayerBits::ProducedByDecayer);
if (particle.has_mothers()) {
mIndicesMother = {particle.mothersIds().front(), particle.mothersIds().back()};
}
Expand All @@ -72,7 +72,6 @@ class OTFParticle
}

// Setters
void setIsPrimary(const bool isPrimary) { mIsPrimary = isPrimary; }
void setCollisionId(const int collisionId) { mCollisionId = collisionId; }
void setPDG(const int pdg) { mPdgCode = pdg; }
void setIndicesMother(const int start, const int stop) { mIndicesMother = {start, stop}; }
Expand Down Expand Up @@ -105,9 +104,9 @@ class OTFParticle
int pdgCode() const { return mPdgCode; }
int globalIndex() const { return mGlobalIndex; }
int collisionId() const { return mCollisionId; }
bool isAlive() const { return mIsAlive; }
bool isPrimary() const { return mIsPrimary; }
bool isFromMcParticles() const { return mIsFromMcParticles; }
bool isAlive() const { return checkBit(DecayerBits::IsAlive); }
bool isPrimary() const { return checkBit(DecayerBits::IsPrimary); }
bool isFromMcParticles() const { return !checkBit(DecayerBits::ProducedByDecayer); }
float weight() const
{
static constexpr float Weight = 1.f;
Expand Down Expand Up @@ -154,8 +153,8 @@ class OTFParticle
int getMotherIndexStop() const { return mIndicesMother[1]; }
int getDaughterIndexStart() const { return mIndicesDaughter[0]; }
int getDaughterIndexStop() const { return mIndicesDaughter[1]; }
std::array<int, 2> getMothers() const { return mIndicesMother; }
std::array<int, 2> getDaughters() const { return mIndicesDaughter; }
const std::array<int, 2>& getMothers() const { return mIndicesMother; }
const std::array<int, 2>& getDaughters() const { return mIndicesDaughter; }
std::span<const int> getMotherSpan() const { return hasMothers() ? std::span<const int>(mIndicesMother.data(), 2) : std::span<const int>(); }

// Checks
Expand All @@ -177,7 +176,7 @@ class OTFParticle
void setBitOn(DecayerBits bit) { mBits.set(static_cast<size_t>(bit), true); }
void setBitOff(DecayerBits bit) { mBits.set(static_cast<size_t>(bit), false); }

std::bitset<8> getBits() const { return mBits; }
const std::bitset<8>& getBits() const { return mBits; }
uint8_t getBitsValue() const { return static_cast<uint8_t>(mBits.to_ulong()); }
void setBits(std::bitset<8> bits) { mBits = bits; }

Expand All @@ -186,8 +185,6 @@ class OTFParticle
int mCollisionId{-1};
float mVx{}, mVy{}, mVz{}, mVt{};
float mPx{}, mPy{}, mPz{}, mE{};
bool mIsAlive{}, mIsFromMcParticles{false};
bool mIsPrimary{};

int mStatusCode{};
uint8_t mFlag{};
Expand Down
1 change: 1 addition & 0 deletions ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -161,9 +161,10 @@

const float trackTimeNS = trackLength / trackVelocity * PicoToNano;
particle.setIndicesDaughter(particlesInDataframe - indexOffset + allParticles.size(), particlesInDataframe - indexOffset + allParticles.size() + (decayStack.size() - 1));
for (auto& daughter : decayStack) {

Check failure on line 164 in ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[const-ref-in-for-loop]

Use constant references for non-modified iterators in range-based for loops.
daughter.setIndicesMother(particlesInDataframe - indexOffset + i, particlesInDataframe - indexOffset + i);
daughter.setCollisionId(particle.collisionId());
daughter.setBitOn(o2::upgrade::DecayerBits::ProducedByDecayer);
daughter.setBitOn(o2::upgrade::DecayerBits::IsAlive);
daughter.setBitOff(o2::upgrade::DecayerBits::IsPrimary);
daughter.setProductionTime(particle.vt() + trackTimeNS);
Expand Down Expand Up @@ -198,7 +199,7 @@
decayParticles(0, allParticles.size());

// Fill output table
for (auto& otfParticle : allParticles) {

Check failure on line 202 in ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[const-ref-in-for-loop]

Use constant references for non-modified iterators in range-based for loops.
otfParticle.setIndexOffset(indexOffset);
if (otfParticle.hasNaN()) {
histos.fill(HIST("hNaNBookkeeping"), 1);
Expand Down
64 changes: 54 additions & 10 deletions ALICE3/TableProducer/OTF/onTheFlyTracker.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@
#include <TGeoGlobalMagField.h>
#include <TH1.h>
#include <TH2.h>
#include <TLorentzVector.h>

Check failure on line 80 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
#include <TMCProcess.h>
#include <TMath.h>
#include <TPDGCode.h>
Expand Down Expand Up @@ -587,6 +587,19 @@
insertHist(histPath + "h2dBRPtResAbs", "h2dPtResAbs;Gen p_{T};#Delta p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
}

if (doprocessDecayer) {
insertHist(histPath + "h2dPrimaryPtRes", "h2dPrimaryPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dSecondaryPtRes", "h2dSecondaryPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dPrimaryElPtRes", "h2dPrimaryElPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dSecondaryElPtRes", "h2dSecondaryElPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dPrimaryPiPtRes", "h2dPrimaryPiPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dSecondaryPiPtRes", "h2dSecondaryPiPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dPrimaryKaPtRes", "h2dPrimaryKaPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dSecondaryKaPtRes", "h2dSecondaryKaPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dPrimaryPrPtRes", "h2dPrimaryPrPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
insertHist(histPath + "h2dSecondaryPrPtRes", "h2dSecondaryPrPtRes;Gen p_{T};#Delta p_{T} / Reco p_{T}", {kTH2D, {{axes.axisMomentum, axes.axisPtRes}}});
}

} // end config loop

// Basic QA
Expand Down Expand Up @@ -763,7 +776,7 @@
/// \param xiDecayVertex the address of the xi decay vertex
/// \param laDecayVertex the address of the la decay vertex
template <typename McParticleType>
void decayCascade(McParticleType particle, o2::track::TrackParCov track, std::vector<TLorentzVector>& decayDaughters, std::vector<double>& xiDecayVertex, std::vector<double>& laDecayVertex)

Check failure on line 779 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
{
const double uXi = rand.Uniform(0, 1);
const double ctauXi = 4.91; // cm
Expand All @@ -786,12 +799,12 @@
xiDecayVertex.push_back(particle.vz() + rxyzXi * (particle.pz() / particle.p()));

std::vector<double> xiDaughters = {o2::constants::physics::MassLambda, o2::constants::physics::MassPionCharged};
TLorentzVector xi(newPx, newPy, particle.pz(), newE);

Check failure on line 802 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
TGenPhaseSpace xiDecay;
xiDecay.SetDecay(xi, 2, xiDaughters.data());
xiDecay.Generate();
decayDaughters.push_back(*xiDecay.GetDecay(1));
TLorentzVector la = *xiDecay.GetDecay(0);

Check failure on line 807 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.

const double uLa = rand.Uniform(0, 1);
const double ctauLa = 7.845; // cm
Expand All @@ -814,7 +827,7 @@
/// \param decayDaughters the address of resulting daughters
/// \param v0DecayVertex the address of the la decay vertex
template <typename McParticleType>
void decayV0Particle(McParticleType particle, std::vector<TLorentzVector>& decayDaughters, std::vector<double>& v0DecayVertex, int pdgCode)

Check failure on line 830 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
{
double u = rand.Uniform(0, 1);
double v0Mass = -1.;
Expand Down Expand Up @@ -848,7 +861,7 @@

const double v0BetaGamma = particle.p() / v0Mass;
const double v0rxyz = (-v0BetaGamma * ctau * std::log(1 - u));
TLorentzVector v0(particle.px(), particle.py(), particle.pz(), particle.e());

Check failure on line 864 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.

v0DecayVertex.push_back(particle.vx() + v0rxyz * (particle.px() / particle.p()));
v0DecayVertex.push_back(particle.vy() + v0rxyz * (particle.py() / particle.p()));
Expand Down Expand Up @@ -923,7 +936,7 @@
o2::upgrade::convertMCParticleToO2Track(mcParticle, trackParCov, pdgDB);
const std::string histPath = "Configuration_" + std::to_string(icfg) + "/";

std::vector<TLorentzVector> cascadeDecayProducts;

Check failure on line 939 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
std::vector<double> xiDecayVertex, laDecayVertex;
static constexpr int kCascProngs = 3;
std::array<o2::track::TrackParCov, kCascProngs> xiDaughterTrackParCovsPerfect;
Expand Down Expand Up @@ -953,7 +966,7 @@
getHist(TH1, histPath + "hXiBuilding")->Fill(0.0f);
}

o2::upgrade::convertTLorentzVectorToO2Track(PDG_t::kPiMinus, cascadeDecayProducts[0], xiDecayVertex, xiDaughterTrackParCovsPerfect[0], pdgDB);

Check failure on line 969 in ALICE3/TableProducer/OTF/onTheFlyTracker.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[root/lorentz-vector]

Do not use the TLorentzVector legacy class. Use std::array with RecoDecay methods or the ROOT::Math::LorentzVector template instead.
xiDaughterTrackParCovsPerfect[0].setPID(pdgCodeToPID(PDG_t::kPiMinus));
o2::upgrade::convertTLorentzVectorToO2Track(PDG_t::kPiMinus, cascadeDecayProducts[1], laDecayVertex, xiDaughterTrackParCovsPerfect[1], pdgDB);
xiDaughterTrackParCovsPerfect[1].setPID(pdgCodeToPID(PDG_t::kPiMinus));
Expand Down Expand Up @@ -2043,6 +2056,7 @@
const bool pdgsToBeHandled = longLivedToBeHandled || (enableNucleiSmearing && nucleiToBeHandled);

o2::upgrade::OTFParticle otfParticle(mcParticle);
otfParticle.setBits(mcParticle.decayerBits_raw());
if (otfParticle.hasNaN()) {
continue;
}
Expand Down Expand Up @@ -2080,17 +2094,16 @@

multiplicityCounter++;
o2::track::TrackParCov trackParCov;
const bool isDecayDaughter = (mcParticle.getProcess() == TMCProcess::kPDecay);
const float time = (eventCollisionTimeNS + gRandom->Gaus(0., timeResolutionNs)) * nsToMus;

bool reconstructed = false;
int nTrkHits = 0;
if (enablePrimarySmearing && otfParticle.checkBit(o2::upgrade::DecayerBits::IsPrimary)) {
const bool isSecondary = !otfParticle.isPrimary() && otfParticle.checkBit(o2::upgrade::DecayerBits::ProducedByDecayer) && otfParticle.isAlive();
if (enablePrimarySmearing && otfParticle.isPrimary()) {
o2::upgrade::convertMCParticleToO2Track(mcParticle, trackParCov, pdgDB);
computeBremsstrahlungLoss(icfg, mcParticle, trackParCov);
reconstructed = mSmearer[icfg]->smearTrack(trackParCov, mcParticle.pdgCode(), dNdEta);
nTrkHits = fastTrackerSettings.minSiliconHits;
} else if (enableSecondarySmearing && !otfParticle.checkBit(o2::upgrade::DecayerBits::IsPrimary) && otfParticle.checkBit(o2::upgrade::DecayerBits::ProducedByDecayer) && otfParticle.checkBit(o2::upgrade::DecayerBits::IsAlive)) {
} else if (enableSecondarySmearing && isSecondary) {
o2::track::TrackParCov perfectTrackParCov;
o2::upgrade::convertMCParticleToO2Track(mcParticle, perfectTrackParCov, pdgDB);
computeBremsstrahlungLoss(icfg, mcParticle, perfectTrackParCov);
Expand All @@ -2114,21 +2127,52 @@

histos.fill(HIST("hNaNBookkeeping"), 0.0f, 1.0f);
if (enablePrimarySmearing) {
const float ptResolution = (trackParCov.getPt() - mcParticle.pt()) / trackParCov.getPt();
getHist(TH1, histPath + "hPtReconstructed")->Fill(trackParCov.getPt());
if (std::abs(mcParticle.pdgCode()) == kElectron)
if (otfParticle.isPrimary()) {
getHist(TH2, histPath + "h2dPrimaryPtRes")->Fill(trackParCov.getPt(), ptResolution);
} else {
getHist(TH2, histPath + "h2dSecondaryPtRes")->Fill(trackParCov.getPt(), ptResolution);
}
if (std::abs(mcParticle.pdgCode()) == kElectron) {
getHist(TH1, histPath + "hPtReconstructedEl")->Fill(trackParCov.getPt());
if (std::abs(mcParticle.pdgCode()) == kPiPlus)
if (otfParticle.isPrimary()) {
getHist(TH2, histPath + "h2dPrimaryElPtRes")->Fill(trackParCov.getPt(), ptResolution);
} else {
getHist(TH2, histPath + "h2dSecondaryElPtRes")->Fill(trackParCov.getPt(), ptResolution);
}
}
if (std::abs(mcParticle.pdgCode()) == kPiPlus) {
getHist(TH1, histPath + "hPtReconstructedPi")->Fill(trackParCov.getPt());
if (std::abs(mcParticle.pdgCode()) == kKPlus)
if (otfParticle.isPrimary()) {
getHist(TH2, histPath + "h2dPrimaryPiPtRes")->Fill(trackParCov.getPt(), ptResolution);
} else {
getHist(TH2, histPath + "h2dSecondaryPiPtRes")->Fill(trackParCov.getPt(), ptResolution);
}
}
if (std::abs(mcParticle.pdgCode()) == kKPlus) {
getHist(TH1, histPath + "hPtReconstructedKa")->Fill(trackParCov.getPt());
if (std::abs(mcParticle.pdgCode()) == kProton)
if (otfParticle.isPrimary()) {
getHist(TH2, histPath + "h2dPrimaryKaPtRes")->Fill(trackParCov.getPt(), ptResolution);
} else {
getHist(TH2, histPath + "h2dSecondaryKaPtRes")->Fill(trackParCov.getPt(), ptResolution);
}
}
if (std::abs(mcParticle.pdgCode()) == kProton) {
getHist(TH1, histPath + "hPtReconstructedPr")->Fill(trackParCov.getPt());
if (otfParticle.isPrimary()) {
getHist(TH2, histPath + "h2dPrimaryPrPtRes")->Fill(trackParCov.getPt(), ptResolution);
} else {
getHist(TH2, histPath + "h2dSecondaryPrPtRes")->Fill(trackParCov.getPt(), ptResolution);
}
}
}

if (reconstructed) {
tracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isDecayDaughter, false, 0, nTrkHits, kRecoPrimary});
tracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isSecondary, false, 0, nTrkHits, kRecoPrimary});
getHist(TH1, histPath + "hPtReconstructedPr")->Fill(trackParCov.getPt());
} else {
ghostTracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isDecayDaughter, false, 0, nTrkHits, kGhostPrimary});
ghostTracksAlice3.push_back(TrackAlice3{trackParCov, mcParticle.globalIndex(), time, timeResolutionUs, isSecondary, false, 0, nTrkHits, kGhostPrimary});
}
}

Expand Down
Loading