diff --git a/PWGLF/CMakeLists.txt b/PWGLF/CMakeLists.txt index 9d96b0b8d33..797ac01500f 100644 --- a/PWGLF/CMakeLists.txt +++ b/PWGLF/CMakeLists.txt @@ -13,4 +13,4 @@ # add_subdirectory(DataModel) add_subdirectory(Tasks) add_subdirectory(TableProducer) - +add_subdirectory(Utils) diff --git a/PWGLF/DataModel/LFStrangenessTables.h b/PWGLF/DataModel/LFStrangenessTables.h index d7e0628d274..efb3bbb51d1 100644 --- a/PWGLF/DataModel/LFStrangenessTables.h +++ b/PWGLF/DataModel/LFStrangenessTables.h @@ -132,6 +132,10 @@ DECLARE_SOA_DYNAMIC_COLUMN(HasTRD, hasTRD, //! Flag to check if track has a TRD [](uint8_t detectorMap) -> bool { return detectorMap & o2::aod::track::TRD; }); DECLARE_SOA_DYNAMIC_COLUMN(HasTOF, hasTOF, //! Flag to check if track has a TOF measurement [](uint8_t detectorMap) -> bool { return detectorMap & o2::aod::track::TOF; }); +DECLARE_SOA_DYNAMIC_COLUMN(HasITSTracker, hasITSTracker, //! Flag to check if track is from ITS tracker + [](uint8_t detectorMap, float itsChi2PerNcl) -> bool { return (detectorMap & o2::aod::track::ITS) ? (itsChi2PerNcl > -1e-3f) : false; }); +DECLARE_SOA_DYNAMIC_COLUMN(HasITSAfterburner, hasITSAfterburner, //! Flag to check if track is from ITS AB + [](uint8_t detectorMap, float itsChi2PerNcl) -> bool { return (detectorMap & o2::aod::track::ITS) ? (itsChi2PerNcl < -1e-3f) : false; }); } // namespace dautrack DECLARE_SOA_TABLE(DauTrackExtras_000, "AOD", "DAUTRACKEXTRA", //! detector properties of decay daughters @@ -144,7 +148,9 @@ DECLARE_SOA_TABLE(DauTrackExtras_000, "AOD", "DAUTRACKEXTRA", //! detector prope dautrack::HasITS, dautrack::HasTPC, dautrack::HasTRD, - dautrack::HasTOF); + dautrack::HasTOF, + dautrack::HasITSTracker, + dautrack::HasITSAfterburner); DECLARE_SOA_TABLE_VERSIONED(DauTrackExtras_001, "AOD", "DAUTRACKEXTRA", 1, //! detector properties of decay daughters dautrack::ITSChi2PerNcl, @@ -157,7 +163,9 @@ DECLARE_SOA_TABLE_VERSIONED(DauTrackExtras_001, "AOD", "DAUTRACKEXTRA", 1, //! d dautrack::HasITS, dautrack::HasTPC, dautrack::HasTRD, - dautrack::HasTOF); + dautrack::HasTOF, + dautrack::HasITSTracker, + dautrack::HasITSAfterburner); DECLARE_SOA_TABLE(DauTrackMCIds, "AOD", "DAUTRACKMCID", // index table when using AO2Ds dautrack::ParticleMCId); diff --git a/PWGLF/TableProducer/Strangeness/lambdakzeromcfinder.cxx b/PWGLF/TableProducer/Strangeness/lambdakzeromcfinder.cxx index 5e0c7494282..c08c7575e44 100644 --- a/PWGLF/TableProducer/Strangeness/lambdakzeromcfinder.cxx +++ b/PWGLF/TableProducer/Strangeness/lambdakzeromcfinder.cxx @@ -183,9 +183,11 @@ struct lambdakzeromcfinder { int positivePdg = 211; int negativePdg = -211; + int relevantProcess = 4; // normal search: decay if (mcParticle.pdgCode() == 22) { positivePdg = -11; negativePdg = +11; + relevantProcess = 5; // look for pair production if photon } if (mcParticle.pdgCode() == 3122) { positivePdg = 2212; @@ -204,7 +206,7 @@ struct lambdakzeromcfinder { auto const& daughters = mcParticle.template daughters_as(); if (daughters.size() >= 2) { for (auto const& daughter : daughters) { // might be better ways of doing this but ok - if (daughter.getProcess() != 4) + if (daughter.getProcess() != relevantProcess) continue; // skip deltarays (if ever), stick to decay products only if (daughter.pdgCode() == positivePdg) { auto const& thisDaughterTracks = daughter.template tracks_as(); diff --git a/PWGLF/TableProducer/converters/CMakeLists.txt b/PWGLF/TableProducer/converters/CMakeLists.txt index 2abb2f6fdbc..7080dc720dd 100644 --- a/PWGLF/TableProducer/converters/CMakeLists.txt +++ b/PWGLF/TableProducer/converters/CMakeLists.txt @@ -14,11 +14,6 @@ o2physics_add_dpl_workflow(stradautrackstofpidconverter PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore COMPONENT_NAME Analysis) -o2physics_add_dpl_workflow(stradautracksextraconverter - SOURCES stradautracksextraconverter.cxx - PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore - COMPONENT_NAME Analysis) - o2physics_add_dpl_workflow(strarawcentsconverter SOURCES strarawcentsconverter.cxx PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore diff --git a/PWGLF/Tasks/QC/CMakeLists.txt b/PWGLF/Tasks/QC/CMakeLists.txt index f2e9bf5646c..5ebd25f4fd5 100755 --- a/PWGLF/Tasks/QC/CMakeLists.txt +++ b/PWGLF/Tasks/QC/CMakeLists.txt @@ -96,7 +96,7 @@ o2physics_add_dpl_workflow(kfstrangenessstudy o2physics_add_dpl_workflow(findable-study SOURCES findableStudy.cxx - PUBLIC_LINK_LIBRARIES O2::DetectorsBase O2Physics::AnalysisCore + PUBLIC_LINK_LIBRARIES O2::DetectorsBase O2Physics::AnalysisCore O2Physics::v0SelectionGroup COMPONENT_NAME Analysis) o2physics_add_dpl_workflow(mcinelgt0 diff --git a/PWGLF/Tasks/QC/findableStudy.cxx b/PWGLF/Tasks/QC/findableStudy.cxx index 4f013ba4857..03d496ef0d1 100644 --- a/PWGLF/Tasks/QC/findableStudy.cxx +++ b/PWGLF/Tasks/QC/findableStudy.cxx @@ -45,6 +45,9 @@ #include "DataFormatsParameters/GRPObject.h" #include "DataFormatsParameters/GRPMagField.h" #include "CCDB/BasicCCDBManager.h" +#include "PWGLF/Utils/v0SelectionBits.h" +#include "PWGLF/Utils/v0SelectionGroup.h" +#include "PWGLF/Utils/v0SelectionTools.h" #include #include @@ -59,55 +62,109 @@ using namespace o2::framework::expressions; using std::array; using recoStraCollisions = soa::Join; -using reconstructedV0s = soa::Join; +using reconstructedV0s = soa::Join; +using reconstructedV0sNoMC = soa::Join; + using dauTracks = soa::Join; -// simple helper -#define bitset(var, nbit) ((var) |= (static_cast(1) << static_cast(nbit))) +// simple checkers, but ensure 64 bit integers +#define bitset(var, nbit) ((var) |= (static_cast(1) << static_cast(nbit))) +#define bitcheck(var, nbit) ((var) & (static_cast(1) << static_cast(nbit))) struct findableStudy { HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject}; // master PDG code selection Configurable pdgCode{"pdgCode", 310, "PDG code to select"}; + Configurable skipITSonly{"skipITSonly", true, "skip reco V0s if an ITS-only (no TPC) prong present"}; ConfigurableAxis axisPt{"axisPt", {VARIABLE_WIDTH, 0.0f, 0.1f, 0.2f, 0.3f, 0.4f, 0.5f, 0.6f, 0.7f, 0.8f, 0.9f, 1.0f, 1.1f, 1.2f, 1.3f, 1.4f, 1.5f, 1.6f, 1.7f, 1.8f, 1.9f, 2.0f, 2.2f, 2.4f, 2.6f, 2.8f, 3.0f, 3.2f, 3.4f, 3.6f, 3.8f, 4.0f, 4.4f, 4.8f, 5.2f, 5.6f, 6.0f, 6.5f, 7.0f, 7.5f, 8.0f, 9.0f, 10.0f, 11.0f, 12.0f, 13.0f, 14.0f, 15.0f, 17.0f, 19.0f, 21.0f, 23.0f, 25.0f, 30.0f, 35.0f, 40.0f, 50.0f}, "pt axis for analysis"}; ConfigurableAxis axisCentrality{"axisCentrality", {VARIABLE_WIDTH, 0.0f, 5.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f}, "Centrality"}; + // +-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+ + // Full wrapper for configurables related to actual analysis + Configurable v0Selections{"v0Selections", {}, "V0 selection criteria for analysis"}; + // +-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+-~-+ + + // pack track quality but separte also afterburner + // dynamic range: 0-31 + enum selection : int { hasTPC = 0, + hasITSTracker, + hasITSAfterburner, + hasTRD, + hasTOF }; + + uint64_t maskTopological; + uint64_t maskTrackProperties; + + uint64_t maskK0ShortSpecific; + uint64_t maskLambdaSpecific; + uint64_t maskAntiLambdaSpecific; + + uint64_t maskSelectionK0Short; + uint64_t maskSelectionLambda; + uint64_t maskSelectionAntiLambda; + void init(InitContext const&) { + v0Selections->PrintSelections(); // for the logs + + v0Selections->provideMasks(maskTopological, maskTrackProperties, maskK0ShortSpecific, maskLambdaSpecific, maskAntiLambdaSpecific); + + // Primary particle selection, central to analysis + maskSelectionK0Short = maskTopological | maskTrackProperties | maskK0ShortSpecific; + maskSelectionLambda = maskTopological | maskTrackProperties | maskLambdaSpecific; + maskSelectionAntiLambda = maskTopological | maskTrackProperties | maskAntiLambdaSpecific; + // Event counting histos.add("hCentrality", "hCentrality", kTH1D, {axisCentrality}); // Duplicate counting histos.add("hNRecoV0s", "hNRecoV0s", kTH1D, {{50, -0.5, 49.5f}}); + histos.add("hNRecoV0sWithTPC", "hNRecoV0sWithTPC", kTH1D, {{50, -0.5, 49.5f}}); histos.add("hNRecoV0sWrongColl", "hNRecoV0sWrongColl", kTH1D, {{50, -0.5, 49.5f}}); - // Findable versus found - histos.add("h2dPtVsCentrality_All_Findable", "hPtVsCentrality_All_Findable", kTH2D, {axisCentrality, axisPt}); - histos.add("h2dPtVsCentrality_All_Found", "hPtVsCentrality_All_Found", kTH2D, {axisCentrality, axisPt}); - - // binary feature study: encode in Y axis - // 0 - positive and negative don't have feature - // 1 - positive has, negative doesn't - // 2 - negative has, positive doesn't - // 3 - both positive and negative have feature - // - // all binary feature presence tests are done such that the "has feature" condition - // is evaluated such that out of any reco V0s, at least one V0 has the desired feature. - // this is required for generality - - const AxisSpec axisBinaryFeature{static_cast(4), -0.5f, +3.5f, ""}; - histos.add("h2dPtVsCentrality_WithTPC_Findable", "hPtVsCentrality_WithTPC_Findable", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithTPC_Found", "hPtVsCentrality_WithTPC_Found", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSTracker_Findable", "hPtVsCentrality_WithITSTracker_Findable", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSTracker_Found", "hPtVsCentrality_WithITSTracker_Found", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSTrackerTPC_Findable", "hPtVsCentrality_WithITSTrackerTPC_Findable", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSTrackerTPC_Found", "hPtVsCentrality_WithITSTrackerTPC_Found", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSABTPC_Findable", "hPtVsCentrality_WithITSABTPC_Findable", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithITSABTPC_Found", "hPtVsCentrality_WithITSABTPC_Found", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithSVertexerOK_Findable", "hPtVsCentrality_WithSVertexerOK_Findable", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); - histos.add("h2dPtVsCentrality_WithSVertexerOK_Found", "hPtVsCentrality_WithSVertexerOK_Found", kTH3D, {axisCentrality, axisPt, axisBinaryFeature}); + // For broad correctness check + histos.add("hFoundVsTracksOK", "hFoundVsTracksOK", kTH2D, {{2, -0.5, 1.5f}, {2, -0.5, 1.5f}}); + + // Global findable + histos.add("h2dPtVsCentrality_Findable", "hPtVsCentrality_Findable", kTH2D, {axisCentrality, axisPt}); + + // Acceptably (for svertexer) tracked + histos.add("h2dPtVsCentrality_AcceptablyTracked", "h2dPtVsCentrality_AcceptablyTracked", kTH2D, {axisCentrality, axisPt}); + + // Found in any capacity, including ITSonly + histos.add("h2dPtVsCentrality_FoundAny", "h2dPtVsCentrality_FoundAny", kTH2D, {axisCentrality, axisPt}); + + // Found with prongs with the TPC (typical analysis setting) + histos.add("h2dPtVsCentrality_Found", "h2dPtVsCentrality_Found", kTH2D, {axisCentrality, axisPt}); + + // Found in loop (may have duplicates, meant as cross-check too) + histos.add("h2dPtVsCentrality_FoundInLoop", "h2dPtVsCentrality_FoundInLoop", kTH2D, {axisCentrality, axisPt}); + + // Passes analysis-level track quality checks + histos.add("h2dPtVsCentrality_PassesTrackQuality", "h2dPtVsCentrality_Analysis_PassesTrackQuality", kTH2D, {axisCentrality, axisPt}); + + // Passes analysis-level topological selection criteria + histos.add("h2dPtVsCentrality_PassesTopological", "h2dPtVsCentrality_PassesTopological", kTH2D, {axisCentrality, axisPt}); + + // Passes analysis-level species-specific + histos.add("h2dPtVsCentrality_PassesThisSpecies", "h2dPtVsCentrality_PassesThisSpecies", kTH2D, {axisCentrality, axisPt}); + + // one-on-one test for topology, operating with track quality and species-specific checks ON + histos.add("h2dPtVsCentrality_V0Radius", "h2dPtVsCentrality_V0Radius", kTH2D, {axisCentrality, axisPt}); + histos.add("h2dPtVsCentrality_V0RadiusMax", "h2dPtVsCentrality_V0RadiusMax", kTH2D, {axisCentrality, axisPt}); + histos.add("h2dPtVsCentrality_V0CosPA", "h2dPtVsCentrality_V0CosPA", kTH2D, {axisCentrality, axisPt}); + histos.add("h2dPtVsCentrality_DcaPosToPV", "h2dPtVsCentrality_DcaPosToPV", kTH2D, {axisCentrality, axisPt}); + histos.add("h2dPtVsCentrality_DcaNegToPV", "h2dPtVsCentrality_DcaNegToPV", kTH2D, {axisCentrality, axisPt}); + histos.add("h2dPtVsCentrality_DcaV0Dau", "h2dPtVsCentrality_DcaV0Dau", kTH2D, {axisCentrality, axisPt}); + + // Track quality tests in steps + histos.add("h2dTrackPropAcceptablyTracked", "h2dTrackPropAcceptablyTracked", kTH3D, {{32, -0.5, 31.5f}, {32, -0.5, 31.5f}, axisCentrality}); + histos.add("h2dTrackPropFound", "h2dTrackPropFound", kTH3D, {{32, -0.5, 31.5f}, {32, -0.5, 31.5f}, axisCentrality}); + histos.add("h2dTrackPropAnalysisTracks", "h2dTrackPropAnalysisTracks", kTH3D, {{32, -0.5, 31.5f}, {32, -0.5, 31.5f}, axisCentrality}); + histos.add("h2dTrackPropAnalysisTopo", "h2dTrackPropAnalysisTopo", kTH3D, {{32, -0.5, 31.5f}, {32, -0.5, 31.5f}, axisCentrality}); + histos.add("h2dTrackPropAnalysisSpecies", "h2dTrackPropAnalysisSpecies", kTH3D, {{32, -0.5, 31.5f}, {32, -0.5, 31.5f}, axisCentrality}); } void processEvents( @@ -117,6 +174,46 @@ struct findableStudy { histos.fill(HIST("hCentrality"), collision.centFT0C()); } + void processDebugCrossCheck( + reconstructedV0sNoMC::iterator const& recV0, // reco V0s for cross-check + dauTracks const& // daughter track extras + ) + { + if (recV0.v0Type() == 1) { + // de-reference daughter track extras + auto pTrack = recV0.posTrackExtra_as(); + auto nTrack = recV0.negTrackExtra_as(); + + // cross-check correctness of new getter + if (pTrack.hasTPC() && !pTrack.hasITS() && !pTrack.hasTRD() && !pTrack.hasTOF()) { + LOGF(info, "X-check: Positive track is TPC only and this is a found V0. Puzzling!"); + } + if (nTrack.hasTPC() && !nTrack.hasITS() && !nTrack.hasTRD() && !nTrack.hasTOF()) { + LOGF(info, "X-check: Negative track is TPC only and this is a found V0. Puzzling!"); + } + } + } + + void processDebugCrossCheckMC( + reconstructedV0s::iterator const& recV0, // reco V0s for cross-check + dauTracks const& // daughter track extras + ) + { + if (recV0.isFound() && recV0.v0Type() == 1) { + // de-reference daughter track extras + auto pTrack = recV0.posTrackExtra_as(); + auto nTrack = recV0.negTrackExtra_as(); + + // cross-check correctness of new getter + if (pTrack.hasTPC() && !pTrack.hasITS() && !pTrack.hasTRD() && !pTrack.hasTOF()) { + LOGF(info, "X-check: Positive track is TPC only and this is a found V0. Puzzling!"); + } + if (nTrack.hasTPC() && !nTrack.hasITS() && !nTrack.hasTRD() && !nTrack.hasTOF()) { + LOGF(info, "X-check: Negative track is TPC only and this is a found V0. Puzzling!"); + } + } + } + void processV0s( aod::V0MCCores::iterator const& v0, // non-duplicated MC V0 decays soa::SmallGroups const& recv0s, // reconstructed versions of the v0 @@ -125,12 +222,30 @@ struct findableStudy { dauTracks const& // daughter track extras ) { - if (v0.pdgCode() != pdgCode || v0.pdgCodePositive() != 211 || v0.pdgCodeNegative() != -211) + int pdgCodePositive = 211; + int pdgCodeNegative = -211; + if (pdgCode == 3122) + pdgCodePositive = 2212; + if (pdgCode == -3122) + pdgCodePositive = -2212; + if (pdgCode == 22) { + pdgCodePositive = -11; + pdgCodeNegative = -11; + } + + if (v0.pdgCode() != pdgCode || v0.pdgCodePositive() != pdgCodePositive || v0.pdgCodeNegative() != pdgCodeNegative) return; if (!v0.isPhysicalPrimary()) return; - float rapidity = RecoDecay::y(std::array{v0.pxPosMC() + v0.pxNegMC(), v0.pyPosMC() + v0.pyNegMC(), v0.pzPosMC() + v0.pzNegMC()}, o2::constants::physics::MassKaonNeutral); + float rapidity = 2.0; + if (pdgCode == 310) + rapidity = RecoDecay::y(std::array{v0.pxPosMC() + v0.pxNegMC(), v0.pyPosMC() + v0.pyNegMC(), v0.pzPosMC() + v0.pzNegMC()}, o2::constants::physics::MassKaonNeutral); + if (pdgCode == 22) + rapidity = RecoDecay::y(std::array{v0.pxPosMC() + v0.pxNegMC(), v0.pyPosMC() + v0.pyNegMC(), v0.pzPosMC() + v0.pzNegMC()}, o2::constants::physics::MassPhoton); + if (pdgCode == 3122 || pdgCode == -3122) + rapidity = RecoDecay::y(std::array{v0.pxPosMC() + v0.pxNegMC(), v0.pyPosMC() + v0.pyNegMC(), v0.pzPosMC() + v0.pzNegMC()}, o2::constants::physics::MassLambda0); + if (std::abs(rapidity) > 0.5f) return; @@ -140,19 +255,31 @@ struct findableStudy { histos.fill(HIST("hNRecoV0s"), recv0s.size()); bool hasWrongCollision = false; float centrality = 100.5f; + bool hasBeenAcceptablyTracked = false; + bool hasBeenFoundAny = false; bool hasBeenFound = false; - - // encode conditionals here - uint32_t withTPC = 0; // if prongs have TPC - uint32_t withITSTracker = 0; // if prongs have been ITS tracked - uint32_t withITSTrackerTPC = 0; // if prongs have TPC and are ITS tracked - uint32_t withITSABTPC = 0; // if prongs have TPC and are ITS afterburned - uint32_t withSVertexerOK = 0; // if prongs have acceptable tracking conditions for svertexer + int nCandidatesWithTPC = 0; for (auto& recv0 : recv0s) { - if (recv0.isFound()) { - hasBeenFound = true; + if (recv0.v0Type() != 1) + continue; // skip anything other than a standard V0 + + // de-reference daughter track extras + auto pTrack = recv0.posTrackExtra_as(); + auto nTrack = recv0.negTrackExtra_as(); + + // skip ITS-only for simplicity + if (skipITSonly) { + if (!pTrack.hasTPC() || !nTrack.hasTPC()) + continue; } + + // define properties for this V0 + bool pTrackOK = false, nTrackOK = false; // tracks are acceptably tracked + + // Detailed analysis level + bool topoV0RadiusOK = false, topoV0RadiusMaxOK = false, topoV0CosPAOK = false, topoDcaPosToPVOK = false, topoDcaNegToPVOK = false, topoDcaV0DauOK = false; + auto coll = recv0.straCollision_as(); int mcCollID_fromCollision = coll.straMCCollisionId(); int mcCollID_fromV0 = recv0.straMCCollisionId(); @@ -164,61 +291,173 @@ struct findableStudy { centrality = coll.centFT0C(); } - // de-reference daughter track extras - auto pTrack = recv0.posTrackExtra_as(); - auto nTrack = recv0.negTrackExtra_as(); - - if (pTrack.hasTPC()) - bitset(withTPC, 0); - if (nTrack.hasTPC()) - bitset(withTPC, 1); + if (recv0.isFound()) { + hasBeenFoundAny = true; // includes also ITS-only, checked before skipITSonly check + } if ( - (pTrack.hasTPC() && pTrack.hasITS()) || - (!pTrack.hasTPC() && pTrack.itsNCls() >= 6)) - bitset(withSVertexerOK, 0); + (pTrack.hasTPC() && pTrack.hasITS()) || // full global track + (pTrack.hasTPC() && pTrack.hasTOF()) || // TPC + TOF is accepted + (pTrack.hasTPC() && pTrack.hasTRD()) || // TPC + TRD is accepted + (!pTrack.hasTPC() && pTrack.itsNCls() >= 6) // long ITS-only + ) { + pTrackOK = true; // for this V0 only + } if ( (nTrack.hasTPC() && nTrack.hasITS()) || - (!nTrack.hasTPC() && nTrack.itsNCls() >= 6)) - bitset(withSVertexerOK, 1); - - if (pTrack.hasITS() && pTrack.itsChi2PerNcl() > -10.0f) - bitset(withITSTracker, 0); - if (nTrack.hasITS() && nTrack.itsChi2PerNcl() > -10.0f) - bitset(withITSTracker, 1); - - if (pTrack.hasTPC() && pTrack.hasITS() && pTrack.itsChi2PerNcl() > -10.0f) - bitset(withITSTrackerTPC, 0); - if (nTrack.hasTPC() && nTrack.hasITS() && nTrack.itsChi2PerNcl() > -10.0f) - bitset(withITSTrackerTPC, 1); - - if (pTrack.hasTPC() && pTrack.hasITS() && pTrack.itsChi2PerNcl() < -10.0f) - bitset(withITSABTPC, 0); - if (nTrack.hasTPC() && nTrack.hasITS() && nTrack.itsChi2PerNcl() < -10.0f) - bitset(withITSABTPC, 1); + (nTrack.hasTPC() && nTrack.hasTOF()) || // TPC + TOF is accepted + (nTrack.hasTPC() && nTrack.hasTRD()) || // TPC + TRD is accepted + (!nTrack.hasTPC() && nTrack.itsNCls() >= 6)) { + nTrackOK = true; // for this V0 only + } + + if (pTrackOK && nTrackOK) + hasBeenAcceptablyTracked = true; + + // cross-check correctness of new getter + if (pTrack.hasITSTracker() && (pTrack.hasITS() && pTrack.itsChi2PerNcl() < -1e-3)) { + LOGF(fatal, "Positive track: inconsistent outcome of ITS tracker getter and explicit check!"); + } + if (nTrack.hasITSTracker() && (pTrack.hasITS() && nTrack.itsChi2PerNcl() < -1e-3)) { + LOGF(fatal, "Negative track: inconsistent outcome of ITS tracker getter and explicit check!"); + } + if (pTrack.hasITSAfterburner() && (pTrack.hasITS() && pTrack.itsChi2PerNcl() > -1e-3)) { + LOGF(fatal, "Positive track: inconsistent outcome of ITS tracker getter and explicit check!"); + } + if (nTrack.hasITSAfterburner() && (pTrack.hasITS() && nTrack.itsChi2PerNcl() > -1e-3)) { + LOGF(fatal, "Negative track: inconsistent outcome of ITS tracker getter and explicit check!"); + } + + // Cross-checking consistency: found should be a subset of nTrack, pTrack OK + histos.fill(HIST("hFoundVsTracksOK"), nTrackOK && pTrackOK, recv0.isFound()); + + // encode conditions of tracks + uint8_t positiveTrackCode = ((uint8_t(pTrack.hasTPC()) << hasTPC) | + (uint8_t(pTrack.hasITSTracker()) << hasITSTracker) | + (uint8_t(pTrack.hasITSAfterburner()) << hasITSAfterburner) | + (uint8_t(pTrack.hasTRD()) << hasTRD) | + (uint8_t(pTrack.hasTOF()) << hasTOF)); + + uint8_t negativeTrackCode = ((uint8_t(nTrack.hasTPC()) << hasTPC) | + (uint8_t(nTrack.hasITSTracker()) << hasITSTracker) | + (uint8_t(nTrack.hasITSAfterburner()) << hasITSAfterburner) | + (uint8_t(nTrack.hasTRD()) << hasTRD) | + (uint8_t(nTrack.hasTOF()) << hasTOF)); + + if (pTrackOK && nTrackOK && ptmc > 1.0 && ptmc < 1.1) { + // this particular V0 reco entry has been acceptably tracked. Do bookkeeping + histos.fill(HIST("h2dTrackPropAcceptablyTracked"), positiveTrackCode, negativeTrackCode, centrality); + } + + // determine if this V0 would go to analysis or not + if (recv0.isFound() && pTrackOK && nTrackOK) { // hack to avoid type check; only interested in found type 1 + // at this stage, this should be REALLY mostly unique (unless you switch skipITSonly to false or so) + // ... but we will cross-check this assumption (hNRecoV0sWithTPC, h2dPtVsCentrality_FoundInLoop) + if (pTrack.hasTPC() && !pTrack.hasITS() && !pTrack.hasTRD() && !pTrack.hasTOF()) { + LOGF(info, "Positive track is TPC only and this is a found V0. Puzzling!"); + } + if (nTrack.hasTPC() && !nTrack.hasITS() && !nTrack.hasTRD() && !nTrack.hasTOF()) { + LOGF(info, "Negative track is TPC only and this is a found V0. Puzzling!"); + } + + nCandidatesWithTPC++; + hasBeenFound = true; + histos.fill(HIST("h2dPtVsCentrality_FoundInLoop"), centrality, ptmc); + if (ptmc > 1.0 && ptmc < 1.1) { + histos.fill(HIST("h2dTrackPropFound"), positiveTrackCode, negativeTrackCode, centrality); + } + + uint64_t selMap = v0data::computeReconstructionBitmap(recv0, pTrack, nTrack, coll, recv0.yLambda(), recv0.yK0Short(), v0Selections); + + // Consider in all cases + selMap = selMap | (uint64_t(1) << v0data::selConsiderK0Short) | (uint64_t(1) << v0data::selConsiderLambda) | (uint64_t(1) << v0data::selConsiderAntiLambda); + + // selection checker: ensure this works on subset of actual svertexer-findable + bool validTrackProperties = v0Selections->verifyMask(selMap, maskTrackProperties) && pTrackOK && nTrackOK; + bool validTopology = v0Selections->verifyMask(selMap, maskTopological); + + uint64_t thisSpeciesMask = maskK0ShortSpecific; + if (pdgCode == 3122) + thisSpeciesMask = maskLambdaSpecific; + if (pdgCode == -3122) + thisSpeciesMask = maskAntiLambdaSpecific; + // add other species masks as necessary + + bool validThisSpecies = v0Selections->verifyMask(selMap, thisSpeciesMask); + + // specific selection (not cumulative) + topoV0RadiusOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selRadius); + topoV0RadiusMaxOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selRadiusMax); + topoV0CosPAOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selCosPA); + topoDcaPosToPVOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selDCAPosToPV); + topoDcaNegToPVOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selDCANegToPV); + topoDcaV0DauOK = v0Selections->verifyMask(selMap, uint64_t(1) << v0data::selDCAV0Dau); + // trackTPCRowsOK = v0Selections->verifyMask(selMap, (uint64_t(1) << v0data::selPosGoodTPCTrack) | (uint64_t(1) << v0data::selNegGoodTPCTrack) ); + + // uint64_t tpcPidMask = (uint64_t(1) << v0data::selTPCPIDPositivePion) | (uint64_t(1) << v0data::selTPCPIDNegativePion); + // if(pdgCode==3122) + // tpcPidMask = (uint64_t(1) << v0data::selTPCPIDPositiveProton) | (uint64_t(1) << v0data::selTPCPIDNegativePion); + // if(pdgCode==-3122) + // tpcPidMask = (uint64_t(1) << v0data::selTPCPIDPositivePion) | (uint64_t(1) << v0data::selTPCPIDNegativeProton); + + // trackTPCPIDOK = v0Selections->verifyMask(selMap, tpcPidMask); + + // Broad level + if (validTrackProperties) { + histos.fill(HIST("h2dPtVsCentrality_PassesTrackQuality"), centrality, ptmc); + if (ptmc > 1.0 && ptmc < 1.1) { + histos.fill(HIST("h2dTrackPropAnalysisTracks"), positiveTrackCode, negativeTrackCode, centrality); + } + } + if (validTrackProperties && validTopology) { + histos.fill(HIST("h2dPtVsCentrality_PassesTopological"), centrality, ptmc); + if (ptmc > 1.0 && ptmc < 1.1) { + histos.fill(HIST("h2dTrackPropAnalysisTopo"), positiveTrackCode, negativeTrackCode, centrality); + } + } + if (validTrackProperties && validTopology && validThisSpecies) { + histos.fill(HIST("h2dPtVsCentrality_PassesThisSpecies"), centrality, ptmc); + if (ptmc > 1.0 && ptmc < 1.1) { + histos.fill(HIST("h2dTrackPropAnalysisSpecies"), positiveTrackCode, negativeTrackCode, centrality); + } + } + + // topological + if (validTrackProperties && validThisSpecies && topoV0RadiusOK) + histos.fill(HIST("h2dPtVsCentrality_V0Radius"), centrality, ptmc); + if (validTrackProperties && validThisSpecies && topoV0RadiusMaxOK) + histos.fill(HIST("h2dPtVsCentrality_V0RadiusMax"), centrality, ptmc); + if (validTrackProperties && validThisSpecies && topoV0CosPAOK) + histos.fill(HIST("h2dPtVsCentrality_V0CosPA"), centrality, ptmc); + if (validTrackProperties && validThisSpecies && topoDcaPosToPVOK) + histos.fill(HIST("h2dPtVsCentrality_DcaPosToPV"), centrality, ptmc); + if (validTrackProperties && validThisSpecies && topoDcaNegToPVOK) + histos.fill(HIST("h2dPtVsCentrality_DcaNegToPV"), centrality, ptmc); + if (validTrackProperties && validThisSpecies && topoDcaV0DauOK) + histos.fill(HIST("h2dPtVsCentrality_DcaV0Dau"), centrality, ptmc); + } } - histos.fill(HIST("h2dPtVsCentrality_All_Findable"), centrality, ptmc); - histos.fill(HIST("h2dPtVsCentrality_WithTPC_Findable"), centrality, ptmc, withTPC); - histos.fill(HIST("h2dPtVsCentrality_WithITSTracker_Findable"), centrality, ptmc, withITSTracker); - histos.fill(HIST("h2dPtVsCentrality_WithITSTrackerTPC_Findable"), centrality, ptmc, withITSTrackerTPC); - histos.fill(HIST("h2dPtVsCentrality_WithITSABTPC_Findable"), centrality, ptmc, withITSABTPC); - histos.fill(HIST("h2dPtVsCentrality_WithSVertexerOK_Findable"), centrality, ptmc, withSVertexerOK); + histos.fill(HIST("hNRecoV0sWithTPC"), nCandidatesWithTPC); + // Major check 1: Findable versus found in some capacity + histos.fill(HIST("h2dPtVsCentrality_Findable"), centrality, ptmc); + if (hasBeenAcceptablyTracked) { + histos.fill(HIST("h2dPtVsCentrality_AcceptablyTracked"), centrality, ptmc); + } + if (hasBeenFoundAny) { + histos.fill(HIST("h2dPtVsCentrality_FoundAny"), centrality, ptmc); + } if (hasBeenFound) { - histos.fill(HIST("h2dPtVsCentrality_All_Found"), centrality, ptmc); - histos.fill(HIST("h2dPtVsCentrality_WithTPC_Found"), centrality, ptmc, withTPC); - histos.fill(HIST("h2dPtVsCentrality_WithITSTracker_Found"), centrality, ptmc, withITSTracker); - histos.fill(HIST("h2dPtVsCentrality_WithITSTrackerTPC_Found"), centrality, ptmc, withITSTrackerTPC); - histos.fill(HIST("h2dPtVsCentrality_WithITSABTPC_Found"), centrality, ptmc, withITSABTPC); - histos.fill(HIST("h2dPtVsCentrality_WithSVertexerOK_Found"), centrality, ptmc, withSVertexerOK); + histos.fill(HIST("h2dPtVsCentrality_Found"), centrality, ptmc); } - if (hasWrongCollision) { histos.fill(HIST("hNRecoV0sWrongColl"), recv0s.size()); } } PROCESS_SWITCH(findableStudy, processEvents, "process collision counters", true); + PROCESS_SWITCH(findableStudy, processDebugCrossCheck, "process debug cross-check of V0 with TPC-only", true); + PROCESS_SWITCH(findableStudy, processDebugCrossCheckMC, "process debug cross-check of V0 with TPC-only, MC version", false); PROCESS_SWITCH(findableStudy, processV0s, "process V0s", true); }; diff --git a/PWGLF/Utils/CMakeLists.txt b/PWGLF/Utils/CMakeLists.txt new file mode 100644 index 00000000000..bc3fd97b824 --- /dev/null +++ b/PWGLF/Utils/CMakeLists.txt @@ -0,0 +1,18 @@ +# Copyright 2019-2020 CERN and copyright holders of ALICE O2. +# See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +# All rights not expressly granted are reserved. +# +# This software is distributed under the terms of the GNU General Public +# License v3 (GPL Version 3), copied verbatim in the file "COPYING". +# +# In applying this license CERN does not waive the privileges and immunities +# granted to it by virtue of its status as an Intergovernmental Organization +# or submit itself to any jurisdiction. + +o2physics_add_library(v0SelectionGroup + SOURCES v0SelectionGroup.cxx + PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore) + +o2physics_target_root_dictionary(v0SelectionGroup + HEADERS v0SelectionGroup.h + LINKDEF v0SelectionGroupLinkDef.h) diff --git a/PWGLF/Utils/v0SelectionBits.h b/PWGLF/Utils/v0SelectionBits.h new file mode 100644 index 00000000000..dd4b48d5a0b --- /dev/null +++ b/PWGLF/Utils/v0SelectionBits.h @@ -0,0 +1,62 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef PWGLF_UTILS_V0SELECTIONBITS_H_ +#define PWGLF_UTILS_V0SELECTIONBITS_H_ + +namespace v0data +{ +// provides simple switches +enum selection : uint64_t { selCosPA = 0, + selRadius, + selRadiusMax, + selDCANegToPV, + selDCAPosToPV, + selDCAV0Dau, + selK0ShortRapidity, + selLambdaRapidity, + selTPCPIDPositivePion, + selTPCPIDNegativePion, + selTPCPIDPositiveProton, + selTPCPIDNegativeProton, + selTOFDeltaTPositiveProtonLambda, + selTOFDeltaTPositivePionLambda, + selTOFDeltaTPositivePionK0Short, + selTOFDeltaTNegativeProtonLambda, + selTOFDeltaTNegativePionLambda, + selTOFDeltaTNegativePionK0Short, + selTOFNSigmaPositiveProtonLambda, // Nsigma + selTOFNSigmaPositivePionLambda, // Nsigma + selTOFNSigmaPositivePionK0Short, // Nsigma + selTOFNSigmaNegativeProtonLambda, // Nsigma + selTOFNSigmaNegativePionLambda, // Nsigma + selTOFNSigmaNegativePionK0Short, // Nsigma + selK0ShortCTau, + selLambdaCTau, + selK0ShortArmenteros, + selPosGoodTPCTrack, // at least min # TPC rows + selNegGoodTPCTrack, // at least min # TPC rows + selPosGoodITSTrack, // at least min # ITS clusters + selNegGoodITSTrack, // at least min # ITS clusters + selPosItsOnly, + selNegItsOnly, + selPosNotTPCOnly, + selNegNotTPCOnly, + selConsiderK0Short, // for mc tagging + selConsiderLambda, // for mc tagging + selConsiderAntiLambda, // for mc tagging + selPhysPrimK0Short, // for mc tagging + selPhysPrimLambda, // for mc tagging + selPhysPrimAntiLambda, // for mc tagging +}; +} // namespace v0data + +#endif // PWGLF_UTILS_V0SELECTIONBITS_H_ diff --git a/PWGLF/Utils/v0SelectionGroup.cxx b/PWGLF/Utils/v0SelectionGroup.cxx new file mode 100644 index 00000000000..6ceffa27482 --- /dev/null +++ b/PWGLF/Utils/v0SelectionGroup.cxx @@ -0,0 +1,107 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// see header for a more detailed description. + +#include "v0SelectionGroup.h" + +// simple checkers, but ensure 64 bit integers +#define bitset(var, nbit) ((var) |= (static_cast(1) << static_cast(nbit))) +#define bitcheck(var, nbit) ((var) & (static_cast(1) << static_cast(nbit))) + +// trivial 64-bit mask checker +bool v0SelectionGroup::verifyMask(uint64_t bitmap, uint64_t mask) const +{ + return (bitmap & mask) == mask; +} + +// provides standard masks given current event selection criteria. +void v0SelectionGroup::provideMasks(uint64_t& maskTopological, uint64_t& maskTrackProperties, uint64_t& maskK0ShortSpecific, uint64_t& maskLambdaSpecific, uint64_t& maskAntiLambdaSpecific) const +{ + maskTopological = (uint64_t(1) << v0data::selCosPA) | (uint64_t(1) << v0data::selRadius) | (uint64_t(1) << v0data::selDCANegToPV) | (uint64_t(1) << v0data::selDCAPosToPV) | (uint64_t(1) << v0data::selDCAV0Dau) | (uint64_t(1) << v0data::selRadiusMax); + maskK0ShortSpecific = (uint64_t(1) << v0data::selK0ShortRapidity) | (uint64_t(1) << v0data::selK0ShortCTau) | (uint64_t(1) << v0data::selK0ShortArmenteros) | (uint64_t(1) << v0data::selConsiderK0Short); + maskLambdaSpecific = (uint64_t(1) << v0data::selLambdaRapidity) | (uint64_t(1) << v0data::selLambdaCTau) | (uint64_t(1) << v0data::selConsiderLambda); + maskAntiLambdaSpecific = (uint64_t(1) << v0data::selLambdaRapidity) | (uint64_t(1) << v0data::selLambdaCTau) | (uint64_t(1) << v0data::selConsiderAntiLambda); + + // ask for specific TPC/TOF PID v0data::selections + maskTrackProperties = 0; + if (requirePosITSonly) { + maskTrackProperties = maskTrackProperties | (uint64_t(1) << v0data::selPosItsOnly) | (uint64_t(1) << v0data::selPosGoodITSTrack); + } else { + maskTrackProperties = maskTrackProperties | (uint64_t(1) << v0data::selPosGoodTPCTrack) | (uint64_t(1) << v0data::selPosGoodITSTrack); + // TPC signal is available: ask for positive track PID + if (TpcPidNsigmaCut < 1e+5) { // safeguard for no cut + maskK0ShortSpecific = maskK0ShortSpecific | (uint64_t(1) << v0data::selTPCPIDPositivePion); + maskLambdaSpecific = maskLambdaSpecific | (uint64_t(1) << v0data::selTPCPIDPositiveProton); + maskAntiLambdaSpecific = maskAntiLambdaSpecific | (uint64_t(1) << v0data::selTPCPIDPositivePion); + } + // TOF PID + if (TofPidNsigmaCutK0Pi < 1e+5) // safeguard for no cut + maskK0ShortSpecific = maskK0ShortSpecific | (uint64_t(1) << v0data::selTOFNSigmaPositivePionK0Short) | (uint64_t(1) << v0data::selTOFDeltaTPositivePionK0Short); + if (TofPidNsigmaCutLaPr < 1e+5) // safeguard for no cut + maskLambdaSpecific = maskLambdaSpecific | (uint64_t(1) << v0data::selTOFNSigmaPositiveProtonLambda) | (uint64_t(1) << v0data::selTOFDeltaTPositiveProtonLambda); + if (TofPidNsigmaCutLaPi < 1e+5) // safeguard for no cut + maskAntiLambdaSpecific = maskAntiLambdaSpecific | (uint64_t(1) << v0data::selTOFNSigmaPositivePionLambda) | (uint64_t(1) << v0data::selTOFDeltaTPositivePionLambda); + } + if (requireNegITSonly) { + maskTrackProperties = maskTrackProperties | (uint64_t(1) << v0data::selNegItsOnly) | (uint64_t(1) << v0data::selNegGoodITSTrack); + } else { + maskTrackProperties = maskTrackProperties | (uint64_t(1) << v0data::selNegGoodTPCTrack) | (uint64_t(1) << v0data::selNegGoodITSTrack); + // TPC signal is available: ask for negative track PID + if (TpcPidNsigmaCut < 1e+5) { // safeguard for no cut + maskK0ShortSpecific = maskK0ShortSpecific | (uint64_t(1) << v0data::selTPCPIDNegativePion); + maskLambdaSpecific = maskLambdaSpecific | (uint64_t(1) << v0data::selTPCPIDNegativePion); + maskAntiLambdaSpecific = maskAntiLambdaSpecific | (uint64_t(1) << v0data::selTPCPIDNegativeProton); + } + // TOF PID + if (TofPidNsigmaCutK0Pi < 1e+5) // safeguard for no cut + maskK0ShortSpecific = maskK0ShortSpecific | (uint64_t(1) << v0data::selTOFNSigmaNegativePionK0Short) | (uint64_t(1) << v0data::selTOFDeltaTNegativePionK0Short); + if (TofPidNsigmaCutLaPr < 1e+5) // safeguard for no cut + maskLambdaSpecific = maskLambdaSpecific | (uint64_t(1) << v0data::selTOFNSigmaNegativePionLambda) | (uint64_t(1) << v0data::selTOFDeltaTNegativePionLambda); + if (TofPidNsigmaCutLaPi < 1e+5) // safeguard for no cut + maskAntiLambdaSpecific = maskAntiLambdaSpecific | (uint64_t(1) << v0data::selTOFNSigmaNegativeProtonLambda) | (uint64_t(1) << v0data::selTOFDeltaTNegativeProtonLambda); + } + + if (skipTPConly) { + maskK0ShortSpecific = maskK0ShortSpecific | (uint64_t(1) << v0data::selPosNotTPCOnly) | (uint64_t(1) << v0data::selNegNotTPCOnly); + maskLambdaSpecific = maskLambdaSpecific | (uint64_t(1) << v0data::selPosNotTPCOnly) | (uint64_t(1) << v0data::selNegNotTPCOnly); + maskAntiLambdaSpecific = maskAntiLambdaSpecific | (uint64_t(1) << v0data::selPosNotTPCOnly) | (uint64_t(1) << v0data::selNegNotTPCOnly); + } +} + +void v0SelectionGroup::PrintSelections() const +{ + LOGF(info, "+++ Phase space selections ++++++++++++++++++++++++"); + LOGF(info, "Rapidity cut ...........................: %.2f", rapidityCut); + LOGF(info, "Eta (daughters) cut ....................: %.2f", daughterEtaCut); + LOGF(info, "+++ Topological selections ++++++++++++++++++++++++"); + LOGF(info, "V0 cosine of PA ........................: %.2f", v0cospa); + LOGF(info, "DCA V0 daughters .......................: %.2f", dcav0dau); + LOGF(info, "Neg track DCA to PV ....................: %.2f", dcanegtopv); + LOGF(info, "Pos track DCA to PV ....................: %.2f", dcapostopv); + LOGF(info, "Minimum radius .........................: %.2f", v0radius); + LOGF(info, "Maximum radius .........................: %.2f", v0radiusMax); + LOGF(info, "+++ Track quality selections ++++++++++++++++++++++"); + LOGF(info, "Minimum TPC rows .......................: %i", minTPCrows); + LOGF(info, "Minimum ITS clusters ...................: %i", minITSclusters); + LOGF(info, "Skip TPC only ..........................: %s", skipTPConly ? "true" : "false"); + LOGF(info, "Require positive ITS only ..............: %s", requirePosITSonly ? "true" : "false"); + LOGF(info, "Require negative ITS only ..............: %s", requireNegITSonly ? "true" : "false"); + LOGF(info, "+++ Particle identification selections ++++++++++++"); + LOGF(info, "TPC PID Nsigma .........................: %.2f", TpcPidNsigmaCut); + LOGF(info, "TOF PID Nsigma LaPr ....................: %.2f", TofPidNsigmaCutLaPr); + LOGF(info, "TOF PID Nsigma LaPi ....................: %.2f", TofPidNsigmaCutLaPi); + LOGF(info, "TOF PID Nsigma K0Pi ....................: %.2f", TofPidNsigmaCutK0Pi); + LOGF(info, "+++ Misc selections ++++++++++++++++++++++++++++++"); + LOGF(info, "K0Short lifetime cut ...................: %.2f", lifetimeCutK0Short); + LOGF(info, "Lambda lifetime cut ....................: %.2f", lifetimeCutLambda); + LOGF(info, "Armenteros podolanski parameter ........: %.2f", armPodCut); +} diff --git a/PWGLF/Utils/v0SelectionGroup.h b/PWGLF/Utils/v0SelectionGroup.h new file mode 100644 index 00000000000..a353a86c742 --- /dev/null +++ b/PWGLF/Utils/v0SelectionGroup.h @@ -0,0 +1,105 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef PWGLF_UTILS_V0SELECTIONGROUP_H_ +#define PWGLF_UTILS_V0SELECTIONGROUP_H_ + +#include +#include +#include +#include "v0SelectionBits.h" +#include "Framework/Logger.h" +#include "CommonConstants/PhysicsConstants.h" + +class v0SelectionGroup +{ + public: + v0SelectionGroup() + : rapidityCut{0.5f}, daughterEtaCut{0.8f}, v0cospa{0.97f}, dcav0dau{1.0f}, dcanegtopv{0.05f}, dcapostopv{0.05f}, v0radius{1.2f}, v0radiusMax{1e+5}, minTPCrows{70}, minITSclusters{-1}, skipTPConly{false}, requirePosITSonly{false}, requireNegITSonly{false}, TpcPidNsigmaCut{5.0f}, TofPidNsigmaCutLaPr{1e+6}, TofPidNsigmaCutLaPi{1e+6}, TofPidNsigmaCutK0Pi{1e+6}, maxDeltaTimeProton{1e+9}, maxDeltaTimePion{1e+9}, lifetimeCutK0Short{20.0f}, lifetimeCutLambda{20.0f}, armPodCut{5.0f} + { + // constructor + } + + void provideMasks(uint64_t& maskTopological, uint64_t& maskTrackProperties, uint64_t& maskK0ShortSpecific, uint64_t& maskLambdaSpecific, uint64_t& maskAntiLambdaSpecific) const; + bool verifyMask(uint64_t bitmap, uint64_t mask) const; + + float getRapidityCut() const { return rapidityCut; } + float getDaughterEtaCut() const { return daughterEtaCut; } + + float getv0cospa() const { return v0cospa; } + float getdcav0dau() const { return dcav0dau; } + float getdcanegtopv() const { return dcanegtopv; } + float getdcapostopv() const { return dcapostopv; } + float getv0radius() const { return v0radius; } + float getv0radiusMax() const { return v0radiusMax; } + + int getminTPCrows() const { return minTPCrows; } + int getminITSclusters() const { return minITSclusters; } + bool getskipTPConly() const { return skipTPConly; } + bool getrequirePosITSonly() const { return requirePosITSonly; } + bool getrequireNegITSonly() const { return requireNegITSonly; } + + float getTpcPidNsigmaCut() const { return TpcPidNsigmaCut; } + float getTofPidNsigmaCutLaPr() const { return TofPidNsigmaCutLaPr; } + float getTofPidNsigmaCutLaPi() const { return TofPidNsigmaCutLaPi; } + float getTofPidNsigmaCutK0Pi() const { return TofPidNsigmaCutK0Pi; } + + float getmaxDeltaTimeProton() const { return maxDeltaTimeProton; } + float getmaxDeltaTimePion() const { return maxDeltaTimePion; } + + float getlifetimeCutK0Short() const { return lifetimeCutK0Short; } + float getlifetimeCutLambda() const { return lifetimeCutLambda; } + + float getarmPodCut() const { return armPodCut; } + + // Helper to print out selections + void PrintSelections() const; + + private: + // Phase space + float rapidityCut; + float daughterEtaCut; + // Topology + float v0cospa; + float dcav0dau; + float dcanegtopv; + float dcapostopv; + float v0radius; + float v0radiusMax; + + // Track quality + int minTPCrows; + int minITSclusters; + bool skipTPConly; + bool requirePosITSonly; + bool requireNegITSonly; + + // PID + float TpcPidNsigmaCut; + float TofPidNsigmaCutLaPr; + float TofPidNsigmaCutLaPi; + float TofPidNsigmaCutK0Pi; + + // PID precursor (compatibility only); + float maxDeltaTimeProton; + float maxDeltaTimePion; + + // Lifetime + float lifetimeCutK0Short; + float lifetimeCutLambda; + + // Armenteros + float armPodCut; + + ClassDefNV(v0SelectionGroup, 1); +}; + +#endif // PWGLF_UTILS_V0SELECTIONGROUP_H_ diff --git a/PWGLF/Utils/v0SelectionGroupLinkDef.h b/PWGLF/Utils/v0SelectionGroupLinkDef.h new file mode 100644 index 00000000000..3c81add2fc2 --- /dev/null +++ b/PWGLF/Utils/v0SelectionGroupLinkDef.h @@ -0,0 +1,25 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef PWGLF_UTILS_V0SELECTIONGROUPLINKDEF_H_ +#define PWGLF_UTILS_V0SELECTIONGROUPLINKDEF_H_ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; +#pragma link C++ nestedclasses; +#pragma link C++ class v0SelectionGroup + ; +#pragma link C++ class std::vector < int> + ; +#pragma link C++ class std::vector < float> + ; +#pragma link C++ class std::vector < double> + ; +#pragma link C++ class std::vector < bool> + ; + +#endif // PWGLF_UTILS_V0SELECTIONGROUPLINKDEF_H_ diff --git a/PWGLF/Utils/v0SelectionTools.h b/PWGLF/Utils/v0SelectionTools.h new file mode 100644 index 00000000000..c88e43eaace --- /dev/null +++ b/PWGLF/Utils/v0SelectionTools.h @@ -0,0 +1,131 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef PWGLF_UTILS_V0SELECTIONTOOLS_H_ +#define PWGLF_UTILS_V0SELECTIONTOOLS_H_ + +#include "v0SelectionBits.h" +#include "v0SelectionGroup.h" + +// simple checkers, but ensure 64 bit integers +#define bitset(var, nbit) ((var) |= (static_cast(1) << static_cast(nbit))) +#define bitcheck(var, nbit) ((var) & (static_cast(1) << static_cast(nbit))) + +namespace v0data +{ + +// utility method to calculate a selection map for the V0s +template +uint64_t computeReconstructionBitmap(TV0 v0, TTrack posTrackExtra, TTrack negTrackExtra, TCollision collision, float rapidityLambda, float rapidityK0Short, const v0SelectionGroup& v0sels) +// precalculate this information so that a check is one mask operation, not many +{ + uint64_t bitMap = 0; + // Base topological variables + if (v0.v0radius() > v0sels.getv0radius()) + bitset(bitMap, v0data::selRadius); + if (v0.v0radius() < v0sels.getv0radiusMax()) + bitset(bitMap, v0data::selRadiusMax); + if (TMath::Abs(v0.dcapostopv()) > v0sels.getdcapostopv()) + bitset(bitMap, v0data::selDCAPosToPV); + if (TMath::Abs(v0.dcanegtopv()) > v0sels.getdcanegtopv()) + bitset(bitMap, v0data::selDCANegToPV); + if (v0.v0cosPA() > v0sels.getv0cospa()) + bitset(bitMap, v0data::selCosPA); + if (v0.dcaV0daughters() < v0sels.getdcav0dau()) + bitset(bitMap, v0data::selDCAV0Dau); + + // rapidity + if (TMath::Abs(rapidityLambda) < v0sels.getRapidityCut()) + bitset(bitMap, v0data::selLambdaRapidity); + if (TMath::Abs(rapidityK0Short) < v0sels.getRapidityCut()) + bitset(bitMap, v0data::selK0ShortRapidity); + + // ITS quality flags + if (posTrackExtra.itsNCls() >= v0sels.getminITSclusters()) + bitset(bitMap, v0data::selPosGoodITSTrack); + if (negTrackExtra.itsNCls() >= v0sels.getminITSclusters()) + bitset(bitMap, v0data::selNegGoodITSTrack); + + // TPC quality flags + if (posTrackExtra.tpcCrossedRows() >= v0sels.getminTPCrows()) + bitset(bitMap, v0data::selPosGoodTPCTrack); + if (negTrackExtra.tpcCrossedRows() >= v0sels.getminTPCrows()) + bitset(bitMap, v0data::selNegGoodTPCTrack); + + // TPC PID + if (fabs(posTrackExtra.tpcNSigmaPi()) < v0sels.getTpcPidNsigmaCut()) + bitset(bitMap, v0data::selTPCPIDPositivePion); + if (fabs(posTrackExtra.tpcNSigmaPr()) < v0sels.getTpcPidNsigmaCut()) + bitset(bitMap, v0data::selTPCPIDPositiveProton); + if (fabs(negTrackExtra.tpcNSigmaPi()) < v0sels.getTpcPidNsigmaCut()) + bitset(bitMap, v0data::selTPCPIDNegativePion); + if (fabs(negTrackExtra.tpcNSigmaPr()) < v0sels.getTpcPidNsigmaCut()) + bitset(bitMap, v0data::selTPCPIDNegativeProton); + + // TOF PID in DeltaT (deprecated, kept for compatibility) + // Positive track + if (fabs(v0.posTOFDeltaTLaPr()) < v0sels.getmaxDeltaTimeProton()) + bitset(bitMap, v0data::selTOFDeltaTPositiveProtonLambda); + if (fabs(v0.posTOFDeltaTLaPi()) < v0sels.getmaxDeltaTimePion()) + bitset(bitMap, v0data::selTOFDeltaTPositivePionLambda); + if (fabs(v0.posTOFDeltaTK0Pi()) < v0sels.getmaxDeltaTimePion()) + bitset(bitMap, v0data::selTOFDeltaTPositivePionK0Short); + // Negative track + if (fabs(v0.negTOFDeltaTLaPr()) < v0sels.getmaxDeltaTimeProton()) + bitset(bitMap, v0data::selTOFDeltaTNegativeProtonLambda); + if (fabs(v0.negTOFDeltaTLaPi()) < v0sels.getmaxDeltaTimePion()) + bitset(bitMap, v0data::selTOFDeltaTNegativePionLambda); + if (fabs(v0.negTOFDeltaTK0Pi()) < v0sels.getmaxDeltaTimePion()) + bitset(bitMap, v0data::selTOFDeltaTNegativePionK0Short); + + // TOF PID in NSigma + // Positive track + if (fabs(v0.tofNSigmaLaPr()) < v0sels.getTofPidNsigmaCutLaPr()) + bitset(bitMap, v0data::selTOFNSigmaPositiveProtonLambda); + if (fabs(v0.tofNSigmaALaPi()) < v0sels.getTofPidNsigmaCutLaPi()) + bitset(bitMap, v0data::selTOFNSigmaPositivePionLambda); + if (fabs(v0.tofNSigmaK0PiPlus()) < v0sels.getTofPidNsigmaCutK0Pi()) + bitset(bitMap, v0data::selTOFNSigmaPositivePionK0Short); + // Negative track + if (fabs(v0.tofNSigmaALaPr()) < v0sels.getTofPidNsigmaCutLaPr()) + bitset(bitMap, v0data::selTOFNSigmaNegativeProtonLambda); + if (fabs(v0.tofNSigmaLaPi()) < v0sels.getTofPidNsigmaCutLaPi()) + bitset(bitMap, v0data::selTOFNSigmaNegativePionLambda); + if (fabs(v0.tofNSigmaK0PiMinus()) < v0sels.getTofPidNsigmaCutK0Pi()) + bitset(bitMap, v0data::selTOFNSigmaNegativePionK0Short); + + // ITS only tag + if (posTrackExtra.tpcCrossedRows() < 1) + bitset(bitMap, v0data::selPosItsOnly); + if (negTrackExtra.tpcCrossedRows() < 1) + bitset(bitMap, v0data::selNegItsOnly); + + // TPC only tag + if (posTrackExtra.detectorMap() != o2::aod::track::TPC) + bitset(bitMap, v0data::selPosNotTPCOnly); + if (negTrackExtra.detectorMap() != o2::aod::track::TPC) + bitset(bitMap, v0data::selNegNotTPCOnly); + + // proper lifetime + if (v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * o2::constants::physics::MassLambda0 < v0sels.getlifetimeCutLambda()) + bitset(bitMap, v0data::selLambdaCTau); + if (v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * o2::constants::physics::MassK0Short < v0sels.getlifetimeCutK0Short()) + bitset(bitMap, v0data::selK0ShortCTau); + + // armenteros + if (v0.qtarm() * v0sels.getarmPodCut() > TMath::Abs(v0.alpha()) || v0sels.getarmPodCut() < 1e-4) + bitset(bitMap, v0data::selK0ShortArmenteros); + + return bitMap; +} +} // namespace v0data + +#endif // PWGLF_UTILS_V0SELECTIONTOOLS_H_