diff --git a/PWGMM/Mult/Tasks/mftReassociationValidation.cxx b/PWGMM/Mult/Tasks/mftReassociationValidation.cxx index 59b9a20e76b..65cb143aa0a 100644 --- a/PWGMM/Mult/Tasks/mftReassociationValidation.cxx +++ b/PWGMM/Mult/Tasks/mftReassociationValidation.cxx @@ -272,6 +272,7 @@ struct MftReassociationValidation { Configurable fillGeneralHistograms{"fillGeneralHistograms", true, "Fill general reassociation efficiencyhistograms"}; Configurable fillSpecificHistograms{"fillSpecificHistograms", true, "Fill specific histograms (DCA, DeltaPos)"}; Configurable fillParametersHistograms{"fillParametersHistograms", true, "Fill parameters histograms"}; + Configurable fillDcaTemplateFit{"fillDcaTemplateFit", true, "Fill DCA template histograms for primary/secondary fraction fit"}; Configurable keepOnlyPhysicalPrimary{"keepOnlyPhysicalPrimary", true, "Keep only physical primary particles"}; Configurable keepOnlySecondaries{"keepOnlySecondaries", false, "Keep only secondary particles"}; Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; @@ -309,8 +310,6 @@ struct MftReassociationValidation { Configurable cutOrphanTracksExplicitly{"cutOrphanTracksExplicitly", false, "if true, cut orphan tracks explicitly"}; Configurable cutOnDcaXY{"cutOnDcaXY", false, "if true, cut on DCA XY"}; Configurable cutOnDcaZ{"cutOnDcaZ", false, "if true, cut on DCA Z"}; - Configurable dcaXYMax{"dcaXYMax", 2.0f, "Maximum value for DCA XY"}; - Configurable dcaZMax{"dcaZMax", 10.0f, "Maximum value for DCA Z"}; Configurable etaMftTrackMax{"etaMftTrackMax", -2.4f, "Maximum value for the eta of MFT tracks when used in cut function"}; Configurable etaMftTrackMin{"etaMftTrackMin", -3.36f, "Minimum value for the eta of MFT tracks when used in cut function"}; Configurable etaMftTrackMaxFilter{"etaMftTrackMaxFilter", -2.0f, "Maximum value for the eta of MFT tracks when used in filter"}; @@ -406,6 +405,8 @@ struct MftReassociationValidation { ConfigurableAxis axisDcaX{"axisDcaX", {500, -0.1, 2.1}, "DCAx binning (cm)"}; ConfigurableAxis axisDcaY{"axisDcaY", {500, -0.1, 2.1}, "DCAy binning (cm)"}; ConfigurableAxis axisDcaZ{"axisDcaZ", {500, -1., 1.}, "DCAz binning (cm)"}; + ConfigurableAxis axisDcaXYfit{"axisDcaXYfit", {500, 0., 2.}, "DCAxy binning (cm) for the template fit"}; + ConfigurableAxis axisDcaZfit{"axisDcaZfit", {500, -1., 1.}, "DCAz binning (cm) for the template fit"}; ConfigurableAxis axisMultiplicity{"axisMultiplicity", {VARIABLE_WIDTH, 0, 5, 10, 20, 30, 40, 50, 100.1}, "multiplicity axis for histograms"}; ConfigurableAxis axisPhi{"axisPhi", {72, 0, TwoPI}, "phi axis for histograms"}; ConfigurableAxis axisPt{"axisPt", {72, 0, 36}, "pt axis for histograms"}; @@ -578,6 +579,11 @@ struct MftReassociationValidation { registry.get(HIST("MC/2D/hIsTrueCollisionAmongCompatibleCollisionsDcaNotReassociated"))->GetXaxis()->SetBinLabel(iBin + 1, labelsMftIsTrueCollisionAmongCompatibleCollisionsStep[iBin].data()); } } + + if (configTask.fillDcaTemplateFit) { + registry.add("MC/2D/DcaFit/hPtVsDcaXY_prim", "Primary DCA_{xy} (2D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("MC/2D/DcaFit/hPtVsDcaXY_sec", "Secondary DCA_{xy} (2D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + } } void addMftMonteCarloHistograms3D() @@ -689,6 +695,20 @@ struct MftReassociationValidation { registry.get(HIST("MC/3D/hIsTrueCollisionAmongCompatibleCollisionsDcaNotReassociated"))->GetXaxis()->SetBinLabel(iBin + 1, labelsMftIsTrueCollisionAmongCompatibleCollisionsStep[iBin].data()); } } + + if (configTask.fillDcaTemplateFit) { + // 3D reassociated + registry.add("MC/3D/DcaFit/hPtVsDcaXY_prim", "Primary DCA_{xy} (3D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("MC/3D/DcaFit/hPtVsDcaXY_sec", "Secondary DCA_{xy} (3D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("MC/3D/DcaFit/hPtVsDcaZ_prim", "Primary DCA_{z} (3D);p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + registry.add("MC/3D/DcaFit/hPtVsDcaZ_sec", "Secondary DCA_{z} (3D);p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + + // Non-ambiguous + registry.add("MC/3D/DcaFit/NonAmb/hPtVsDcaXY_prim", "NonAmb Primary DCA_{xy};p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("MC/3D/DcaFit/NonAmb/hPtVsDcaXY_sec", "NonAmb Secondary DCA_{xy};p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("MC/3D/DcaFit/NonAmb/hPtVsDcaZ_prim", "NonAmb Primary DCA_{z};p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + registry.add("MC/3D/DcaFit/NonAmb/hPtVsDcaZ_sec", "NonAmb Secondary DCA_{z};p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + } } template @@ -833,6 +853,10 @@ struct MftReassociationValidation { if (doprocessData) { addMftHistograms(); + + if (configTask.fillDcaTemplateFit) { + registry.add("Data/2D/DcaFit/hPtVsDcaXY", "Data DCA_{xy} (2D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + } } if (doprocessMcReassociated2d) { @@ -857,6 +881,18 @@ struct MftReassociationValidation { addMftMonteCarloHistograms3D(); } + if (doprocessData3d) { + addMftHistograms(); + + if (configTask.fillDcaTemplateFit) { + registry.add("Data/3D/DcaFit/hPtVsDcaXY", "Data DCA_{xy} (3D);p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("Data/3D/DcaFit/hPtVsDcaZ", "Data DCA_{z} (3D);p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + + registry.add("Data/3D/DcaFit/NonAmb/hPtVsDcaXY", "NonAmb Data DCA_{xy};p_{T} (GeV/c);DCA_{xy} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaXYfit}); + registry.add("Data/3D/DcaFit/NonAmb/hPtVsDcaZ", "NonAmb Data DCA_{z};p_{T} (GeV/c);DCA_{z} (cm)", kTH2F, {configAxis.axisPt, configAxis.axisDcaZfit}); + } + } + } // End of init() function // ========================= @@ -1108,11 +1144,11 @@ struct MftReassociationValidation { } } - if (configMft.cutOnDcaXY && std::abs(dcaXY) > configMft.dcaXYMax) { + if (configMft.cutOnDcaXY && std::abs(dcaXY) > configMft.mftMaxDCAxy) { return false; } registry.fill(HIST("MC/") + HIST(WhatReassociationMethod[ReassociationMethod]) + HIST("hPreciseTrackSelectionCounter"), SpecificTrackSelectionStep::AfterDcaXYCut); - if (configMft.cutOnDcaZ && std::abs(dcaZ) > configMft.dcaZMax) { + if (configMft.cutOnDcaZ && std::abs(dcaZ) > configMft.mftMaxDCAz) { return false; } registry.fill(HIST("MC/") + HIST(WhatReassociationMethod[ReassociationMethod]) + HIST("hPreciseTrackSelectionCounter"), SpecificTrackSelectionStep::AfterDcaZCut); @@ -1199,31 +1235,77 @@ struct MftReassociationValidation { for (const auto& reassociated2dMftTrack : reassociated2dMftTracks) { - registry.fill(HIST("Data/hAmbiguityOfMftTracks"), MftTrackAmbiguityStep::AllMftTracks); auto templatedMftTrack = reassociated2dMftTrack.template mfttrack_as(); o2::track::TrackParCovFwd trackPar = o2::aod::fwdtrackutils::getTrackParCovFwdShift(templatedMftTrack, mZShift); std::array dcaInfOrig; trackPar.propagateToDCAhelix(bZ, {collision.posX(), collision.posY(), collision.posZ()}, dcaInfOrig); - auto dcaXYoriginal = 999.f; - dcaXYoriginal = std::sqrt(dcaInfOrig[0] * dcaInfOrig[0] + dcaInfOrig[1] * dcaInfOrig[1]); + // auto dcaXYoriginal = 999.f; + // dcaXYoriginal = std::sqrt(dcaInfOrig[0] * dcaInfOrig[0] + dcaInfOrig[1] * dcaInfOrig[1]); - if (!isAcceptedMftTrack(templatedMftTrack, false, true, dcaXYoriginal, dcaInfOrig[2])) { + if (!isAcceptedMftTrack(templatedMftTrack, false, true, reassociated2dMftTrack.bestDCAXY(), dcaInfOrig[2])) { continue; } - registry.fill(HIST("Data/hAmbiguityOfMftTracks"), MftTrackAmbiguityStep::AfterTrackSelection); - - if (isAmbiguousMftTrack(reassociated2dMftTrack, true)) { - registry.fill(HIST("Data/hReassociation2dMftTracks"), Reassociation2dMftTracks::NotReassociated2dMftTracks); + if (configTask.fillDcaTemplateFit) { + registry.fill(HIST("Data/2D/DcaFit/hPtVsDcaXY"), templatedMftTrack.pt(), reassociated2dMftTrack.bestDCAXY()); } } // end of loop over reassociated MFT tracks - - // TO-DO the same for reassociated3d (change the histograms) } PROCESS_SWITCH(MftReassociationValidation, processData, "Process MFT reassociation validation for DATA", false); + void processData3d(FilteredCollisionsWSelMult::iterator const& collision, + FilteredMftTracks const& /*mftTracks*/, + soa::SmallGroups const& reassociated3dMftTracks, + aod::BCsWithTimestamps const&) + { + auto bc = collision.template bc_as(); + if (areCorrectionsLoaded == false) { + loadZVertexShiftCorrection(bc); + areCorrectionsLoaded = true; + } + + if (!isAcceptedCollision(collision, false)) { + return; + } + + for (const auto& reassociated3dMftTrack : reassociated3dMftTracks) { + auto templatedMftTrack = reassociated3dMftTrack.template mfttrack_as(); + + // // Compute DCA at original collision for isAcceptedMftTrack + // o2::track::TrackParCovFwd trackPar = + // o2::aod::fwdtrackutils::getTrackParCovFwdShift(templatedMftTrack, mZShift); + // std::array dcaInfOrig; + // trackPar.propagateToDCAhelix( + // bZ, + // {collision.posX(), collision.posY(), collision.posZ()}, + // dcaInfOrig); + // const float dcaXYoriginal = + // std::sqrt(dcaInfOrig[0] * dcaInfOrig[0] + + // dcaInfOrig[1] * dcaInfOrig[1]); + + // Use existing track selection function (unchanged) + if (!isAcceptedMftTrack(templatedMftTrack, false, true, reassociated3dMftTrack.bestDCAXY(), reassociated3dMftTrack.bestDCAZ())) { + continue; + } + + if (configTask.fillDcaTemplateFit) { + registry.fill(HIST("Data/3D/DcaFit/hPtVsDcaXY"), templatedMftTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("Data/3D/DcaFit/hPtVsDcaZ"), templatedMftTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + } + + // The same for nonambi only + if (reassociated3dMftTrack.ambDegree() == 1) { + if (configTask.fillDcaTemplateFit) { + registry.fill(HIST("Data/3D/DcaFit/NonAmb/hPtVsDcaXY"), templatedMftTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("Data/3D/DcaFit/NonAmb/hPtVsDcaZ"), templatedMftTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + } + } + } + } + PROCESS_SWITCH(MftReassociationValidation, processData3d, "Process DCA fit histograms for DATA 3D reassociation", false); + void processCreateLookupTable(FilteredCollisionsWSelMultMcLabels const& collisions, soa::Join const& mcCollisions) { recoVtxX.clear(); @@ -1347,6 +1429,15 @@ struct MftReassociationValidation { registry.fill(HIST("MC/2D/hPreciseTrackSelectionCounter"), SpecificTrackSelectionStep::HasMcParticle); auto particle = templatedTrack.template mcParticle_as(); + + if (configTask.fillDcaTemplateFit) { + if (particle.isPhysicalPrimary()) { + registry.fill(HIST("MC/2D/DcaFit/hPtVsDcaXY_prim"), templatedTrack.pt(), reassociated2dMftTrack.bestDCAXY()); + } else { + registry.fill(HIST("MC/2D/DcaFit/hPtVsDcaXY_sec"), templatedTrack.pt(), reassociated2dMftTrack.bestDCAXY()); + } + } + float deltaX = -999.f; float deltaY = -999.f; float deltaZ = -999.f; @@ -1734,6 +1825,25 @@ struct MftReassociationValidation { registry.fill(HIST("MC/3D/hPreciseTrackSelectionCounter"), SpecificTrackSelectionStep::HasMcParticle); auto particle = templatedTrack.template mcParticle_as(); + + if (configTask.fillDcaTemplateFit) { + if (particle.isPhysicalPrimary()) { + registry.fill(HIST("MC/3D/DcaFit/hPtVsDcaXY_prim"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("MC/3D/DcaFit/hPtVsDcaZ_prim"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + if (reassociated3dMftTrack.ambDegree() == 1) { + registry.fill(HIST("MC/3D/DcaFit/NonAmb/hPtVsDcaXY_prim"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("MC/3D/DcaFit/NonAmb/hPtVsDcaZ_prim"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + } + } else { + registry.fill(HIST("MC/3D/DcaFit/hPtVsDcaXY_sec"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("MC/3D/DcaFit/hPtVsDcaZ_sec"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + if (reassociated3dMftTrack.ambDegree() == 1) { + registry.fill(HIST("MC/3D/DcaFit/NonAmb/hPtVsDcaXY_sec"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAXY()); + registry.fill(HIST("MC/3D/DcaFit/NonAmb/hPtVsDcaZ_sec"), templatedTrack.pt(), reassociated3dMftTrack.bestDCAZ()); + } + } + } + float deltaX = -999.f; float deltaY = -999.f; float deltaZ = -999.f;