diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidTableProducer.cxx b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidTableProducer.cxx index 7e4cb37c9ef..92e25c9fb32 100644 --- a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidTableProducer.cxx +++ b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidTableProducer.cxx @@ -24,13 +24,8 @@ #include #include #include -#include #include -#include #include -#include -#include -#include #include #include @@ -51,6 +46,7 @@ using namespace o2::constants::physics; struct HmpidTableProducer { Produces hmpidAnalysis; + Produces hmpidAnalysisMC; HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; @@ -80,6 +76,16 @@ struct HmpidTableProducer { aod::pidTOFFullPi, aod::pidTOFFullKa, aod::pidTOFFullPr, aod::pidTOFFullDe>; + using TrackCandidatesMC = soa::Join; + + std::unordered_set mCollisionsWithHmpid; + void init(o2::framework::InitContext&) { ccdb->setURL(ccdbConfig.ccdbUrl); @@ -294,42 +300,37 @@ struct HmpidTableProducer { } PROCESS_SWITCH(HmpidTableProducer, processEvent, "Process event level", true); - void processHmpid( + template + void runHmpidAnalysis( aod::HMPIDs const& hmpids, - TrackCandidates const&, + TTrackTable const&, CollisionCandidates const&, aod::BCsWithTimestamps const&) { - static std::unordered_set collisionsWithHmpid; - for (auto const& t : hmpids) { - const auto& globalTrack = t.track_as(); + const auto& globalTrack = t.template track_as(); if (!globalTrack.has_collision()) continue; - const auto& col = globalTrack.collision_as(); - initCCDB(col.bc_as()); + const auto& col = globalTrack.template collision_as(); + initCCDB(col.template bc_as()); uint32_t collId = col.globalIndex(); - // Track quality selection if ((requireITS && !globalTrack.hasITS()) || (requireTPC && !globalTrack.hasTPC()) || - (requireTOF && !globalTrack.hasTOF())) { + (requireTOF && !globalTrack.hasTOF())) continue; - } - if (collisionsWithHmpid.insert(collId).second) { + if (mCollisionsWithHmpid.insert(collId).second) histos.fill(HIST("eventsHmpid"), 0.5); - } // clusSize diagnostics histos.fill(HIST("hClusSize"), t.hmpidClusSize()); bool isCorrupt = (t.hmpidClusSize() <= 0); - if (isCorrupt) { + if (isCorrupt) histos.fill(HIST("hClusSizeCorrupt"), t.hmpidClusSize()); - } // --- M2: clusSize encoding --- int chamberM2 = t.hmpidClusSize() / 1000000; @@ -354,7 +355,7 @@ struct HmpidTableProducer { int16_t charge = globalTrack.sign(); auto prop = o2::base::Propagator::Instance(); - double bz = static_cast(prop->getNominalBz()); // positive sign + double bz = static_cast(prop->getNominalBz()); int chamberM1 = getHmpidChamber(x, p, bz, charge); @@ -374,35 +375,30 @@ struct HmpidTableProducer { // Legend: // bin 0 = clusSize > 0, chamber found (M2 ok) - // bin 1 = clusSize > 0, chamber not found (non dovrebbe mai accadere) + // bin 1 = clusSize > 0, chamber not found // bin 2 = clusSize <= 0, M1 recovery (corrupt, M1 ok) // bin 3 = clusSize <= 0, M1 fails (corrupt, skipped) - if (!isCorrupt && chamberM3 >= 0) { + if (!isCorrupt && chamberM3 >= 0) histos.fill(HIST("hChamberAssignment"), 0.); - } else if (!isCorrupt && chamberM3 < 0) { + else if (!isCorrupt && chamberM3 < 0) histos.fill(HIST("hChamberAssignment"), 1.); - } else if (isCorrupt && chamberM3 >= 0) { + else if (isCorrupt && chamberM3 >= 0) histos.fill(HIST("hChamberAssignment"), 2.); - } else { + else histos.fill(HIST("hChamberAssignment"), 3.); - } histos.fill(HIST("hChamberM3"), chamberM3); histos.fill(HIST("hChamberM3vsM2"), chamberM2, chamberM3); - // Skip track if chamber undetermined - if (chamberM3 < 0) { + if (chamberM3 < 0) continue; - } - // Fill photon charges float hmpidPhotsCharge2[o2::aod::kDimPhotonsCharge]; - for (int i = 0; i < o2::aod::kDimPhotonsCharge; i++) { + for (int i = 0; i < o2::aod::kDimPhotonsCharge; i++) hmpidPhotsCharge2[i] = t.hmpidPhotsCharge()[i]; - } - // Fill output table + // fill hmpid table hmpidAnalysis( t.hmpidSignal(), t.hmpidMom(), globalTrack.p(), t.hmpidXTrack(), t.hmpidYTrack(), @@ -422,9 +418,40 @@ struct HmpidTableProducer { globalTrack.tpcNSigmaDe(), globalTrack.tofNSigmaDe(), col.centFV0A()); + // fill hmpid table for mc tracks if running on MC + if constexpr (isMC) { + if (globalTrack.has_mcParticle()) { + const auto& mc = globalTrack.mcParticle(); + hmpidAnalysisMC(mc.pdgCode(), mc.vx(), mc.vy(), mc.vz(), mc.isPhysicalPrimary(), mc.getProcess()); + } else { + hmpidAnalysisMC(-1, 0.f, 0.f, 0.f, false, -100); + } + } + } // end HMPID loop } + + // process real data + void processHmpid(aod::HMPIDs const& hmpids, + TrackCandidates const& tracks, + CollisionCandidates const& cols, + aod::BCsWithTimestamps const& bcs) + { + // isMC False + runHmpidAnalysis(hmpids, tracks, cols, bcs); + } PROCESS_SWITCH(HmpidTableProducer, processHmpid, "Process HMPID entries", true); + + // process MC + void processHmpidMC(aod::HMPIDs const& hmpids, + TrackCandidatesMC const& tracks, + CollisionCandidates const& cols, + aod::BCsWithTimestamps const& bcs, + aod::McParticles const&) + { + runHmpidAnalysis(hmpids, tracks, cols, bcs); + } + PROCESS_SWITCH(HmpidTableProducer, processHmpidMC, "Process HMPID MC entries", true); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfg) diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/tableHMPID.h b/DPG/Tasks/AOTTrack/PID/HMPID/tableHMPID.h index 6844da6c519..1fcf9980763 100644 --- a/DPG/Tasks/AOTTrack/PID/HMPID/tableHMPID.h +++ b/DPG/Tasks/AOTTrack/PID/HMPID/tableHMPID.h @@ -34,22 +34,18 @@ DECLARE_SOA_COLUMN(ChargeMip, chargeMip, float); DECLARE_SOA_COLUMN(ClusterSize, clusterSize, float); DECLARE_SOA_COLUMN(Chamber, chamber, float); DECLARE_SOA_COLUMN(PhotonsCharge, photonsCharge, float[kDimPhotonsCharge]); - DECLARE_SOA_COLUMN(EtaTrack, etaTrack, float); DECLARE_SOA_COLUMN(PhiTrack, phiTrack, float); DECLARE_SOA_COLUMN(Px, px, float); DECLARE_SOA_COLUMN(Py, py, float); DECLARE_SOA_COLUMN(Pz, pz, float); - DECLARE_SOA_COLUMN(ItsNCluster, itsNCluster, float); DECLARE_SOA_COLUMN(TpcNCluster, tpcNCluster, float); DECLARE_SOA_COLUMN(TpcNClsCrossedRows, tpcNClsCrossedRows, float); DECLARE_SOA_COLUMN(TpcChi2, tpcChi2, float); DECLARE_SOA_COLUMN(ItsChi2, itsChi2, float); - DECLARE_SOA_COLUMN(DcaXY, dcaXY, float); DECLARE_SOA_COLUMN(DcaZ, dcaZ, float); - DECLARE_SOA_COLUMN(TpcNSigmaPi, tpcNSigmaPi, float); DECLARE_SOA_COLUMN(TofNSigmaPi, tofNSigmaPi, float); DECLARE_SOA_COLUMN(TpcNSigmaKa, tpcNSigmaKa, float); @@ -59,7 +55,6 @@ DECLARE_SOA_COLUMN(TofNSigmaPr, tofNSigmaPr, float); DECLARE_SOA_COLUMN(TpcNSigmaDe, tpcNSigmaDe, float); DECLARE_SOA_COLUMN(TofNSigmaDe, tofNSigmaDe, float); DECLARE_SOA_COLUMN(Centrality, centrality, float); - } // namespace variables_table DECLARE_SOA_TABLE(HmpidAnalysis, "AOD", "HMPIDANALYSIS", @@ -97,6 +92,28 @@ DECLARE_SOA_TABLE(HmpidAnalysis, "AOD", "HMPIDANALYSIS", variables_table::TofNSigmaDe, variables_table::Centrality); +// ----------------------------------------------------------------------- +// MC truth table +// ----------------------------------------------------------------------- +namespace hmpid_mc +{ +DECLARE_SOA_COLUMN(PdgCode, pdgCode, int); +DECLARE_SOA_COLUMN(McVx, mcVx, float); +DECLARE_SOA_COLUMN(McVy, mcVy, float); +DECLARE_SOA_COLUMN(McVz, mcVz, float); +DECLARE_SOA_COLUMN(IsPhysPrimary, isPhysPrimary, bool); +DECLARE_SOA_COLUMN(ProcessCode, processCode, int); + +} // namespace hmpid_mc + +DECLARE_SOA_TABLE(HmpidAnalysisMC, "AOD", "HMPIDANALYSISMC", + hmpid_mc::PdgCode, + hmpid_mc::McVx, + hmpid_mc::McVy, + hmpid_mc::McVz, + hmpid_mc::IsPhysPrimary, + hmpid_mc::ProcessCode); + } // namespace o2::aod #endif // DPG_TASKS_AOTTRACK_PID_HMPID_TABLEHMPID_H_