From 11a32367cc82a6af6c987e789be0d3867a1aa3f0 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Wed, 22 Jan 2025 19:03:37 +0100 Subject: [PATCH 01/11] Add ITS PID for f1 and double-phi table producer --- PWGLF/TableProducer/Resonances/doublephitable.cxx | 8 ++++++++ PWGLF/TableProducer/Resonances/f1protonreducedtable.cxx | 8 +++++--- 2 files changed, 13 insertions(+), 3 deletions(-) diff --git a/PWGLF/TableProducer/Resonances/doublephitable.cxx b/PWGLF/TableProducer/Resonances/doublephitable.cxx index ae4e75a26a5..d02ab3d3734 100644 --- a/PWGLF/TableProducer/Resonances/doublephitable.cxx +++ b/PWGLF/TableProducer/Resonances/doublephitable.cxx @@ -45,6 +45,7 @@ #include "DataFormatsTPC/BetheBlochAleph.h" #include "CCDB/BasicCCDBManager.h" #include "CCDB/CcdbApi.h" +#include "Common/DataModel/PIDResponseITS.h" using namespace o2; using namespace o2::framework; @@ -143,6 +144,7 @@ struct doublephitable { ROOT::Math::PxPyPzMVector KaonPlus, KaonMinus, PhiMesonMother, PhiVectorDummy, Phid1dummy, Phid2dummy; void processPhiReducedTable(EventCandidates::iterator const& collision, TrackCandidates const&, aod::BCsWithTimestamps const&) { + o2::aod::ITSResponse itsResponse; bool keepEventDoublePhi = false; int numberPhi = 0; auto currentRunNumber = collision.bc_as().runNumber(); @@ -179,6 +181,9 @@ struct doublephitable { if (!selectionPID(track1)) { continue; } + if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + continue; + } Npostrack = Npostrack + 1; qaRegistry.fill(HIST("hNsigmaPtkaonTPC"), track1.tpcNSigmaKa(), track1.pt()); if (track1.hasTOF()) { @@ -201,6 +206,9 @@ struct doublephitable { if (track2ID == track1ID) { continue; } + if (!(itsResponse.nSigmaITS(track2) > -3.0 && itsResponse.nSigmaITS(track2) < 3.0)) { + continue; + } if (!selectionPair(track1, track2)) { continue; } diff --git a/PWGLF/TableProducer/Resonances/f1protonreducedtable.cxx b/PWGLF/TableProducer/Resonances/f1protonreducedtable.cxx index 7e105ab4a4d..ce0cdcf3af5 100644 --- a/PWGLF/TableProducer/Resonances/f1protonreducedtable.cxx +++ b/PWGLF/TableProducer/Resonances/f1protonreducedtable.cxx @@ -46,6 +46,7 @@ #include "DataFormatsTPC/BetheBlochAleph.h" #include "CCDB/BasicCCDBManager.h" #include "CCDB/CcdbApi.h" +#include "Common/DataModel/PIDResponseITS.h" using namespace o2; using namespace o2::framework; @@ -502,6 +503,7 @@ struct f1protonreducedtable { void processF1ProtonReducedTable(EventCandidates::iterator const& collision, aod::BCsWithTimestamps const&, PrimaryTrackCandidates const& tracks, ResoV0s const& V0s) { + o2::aod::ITSResponse itsResponse; bool keepEventF1Proton = false; int numberF1 = 0; // keep track of indices @@ -629,7 +631,7 @@ struct f1protonreducedtable { nTPCSigmaN[0] = updatePID(track, bgScalingPion, BBAntipion); } - if ((track.sign() > 0 && SelectionPID(track, strategyPIDPion, 0, nTPCSigmaP[0])) || (track.sign() < 0 && SelectionPID(track, strategyPIDPion, 0, nTPCSigmaN[0]))) { + if ((track.sign() > 0 && SelectionPID(track, strategyPIDPion, 0, nTPCSigmaP[0]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0)) || (track.sign() < 0 && SelectionPID(track, strategyPIDPion, 0, nTPCSigmaN[0]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0))) { ROOT::Math::PtEtaPhiMVector temp(track.pt(), track.eta(), track.phi(), massPi); pions.push_back(temp); PionIndex.push_back(track.globalIndex()); @@ -650,7 +652,7 @@ struct f1protonreducedtable { PionTOFHit.push_back(PionTOF); } - if ((track.pt() > cMinKaonPt && track.sign() > 0 && SelectionPID(track, strategyPIDKaon, 1, nTPCSigmaP[1])) || (track.pt() > cMinKaonPt && track.sign() < 0 && SelectionPID(track, strategyPIDKaon, 1, nTPCSigmaN[1]))) { + if ((track.pt() > cMinKaonPt && track.sign() > 0 && SelectionPID(track, strategyPIDKaon, 1, nTPCSigmaP[1]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0)) || (track.pt() > cMinKaonPt && track.sign() < 0 && SelectionPID(track, strategyPIDKaon, 1, nTPCSigmaN[1]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0))) { ROOT::Math::PtEtaPhiMVector temp(track.pt(), track.eta(), track.phi(), massKa); kaons.push_back(temp); KaonIndex.push_back(track.globalIndex()); @@ -673,7 +675,7 @@ struct f1protonreducedtable { KaonTOFHit.push_back(KaonTOF); } - if ((track.pt() < cMaxProtonPt && track.sign() > 0 && SelectionPID(track, strategyPIDProton, 2, nTPCSigmaP[2])) || (track.pt() < cMaxProtonPt && track.sign() < 0 && SelectionPID(track, strategyPIDProton, 2, nTPCSigmaN[2]))) { + if ((track.pt() < cMaxProtonPt && track.sign() > 0 && SelectionPID(track, strategyPIDProton, 2, nTPCSigmaP[2]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0)) || (track.pt() < cMaxProtonPt && track.sign() < 0 && SelectionPID(track, strategyPIDProton, 2, nTPCSigmaN[2]) && (itsResponse.nSigmaITS(track) > -3.0 && itsResponse.nSigmaITS(track) < 3.0))) { ROOT::Math::PtEtaPhiMVector temp(track.pt(), track.eta(), track.phi(), massPr); protons.push_back(temp); ProtonIndex.push_back(track.globalIndex()); From 28a8f2a54bfd7172830141495556bb8216843b3a Mon Sep 17 00:00:00 2001 From: skundu692 Date: Sat, 25 Jan 2025 17:12:19 +0100 Subject: [PATCH 02/11] Code clean up and new process function --- PWGLF/Tasks/Resonances/CMakeLists.txt | 4 +- PWGLF/Tasks/Resonances/highmasslambda.cxx | 663 +++++++++++----------- 2 files changed, 323 insertions(+), 344 deletions(-) diff --git a/PWGLF/Tasks/Resonances/CMakeLists.txt b/PWGLF/Tasks/Resonances/CMakeLists.txt index 3982955c885..4f00934dba0 100644 --- a/PWGLF/Tasks/Resonances/CMakeLists.txt +++ b/PWGLF/Tasks/Resonances/CMakeLists.txt @@ -151,7 +151,7 @@ o2physics_add_dpl_workflow(kaonkaonanalysis o2physics_add_dpl_workflow(highmasslambda SOURCES highmasslambda.cxx - PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore + PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DCAFitter COMPONENT_NAME Analysis) o2physics_add_dpl_workflow(lambdav2 @@ -187,4 +187,4 @@ o2physics_add_dpl_workflow(kshortlambda o2physics_add_dpl_workflow(rho770analysis SOURCES rho770analysis.cxx PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore - COMPONENT_NAME Analysis) \ No newline at end of file + COMPONENT_NAME Analysis) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 73be952f278..4a6fec731e4 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -55,57 +55,52 @@ #include "CCDB/BasicCCDBManager.h" #include "PWGLF/DataModel/LFStrangenessTables.h" #include "Common/DataModel/PIDResponseITS.h" +// #include "PWGHF/Utils/utilsBfieldCCDB.h" +#include "PWGHF/Utils/utilsEvSelHf.h" +#include "PWGHF/Utils/utilsTrkCandHf.h" +#include "ReconstructionDataFormats/DCA.h" +#include "ReconstructionDataFormats/V0.h" +#include "DCAFitter/DCAFitterN.h" using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; using std::array; struct highmasslambda { - - int mRunNumber; int multEstimator; - float d_bz; Service ccdb; Service pdg; - // CCDB options - // Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; - // Configurable grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"}; - // Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; - // Configurable lutPath{"lutPath", "GLO/Param/MatLUT", "Path of the Lut parametrization"}; - // Configurable geoPath{"geoPath", "GLO/Config/GeometryAligned", "Path of the geometry file"}; + o2::base::Propagator::MatCorrType noMatCorr = o2::base::Propagator::MatCorrType::USEMatCorrNONE; + o2::vertexing::DCAFitterN<2> df; + int runNumber{0}; + double bz{0.}; + + Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; + Configurable ccdbPathLut{"ccdbPathLut", "GLO/Param/MatLUT", "Path for LUT parametrization"}; + Configurable ccdbPathGrp{"ccdbPathGrp", "GLO/GRP/GRP", "Path of the grp file (Run 2)"}; + Configurable ccdbPathGrpMag{"ccdbPathGrpMag", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object (Run 3)"}; + Configurable isRun2{"isRun2", false, "enable Run 2 or Run 3 GRP objects for magnetic field"}; + Configurable cnfabsdca{"cnfabsdca", false, "Use Abs DCA for secondary vertex fitting"}; // fill output - Configurable fillQA{"fillQA", false, "fillQA"}; - Configurable useSP{"useSP", false, "useSP"}; - // Configurable additionalEvselITS{"additionalEvselITS", true, "additionalEvselITS"}; - Configurable useSignDCAV0{"useSignDCAV0", true, "useSignDCAV0"}; - Configurable fillDefault{"fillDefault", false, "fill Occupancy"}; Configurable cfgOccupancyCut{"cfgOccupancyCut", 2500, "Occupancy cut"}; - Configurable fillDecayLength{"fillDecayLength", true, "fill decay length"}; - Configurable fillPolarization{"fillPolarization", false, "fill polarization"}; Configurable fillRotation{"fillRotation", false, "fill rotation"}; // events Configurable cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"}; Configurable cfgCutCentralityMax{"cfgCutCentralityMax", 50.0f, "Accepted maximum Centrality"}; Configurable cfgCutCentralityMin{"cfgCutCentralityMin", 30.0f, "Accepted minimum Centrality"}; + Configurable additionalEvSel{"additionalEvSel", true, "additionalEvSel"}; // proton track cut - Configurable rejectPID{"rejectPID", true, "pion, kaon, electron rejection"}; - Configurable cfgCutTOFBeta{"cfgCutTOFBeta", 0.0, "cut TOF beta"}; - Configurable ispTdifferentialDCA{"ispTdifferentialDCA", true, "is pT differential DCA"}; - Configurable useGlobalTrack{"useGlobalTrack", true, "use global track"}; Configurable confMinRot{"confMinRot", 5.0 * TMath::Pi() / 6.0, "Minimum of rotation"}; Configurable confMaxRot{"confMaxRot", 7.0 * TMath::Pi() / 6.0, "Maximum of rotation"}; Configurable confRapidity{"confRapidity", 0.8, "cut on Rapidity"}; - Configurable cfgCutPT{"cfgCutPT", 0.3, "PT cut on daughter track"}; + Configurable cfgCutPT{"cfgCutPT", 0.4, "PT cut on daughter track"}; Configurable cfgCutEta{"cfgCutEta", 0.8, "Eta cut on daughter track"}; Configurable cfgCutDCAxymin1{"cfgCutDCAxymin1", 0.005f, "Minimum DCAxy range for tracks pt 0 to 0.5"}; Configurable cfgCutDCAxymin2{"cfgCutDCAxymin2", 0.003f, "Minimum DCAxy range for tracks pt 0.5 to 1"}; Configurable cfgCutDCAxymin3{"cfgCutDCAxymin3", 0.003f, "Minimum DCAxy range for tracks pt 1.0 to 1.5"}; Configurable cfgCutDCAxymin4{"cfgCutDCAxymin4", 0.002f, "Minimum DCAxy range for tracks pt 1.5 to 2.0"}; - Configurable cfgCutDCAxymin5{"cfgCutDCAxymin5", 0.001f, "Minimum DCAxy range for tracks pt 2.0 to 2.5"}; - Configurable cfgCutDCAxymin6{"cfgCutDCAxymin6", 0.0003f, "Minimum DCAxy range for tracks pt 2.5 to 3.0"}; - Configurable cfgCutDCAxymin7{"cfgCutDCAxymin7", 0.0003f, "Minimum DCAxy range for tracks pt 3.0 to 4.0"}; - Configurable cfgCutDCAxymin8{"cfgCutDCAxymin8", 0.0003f, "Minimum DCAxy range for tracks pt 4.0 to 10.0"}; + Configurable cfgCutDCAxymin5{"cfgCutDCAxymin5", 0.001f, "Minimum DCAxy range for tracks pt 2.0 to 1000.5"}; Configurable cfgCutDCAxy{"cfgCutDCAxy", 0.1f, "DCAxy range for tracks"}; Configurable cfgCutDCAz{"cfgCutDCAz", 1.0f, "DCAz range for tracks"}; Configurable cfgITScluster{"cfgITScluster", 5, "Number of ITS cluster"}; @@ -113,9 +108,6 @@ struct highmasslambda { Configurable PIDstrategy{"PIDstrategy", 0, "0: TOF Veto, 1: TOF Veto opti, 2: TOF, 3: TOF loose 1, 4: TOF loose 2, 5: old pt dep"}; Configurable nsigmaCutTPC{"nsigmacutTPC", 3.0, "Value of the TPC Nsigma cut"}; Configurable nsigmaCutTOF{"nsigmaCutTOF", 3.0, "TOF PID"}; - Configurable nsigmaCutTPCPre{"nsigmacutTPCPre", 3.0, "Value of the TPC Nsigma cut Pre filter"}; - Configurable minnsigmaCutTPCPre{"minnsigmacutTPCPre", -2.0, "Minimum Value of the TPC Nsigma cut Pre filter"}; - Configurable kaonrejpar{"kaonrejpar", 1.0, "Kaon rej. par"}; // Configs for V0 Configurable ConfV0PtMin{"ConfV0PtMin", 0.f, "Minimum transverse momentum of V0"}; Configurable ConfV0DCADaughMax{"ConfV0DCADaughMax", 0.2f, "Maximum DCA between the V0 daughters"}; @@ -128,136 +120,103 @@ struct highmasslambda { Configurable cMinLambdaMass{"cMinLambdaMass", 2.18, "Minimum lambda mass"}; Configurable cMaxLambdaMass{"cMaxLambdaMass", 2.42, "Maximum lambda mass"}; // config for V0 daughters - Configurable ConfDaughEta{"ConfDaughEta", 0.8f, "V0 Daugh sel: max eta"}; - Configurable ConfDaughPt{"ConfDaughPt", 0.1f, "V0 Daugh sel: min pt"}; - Configurable ConfDaughTPCnclsMin{"ConfDaughTPCnclsMin", 50.f, "V0 Daugh sel: Min. nCls TPC"}; Configurable ConfDaughDCAMin{"ConfDaughDCAMin", 0.08f, "V0 Daugh sel: Max. DCA Daugh to PV (cm)"}; Configurable ConfDaughPIDCuts{"ConfDaughPIDCuts", 3, "PID selections for KS0 daughters"}; - // Fill strategy - // Configurable cfgSelectDaughterTopology{"cfgSelectDaughterTopology", 2, "Select daughter for topology"}; + // config SVx + Configurable ConfMaxDecayLength{"ConfMaxDecayLength", 0.1f, "Maximum decay length (cm)"}; + Configurable ConfMinCPA{"ConfMinCPA", 0.9f, "Minimum CPA"}; // Mixed event Configurable cfgNoMixedEvents{"cfgNoMixedEvents", 1, "Number of mixed events per event"}; /// activate rotational background Configurable nBkgRotations{"nBkgRotations", 9, "Number of rotated copies (background) per each original candidate"}; - // THnsparse bining ConfigurableAxis configThnAxisInvMass{"configThnAxisInvMass", {60, 2.15, 2.45}, "#it{M} (GeV/#it{c}^{2})"}; - ConfigurableAxis configThnAxisPt{"configThnAxisPt", {24, 1.0, 25.}, "#it{p}_{T} (GeV/#it{c})"}; - ConfigurableAxis configThnAxisCosThetaStar{"configThnAxisCosThetaStar", {10, -1.0, 1.}, "cos(#vartheta)"}; - ConfigurableAxis configThnAxisCentrality{"configThnAxisCentrality", {1, 30., 50}, "Centrality"}; - ConfigurableAxis configThnAxisPhiminusPsi{"configThnAxisPhiminusPsi", {6, 0.0, TMath::Pi()}, "#phi - #psi"}; + ConfigurableAxis configThnAxisPt{"configThnAxisPt", {5, 1.0, 6.}, "#it{p}_{T} (GeV/#it{c})"}; ConfigurableAxis configThnAxisV2{"configThnAxisV2", {80, -1, 1}, "V2"}; - ConfigurableAxis configThnAxisSA{"configThnAxisSA", {100, -1, 1}, "SA"}; - ConfigurableAxis cnfigThnAxisDecayLength{"cnfigThnAxisDecayLength", {150, 0.0, 0.3}, "SA"}; - ConfigurableAxis cnfigThnAxisDCASum{"cnfigThnAxisDCASum", {150, -0.3, 0.3}, "SA"}; + ConfigurableAxis cnfigThnAxisDCA{"cnfigThnAxisDCA", {100, 0.0, 0.1}, "DCA"}; + ConfigurableAxis cnfigThnAxisDecayLength{"cnfigThnAxisDecayLength", {150, 0.0, 0.3}, "decay length"}; ConfigurableAxis cnfigThnAxisPtProton{"cnfigThnAxisPtProton", {16, 0.0, 8.0}, "pT"}; - ConfigurableAxis configThnAxisSP{"configThnAxisSP", {400, -12, 12}, "SP"}; + ConfigurableAxis cnfigThnAxisCPA{"cnfigThnAxisCPA", {300, 0.8, 1.1}, "CPA"}; + // ConfigurableAxis configThnAxisCosThetaStar{"configThnAxisCosThetaStar", {10, -1.0, 1.}, "cos(#vartheta)"}; + // ConfigurableAxis configThnAxisSA{"configThnAxisSA", {100, -1, 1}, "SA"}; Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex; Filter centralityFilter = (nabs(aod::cent::centFT0C) < cfgCutCentralityMax && nabs(aod::cent::centFT0C) > cfgCutCentralityMin); Filter acceptanceFilter = (nabs(aod::track::eta) < cfgCutEta && nabs(aod::track::pt) > cfgCutPT); Filter dcaCutFilter = (nabs(aod::track::dcaXY) < cfgCutDCAxy) && (nabs(aod::track::dcaZ) < cfgCutDCAz); - Filter pidFilter = aod::pidtpc::tpcNSigmaPr > minnsigmaCutTPCPre&& aod::pidtpc::tpcNSigmaPr < nsigmaCutTPCPre; + Filter pidFilter = nabs(aod::pidtpc::tpcNSigmaPr) < nsigmaCutTPC; using EventCandidates = soa::Filtered>; - using TrackCandidates = soa::Filtered>; + + using TrackCandidates = soa::Filtered>; using AllTrackCandidates = soa::Join; using ResoV0s = aod::V0Datas; + using TrackCandidatesSvx = soa::Filtered>; + using AllTrackCandidatesSvx = soa::Join; + using ResoV0sSvx = soa::Join; + SliceCache cache; - // Partition posTracks = aod::track::signed1Pt > cfgCutCharge; - // Partition negTracks = aod::track::signed1Pt < cfgCutCharge; HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; void init(o2::framework::InitContext&) { - std::vector occupancyBinning = {-0.5, 500.0, 1000.0, 1500.0, 2000.0, 3000.0, 4000.0, 5000.0, 50000.0}; + // std::vector dcaBinning = {0.0, 0.0005, 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009, 0.01, 0.012, 0.014, 0.016, 0.02, 0.03, 0.05, 0.1, 0.5, 1.0}; - std::vector dcaBinning = {0.0, 0.0005, 0.001, 0.0012, 0.0014, 0.0016, 0.002, 0.0025, 0.003, 0.004, 0.005, 0.006, 0.008, 0.01, 0.015, 0.02, 0.04, 0.05, 0.06, 0.08, 0.1, 0.3, 1.0}; - std::vector ptProtonBinning = {0.2, 0.3, 0.5, 0.6, 0.8, 1.2, 1.4, 1.6, 2.0, 3.0, 4.0, 6.0}; - std::vector ptLambdaBinning = {2.0, 3.0, 4.0, 5.0, 6.0}; + // std::vector dcaBinning = {0.0, 0.0005, 0.001, 0.0012, 0.0014, 0.0016, 0.002, 0.0025, 0.003, 0.004, 0.005, 0.006, 0.008, 0.01, 0.015, 0.02, 0.04, 0.05, 0.06, 0.08, 0.1, 0.3, 1.0}; + // std::vector ptProtonBinning = {0.2, 0.3, 0.5, 0.6, 0.8, 1.2, 1.4, 1.6, 2.0, 3.0, 4.0, 6.0}; + // std::vector ptLambdaBinning = {2.0, 3.0, 4.0, 5.0, 6.0}; + std::vector occupancyBinning = {-0.5, 500.0, 1000.0, 1500.0, 2000.0, 3000.0, 4000.0, 5000.0, 50000.0}; + AxisSpec resAxis = {1600, -16, 16, "Res"}; + AxisSpec phiAxis = {500, -6.28, 6.28, "phi"}; + AxisSpec centAxis = {8, 0, 80, "V0M (%)"}; const AxisSpec thnAxisInvMass{configThnAxisInvMass, "#it{M} (GeV/#it{c}^{2})"}; const AxisSpec thnAxisPt{configThnAxisPt, "#it{p}_{T} (GeV/#it{c})"}; - const AxisSpec thnAxisCosThetaStar{configThnAxisCosThetaStar, "cos(#vartheta)"}; - const AxisSpec thnAxisPhiminusPsi{configThnAxisPhiminusPsi, "#phi - #psi"}; - const AxisSpec thnAxisCentrality{configThnAxisCentrality, "Centrality (%)"}; const AxisSpec thnAxisV2{configThnAxisV2, "V2"}; - const AxisSpec thnAxisSA{configThnAxisSA, "SA"}; + const AxisSpec thnAxisDCA{cnfigThnAxisDCA, "DCAxy"}; const AxisSpec thnAxisDecayLength{cnfigThnAxisDecayLength, "Decay Length"}; - const AxisSpec thnAxisDCASum{cnfigThnAxisDCASum, "DCA Sum"}; const AxisSpec thnAxisPtProton{cnfigThnAxisPtProton, "Proton Pt"}; - const AxisSpec thnAxisSP{configThnAxisSP, "SP"}; - - AxisSpec phiAxis = {500, -6.28, 6.28, "phi"}; - AxisSpec resAxis = {1000, -10, 10, "Res"}; - AxisSpec centAxis = {8, 0, 80, "V0M (%)"}; - // AxisSpec ptProtonAxis = {16, 0.0, 8, "V0M (%)"}; - AxisSpec dcaAxis = {dcaBinning, "DCAxy"}; - AxisSpec dcatoPVAxis = {50, 0.0, 0.5, "V0M (%)"}; + const AxisSpec thnAxisCPA{cnfigThnAxisCPA, "CPA"}; AxisSpec occupancyAxis = {occupancyBinning, "occupancy"}; - AxisSpec ptProtonAxis = {ptProtonBinning, "pT proton"}; - histos.add("hMomCorr", "hMomCorr", kTH3F, {{200, -10.0f, 10.0f}, {200, -10.0f, 10.0f}, {8, 0.0f, 80.0f}}); + // const AxisSpec thnAxisCosThetaStar{configThnAxisCosThetaStar, "cos(#vartheta)"}; + // const AxisSpec thnAxisPhiminusPsi{configThnAxisPhiminusPsi, "#phi - #psi"}; + // const AxisSpec thnAxisCentrality{configThnAxisCentrality, "Centrality (%)"}; + // const AxisSpec thnAxisSA{configThnAxisSA, "SA"}; + + histos.add("hMomCorr", "hMomCorr", kTH2F, {{200, -10.0f, 10.0f}, {200, -10.0f, 10.0f}}); histos.add("hInvMassKs0", "hInvMassKs0", kTH1F, {{200, 0.4f, 0.6f}}); histos.add("hV0Dca", "hV0Dca", kTH1F, {{2000, -1.0f, 1.0f}}); histos.add("hpTvsRapidity", "pT vs Rapidity", kTH2F, {{100, 0.0f, 10.0f}, {300, -1.5f, 1.5f}}); - histos.add("hFTOCvsTPCNoCut", "Mult correlation FT0C vs. TPC without any cut", kTH2F, {{80, 0.0f, 80.0f}, {100, -0.5f, 5999.5f}}); - histos.add("hFTOCvsTPC", "Mult correlation FT0C vs. TPC", kTH2F, {{80, 0.0f, 80.0f}, {100, -0.5f, 5999.5f}}); - histos.add("hFTOCvsTPCSelected", "Mult correlation FT0C vs. TPC after selection", kTH2F, {{80, 0.0f, 80.0f}, {100, -0.5f, 5999.5f}}); histos.add("hCentrality", "Centrality distribution", kTH1F, {{200, 0.0, 200.0}}); histos.add("hVtxZ", "Vertex distribution in Z;Z (cm)", kTH1F, {{400, -20.0, 20.0}}); - histos.add("hOccupancy", "Occupancy", kTH1F, {{5000, 0.0, 50000.0}}); + histos.add("hOccupancy", "Occupancy", kTH1F, {{5000, -0.5, 50000.5}}); histos.add("hRotation", "hRotation", kTH1F, {{360, 0.0, 2.0 * TMath::Pi()}}); histos.add("hEta", "Eta distribution", kTH1F, {{200, -1.0f, 1.0f}}); histos.add("hDcaxy", "Dcaxy distribution", kTH1F, {{1000, -0.5f, 0.5f}}); histos.add("hDcaz", "Dcaz distribution", kTH1F, {{1000, -0.5f, 0.5f}}); - histos.add("hNsigmaProtonTPCDiff", "Difference NsigmaProton NsigmaKaon TPC distribution", kTH3F, {{100, -5.0f, 5.0f}, {100, -5.0f, 5.0f}, {80, 0.0f, 8.0f}}); - - histos.add("hNsigmaProtonElectronTPC", "NsigmaProton-Electron TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonPionTPC", "NsigmaProton-Pion TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonKaonTPC", "NsigmaProton-Kaon TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - - histos.add("hNsigmaProtonElectronTPC_afterKa", "NsigmaProton-Electron TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonPionTPC_afterKa", "NsigmaProton-Pion TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonKaonTPC_afterKa", "NsigmaProton-Kaon TPC distribution", kTH3F, {{60, -3.0f, 3.0f}, {200, -10.0f, 10.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonTPC", "NsigmaProton TPC distribution", kTH2F, {{100, -5.0f, 5.0f}, {60, 0.0f, 6.0f}}); histos.add("hNsigmaProtonTOF", "NsigmaProton TOF distribution", kTH2F, {{1000, -50.0f, 50.0f}, {60, 0.0f, 6.0f}}); - histos.add("hNsigmaProtonTOFPre", "NsigmaProton TOF distribution Pre sel", kTH2F, {{1000, -50.0f, 50.0f}, {60, 0.0f, 6.0f}}); - histos.add("hMassvsDecaySum", "hMassvsDecaySum", kTH2F, {thnAxisInvMass, thnAxisDCASum}); + histos.add("hNsigmaProtonTPCPre", "NsigmaProton TPC distribution Pre sel", kTH2F, {{1000, -50.0f, 50.0f}, {60, 0.0f, 6.0f}}); histos.add("hPsiFT0C", "PsiFT0C", kTH3F, {centAxis, phiAxis, occupancyAxis}); histos.add("hPsiFT0A", "PsiFT0A", kTH3F, {centAxis, phiAxis, occupancyAxis}); histos.add("hPsiTPC", "PsiTPC", kTH3F, {centAxis, phiAxis, occupancyAxis}); - histos.add("hPsiTPCR", "PsiTPCR", kTH3F, {centAxis, phiAxis, occupancyAxis}); - histos.add("hPsiTPCL", "PsiTPCL", kTH3F, {centAxis, phiAxis, occupancyAxis}); - - if (fillDefault) { - histos.add("hSparseV2SASameEvent_V2", "hSparseV2SASameEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcaAxis, ptProtonAxis}); - histos.add("hSparseV2SASameEventRotational_V2", "hSparseV2SASameEventRotational", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcaAxis, ptProtonAxis}); - histos.add("hSparseV2SAMixedEvent_V2", "hSparseV2SAMixedEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcaAxis, ptProtonAxis}); - } - if (fillDecayLength) { - histos.add("hSparseV2SASameEvent_V2_dcatopv", "hSparseV2SASameEvent_V2_dcatopv", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcatoPVAxis, ptProtonAxis}); - histos.add("hSparseV2SASameEventRotational_V2_dcatopv", "hSparseV2SASameEventRotational_V2_dcatopv", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcatoPVAxis, ptProtonAxis}); - histos.add("hSparseV2SAMixedEvent_V2_dcatopv", "hSparseV2SAMixedEvent_V2_dcatopv", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSP, dcatoPVAxis, ptProtonAxis}); - } - if (fillPolarization) { - histos.add("hSparseV2SASameEventplus_SA", "hSparseV2SASameEventplus_SA", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SASameEventplus_SA_A0", "hSparseV2SASameEventplus_SA_A0", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SASameEventplus_SA_azimuth", "hSparseV2SASameEventplus_SA_azimuth", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisCentrality}); - histos.add("hSparseV2SASameEventminus_SA", "hSparseV2SASameEventminus_SA", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SASameEventminus_SA_A0", "hSparseV2SASameEventminus_SA_A0", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SASameEventminus_SA_azimuth", "hSparseV2SASameEventminus_SA_azimuth", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisCentrality}); - - histos.add("hSparseV2SAMixedEventplus_SA", "hSparseV2SAMixedEventplus_SA", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SAMixedEventplus_SA_A0", "hSparseV2SAMixedEventplus_SA_A0", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SAMixedEventplus_SA_azimuth", "hSparseV2SAMixedEventplus_SA_azimuth", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisCentrality}); - histos.add("hSparseV2SAMixedEventminus_SA", "hSparseV2SAMixedEventminus_SA", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SAMixedEventminus_SA_A0", "hSparseV2SAMixedEventminus_SA_A0", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisPhiminusPsi, thnAxisCentrality}); - histos.add("hSparseV2SAMixedEventminus_SA_azimuth", "hSparseV2SAMixedEventminus_SA_azimuth", HistType::kTHnSparseF, {thnAxisInvMass, ptLambdaBinning, thnAxisSA, thnAxisCentrality}); - } + + // SVX histo + histos.add("hDecayLengthxy", "Decay length xy", kTH1F, {{500, 0.0f, 0.1f}}); + histos.add("hDecayLength", "Decay length", kTH1F, {{500, 0.0f, 0.1f}}); + histos.add("hImpactPar0", "hImpactPar0", kTH1F, {{500, 0.0f, 0.1f}}); + histos.add("hImpactPar1", "hImpactPar1", kTH1F, {{500, 0.0f, 0.1f}}); + histos.add("hCPA", "hCPA", kTH1F, {{220, -1.1f, 1.1f}}); + + histos.add("hSparseV2SASameEvent_V2", "hSparseV2SASameEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); + histos.add("hSparseV2SASameEventRotational_V2", "hSparseV2SASameEventRotational", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); + histos.add("hSparseV2SAMixedEvent_V2", "hSparseV2SAMixedEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); + + histos.add("hSparseV2SASameEvent_V2_SVX", "hSparseV2SASameEvent_V2_SVX", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDecayLength, thnAxisCPA}); + histos.add("hSparseV2SASameEventRotational_V2_SVX", "hSparseV2SASameEventRotational_SVX", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDecayLength, thnAxisCPA}); // histogram for resolution histos.add("ResFT0CTPC", "ResFT0CTPC", kTH2F, {centAxis, resAxis}); @@ -266,121 +225,41 @@ struct highmasslambda { histos.add("ResTPCRTPCL", "ResTPCRTPCL", kTH2F, {centAxis, resAxis}); histos.add("ResFT0CFT0A", "ResFT0CFT0A", kTH2F, {centAxis, resAxis}); histos.add("ResFT0ATPC", "ResFT0ATPC", kTH2F, {centAxis, resAxis}); - } - template - bool selectionTrack(const T& candidate) - { - if (useGlobalTrack && !(candidate.isGlobalTrackWoDCA() && candidate.itsNCls() > cfgITScluster && candidate.tpcNClsFound() > cfgTPCcluster)) { - return false; - } - if (!useGlobalTrack && !(candidate.itsNCls() > cfgITScluster && candidate.tpcNClsFound() > cfgTPCcluster)) { - return false; - } - if (ispTdifferentialDCA) { - if (candidate.pt() > 0.0 && candidate.pt() < 0.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin1) { - return false; - } - if (candidate.pt() >= 0.5 && candidate.pt() < 1.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin2) { - return false; - } - if (candidate.pt() >= 1.0 && candidate.pt() < 1.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin3) { - return false; - } - if (candidate.pt() >= 1.5 && candidate.pt() < 2.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin4) { - return false; - } - if (candidate.pt() >= 2.0 && candidate.pt() < 2.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin5) { - return false; - } - if (candidate.pt() >= 2.5 && candidate.pt() < 3.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin6) { - return false; - } - if (candidate.pt() >= 3.0 && candidate.pt() < 4.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin7) { - return false; - } - if (candidate.pt() >= 4.0 && candidate.pt() < 10.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin8) { - return false; - } - } - - if (!ispTdifferentialDCA) { - if (candidate.pt() > 0.0 && candidate.pt() < 0.5 && TMath::Abs(candidate.dcaXY()) < 0.004) { - return false; - } - if (candidate.pt() >= 0.5 && candidate.pt() < 1.0 && TMath::Abs(candidate.dcaXY()) < 0.003) { - return false; - } - if (candidate.pt() >= 1.0 && candidate.pt() < 2.0 && TMath::Abs(candidate.dcaXY()) < 0.002) { - return false; - } - } - return true; + df.setPropagateToPCA(true); + df.setMaxR(200); + df.setMaxDZIni(4); + df.setMinParamChange(1.e-3); + df.setMinRelChi2Change(0.9); + df.setUseAbsDCA(cnfabsdca); + df.setWeightedFinalPCA(cnfabsdca); + df.setMatCorrType(noMatCorr); + + ccdb->setURL(ccdburl); + ccdb->setCaching(true); + ccdb->setLocalObjectValidityChecking(); + runNumber = 0; + bz = 0; } - template - bool rejectPi(const T& candidate) - { - if (candidate.tpcInnerParam() > 0.9 && candidate.tpcInnerParam() < 1.0 && candidate.tpcNSigmaPi() < 6.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.0 && candidate.tpcInnerParam() < 1.1 && candidate.tpcNSigmaPi() < 4.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.1 && candidate.tpcInnerParam() < 1.2 && candidate.tpcNSigmaPi() < 3.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.2 && candidate.tpcInnerParam() < 1.4 && candidate.tpcNSigmaPi() < 1.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.4 && candidate.tpcInnerParam() < 1.5 && candidate.tpcNSigmaPi() < 0.5) { - return false; - } - return true; - } - - template - bool rejectEl(const T& candidate) - { - - if (candidate.tpcInnerParam() > 0.7 && candidate.tpcInnerParam() < 0.8 && candidate.tpcNSigmaEl() < 2.0) { - return false; - } - if (candidate.tpcInnerParam() > 0.8 && candidate.tpcInnerParam() < 0.9 && candidate.tpcNSigmaEl() < 0.0) { - return false; - } - if (candidate.tpcInnerParam() > 0.9 && candidate.tpcInnerParam() < 1.0 && candidate.tpcNSigmaEl() < -1.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.0 && candidate.tpcInnerParam() < 1.1 && candidate.tpcNSigmaEl() < -2.0) { - return false; - } - if (candidate.tpcInnerParam() > 1.1 && candidate.tpcInnerParam() < 1.2 && candidate.tpcNSigmaEl() < -3.0) { - return false; - } - - return true; - } - - template - bool rejectKa(const T& candidate) + bool selectionTrack(const T& candidate) { - if (candidate.tpcInnerParam() > 0.7 && candidate.tpcInnerParam() < 0.8 && candidate.tpcNSigmaKa() < 7.5) { + if (!(candidate.isGlobalTrackWoDCA() && candidate.itsNCls() > cfgITScluster && candidate.tpcNClsFound() > cfgTPCcluster)) { return false; } - if (candidate.tpcInnerParam() > 0.8 && candidate.tpcInnerParam() < 0.9 && candidate.tpcNSigmaKa() < 6.0) { + if (candidate.pt() > 0.0 && candidate.pt() < 0.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin1) { return false; } - if (candidate.tpcInnerParam() > 0.9 && candidate.tpcInnerParam() < 1.1 && candidate.tpcNSigmaKa() < 5.0) { + if (candidate.pt() >= 0.5 && candidate.pt() < 1.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin2) { return false; } - if (candidate.tpcInnerParam() > 1.1 && candidate.tpcInnerParam() < 1.2 && candidate.tpcNSigmaKa() < 3.5) { + if (candidate.pt() >= 1.0 && candidate.pt() < 1.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin3) { return false; } - if (candidate.tpcInnerParam() > 1.2 && candidate.tpcInnerParam() < 1.4 && candidate.tpcNSigmaKa() < 3.0) { + if (candidate.pt() >= 1.5 && candidate.pt() < 2.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin4) { return false; } - if (candidate.tpcInnerParam() > 1.4 && candidate.tpcInnerParam() < 1.5 && candidate.tpcNSigmaKa() < 2.5) { + if (candidate.pt() >= 2.0 && candidate.pt() < 10000000.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin5) { return false; } return true; @@ -391,12 +270,15 @@ struct highmasslambda { bool selectionPID1(const T& candidate) { if (candidate.hasTOF()) { - if (candidate.beta() > cfgCutTOFBeta && candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.p() < 0.5 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + return true; + } + if (candidate.p() >= 0.5 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } } if (!candidate.hasTOF()) { - if (candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + if (TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } } @@ -407,10 +289,10 @@ struct highmasslambda { template bool selectionPID2(const T& candidate) { - if (candidate.tpcInnerParam() < 1.0 && candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + if (candidate.tpcInnerParam() < 1.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } - if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 2.0 && candidate.beta() > cfgCutTOFBeta && candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 2.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } if (candidate.tpcInnerParam() >= 2.0 && candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC) { @@ -469,13 +351,13 @@ struct highmasslambda { if (charge > 0 && sign < 0) { return false; } - if (TMath::Abs(eta) > ConfDaughEta) { + if (TMath::Abs(eta) > 0.8) { return false; } - if (TMath::Abs(pt) < ConfDaughPt) { + if (TMath::Abs(pt) < 0.15) { return false; } - if (tpcNClsF < ConfDaughTPCnclsMin) { + if (tpcNClsF < 50) { return false; } if (TMath::Abs(dcaXY) < ConfDaughDCAMin) { @@ -506,7 +388,6 @@ struct highmasslambda { using BinningTypeVertexContributor = ColumnBinningPolicy; ROOT::Math::PxPyPzMVector Lambdac, Proton, Kshort, LambdacRot, KshortRot, fourVecDauCM; ROOT::Math::XYZVector threeVecDauCM, threeVecDauCMXY, eventplaneVec, eventplaneVecNorm, beamvector; - double massPi = TDatabasePDG::Instance()->GetParticle(kPiPlus)->Mass(); // FIXME: Get from the common header double massPr = TDatabasePDG::Instance()->GetParticle(kProton)->Mass(); // FIXME: Get from the common header double massK0s = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(); // FIXME: Get from the common header double v2, v2Rot; @@ -515,13 +396,12 @@ struct highmasslambda { if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { return; } - if (!collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + if (additionalEvSel && !collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { return; } o2::aod::ITSResponse itsResponse; auto centrality = collision.centFT0C(); auto multTPC = collision.multNTracksPV(); - histos.fill(HIST("hFTOCvsTPCNoCut"), centrality, multTPC); int occupancy = collision.trackOccupancyInTimeRange(); if (occupancy > cfgOccupancyCut) { return; @@ -537,14 +417,9 @@ struct highmasslambda { auto QTPC = collision.qTPC(); auto QTPCR = collision.qTPCR(); auto QTPCL = collision.qTPCL(); - - histos.fill(HIST("hFTOCvsTPC"), centrality, multTPC); - histos.fill(HIST("hFTOCvsTPCSelected"), centrality, multTPC); histos.fill(HIST("hPsiFT0C"), centrality, psiFT0C, occupancy); histos.fill(HIST("hPsiFT0A"), centrality, psiFT0A, occupancy); histos.fill(HIST("hPsiTPC"), centrality, psiTPC, occupancy); - histos.fill(HIST("hPsiTPCR"), centrality, psiTPCR, occupancy); - histos.fill(HIST("hPsiTPCL"), centrality, psiTPCL, occupancy); histos.fill(HIST("ResFT0CTPC"), centrality, QFT0C * QTPC * TMath::Cos(2.0 * (psiFT0C - psiTPC))); histos.fill(HIST("ResFT0CTPCR"), centrality, QFT0C * QTPCR * TMath::Cos(2.0 * (psiFT0C - psiTPCR))); histos.fill(HIST("ResFT0CTPCL"), centrality, QFT0C * QTPCL * TMath::Cos(2.0 * (psiFT0C - psiTPCL))); @@ -559,6 +434,7 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } + histos.fill(HIST("hNsigmaProtonTPCPre"), track1.tpcNSigmaPr(), track1.pt()); // PID check if (PIDstrategy == 0 && !selectionPID1(track1)) { continue; @@ -569,29 +445,9 @@ struct highmasslambda { if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { continue; } - if (track1.hasTOF()) { - histos.fill(HIST("hNsigmaProtonTOFPre"), track1.tofNSigmaPr(), track1.pt()); - } - if (!track1.hasTOF()) { - if (fillQA) { - histos.fill(HIST("hNsigmaProtonElectronTPC"), track1.tpcNSigmaPr(), track1.tpcNSigmaEl(), track1.tpcInnerParam()); - histos.fill(HIST("hNsigmaProtonPionTPC"), track1.tpcNSigmaPr(), track1.tpcNSigmaPi(), track1.tpcInnerParam()); - histos.fill(HIST("hNsigmaProtonKaonTPC"), track1.tpcNSigmaPr(), track1.tpcNSigmaKa(), track1.tpcInnerParam()); - } - if (rejectPID && !rejectKa(track1)) { - continue; - } - if (fillQA) { - histos.fill(HIST("hNsigmaProtonElectronTPC_afterKa"), track1.tpcNSigmaPr(), track1.tpcNSigmaEl(), track1.tpcInnerParam()); - histos.fill(HIST("hNsigmaProtonPionTPC_afterKa"), track1.tpcNSigmaPr(), track1.tpcNSigmaPi(), track1.tpcInnerParam()); - histos.fill(HIST("hNsigmaProtonKaonTPC_afterKa"), track1.tpcNSigmaPr(), track1.tpcNSigmaKa(), track1.tpcInnerParam()); - } - } - histos.fill(HIST("hMomCorr"), track1.p() / track1.sign(), track1.p() - track1.tpcInnerParam(), centrality); histos.fill(HIST("hEta"), track1.eta()); histos.fill(HIST("hDcaxy"), track1.dcaXY()); histos.fill(HIST("hDcaz"), track1.dcaZ()); - histos.fill(HIST("hNsigmaProtonTPCDiff"), track1.tpcNSigmaPr(), track1.tpcNSigmaKa(), track1.pt()); histos.fill(HIST("hNsigmaProtonTPC"), track1.tpcNSigmaPr(), track1.pt()); if (track1.hasTOF()) { histos.fill(HIST("hNsigmaProtonTOF"), track1.tofNSigmaPr(), track1.pt()); @@ -601,10 +457,8 @@ struct highmasslambda { if (!SelectionV0(collision, v0)) { continue; } - auto anglesign = (v0.x() - collision.posX()) * v0.px() + (v0.y() - collision.posY()) * v0.py() + (v0.z() - collision.posZ()) * v0.pz(); - anglesign = anglesign / TMath::Abs(anglesign); if (firstprimarytrack == 0) { - histos.fill(HIST("hV0Dca"), anglesign * v0.dcav0topv()); + histos.fill(HIST("hV0Dca"), v0.dcav0topv()); } auto postrack = v0.template posTrack_as(); auto negtrack = v0.template negTrack_as(); @@ -628,28 +482,9 @@ struct highmasslambda { Kshort = ROOT::Math::PxPyPzMVector(v0.px(), v0.py(), v0.pz(), v0.mK0Short()); Lambdac = Proton + Kshort; auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); - - if (useSP) { - v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - } - if (!useSP) { - v2 = TMath::Cos(2.0 * phiminuspsi); - } - auto dcasum = 0.0; - if (useSignDCAV0) { - dcasum = anglesign * (v0.dcav0topv()) - track1.dcaXY(); - } - if (!useSignDCAV0) { - dcasum = v0.dcav0topv() - track1.dcaXY(); - } - histos.fill(HIST("hMassvsDecaySum"), Lambdac.M(), dcasum); + v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - if (fillDefault) { - histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); - } - if (fillDecayLength) { - histos.fill(HIST("hSparseV2SASameEvent_V2_dcatopv"), Lambdac.M(), Lambdac.Pt(), v2, v0.dcav0topv(), Proton.Pt()); - } + histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY())); } if (fillRotation) { for (int nrotbkg = 0; nrotbkg < nBkgRotations; nrotbkg++) { @@ -663,43 +498,12 @@ struct highmasslambda { KshortRot = ROOT::Math::PxPyPzMVector(rotKaonPx, rotKaonPy, Kshort.pz(), massK0s); LambdacRot = Proton + KshortRot; auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); - if (useSP) { - v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; - } - if (!useSP) { - v2Rot = TMath::Cos(2.0 * phiminuspsiRot); - } - + v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && TMath::Abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { - if (fillDefault) { - histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, TMath::Abs(track1.dcaXY()), Proton.Pt()); - } - if (fillDecayLength) { - histos.fill(HIST("hSparseV2SASameEventRotational_V2_dcatopv"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, v0.dcav0topv(), Proton.Pt()); - } + histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, TMath::Abs(track1.dcaXY())); } } } - if (fillPolarization) { - ROOT::Math::Boost boost{Lambdac.BoostToCM()}; - fourVecDauCM = boost(Kshort); - threeVecDauCM = fourVecDauCM.Vect(); - threeVecDauCMXY = ROOT::Math::XYZVector(threeVecDauCM.X(), threeVecDauCM.Y(), 0.); - beamvector = ROOT::Math::XYZVector(0, 0, 1); - auto cosThetaStar = beamvector.Dot(threeVecDauCM) / std::sqrt(threeVecDauCM.Mag2()) / std::sqrt(beamvector.Mag2()); - auto SA = cosThetaStar * TMath::Sin(2.0 * phiminuspsi); - - if (track1.sign() > 0) { - histos.fill(HIST("hSparseV2SASameEventplus_SA"), Lambdac.M(), Lambdac.Pt(), cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SASameEventplus_SA_A0"), Lambdac.M(), Lambdac.Pt(), cosThetaStar * cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SASameEventplus_SA_azimuth"), Lambdac.M(), Lambdac.Pt(), SA, centrality); - } - if (track1.sign() < 0) { - histos.fill(HIST("hSparseV2SASameEventminus_SA"), Lambdac.M(), Lambdac.Pt(), cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SASameEventminus_SA_A0"), Lambdac.M(), Lambdac.Pt(), cosThetaStar * cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SASameEventminus_SA_azimuth"), Lambdac.M(), Lambdac.Pt(), SA, centrality); - } - } } } } @@ -719,14 +523,14 @@ struct highmasslambda { if (collision1.bcId() == collision2.bcId()) { continue; } - if (!collision1.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + if (additionalEvSel && !collision1.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { continue; } - if (!collision2.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + if (additionalEvSel && !collision2.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { continue; } o2::aod::ITSResponse itsResponse; - auto centrality = collision1.centFT0C(); + auto psiFT0C = collision1.psiFT0C(); auto QFT0C = collision1.qFT0C(); int occupancy1 = collision1.trackOccupancyInTimeRange(); @@ -753,11 +557,6 @@ struct highmasslambda { continue; } - if (!track1.hasTOF()) { - if (rejectPID && !rejectEl(track1)) { - continue; - } - } if (!SelectionV0(collision2, v0)) { continue; } @@ -779,47 +578,227 @@ struct highmasslambda { continue; } auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); - if (useSP) { - v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - } - if (!useSP) { - v2 = TMath::Cos(2.0 * phiminuspsi); + v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; + if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY())); } - auto anglesign = (v0.x() - collision1.posX()) * v0.px() + (v0.y() - collision1.posY()) * v0.py() + (v0.z() - collision1.posZ()) * v0.pz(); - anglesign = anglesign / TMath::Abs(anglesign); + } + } + } + PROCESS_SWITCH(highmasslambda, processMixedEventOpti, "Process Mixed event new", false); + void processSameEventSvx(EventCandidates::iterator const& collision, TrackCandidatesSvx const& tracks, AllTrackCandidatesSvx const&, ResoV0sSvx const& V0s, aod::BCsWithTimestamps const&) + { + if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { + return; + } + if (additionalEvSel && !collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + return; + } - if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - if (fillDefault) { - histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); - } - if (fillDecayLength) { - histos.fill(HIST("hSparseV2SAMixedEvent_V2_dcatopv"), Lambdac.M(), Lambdac.Pt(), v2, v0.dcav0topv(), Proton.Pt()); - } + /// Set the magnetic field from ccdb. + auto bc = collision.template bc_as(); + if (runNumber != bc.runNumber()) { + o2::parameters::GRPMagField* grpo = ccdb->getForTimeStamp(ccdbPathGrpMag, bc.timestamp()); + if (grpo == nullptr) { + LOGF(fatal, "Run 3 GRP object (type o2::parameters::GRPMagField) is not available in CCDB for run=%d at timestamp=%llu", bc.runNumber(), bc.timestamp()); + } + o2::base::Propagator::initFieldFromGRP(grpo); + bz = o2::base::Propagator::Instance()->getNominalBz(); + runNumber = bc.runNumber(); + } + df.setBz(bz); + int occupancy = collision.trackOccupancyInTimeRange(); + if (occupancy > cfgOccupancyCut) { + return; + } + o2::aod::ITSResponse itsResponse; + auto centrality = collision.centFT0C(); + auto multTPC = collision.multNTracksPV(); + auto psiFT0C = collision.psiFT0C(); + auto psiFT0A = collision.psiFT0A(); + auto psiTPC = collision.psiTPC(); + auto psiTPCR = collision.psiTPCR(); + auto psiTPCL = collision.psiTPCL(); + auto QFT0C = collision.qFT0C(); + auto QFT0A = collision.qFT0A(); + auto QTPC = collision.qTPC(); + auto QTPCR = collision.qTPCR(); + auto QTPCL = collision.qTPCL(); + histos.fill(HIST("hPsiFT0C"), centrality, psiFT0C, occupancy); + histos.fill(HIST("hPsiFT0A"), centrality, psiFT0A, occupancy); + histos.fill(HIST("hPsiTPC"), centrality, psiTPC, occupancy); + histos.fill(HIST("ResFT0CTPC"), centrality, QFT0C * QTPC * TMath::Cos(2.0 * (psiFT0C - psiTPC))); + histos.fill(HIST("ResFT0CTPCR"), centrality, QFT0C * QTPCR * TMath::Cos(2.0 * (psiFT0C - psiTPCR))); + histos.fill(HIST("ResFT0CTPCL"), centrality, QFT0C * QTPCL * TMath::Cos(2.0 * (psiFT0C - psiTPCL))); + histos.fill(HIST("ResTPCRTPCL"), centrality, QTPCR * QTPCL * TMath::Cos(2.0 * (psiTPCR - psiTPCL))); + histos.fill(HIST("ResFT0CFT0A"), centrality, QFT0C * QFT0A * TMath::Cos(2.0 * (psiFT0C - psiFT0A))); + histos.fill(HIST("ResFT0ATPC"), centrality, QTPC * QFT0A * TMath::Cos(2.0 * (psiTPC - psiFT0A))); + histos.fill(HIST("hCentrality"), centrality); + histos.fill(HIST("hVtxZ"), collision.posZ()); + histos.fill(HIST("hOccupancy"), occupancy); + auto firstprimarytrack = 0; + for (auto track1 : tracks) { + if (!selectionTrack(track1)) { + continue; + } + histos.fill(HIST("hNsigmaProtonTPCPre"), track1.tpcNSigmaPr(), track1.pt()); + // PID check + if (PIDstrategy == 0 && !selectionPID1(track1)) { + continue; + } + if (PIDstrategy == 1 && !selectionPID2(track1)) { + continue; + } + if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { + continue; + } + histos.fill(HIST("hEta"), track1.eta()); + histos.fill(HIST("hDcaxy"), track1.dcaXY()); + histos.fill(HIST("hDcaz"), track1.dcaZ()); + histos.fill(HIST("hNsigmaProtonTPC"), track1.tpcNSigmaPr(), track1.pt()); + if (track1.hasTOF()) { + histos.fill(HIST("hNsigmaProtonTOF"), track1.tofNSigmaPr(), track1.pt()); + } + auto track1ID = track1.globalIndex(); + auto trackParCovBach = getTrackParCov(track1); + + for (auto v0 : V0s) { + if (!SelectionV0(collision, v0)) { + continue; + } + if (firstprimarytrack == 0) { + histos.fill(HIST("hV0Dca"), v0.dcav0topv()); + } + auto postrack = v0.template posTrack_as(); + auto negtrack = v0.template negTrack_as(); + if (!isSelectedV0Daughter(postrack, 1)) { + continue; + } + if (!isSelectedV0Daughter(negtrack, -1)) { + continue; + } + if (track1ID == postrack.globalIndex()) { + continue; } - if (fillPolarization) { - ROOT::Math::Boost boost{Lambdac.BoostToCM()}; - fourVecDauCM = boost(Kshort); - threeVecDauCM = fourVecDauCM.Vect(); - threeVecDauCMXY = ROOT::Math::XYZVector(threeVecDauCM.X(), threeVecDauCM.Y(), 0.); - beamvector = ROOT::Math::XYZVector(0, 0, 1); - auto cosThetaStar = beamvector.Dot(threeVecDauCM) / std::sqrt(threeVecDauCM.Mag2()) / std::sqrt(beamvector.Mag2()); - auto SA = cosThetaStar * TMath::Sin(2.0 * phiminuspsi); - - if (track1.sign() > 0) { - histos.fill(HIST("hSparseV2SAMixedEventplus_SA"), Lambdac.M(), Lambdac.Pt(), cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SAMixedEventplus_SA_A0"), Lambdac.M(), Lambdac.Pt(), cosThetaStar * cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SAMixedEventplus_SA_azimuth"), Lambdac.M(), Lambdac.Pt(), SA, centrality); + if (track1ID == negtrack.globalIndex()) { + continue; + } + if (firstprimarytrack == 0) { + histos.fill(HIST("hInvMassKs0"), v0.mK0Short()); + } + firstprimarytrack = firstprimarytrack + 1; + float v0x, v0y, v0z, v0px, v0py, v0pz; + float posTrackX, negTrackX; + o2::track::TrackParCov trackParCovV0DaughPos; + o2::track::TrackParCov trackParCovV0DaughNeg; + trackParCovV0DaughPos = getTrackParCov(postrack); // check that aod::TracksWCov does not need TracksDCA! + trackParCovV0DaughNeg = getTrackParCov(negtrack); // check that aod::TracksWCov does not need TracksDCA! + posTrackX = v0.posX(); + negTrackX = v0.negX(); + v0x = v0.x(); + v0y = v0.y(); + v0z = v0.z(); + const std::array vertexV0 = {v0x, v0y, v0z}; + + v0px = v0.px(); + v0py = v0.py(); + v0pz = v0.pz(); + const std::array momentumV0 = {v0px, v0py, v0pz}; + + std::array covV0Pos = {0.}; + for (int i = 0; i < 6; i++) { + covV0Pos[i] = v0.positionCovMat()[i]; + } + trackParCovV0DaughPos.propagateTo(posTrackX, bz); // propagate the track to the X closest to the V0 vertex + trackParCovV0DaughNeg.propagateTo(negTrackX, bz); // propagate the track to the X closest to the V0 vertex + + // we build the neutral track to then build the cascade + // auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, {0, 0, 0, 0, 0, 0}, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) + auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, covV0Pos, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) + std::array pVecV0 = {0., 0., 0.}; + std::array pVecBach = {0., 0., 0.}; + std::array pVecCand = {0., 0., 0.}; + try { + if (df.process(trackV0, trackParCovBach) == 0) { + continue; + } else { } - if (track1.sign() < 0) { - histos.fill(HIST("hSparseV2SAMixedEventminus_SA"), Lambdac.M(), Lambdac.Pt(), cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SAMixedEventminus_SA_A0"), Lambdac.M(), Lambdac.Pt(), cosThetaStar * cosThetaStar, phiminuspsi, centrality); - histos.fill(HIST("hSparseV2SAMixedEventminus_SA_azimuth"), Lambdac.M(), Lambdac.Pt(), SA, centrality); + } catch (const std::runtime_error& error) { + continue; + } + df.propagateTracksToVertex(); // propagate the bachelor and V0 to the Lambdac vertex + trackV0.getPxPyPzGlo(pVecV0); // momentum of D0 at the Lambdac vertex + trackParCovBach.getPxPyPzGlo(pVecBach); // momentum of proton at the Lambdac vertex + pVecCand = RecoDecay::pVec(pVecV0, pVecBach); + const auto& secondaryVertex = df.getPCACandidate(); + + // get track impact parameters + // This modifies track momenta! + auto primaryVertex = getPrimaryVertex(collision); + o2::dataformats::DCA impactParameter0; + o2::dataformats::DCA impactParameter1; + trackV0.propagateToDCA(primaryVertex, bz, &impactParameter0); + trackParCovBach.propagateToDCA(primaryVertex, bz, &impactParameter1); + double phi, theta; + getPointDirection(std::array{collision.posX(), collision.posY(), collision.posZ()}, secondaryVertex, phi, theta); + // auto errorDecayLength = std::sqrt(getRotatedCovMatrixXX(covMatrixPV, phi, theta) + getRotatedCovMatrixXX(covMatrixPCA, phi, theta)); + // auto errorDecayLengthXY = std::sqrt(getRotatedCovMatrixXX(covMatrixPV, phi, 0.) + getRotatedCovMatrixXX(covMatrixPCA, phi, 0.)); + + Kshort = ROOT::Math::PxPyPzMVector(pVecV0[0], pVecV0[1], pVecV0[2], massK0s); + Proton = ROOT::Math::PxPyPzMVector(pVecBach[0], pVecBach[1], pVecBach[2], massPr); + Lambdac = Proton + Kshort; + auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); + v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; + + double protonimpactparameter = impactParameter1.getY(); + double kshortimpactparameter = impactParameter0.getY(); + + double decaylengthx = secondaryVertex[0] - collision.posX(); + double decaylengthy = secondaryVertex[1] - collision.posY(); + double decaylengthz = secondaryVertex[2] - collision.posZ(); + double decaylength = TMath::Sqrt(decaylengthx * decaylengthx + decaylengthy * decaylengthy + decaylengthz * decaylengthz); + double decaylengthxy = TMath::Sqrt(decaylengthx * decaylengthx + decaylengthy * decaylengthy); + double anglesign = decaylengthx * Lambdac.Px() + decaylengthy * Lambdac.Py() + decaylengthz * Lambdac.Pz(); + double CPAlambdac = anglesign / (decaylength * Lambdac.P()); + + histos.fill(HIST("hDecayLengthxy"), decaylengthxy); + histos.fill(HIST("hDecayLength"), decaylength); + histos.fill(HIST("hImpactPar0"), protonimpactparameter); + histos.fill(HIST("hImpactPar1"), kshortimpactparameter); + histos.fill(HIST("hCPA"), CPAlambdac); + histos.fill(HIST("hMomCorr"), Proton.P() - track1.p(), v0.p() - Kshort.P()); + if (decaylength > ConfMaxDecayLength) { + continue; + } + if (CPAlambdac < ConfMinCPA) { + continue; + } + + if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEvent_V2_SVX"), Lambdac.M(), Lambdac.Pt(), v2, decaylengthxy, CPAlambdac); + } + if (fillRotation) { + for (int nrotbkg = 0; nrotbkg < nBkgRotations; nrotbkg++) { + auto anglestart = confMinRot; + auto angleend = confMaxRot; + auto anglestep = (angleend - anglestart) / (1.0 * (nBkgRotations - 1)); + auto rotangle = anglestart + nrotbkg * anglestep; + histos.fill(HIST("hRotation"), rotangle); + auto rotKaonPx = Kshort.px() * std::cos(rotangle) - Kshort.py() * std::sin(rotangle); + auto rotKaonPy = Kshort.px() * std::sin(rotangle) + Kshort.py() * std::cos(rotangle); + KshortRot = ROOT::Math::PxPyPzMVector(rotKaonPx, rotKaonPy, Kshort.pz(), massK0s); + LambdacRot = Proton + KshortRot; + auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); + v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; + if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && TMath::Abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEventRotational_V2_SVX"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, decaylengthxy, CPAlambdac); + } } } } } } - PROCESS_SWITCH(highmasslambda, processMixedEventOpti, "Process Mixed event new", true); + PROCESS_SWITCH(highmasslambda, processSameEventSvx, "Process Same event SVX", true); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { From 2ce7e919cd6285a39b59b00f6a56d0ebe3dd7e26 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Sat, 25 Jan 2025 17:22:05 +0100 Subject: [PATCH 03/11] fix megalinter error --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 1 + 1 file changed, 1 insertion(+) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 4a6fec731e4..c68c57d06d0 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -26,6 +26,7 @@ #include #include #include +#include #include "TRandom3.h" #include "Math/Vector3D.h" From a50014645a4e3fd70db88461a1f653f1cf719f6c Mon Sep 17 00:00:00 2001 From: skundu692 Date: Sat, 25 Jan 2025 18:18:32 +0100 Subject: [PATCH 04/11] Fix PID bug --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index c68c57d06d0..91ade24c866 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -293,11 +293,16 @@ struct highmasslambda { if (candidate.tpcInnerParam() < 1.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } - if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 2.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 2.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && candidate.hasTOF() && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } - if (candidate.tpcInnerParam() >= 2.0 && candidate.tpcNSigmaPr() > -nsigmaCutTPC && candidate.tpcNSigmaPr() < nsigmaCutTPC) { - return true; + if (candidate.tpcInnerParam() >= 2.0) { + if (candidate.hasTOF() && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + return true; + } + if (!candidate.hasTOF() && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + return true; + } } return false; } From 058af4511988a24b5ed1f4636253ba6f0f3c9a31 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Sat, 25 Jan 2025 18:31:24 +0100 Subject: [PATCH 05/11] Remove unused variable --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 91ade24c866..07fc3ed396e 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -407,7 +407,7 @@ struct highmasslambda { } o2::aod::ITSResponse itsResponse; auto centrality = collision.centFT0C(); - auto multTPC = collision.multNTracksPV(); + // auto multTPC = collision.multNTracksPV(); int occupancy = collision.trackOccupancyInTimeRange(); if (occupancy > cfgOccupancyCut) { return; @@ -619,7 +619,7 @@ struct highmasslambda { } o2::aod::ITSResponse itsResponse; auto centrality = collision.centFT0C(); - auto multTPC = collision.multNTracksPV(); + // auto multTPC = collision.multNTracksPV(); auto psiFT0C = collision.psiFT0C(); auto psiFT0A = collision.psiFT0A(); auto psiTPC = collision.psiTPC(); From 31fde7703a92d843dad13f87e4be784169d1af15 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Sun, 26 Jan 2025 09:20:04 +0100 Subject: [PATCH 06/11] Fix histogram range --- PWGLF/Tasks/Resonances/phipbpb.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGLF/Tasks/Resonances/phipbpb.cxx b/PWGLF/Tasks/Resonances/phipbpb.cxx index 5d1cb842fe0..cd9fe9861cd 100644 --- a/PWGLF/Tasks/Resonances/phipbpb.cxx +++ b/PWGLF/Tasks/Resonances/phipbpb.cxx @@ -156,7 +156,7 @@ struct phipbpb { const AxisSpec thnAxisSA{configThnAxisSA, "SA"}; AxisSpec cumulantAxis = {200, -1, 1, "phi"}; AxisSpec phiAxis = {500, -6.28, 6.28, "phi"}; - AxisSpec resAxis = {2000, -10, 10, "Res"}; + AxisSpec resAxis = {6000, -30, 30, "Res"}; AxisSpec centAxis = {8, 0, 80, "V0M (%)"}; AxisSpec occupancyAxis = {occupancyBinning, "Occupancy"}; From 333601706860d60545c7514de06a1ad03c767ae3 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Tue, 28 Jan 2025 19:50:54 +0100 Subject: [PATCH 07/11] Improved DCA fitter for secondary reconstruction and add new PID selection --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 156 +++++++++++++++------- 1 file changed, 109 insertions(+), 47 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 07fc3ed396e..efdb3418447 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -170,7 +170,7 @@ struct highmasslambda { // std::vector ptLambdaBinning = {2.0, 3.0, 4.0, 5.0, 6.0}; std::vector occupancyBinning = {-0.5, 500.0, 1000.0, 1500.0, 2000.0, 3000.0, 4000.0, 5000.0, 50000.0}; - AxisSpec resAxis = {1600, -16, 16, "Res"}; + AxisSpec resAxis = {1600, -30, 30, "Res"}; AxisSpec phiAxis = {500, -6.28, 6.28, "phi"}; AxisSpec centAxis = {8, 0, 80, "V0M (%)"}; const AxisSpec thnAxisInvMass{configThnAxisInvMass, "#it{M} (GeV/#it{c}^{2})"}; @@ -212,9 +212,9 @@ struct highmasslambda { histos.add("hImpactPar1", "hImpactPar1", kTH1F, {{500, 0.0f, 0.1f}}); histos.add("hCPA", "hCPA", kTH1F, {{220, -1.1f, 1.1f}}); - histos.add("hSparseV2SASameEvent_V2", "hSparseV2SASameEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); - histos.add("hSparseV2SASameEventRotational_V2", "hSparseV2SASameEventRotational", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); - histos.add("hSparseV2SAMixedEvent_V2", "hSparseV2SAMixedEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA}); + histos.add("hSparseV2SASameEvent_V2", "hSparseV2SASameEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA, thnAxisPtProton}); + histos.add("hSparseV2SASameEventRotational_V2", "hSparseV2SASameEventRotational", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA, thnAxisPtProton}); + histos.add("hSparseV2SAMixedEvent_V2", "hSparseV2SAMixedEvent_V2", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDCA, thnAxisPtProton}); histos.add("hSparseV2SASameEvent_V2_SVX", "hSparseV2SASameEvent_V2_SVX", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDecayLength, thnAxisCPA}); histos.add("hSparseV2SASameEventRotational_V2_SVX", "hSparseV2SASameEventRotational_SVX", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisV2, thnAxisDecayLength, thnAxisCPA}); @@ -307,6 +307,27 @@ struct highmasslambda { return false; } + // TPC TOF + template + bool selectionPID3(const T& candidate) + { + if (candidate.hasTOF() && candidate.tpcInnerParam() > 0.6 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + return true; + } + if (!candidate.hasTOF()) { + if (candidate.tpcInnerParam() < 1.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + return true; + } + if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 1.6 && candidate.tpcNSigmaPr() < nsigmaCutTPC && candidate.tpcNSigmaPr() > -1.5) { + return true; + } + if (candidate.tpcInnerParam() >= 1.6 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + return true; + } + } + return false; + } + template bool SelectionV0(Collision const& collision, V0 const& candidate) { @@ -448,6 +469,9 @@ struct highmasslambda { if (PIDstrategy == 1 && !selectionPID2(track1)) { continue; } + if (PIDstrategy == 2 && !selectionPID3(track1)) { + continue; + } if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { continue; } @@ -490,7 +514,7 @@ struct highmasslambda { auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY())); + histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); } if (fillRotation) { for (int nrotbkg = 0; nrotbkg < nBkgRotations; nrotbkg++) { @@ -506,7 +530,7 @@ struct highmasslambda { auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && TMath::Abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, TMath::Abs(track1.dcaXY())); + histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, TMath::Abs(track1.dcaXY()), Proton.Pt()); } } } @@ -559,6 +583,9 @@ struct highmasslambda { if (PIDstrategy == 1 && !selectionPID2(track1)) { continue; } + if (PIDstrategy == 2 && !selectionPID3(track1)) { + continue; + } if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { continue; } @@ -586,7 +613,7 @@ struct highmasslambda { auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY())); + histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); } } } @@ -655,6 +682,9 @@ struct highmasslambda { if (PIDstrategy == 1 && !selectionPID2(track1)) { continue; } + if (PIDstrategy == 2 && !selectionPID3(track1)) { + continue; + } if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { continue; } @@ -667,6 +697,7 @@ struct highmasslambda { } auto track1ID = track1.globalIndex(); auto trackParCovBach = getTrackParCov(track1); + // auto trackParCovBach = getTrackParCov(bach); for (auto v0 : V0s) { if (!SelectionV0(collision, v0)) { @@ -693,53 +724,84 @@ struct highmasslambda { histos.fill(HIST("hInvMassKs0"), v0.mK0Short()); } firstprimarytrack = firstprimarytrack + 1; - float v0x, v0y, v0z, v0px, v0py, v0pz; - float posTrackX, negTrackX; - o2::track::TrackParCov trackParCovV0DaughPos; - o2::track::TrackParCov trackParCovV0DaughNeg; - trackParCovV0DaughPos = getTrackParCov(postrack); // check that aod::TracksWCov does not need TracksDCA! - trackParCovV0DaughNeg = getTrackParCov(negtrack); // check that aod::TracksWCov does not need TracksDCA! - posTrackX = v0.posX(); - negTrackX = v0.negX(); - v0x = v0.x(); - v0y = v0.y(); - v0z = v0.z(); - const std::array vertexV0 = {v0x, v0y, v0z}; - - v0px = v0.px(); - v0py = v0.py(); - v0pz = v0.pz(); - const std::array momentumV0 = {v0px, v0py, v0pz}; - - std::array covV0Pos = {0.}; - for (int i = 0; i < 6; i++) { - covV0Pos[i] = v0.positionCovMat()[i]; - } - trackParCovV0DaughPos.propagateTo(posTrackX, bz); // propagate the track to the X closest to the V0 vertex - trackParCovV0DaughNeg.propagateTo(negTrackX, bz); // propagate the track to the X closest to the V0 vertex - - // we build the neutral track to then build the cascade - // auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, {0, 0, 0, 0, 0, 0}, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) - auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, covV0Pos, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) + // LOGF(info, "Before dca fitter"); std::array pVecV0 = {0., 0., 0.}; std::array pVecBach = {0., 0., 0.}; std::array pVecCand = {0., 0., 0.}; + const std::array vertexV0 = {v0.x(), v0.y(), v0.z()}; + const std::array momentumV0 = {v0.px(), v0.py(), v0.pz()}; + // we build the neutral track to then build the cascade + std::array covV = {0.}; + constexpr int MomInd[6] = {9, 13, 14, 18, 19, 20}; // cov matrix elements for momentum component + for (int i = 0; i < 6; i++) { + covV[MomInd[i]] = v0.momentumCovMat()[i]; + covV[i] = v0.positionCovMat()[i]; + } + auto trackV0 = o2::track::TrackParCov(vertexV0, momentumV0, covV, 0, true); + trackV0.setAbsCharge(0); + trackV0.setPID(o2::track::PID::K0); + + int nCand2 = 0; try { - if (df.process(trackV0, trackParCovBach) == 0) { - continue; - } else { - } - } catch (const std::runtime_error& error) { + nCand2 = df.process(trackV0, trackParCovBach); + } catch (...) { continue; } - df.propagateTracksToVertex(); // propagate the bachelor and V0 to the Lambdac vertex - trackV0.getPxPyPzGlo(pVecV0); // momentum of D0 at the Lambdac vertex - trackParCovBach.getPxPyPzGlo(pVecBach); // momentum of proton at the Lambdac vertex - pVecCand = RecoDecay::pVec(pVecV0, pVecBach); - const auto& secondaryVertex = df.getPCACandidate(); - // get track impact parameters - // This modifies track momenta! + if (nCand2 == 0) { + continue; + } + df.propagateTracksToVertex(); // propagate the bach and V0 to the Lc vertex + df.getTrack(0).getPxPyPzGlo(pVecV0); // take the momentum at the Lc vertex + df.getTrack(1).getPxPyPzGlo(pVecBach); + // LOGF(info, "after dca fitter"); + + /* + float v0x, v0y, v0z, v0px, v0py, v0pz; + float posTrackX, negTrackX; + o2::track::TrackParCov trackParCovV0DaughPos; + o2::track::TrackParCov trackParCovV0DaughNeg; + trackParCovV0DaughPos = getTrackParCov(postrack); // check that aod::TracksWCov does not need TracksDCA! + trackParCovV0DaughNeg = getTrackParCov(negtrack); // check that aod::TracksWCov does not need TracksDCA! + posTrackX = v0.posX(); + negTrackX = v0.negX(); + v0x = v0.x(); + v0y = v0.y(); + v0z = v0.z(); + const std::array vertexV0 = {v0x, v0y, v0z}; + + v0px = v0.px(); + v0py = v0.py(); + v0pz = v0.pz(); + const std::array momentumV0 = {v0px, v0py, v0pz}; + + std::array covV0Pos = {0.}; + for (int i = 0; i < 6; i++) { + covV0Pos[i] = v0.positionCovMat()[i]; + } + trackParCovV0DaughPos.propagateTo(posTrackX, bz); // propagate the track to the X closest to the V0 vertex + trackParCovV0DaughNeg.propagateTo(negTrackX, bz); // propagate the track to the X closest to the V0 vertex + + // we build the neutral track to then build the cascade + // auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, {0, 0, 0, 0, 0, 0}, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) + auto trackV0 = o2::dataformats::V0(vertexV0, momentumV0, covV0Pos, trackParCovV0DaughPos, trackParCovV0DaughNeg); // build the V0 track (indices for v0 daughters set to 0 for now) + std::array pVecV0 = {0., 0., 0.}; + std::array pVecBach = {0., 0., 0.}; + std::array pVecCand = {0., 0., 0.}; + try { + if (df.process(trackV0, trackParCovBach) == 0) { + continue; + } else { + } + } catch (const std::runtime_error& error) { + continue; + } + df.propagateTracksToVertex(); // propagate the bachelor and V0 to the Lambdac vertex + trackV0.getPxPyPzGlo(pVecV0); // momentum of D0 at the Lambdac vertex + trackParCovBach.getPxPyPzGlo(pVecBach); // momentum of proton at the Lambdac vertex + */ + + const auto& secondaryVertex = df.getPCACandidate(); auto primaryVertex = getPrimaryVertex(collision); o2::dataformats::DCA impactParameter0; o2::dataformats::DCA impactParameter1; From bf8e141acf4bfbab9dc1731f13f17ac6000a2606 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Wed, 29 Jan 2025 14:29:36 +0100 Subject: [PATCH 08/11] Improve PID for proton selection --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 119 ++++++++++++++-------- 1 file changed, 75 insertions(+), 44 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 2df58a81c28..e8c6b0231f9 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -248,19 +248,19 @@ struct highmasslambda { if (!(candidate.isGlobalTrackWoDCA() && candidate.itsNCls() > cfgITScluster && candidate.tpcNClsFound() > cfgTPCcluster)) { return false; } - if (candidate.pt() > 0.0 && candidate.pt() < 0.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin1) { + if (candidate.pt() > 0.0 && candidate.pt() < 0.5 && std::abs(candidate.dcaXY()) < cfgCutDCAxymin1) { return false; } - if (candidate.pt() >= 0.5 && candidate.pt() < 1.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin2) { + if (candidate.pt() >= 0.5 && candidate.pt() < 1.0 && std::abs(candidate.dcaXY()) < cfgCutDCAxymin2) { return false; } - if (candidate.pt() >= 1.0 && candidate.pt() < 1.5 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin3) { + if (candidate.pt() >= 1.0 && candidate.pt() < 1.5 && std::abs(candidate.dcaXY()) < cfgCutDCAxymin3) { return false; } - if (candidate.pt() >= 1.5 && candidate.pt() < 2.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin4) { + if (candidate.pt() >= 1.5 && candidate.pt() < 2.0 && std::abs(candidate.dcaXY()) < cfgCutDCAxymin4) { return false; } - if (candidate.pt() >= 2.0 && candidate.pt() < 10000000.0 && TMath::Abs(candidate.dcaXY()) < cfgCutDCAxymin5) { + if (candidate.pt() >= 2.0 && candidate.pt() < 10000000.0 && std::abs(candidate.dcaXY()) < cfgCutDCAxymin5) { return false; } return true; @@ -271,15 +271,15 @@ struct highmasslambda { bool selectionPID1(const T& candidate) { if (candidate.hasTOF()) { - if (candidate.p() < 0.5 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (candidate.p() < 0.5 && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } - if (candidate.p() >= 0.5 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.p() >= 0.5 && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && std::abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } } if (!candidate.hasTOF()) { - if (TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } } @@ -290,17 +290,33 @@ struct highmasslambda { template bool selectionPID2(const T& candidate) { - if (candidate.tpcInnerParam() < 1.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (candidate.pt() < 0.7 && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } - if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 2.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && candidate.hasTOF() && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { - return true; + if (candidate.pt() >= 0.7 && candidate.pt() < 1.4) { + if (candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && std::abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + return true; + } + if (!candidate.hasTOF()) { + if (candidate.pt() >= 0.7 && candidate.pt() < 0.8 && candidate.tpcNSigmaPr() > -1.8 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + if (candidate.pt() >= 0.8 && candidate.pt() < 0.9 && candidate.tpcNSigmaPr() > -1.7 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + if (candidate.pt() >= 0.9 && candidate.pt() < 1.0 && candidate.tpcNSigmaPr() > -1.6 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + if (candidate.pt() >= 1.0 && candidate.pt() < 1.4 && candidate.tpcNSigmaPr() > -1.0 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + } } - if (candidate.tpcInnerParam() >= 2.0) { - if (candidate.hasTOF() && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.pt() >= 1.4) { + if (candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && std::abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } - if (!candidate.hasTOF() && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (!candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } } @@ -311,17 +327,30 @@ struct highmasslambda { template bool selectionPID3(const T& candidate) { - if (candidate.hasTOF() && candidate.tpcInnerParam() > 0.6 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { + if (candidate.pt() < 0.7 && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } - if (!candidate.hasTOF()) { - if (candidate.tpcInnerParam() < 1.0 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (candidate.pt() >= 0.7 && candidate.pt() < 1.8) { + if (candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && std::abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } - if (candidate.tpcInnerParam() >= 1.0 && candidate.tpcInnerParam() < 1.6 && candidate.tpcNSigmaPr() < nsigmaCutTPC && candidate.tpcNSigmaPr() > -1.5) { + if (!candidate.hasTOF()) { + if (candidate.pt() >= 0.7 && candidate.pt() < 0.8 && candidate.tpcNSigmaPr() > -1.8 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + if (candidate.pt() >= 0.8 && candidate.pt() < 0.9 && candidate.tpcNSigmaPr() > -1.7 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + if (candidate.pt() >= 0.9 && candidate.pt() < 1.0 && candidate.tpcNSigmaPr() > -1.6 && candidate.tpcNSigmaPr() < nsigmaCutTPC) { + return true; + } + } + } + if (candidate.pt() >= 1.8) { + if (candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC && std::abs(candidate.tofNSigmaPr()) < nsigmaCutTOF) { return true; } - if (candidate.tpcInnerParam() >= 1.6 && TMath::Abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { + if (!candidate.hasTOF() && std::abs(candidate.tpcNSigmaPr()) < nsigmaCutTPC) { return true; } } @@ -331,13 +360,13 @@ struct highmasslambda { template bool SelectionV0(Collision const& collision, V0 const& candidate) { - if (TMath::Abs(candidate.dcav0topv()) > cMaxV0DCA) { + if (std::abs(candidate.dcav0topv()) > cMaxV0DCA) { return false; } const float pT = candidate.pt(); const std::vector decVtx = {candidate.x(), candidate.y(), candidate.z()}; const float tranRad = candidate.v0radius(); - const double dcaDaughv0 = TMath::Abs(candidate.dcaV0daughters()); + const double dcaDaughv0 = std::abs(candidate.dcaV0daughters()); const double cpav0 = candidate.v0cosPA(); float CtauK0s = candidate.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(); // FIXME: Get from the common header @@ -359,7 +388,7 @@ struct highmasslambda { if (tranRad > ConfV0TranRadV0Max) { return false; } - if (TMath::Abs(CtauK0s) > cMaxV0LifeTime || candidate.mK0Short() < lowmasscutks0 || candidate.mK0Short() > highmasscutks0) { + if (std::abs(CtauK0s) > cMaxV0LifeTime || candidate.mK0Short() < lowmasscutks0 || candidate.mK0Short() > highmasscutks0) { return false; } return true; @@ -378,19 +407,19 @@ struct highmasslambda { if (charge > 0 && sign < 0) { return false; } - if (TMath::Abs(eta) > 0.8) { + if (std::abs(eta) > 0.8) { return false; } - if (TMath::Abs(pt) < 0.15) { + if (std::abs(pt) < 0.15) { return false; } if (tpcNClsF < 50) { return false; } - if (TMath::Abs(dcaXY) < ConfDaughDCAMin) { + if (std::abs(dcaXY) < ConfDaughDCAMin) { return false; } - if (TMath::Abs(track.tpcNSigmaPi()) > ConfDaughPIDCuts) { + if (std::abs(track.tpcNSigmaPi()) > ConfDaughPIDCuts) { return false; } return true; @@ -461,6 +490,9 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } + if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + continue; + } histos.fill(HIST("hNsigmaProtonTPCPre"), track1.tpcNSigmaPr(), track1.pt()); // PID check if (PIDstrategy == 0 && !selectionPID1(track1)) { @@ -472,9 +504,6 @@ struct highmasslambda { if (PIDstrategy == 2 && !selectionPID3(track1)) { continue; } - if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { - continue; - } histos.fill(HIST("hEta"), track1.eta()); histos.fill(HIST("hDcaxy"), track1.dcaXY()); histos.fill(HIST("hDcaz"), track1.dcaZ()); @@ -513,8 +542,8 @@ struct highmasslambda { Lambdac = Proton + Kshort; auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); + if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, std::abs(track1.dcaXY()), Proton.Pt()); } if (fillRotation) { for (int nrotbkg = 0; nrotbkg < nBkgRotations; nrotbkg++) { @@ -529,8 +558,8 @@ struct highmasslambda { LambdacRot = Proton + KshortRot; auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; - if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && TMath::Abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, TMath::Abs(track1.dcaXY()), Proton.Pt()); + if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, std::abs(track1.dcaXY()), Proton.Pt()); } } } @@ -576,6 +605,10 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } + if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + continue; + } + // PID check if (PIDstrategy == 0 && !selectionPID1(track1)) { continue; @@ -586,9 +619,6 @@ struct highmasslambda { if (PIDstrategy == 2 && !selectionPID3(track1)) { continue; } - if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { - continue; - } if (!SelectionV0(collision2, v0)) { continue; @@ -607,13 +637,13 @@ struct highmasslambda { if (Lambdac.Pt() > 6.0 || Lambdac.Pt() < 2.0) { continue; } - if (TMath::Abs(Lambdac.Rapidity()) > confRapidity) { + if (std::abs(Lambdac.Rapidity()) > confRapidity) { continue; } auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, TMath::Abs(track1.dcaXY()), Proton.Pt()); + if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, std::abs(track1.dcaXY()), Proton.Pt()); } } } @@ -674,6 +704,10 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } + if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + continue; + } + histos.fill(HIST("hNsigmaProtonTPCPre"), track1.tpcNSigmaPr(), track1.pt()); // PID check if (PIDstrategy == 0 && !selectionPID1(track1)) { @@ -685,9 +719,6 @@ struct highmasslambda { if (PIDstrategy == 2 && !selectionPID3(track1)) { continue; } - if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -2.5 && itsResponse.nSigmaITS(track1) < 2.5)) { - continue; - } histos.fill(HIST("hEta"), track1.eta()); histos.fill(HIST("hDcaxy"), track1.dcaXY()); histos.fill(HIST("hDcaz"), track1.dcaZ()); @@ -842,7 +873,7 @@ struct highmasslambda { continue; } - if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && TMath::Abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { histos.fill(HIST("hSparseV2SASameEvent_V2_SVX"), Lambdac.M(), Lambdac.Pt(), v2, decaylengthxy, CPAlambdac); } if (fillRotation) { @@ -858,7 +889,7 @@ struct highmasslambda { LambdacRot = Proton + KshortRot; auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; - if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && TMath::Abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { + if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { histos.fill(HIST("hSparseV2SASameEventRotational_V2_SVX"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, decaylengthxy, CPAlambdac); } } From 52a421f3690787b25b69ebc550848ab9c149b246 Mon Sep 17 00:00:00 2001 From: skundu692 Date: Wed, 29 Jan 2025 17:30:55 +0100 Subject: [PATCH 09/11] Separate additional event selection --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 23 ++++++++++++++--------- 1 file changed, 14 insertions(+), 9 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index e8c6b0231f9..072eb11d315 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -449,10 +449,10 @@ struct highmasslambda { double v2, v2Rot; void processSameEvent(EventCandidates::iterator const& collision, TrackCandidates const& tracks, AllTrackCandidates const&, ResoV0s const& V0s, aod::BCs const&) { - if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { + if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { return; } - if (additionalEvSel && !collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + if (additionalEvSel && (!collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) { return; } o2::aod::ITSResponse itsResponse; @@ -573,11 +573,17 @@ struct highmasslambda { BinningTypeVertexContributor binningOnPositions{{axisVertex, axisMultiplicityClass, axisEPAngle}, true}; Pair pairs{binningOnPositions, cfgNoMixedEvents, -1, collisions, tracksV0sTuple, &cache}; // -1 is the number of the bin to skip for (auto& [collision1, tracks1, collision2, tracks2] : pairs) { - if (!collision1.sel8() || !collision1.triggereventep() || !collision1.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision1.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision1.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision1.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision1.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { - continue; + if (!collision1.sel8() || !collision1.triggereventep() || !collision1.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision1.selection_bit(aod::evsel::kNoITSROFrameBorder)) { + return; } - if (!collision2.sel8() || !collision2.triggereventep() || !collision2.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision2.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision2.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision2.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision2.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { - continue; + if (additionalEvSel && (!collision1.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll) || !collision1.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision1.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision1.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) { + return; + } + if (!collision2.sel8() || !collision2.triggereventep() || !collision2.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision2.selection_bit(aod::evsel::kNoITSROFrameBorder)) { + return; + } + if (additionalEvSel && (!collision2.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll) || !collision2.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision2.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision2.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) { + return; } if (collision1.bcId() == collision2.bcId()) { continue; @@ -651,13 +657,12 @@ struct highmasslambda { PROCESS_SWITCH(highmasslambda, processMixedEventOpti, "Process Mixed event new", false); void processSameEventSvx(EventCandidates::iterator const& collision, TrackCandidatesSvx const& tracks, AllTrackCandidatesSvx const&, ResoV0sSvx const& V0s, aod::BCsWithTimestamps const&) { - if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)) { + if (!collision.sel8() || !collision.triggereventep() || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { return; } - if (additionalEvSel && !collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) { + if (additionalEvSel && (!collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll) || !collision.selection_bit(aod::evsel::kNoSameBunchPileup) || !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV) || !collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) { return; } - /// Set the magnetic field from ccdb. auto bc = collision.template bc_as(); if (runNumber != bc.runNumber()) { From 1a7ec0da8af2687fe9d51efaf23bb97869d85cae Mon Sep 17 00:00:00 2001 From: skundu692 Date: Thu, 30 Jan 2025 19:46:52 +0100 Subject: [PATCH 10/11] Improved rotational background --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 72 +++++++++++++++++++---- 1 file changed, 60 insertions(+), 12 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 072eb11d315..22f83409e44 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -86,6 +86,7 @@ struct highmasslambda { // fill output Configurable cfgOccupancyCut{"cfgOccupancyCut", 2500, "Occupancy cut"}; Configurable fillRotation{"fillRotation", false, "fill rotation"}; + Configurable useSP{"useSP", false, "useSP"}; // events Configurable cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"}; Configurable cfgCutCentralityMax{"cfgCutCentralityMax", 50.0f, "Accepted maximum Centrality"}; @@ -109,6 +110,8 @@ struct highmasslambda { Configurable PIDstrategy{"PIDstrategy", 0, "0: TOF Veto, 1: TOF Veto opti, 2: TOF, 3: TOF loose 1, 4: TOF loose 2, 5: old pt dep"}; Configurable nsigmaCutTPC{"nsigmacutTPC", 3.0, "Value of the TPC Nsigma cut"}; Configurable nsigmaCutTOF{"nsigmaCutTOF", 3.0, "TOF PID"}; + Configurable nsigmaCutITS{"nsigmaCutITS", 3.0, "Value of the ITS Nsigma cut"}; + // Configs for V0 Configurable ConfV0PtMin{"ConfV0PtMin", 0.f, "Minimum transverse momentum of V0"}; Configurable ConfV0DCADaughMax{"ConfV0DCADaughMax", 0.2f, "Maximum DCA between the V0 daughters"}; @@ -490,7 +493,7 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } - if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -nsigmaCutITS && itsResponse.nSigmaITS(track1) < nsigmaCutITS)) { continue; } histos.fill(HIST("hNsigmaProtonTPCPre"), track1.tpcNSigmaPr(), track1.pt()); @@ -542,7 +545,10 @@ struct highmasslambda { Lambdac = Proton + Kshort; auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + if (useSP) { + v2 = TMath::Cos(2.0 * phiminuspsi); + } + if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 1.0 && Lambdac.Pt() <= 6.0) { histos.fill(HIST("hSparseV2SASameEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, std::abs(track1.dcaXY()), Proton.Pt()); } if (fillRotation) { @@ -558,7 +564,10 @@ struct highmasslambda { LambdacRot = Proton + KshortRot; auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; - if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { + if (useSP) { + v2Rot = TMath::Cos(2.0 * phiminuspsiRot); + } + if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 1.0 && LambdacRot.Pt() <= 6.0) { histos.fill(HIST("hSparseV2SASameEventRotational_V2"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, std::abs(track1.dcaXY()), Proton.Pt()); } } @@ -611,7 +620,7 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } - if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -nsigmaCutITS && itsResponse.nSigmaITS(track1) < nsigmaCutITS)) { continue; } @@ -648,7 +657,10 @@ struct highmasslambda { } auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { + if (useSP) { + v2 = TMath::Cos(2.0 * phiminuspsi); + } + if (occupancy1 < cfgOccupancyCut && occupancy2 < cfgOccupancyCut && Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 1.0 && Lambdac.Pt() <= 6.0) { histos.fill(HIST("hSparseV2SAMixedEvent_V2"), Lambdac.M(), Lambdac.Pt(), v2, std::abs(track1.dcaXY()), Proton.Pt()); } } @@ -709,7 +721,7 @@ struct highmasslambda { if (!selectionTrack(track1)) { continue; } - if (!(itsResponse.nSigmaITS(track1) > -3.0 && itsResponse.nSigmaITS(track1) < 3.0)) { + if (track1.p() < 1.0 && !(itsResponse.nSigmaITS(track1) > -nsigmaCutITS && itsResponse.nSigmaITS(track1) < nsigmaCutITS)) { continue; } @@ -853,7 +865,9 @@ struct highmasslambda { Lambdac = Proton + Kshort; auto phiminuspsi = GetPhiInRange(Lambdac.Phi() - psiFT0C); v2 = TMath::Cos(2.0 * phiminuspsi) * QFT0C; - + if (useSP) { + v2 = TMath::Cos(2.0 * phiminuspsi); + } double protonimpactparameter = impactParameter1.getY(); double kshortimpactparameter = impactParameter0.getY(); @@ -878,8 +892,8 @@ struct highmasslambda { continue; } - if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 2.0 && Lambdac.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEvent_V2_SVX"), Lambdac.M(), Lambdac.Pt(), v2, decaylengthxy, CPAlambdac); + if (Lambdac.M() > cMinLambdaMass && Lambdac.M() <= cMaxLambdaMass && std::abs(Lambdac.Rapidity()) < confRapidity && Lambdac.Pt() > 1.0 && Lambdac.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEvent_V2_SVX"), Lambdac.M(), Lambdac.Pt(), v2, decaylength, CPAlambdac); } if (fillRotation) { for (int nrotbkg = 0; nrotbkg < nBkgRotations; nrotbkg++) { @@ -890,12 +904,46 @@ struct highmasslambda { histos.fill(HIST("hRotation"), rotangle); auto rotKaonPx = Kshort.px() * std::cos(rotangle) - Kshort.py() * std::sin(rotangle); auto rotKaonPy = Kshort.px() * std::sin(rotangle) + Kshort.py() * std::cos(rotangle); - KshortRot = ROOT::Math::PxPyPzMVector(rotKaonPx, rotKaonPy, Kshort.pz(), massK0s); + ////////// DCA fitter //////////////// + // LOGF(info, "Before dca fitter"); + std::array pVecV0rot = {0., 0., 0.}; + std::array pVecBachrot = {0., 0., 0.}; + const std::array momentumV0rot = {rotKaonPx, rotKaonPy, v0.pz()}; + auto trackV0rot = o2::track::TrackParCov(vertexV0, momentumV0rot, covV, 0, true); + trackV0rot.setAbsCharge(0); + trackV0rot.setPID(o2::track::PID::K0); + int nCand2rot = 0; + try { + nCand2rot = df.process(trackV0rot, trackParCovBach); + } catch (...) { + continue; + } + if (nCand2rot == 0) { + continue; + } + df.propagateTracksToVertex(); // propagate the bach and V0 to the Lc vertex + df.getTrack(0).getPxPyPzGlo(pVecV0rot); // take the momentum at the Lc vertex + df.getTrack(1).getPxPyPzGlo(pVecBachrot); + const auto& secondaryVertexrot = df.getPCACandidate(); + double phirot, thetarot; + getPointDirection(std::array{collision.posX(), collision.posY(), collision.posZ()}, secondaryVertexrot, phirot, thetarot); + KshortRot = ROOT::Math::PxPyPzMVector(pVecV0rot[0], pVecV0rot[1], pVecV0rot[2], massK0s); + Proton = ROOT::Math::PxPyPzMVector(pVecBachrot[0], pVecBachrot[1], pVecBachrot[2], massPr); LambdacRot = Proton + KshortRot; auto phiminuspsiRot = GetPhiInRange(LambdacRot.Phi() - psiFT0C); v2Rot = TMath::Cos(2.0 * phiminuspsiRot) * QFT0C; - if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 2.0 && LambdacRot.Pt() <= 6.0) { - histos.fill(HIST("hSparseV2SASameEventRotational_V2_SVX"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, decaylengthxy, CPAlambdac); + if (useSP) { + v2Rot = TMath::Cos(2.0 * phiminuspsiRot); + } + double decaylengthxrot = secondaryVertexrot[0] - collision.posX(); + double decaylengthyrot = secondaryVertexrot[1] - collision.posY(); + double decaylengthzrot = secondaryVertexrot[2] - collision.posZ(); + double decaylengthrot = TMath::Sqrt(decaylengthxrot * decaylengthxrot + decaylengthyrot * decaylengthyrot + decaylengthzrot * decaylengthzrot); + double decaylengthxyrot = TMath::Sqrt(decaylengthxrot * decaylengthxrot + decaylengthyrot * decaylengthyrot); + double anglesignrot = decaylengthxrot * LambdacRot.Px() + decaylengthyrot * LambdacRot.Py() + decaylengthzrot * LambdacRot.Pz(); + double CPAlambdacrot = anglesignrot / (decaylengthrot * LambdacRot.P()); + if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 1.0 && LambdacRot.Pt() <= 6.0) { + histos.fill(HIST("hSparseV2SASameEventRotational_V2_SVX"), LambdacRot.M(), LambdacRot.Pt(), v2Rot, decaylengthrot, CPAlambdacrot); } } } From 70746bfe60b516a2079bd807cd24d39853e25091 Mon Sep 17 00:00:00 2001 From: skundu692 <86804743+skundu692@users.noreply.github.com> Date: Thu, 30 Jan 2025 20:42:48 +0100 Subject: [PATCH 11/11] Remove unused variable --- PWGLF/Tasks/Resonances/highmasslambda.cxx | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/PWGLF/Tasks/Resonances/highmasslambda.cxx b/PWGLF/Tasks/Resonances/highmasslambda.cxx index 22f83409e44..85dda8d4fa7 100644 --- a/PWGLF/Tasks/Resonances/highmasslambda.cxx +++ b/PWGLF/Tasks/Resonances/highmasslambda.cxx @@ -902,8 +902,8 @@ struct highmasslambda { auto anglestep = (angleend - anglestart) / (1.0 * (nBkgRotations - 1)); auto rotangle = anglestart + nrotbkg * anglestep; histos.fill(HIST("hRotation"), rotangle); - auto rotKaonPx = Kshort.px() * std::cos(rotangle) - Kshort.py() * std::sin(rotangle); - auto rotKaonPy = Kshort.px() * std::sin(rotangle) + Kshort.py() * std::cos(rotangle); + float rotKaonPx = Kshort.px() * std::cos(rotangle) - Kshort.py() * std::sin(rotangle); + float rotKaonPy = Kshort.px() * std::sin(rotangle) + Kshort.py() * std::cos(rotangle); ////////// DCA fitter //////////////// // LOGF(info, "Before dca fitter"); std::array pVecV0rot = {0., 0., 0.}; @@ -939,7 +939,7 @@ struct highmasslambda { double decaylengthyrot = secondaryVertexrot[1] - collision.posY(); double decaylengthzrot = secondaryVertexrot[2] - collision.posZ(); double decaylengthrot = TMath::Sqrt(decaylengthxrot * decaylengthxrot + decaylengthyrot * decaylengthyrot + decaylengthzrot * decaylengthzrot); - double decaylengthxyrot = TMath::Sqrt(decaylengthxrot * decaylengthxrot + decaylengthyrot * decaylengthyrot); + // double decaylengthxyrot = TMath::Sqrt(decaylengthxrot * decaylengthxrot + decaylengthyrot * decaylengthyrot); double anglesignrot = decaylengthxrot * LambdacRot.Px() + decaylengthyrot * LambdacRot.Py() + decaylengthzrot * LambdacRot.Pz(); double CPAlambdacrot = anglesignrot / (decaylengthrot * LambdacRot.P()); if (LambdacRot.M() > cMinLambdaMass && LambdacRot.M() <= cMaxLambdaMass && std::abs(LambdacRot.Rapidity()) < confRapidity && LambdacRot.Pt() > 1.0 && LambdacRot.Pt() <= 6.0) {