1616
1717// #include <TDatabasePDG.h>
1818#include " PWGLF/DataModel/LFStrangenessTables.h"
19+ #include " PWGLF/DataModel/mcCentrality.h"
20+ #include " PWGLF/Utils/inelGt.h"
1921
2022#include " Common/Core/TrackSelection.h"
2123#include " Common/Core/trackUtilities.h"
@@ -105,6 +107,7 @@ struct Kstarqa {
105107 Configurable<bool > cfgUseITSTPCRefit{" cfgUseITSTPCRefit" , false , " Require ITS Refit" };
106108 Configurable<bool > isNoCollInTimeRangeStandard{" isNoCollInTimeRangeStandard" , false , " No collision in time range standard" };
107109 Configurable<bool > isApplyPtDepDCAxyCut{" isApplyPtDepDCAxyCut" , false , " Apply pT dependent DCAxy cut" };
110+ Configurable<bool > isGoldenChi2{" isGoldenChi2" , false , " Apply golden chi2 cut" };
108111
109112 // cuts on mother
110113 Configurable<bool > isApplyCutsOnMother{" isApplyCutsOnMother" , false , " Enable additional cuts on Kstar mother" };
@@ -427,6 +430,8 @@ struct Kstarqa {
427430 if (std::abs (candidate.dcaXY ()) > (0.0105 + 0.035 / std::pow (candidate.pt (), 1.1 )))
428431 return false ;
429432 }
433+ if (selectionConfig.isGoldenChi2 && candidate.passedGoldenChi2 ())
434+ return false ;
430435 if (std::abs (candidate.dcaZ ()) > selectionConfig.cfgCutDCAz )
431436 return false ;
432437 if (candidate.tpcCrossedRowsOverFindableCls () < selectionConfig.cfgRCRFC )
@@ -660,12 +665,13 @@ struct Kstarqa {
660665 Filter acceptanceFilter = (nabs(aod::track::eta) < selectionConfig.cfgCutEta && nabs(aod::track::pt) > selectionConfig.cfgCutPT);
661666 Filter fDCAcutFilter = (nabs(aod::track::dcaXY) < selectionConfig.cfgCutDCAxy) && (nabs(aod::track::dcaZ) < selectionConfig.cfgCutDCAz);
662667
663- using EventCandidates = soa::Join<aod::Collisions, aod::EvSels, aod::FT0Mults, aod::MultZeqs, aod::CentFT0Ms, aod::CentFT0Cs, aod::CentFT0As, aod::CentFV0As, aod::PVMults>; // aod::CentNGlobals, aod::CentNTPVs, aod::CentMFTs
664- using EventCandidatesMix = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::FT0Mults, aod::MultZeqs, aod::CentFT0Ms, aod::CentFT0Cs, aod::CentFT0As, aod::CentFV0As>>; // aod::CentNGlobals, aod::CentNTPVs, aod::CentMFTs
665- using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullPi, aod::pidTOFFullPi, aod::pidTPCFullKa, aod::pidTOFFullKa, aod::pidTOFbeta>>;
666- using EventCandidatesMC = soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels, aod::CentFT0Ms, aod::CentFT0As, aod::CentFT0Cs, aod::FT0Mults, aod::PVMults>;
668+ using EventCandidates = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::FT0Mults, aod::MultZeqs, aod::CentFT0Ms, aod::CentFT0Cs, aod::CentFT0As, aod::CentFV0As, aod::PVMults>>; // aod::CentNGlobals, aod::CentNTPVs, aod::CentMFTs
669+ using EventCandidatesMix = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::FT0Mults, aod::MultZeqs, aod::CentFT0Ms, aod::CentFT0Cs, aod::CentFT0As, aod::CentFV0As>>; // aod::CentNGlobals, aod::CentNTPVs, aod::CentMFTs
670+ using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullPi, aod::pidTOFFullPi, aod::pidTPCFullKa, aod::pidTOFFullKa, aod::pidTOFbeta, aod::TrackSelectionExtension>>;
671+ using EventCandidatesMC = soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels, aod::CentFT0Ms, aod::CentFT0As, aod::CentFT0Cs, aod::MultZeqs, aod::FT0Mults, aod::PVMults, aod::CentFV0As>;
672+ // using EventCandidatesMC = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels, aod::CentFT0Ms, aod::CentFT0As, aod::CentFT0Cs, aod::FT0Mults, aod::PVMults, aod::CentFV0As>>;
667673
668- using TrackCandidatesMC = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullKa, aod::pidTOFFullKa, aod::pidTPCFullPi, aod::pidTOFFullPi, aod::McTrackLabels, aod::pidTOFbeta>>;
674+ using TrackCandidatesMC = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullKa, aod::pidTOFFullKa, aod::pidTPCFullPi, aod::pidTOFFullPi, aod::McTrackLabels, aod::pidTOFbeta, aod::TrackSelectionExtension >>;
669675 using EventMCGenerated = soa::Join<aod::McCollisions, aod::MultsExtraMC>; // aod::CentNGlobals, aod::CentNTPVs, aod::CentMFTs
670676
671677 // *********Varibles declaration***************
@@ -1005,17 +1011,21 @@ struct Kstarqa {
10051011 using BinningTypeVertexContributor = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0M>;
10061012 using BinningTypeFT0A = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0A>;
10071013 using BinningTypeFV0A = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFV0A>;
1014+ using BinningTypeMC = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0M>;
10081015
10091016 BinningTypeVertexContributor binningOnPositions{{axisVertex, axisMultiplicity}, true };
10101017 BinningTypeCentralityM binningOnCentrality{{axisVertex, axisMultiplicity}, true };
10111018 BinningTypeFT0A binningOnFT0A{{axisVertex, axisMultiplicity}, true };
10121019 BinningTypeFV0A binningOnFV0A{{axisVertex, axisMultiplicity}, true };
1020+ BinningTypeMC binningOnMC{{axisVertex, axisMultiplicity}, true };
10131021
10141022 SameKindPair<EventCandidates, TrackCandidates, BinningTypeVertexContributor> pair1{binningOnPositions, selectionConfig.cfgNoMixedEvents , -1 , &cache};
10151023 SameKindPair<EventCandidates, TrackCandidates, BinningTypeCentralityM> pair2{binningOnCentrality, selectionConfig.cfgNoMixedEvents , -1 , &cache};
10161024 SameKindPair<EventCandidates, TrackCandidates, BinningTypeFT0A> pair3{binningOnFT0A, selectionConfig.cfgNoMixedEvents , -1 , &cache};
10171025 SameKindPair<EventCandidates, TrackCandidates, BinningTypeFV0A> pair4{binningOnFV0A, selectionConfig.cfgNoMixedEvents , -1 , &cache};
10181026
1027+ SameKindPair<EventCandidatesMC, TrackCandidatesMC, BinningTypeMC> pairmc{binningOnMC, selectionConfig.cfgNoMixedEvents , -1 , &cache};
1028+
10191029 void processME (EventCandidatesMix const &, TrackCandidates const &)
10201030 {
10211031 // Map estimator to pair and multiplicity accessor
@@ -1061,9 +1071,193 @@ struct Kstarqa {
10611071 runMixing (pair4, [](const auto & c) { return c.centFV0A (); });
10621072 }
10631073 }
1064-
10651074 PROCESS_SWITCH (Kstarqa, processME, " Process Mixed event" , true );
10661075
1076+ void processSEMC (EventCandidatesMC::iterator const & collision, TrackCandidatesMC const & tracks, aod::McParticles const &, aod::McCollisions const & /* mcCollisions*/ )
1077+ {
1078+ int occupancy = collision.trackOccupancyInTimeRange ();
1079+ rEventSelection.fill (HIST (" hOccupancy" ), occupancy);
1080+
1081+ if (!selectionEvent (collision, true )) {
1082+ return ;
1083+ }
1084+
1085+ multiplicity = -1 ;
1086+
1087+ if (cSelectMultEstimator == kFT0M ) {
1088+ multiplicity = collision.centFT0M ();
1089+ } else if (cSelectMultEstimator == kFT0A ) {
1090+ multiplicity = collision.centFT0A ();
1091+ } else if (cSelectMultEstimator == kFT0C ) {
1092+ multiplicity = collision.centFT0C ();
1093+ } else if (cSelectMultEstimator == kFV0A ) {
1094+ multiplicity = collision.centFV0A ();
1095+ } else {
1096+ multiplicity = collision.centFT0M (); // default
1097+ }
1098+
1099+ // Fill the event counter
1100+ if (cQAevents) {
1101+ rEventSelection.fill (HIST (" hVertexZRec" ), collision.posZ ());
1102+ rEventSelection.fill (HIST (" hMultiplicity" ), multiplicity);
1103+ rEventSelection.fill (HIST (" multdist_FT0M" ), collision.multFT0M ());
1104+ // rEventSelection.fill(HIST("multdist_FT0A"), collision.multFT0A());
1105+ // rEventSelection.fill(HIST("multdist_FT0C"), collision.multFT0C());
1106+ // rEventSelection.fill(HIST("hNcontributor"), collision.numContrib());
1107+ }
1108+
1109+ for (const auto & [track1, track2] : combinations (CombinationsFullIndexPolicy (tracks, tracks))) {
1110+ rEventSelection.fill (HIST (" tracksCheckData" ), 0.5 );
1111+ if (!selectionTrack (track1)) {
1112+ continue ;
1113+ }
1114+ if (!selectionTrack (track2)) {
1115+ continue ;
1116+ }
1117+ rEventSelection.fill (HIST (" tracksCheckData" ), 1.5 );
1118+
1119+ if (cQAplots) {
1120+ hPID.fill (HIST (" Before/hNsigmaTPC_Ka_before" ), track1.pt (), track1.tpcNSigmaKa ());
1121+ hPID.fill (HIST (" Before/hNsigmaTOF_Ka_before" ), track1.pt (), track1.tofNSigmaKa ());
1122+ hPID.fill (HIST (" Before/hNsigmaTPC_Pi_before" ), track2.pt (), track2.tpcNSigmaPi ());
1123+ hPID.fill (HIST (" Before/hNsigmaTOF_Pi_before" ), track2.pt (), track2.tofNSigmaPi ());
1124+ hPID.fill (HIST (" Before/hNsigma_TPC_TOF_Ka_before" ), track1.tpcNSigmaKa (), track1.tofNSigmaKa ());
1125+ hPID.fill (HIST (" Before/hNsigma_TPC_TOF_Pi_before" ), track2.tpcNSigmaPi (), track2.tofNSigmaPi ());
1126+
1127+ hPID.fill (HIST (" Before/hTPCnsigKa_mult_pt" ), track1.tpcNSigmaKa (), multiplicity, track1.pt ());
1128+ hPID.fill (HIST (" Before/hTPCnsigPi_mult_pt" ), track2.tpcNSigmaPi (), multiplicity, track2.pt ());
1129+ hPID.fill (HIST (" Before/hTOFnsigKa_mult_pt" ), track1.tofNSigmaKa (), multiplicity, track1.pt ());
1130+ hPID.fill (HIST (" Before/hTOFnsigPi_mult_pt" ), track2.tofNSigmaKa (), multiplicity, track2.pt ());
1131+
1132+ hOthers.fill (HIST (" hCRFC_before" ), track1.tpcCrossedRowsOverFindableCls ());
1133+ hOthers.fill (HIST (" dE_by_dx_TPC" ), track1.p (), track1.tpcSignal ());
1134+ hOthers.fill (HIST (" hphi" ), track1.phi ());
1135+
1136+ if (track1.sign () < 0 ) {
1137+ hPID.fill (HIST (" Before/h1PID_TPC_neg_kaon" ), track1.tpcNSigmaKa ());
1138+ hPID.fill (HIST (" Before/h1PID_TPC_neg_pion" ), track2.tpcNSigmaPi ());
1139+ hPID.fill (HIST (" Before/h1PID_TOF_neg_kaon" ), track1.tofNSigmaKa ());
1140+ hPID.fill (HIST (" Before/h1PID_TOF_neg_pion" ), track2.tofNSigmaPi ());
1141+ } else {
1142+ hPID.fill (HIST (" Before/h1PID_TPC_pos_kaon" ), track1.tpcNSigmaKa ());
1143+ hPID.fill (HIST (" Before/h1PID_TPC_pos_pion" ), track2.tpcNSigmaPi ());
1144+ hPID.fill (HIST (" Before/h1PID_TOF_pos_kaon" ), track1.tofNSigmaKa ());
1145+ hPID.fill (HIST (" Before/h1PID_TOF_pos_pion" ), track2.tofNSigmaPi ());
1146+ }
1147+ }
1148+
1149+ if (cQAevents) {
1150+ rEventSelection.fill (HIST (" hDcaxy" ), track1.dcaXY ());
1151+ rEventSelection.fill (HIST (" hDcaz" ), track1.dcaZ ());
1152+ }
1153+
1154+ // since we are using combinations full index policy, so repeated pairs are allowed, so we can check one with Kaon and other with pion
1155+ if (!applypTdepPID && !selectionPID (track1, 1 )) // Track 1 is checked with Kaon
1156+ continue ;
1157+ if (!applypTdepPID && !selectionPID (track2, 0 )) // Track 2 is checked with Pion
1158+ continue ;
1159+
1160+ if (applypTdepPID && !selectionPIDNew (track1, 1 )) // Track 1 is checked with Kaon
1161+ continue ;
1162+ if (applypTdepPID && !selectionPIDNew (track2, 0 )) // Track 2 is checked with Pion
1163+ continue ;
1164+
1165+ rEventSelection.fill (HIST (" tracksCheckData" ), 2.5 );
1166+
1167+ if (cFakeTrack && isFakeTrack (track1, 1 )) // Kaon
1168+ continue ;
1169+ if (cFakeTrack && isFakeTrack (track2, 0 )) // Pion
1170+ continue ;
1171+
1172+ // if (cMID) {
1173+ // if (cMIDselectionPID(track1, 0)) // Kaon misidentified as pion
1174+ // continue;
1175+ // if (cMIDselectionPID(track1, 2)) // Kaon misidentified as proton
1176+ // continue;
1177+ // if (cMIDselectionPID(track2, 1)) // Pion misidentified as kaon
1178+ // continue;
1179+ // }
1180+
1181+ rEventSelection.fill (HIST (" tracksCheckData" ), 3.5 );
1182+
1183+ if (cQAplots) {
1184+ hPID.fill (HIST (" After/hDcaxyPi" ), track2.dcaXY ());
1185+ hPID.fill (HIST (" After/hDcaxyKa" ), track1.dcaXY ());
1186+ hPID.fill (HIST (" After/hDcazPi" ), track2.dcaZ ());
1187+ hPID.fill (HIST (" After/hDcazKa" ), track1.dcaZ ());
1188+
1189+ hPID.fill (HIST (" After/hTPCnsigKa_mult_pt" ), track1.tpcNSigmaKa (), multiplicity, track1.pt ());
1190+ hPID.fill (HIST (" After/hTPCnsigPi_mult_pt" ), track2.tpcNSigmaPi (), multiplicity, track2.pt ());
1191+ hPID.fill (HIST (" After/hTOFnsigKa_mult_pt" ), track1.tofNSigmaKa (), multiplicity, track1.pt ());
1192+ hPID.fill (HIST (" After/hTOFnsigPi_mult_pt" ), track2.tofNSigmaKa (), multiplicity, track2.pt ());
1193+ hOthers.fill (HIST (" hEta_after" ), track1.eta ());
1194+ hOthers.fill (HIST (" hCRFC_after" ), track1.tpcCrossedRowsOverFindableCls ());
1195+ hPID.fill (HIST (" After/hNsigmaKaonTPC_after" ), track1.pt (), track1.tpcNSigmaKa ());
1196+ hPID.fill (HIST (" After/hNsigmaKaonTOF_after" ), track1.pt (), track1.tofNSigmaKa ());
1197+ hPID.fill (HIST (" After/hNsigmaPionTPC_after" ), track2.pt (), track2.tpcNSigmaPi ());
1198+ hPID.fill (HIST (" After/hNsigmaPionTOF_after" ), track2.pt (), track2.tofNSigmaPi ());
1199+ hPID.fill (HIST (" After/hNsigma_TPC_TOF_Ka_after" ), track1.tpcNSigmaKa (), track1.tofNSigmaKa ());
1200+ hPID.fill (HIST (" After/hNsigma_TPC_TOF_Pi_after" ), track2.tpcNSigmaPi (), track2.tofNSigmaPi ());
1201+ }
1202+
1203+ if (track1.globalIndex () == track2.globalIndex ())
1204+ continue ;
1205+
1206+ rEventSelection.fill (HIST (" tracksCheckData" ), 4.5 );
1207+
1208+ daughter1 = ROOT::Math::PxPyPzMVector (track1.px (), track1.py (), track1.pz (), massKa);
1209+ daughter2 = ROOT::Math::PxPyPzMVector (track2.px (), track2.py (), track2.pz (), massPi);
1210+ mother = daughter1 + daughter2; // Kstar meson
1211+
1212+ if (selectionConfig.isApplyCutsOnMother ) {
1213+ if (mother.Pt () >= selectionConfig.cMaxPtMotherCut ) // excluding candidates in overflow
1214+ continue ;
1215+ if (mother.M () >= selectionConfig.cMaxMinvMotherCut ) // excluding candidates in overflow
1216+ continue ;
1217+ }
1218+
1219+ hOthers.fill (HIST (" hKstar_Rap" ), mother.Rapidity ());
1220+ hOthers.fill (HIST (" hKstar_Eta" ), mother.Eta ());
1221+
1222+ isMix = false ;
1223+ fillInvMass (daughter1, daughter2, mother, multiplicity, isMix, track1, track2);
1224+ }
1225+ }
1226+ PROCESS_SWITCH (Kstarqa, processSEMC, " Process same event in MC" , true );
1227+
1228+ void processMEMC (EventCandidatesMC const &, TrackCandidatesMC const &, aod::McParticles const &, aod::McCollisions const &)
1229+ {
1230+ for (const auto & [c1, tracks1, c2, tracks2] : pairmc) {
1231+
1232+ if (!selectionEvent (c1, false ) || !selectionEvent (c2, false )) {
1233+ continue ;
1234+ }
1235+
1236+ // multiplicity = multiplicityGetter(c1);
1237+ multiplicity = c1.centFT0M (); // default, can be changed later
1238+
1239+ for (const auto & [t1, t2] : o2::soa::combinations (o2::soa::CombinationsFullIndexPolicy (tracks1, tracks2))) {
1240+ if (!selectionTrack (t1) || !selectionTrack (t2))
1241+ continue ;
1242+ if (!selectionPID (t1, 1 ) || !selectionPID (t2, 0 ))
1243+ continue ;
1244+
1245+ daughter1 = ROOT::Math::PxPyPzMVector (t1.px (), t1.py (), t1.pz (), massKa);
1246+ daughter2 = ROOT::Math::PxPyPzMVector (t2.px (), t2.py (), t2.pz (), massPi);
1247+ mother = daughter1 + daughter2;
1248+
1249+ isMix = true ;
1250+
1251+ if (std::abs (mother.Rapidity ()) < selectionConfig.rapidityMotherData ) {
1252+ fillInvMass (daughter1, daughter2, mother, multiplicity, isMix, t1, t2);
1253+ }
1254+ }
1255+ }
1256+ }
1257+ PROCESS_SWITCH (Kstarqa, processMEMC, " Process mixed-event in MC" , true );
1258+
1259+ Service<o2::framework::O2DatabasePDG> pdgDB;
1260+
10671261 // void processGen(EventMCGenerated::iterator const& mcCollision, aod::McParticles const& mcParticles, const soa::SmallGroups<EventCandidatesMC>& collisions)
10681262 void processGen (aod::McCollision const & mcCollision, aod::McParticles const & mcParticles, const soa::SmallGroups<EventCandidatesMC>& collisions)
10691263 {
@@ -1086,6 +1280,15 @@ struct Kstarqa {
10861280 multiplicity = -1.0 ;
10871281 // float impactParameter = mcCollision.impactParameter();
10881282
1283+ bool isINELgt0true = false ;
1284+
1285+ if (pwglf::isINELgtNmc (mcParticles, 0 , pdgDB)) {
1286+ isINELgt0true = true ;
1287+ }
1288+ if (selectionConfig.isINELgt0 && !isINELgt0true) {
1289+ return ;
1290+ }
1291+
10891292 // if (selectionConfig.isINELgt0 && !mcCollision.isInelGt0()) {
10901293 // return;
10911294 // }
@@ -1129,16 +1332,12 @@ struct Kstarqa {
11291332
11301333 for (const auto & mcParticle : mcParticles) {
11311334 if (std::abs (mcParticle.y ()) < selectionConfig.rapidityMotherData && std::abs (mcParticle.pdgCode ()) == o2::constants::physics::kK0Star892 ) {
1132- // if (inelgt0MCgen) {
11331335 hInvMass.fill (HIST (" hAllKstarGenCollisisons" ), multiplicity, mcParticle.pt ());
1134- // }
11351336 }
11361337 }
11371338
11381339 const auto evtReconstructedAndSelected = std::find (selectedEvents.begin (), selectedEvents.end (), mcCollision.globalIndex ()) != selectedEvents.end ();
1139- // if (inelgt0MCgen) {
11401340 hInvMass.fill (HIST (" hAllGenCollisions" ), multiplicity);
1141- // }
11421341 if (!cAllGenCollisions && !evtReconstructedAndSelected) { // Check that the event is reconstructed and that the reconstructed events pass the selection
11431342 return ;
11441343 }
@@ -1201,6 +1400,15 @@ struct Kstarqa {
12011400 // return;
12021401 // }
12031402
1403+ bool isINELgt0true = false ;
1404+
1405+ if (pwglf::isINELgtNmc (mcParticles, 0 , pdgDB)) {
1406+ isINELgt0true = true ;
1407+ }
1408+ if (selectionConfig.isINELgt0 && !isINELgt0true) {
1409+ return ;
1410+ }
1411+
12041412 auto impactPar = mcCollision.impactParameter ();
12051413 hInvMass.fill (HIST (" MCcorrections/hImpactParameterGen" ), impactPar);
12061414
0 commit comments