@@ -78,6 +78,8 @@ struct Derivedupcanalysis {
7878 Configurable<bool > analyseOmega{" analyseOmega" , true , " process Omega-like candidates" };
7979 Configurable<bool > analyseAntiOmega{" analyseAntiOmega" , true , " process AntiOmega-like candidates" };
8080
81+ Configurable<std::vector<int >> generatorIds{" generatorIds" , std::vector<int >{-1 }, " MC generatorIds to process" };
82+
8183 // Event selections
8284 Configurable<bool > rejectITSROFBorder{" rejectITSROFBorder" , true , " reject events at ITS ROF border" };
8385 Configurable<bool > rejectTFBorder{" rejectTFBorder" , true , " reject events at TF border" };
@@ -861,6 +863,7 @@ struct Derivedupcanalysis {
861863 histos.add (" eventQA/hTracksPVeta1VsTracksGlobal" , " hTracksPVeta1VsTracksGlobal" , kTH3F , {axisNTracksPVeta1, axisNTracksGlobal, axisSelGap});
862864 histos.add (" eventQA/hCentralityVsTracksGlobal" , " hCentralityVsTracksGlobal" , kTH3F , {axisFT0Cqa, axisNTracksGlobal, axisSelGap});
863865 histos.add (" eventQA/hGapSide" , " Gap side; Entries" , kTH1F , {{5 , -0.5 , 4.5 }});
866+ histos.add (" eventQA/hRawGapSide" , " Raw Gap side; Entries" , kTH1F , {{6 , -1.5 , 4.5 }});
864867 histos.add (" eventQA/hSelGapSide" , " Selected gap side; Entries" , kTH1F , {axisSelGap});
865868 histos.add (" eventQA/hPosX" , " Vertex position in x" , kTH2F , {{100 , -0.1 , 0.1 }, axisSelGap});
866869 histos.add (" eventQA/hPosY" , " Vertex position in y" , kTH2F , {{100 , -0.1 , 0.1 }, axisSelGap});
@@ -977,14 +980,16 @@ struct Derivedupcanalysis {
977980 }
978981
979982 template <typename TCollision>
980- int getGapSide (TCollision const & collision)
983+ int getGapSide (TCollision const & collision, bool fillQA )
981984 {
982985 int selGapSide = sgSelector.trueGap (collision, upcCuts.fv0a , upcCuts.ft0a , upcCuts.ft0c , upcCuts.zdc );
983- histos.fill (HIST (" eventQA/hGapSide" ), collision.gapSide ());
984- histos.fill (HIST (" eventQA/hSelGapSide" ), selGapSide);
985- histos.fill (HIST (" eventQA/hFT0" ), collision.totalFT0AmplitudeA (), collision.totalFT0AmplitudeC (), selGapSide);
986- histos.fill (HIST (" eventQA/hFDD" ), collision.totalFDDAmplitudeA (), collision.totalFDDAmplitudeC (), selGapSide);
987- histos.fill (HIST (" eventQA/hZN" ), collision.energyCommonZNA (), collision.energyCommonZNC (), selGapSide);
986+ if (fillQA) {
987+ histos.fill (HIST (" eventQA/hGapSide" ), collision.gapSide ());
988+ histos.fill (HIST (" eventQA/hSelGapSide" ), selGapSide);
989+ histos.fill (HIST (" eventQA/hFT0" ), collision.totalFT0AmplitudeA (), collision.totalFT0AmplitudeC (), selGapSide);
990+ histos.fill (HIST (" eventQA/hFDD" ), collision.totalFDDAmplitudeA (), collision.totalFDDAmplitudeC (), selGapSide);
991+ histos.fill (HIST (" eventQA/hZN" ), collision.energyCommonZNA (), collision.energyCommonZNC (), selGapSide);
992+ }
988993 return selGapSide;
989994 }
990995
@@ -1014,7 +1019,7 @@ struct Derivedupcanalysis {
10141019 float qaBin;
10151020 };
10161021
1017- const std::array<SelectionCheck, 15 > checks = {{
1022+ const std::array<SelectionCheck, 14 > checks = {{
10181023 {true , true , 0.0 }, // All collisions
10191024 {requireIsTriggerTVX, collision.selection_bit (aod::evsel::kIsTriggerTVX ), 1.0 }, // Triggered by FT0M
10201025 {true , std::fabs (collision.posZ ()) < maxZVtxPosition, 2.0 }, // Vertex-Z selected
@@ -1029,7 +1034,6 @@ struct Derivedupcanalysis {
10291034 {requireNoCollInTimeRangeNarrow, collision.selection_bit (o2::aod::evsel::kNoCollInTimeRangeNarrow ), 11.0 }, // No collision within +-4 µs
10301035 {minOccupancy >= 0 , collision.trackOccupancyInTimeRange () >= minOccupancy, 12.0 }, // Above min occupancy
10311036 {maxOccupancy > 0 , collision.trackOccupancyInTimeRange () < maxOccupancy, 13.0 }, // Below max occupancy
1032- {studyUPConly, collision.isUPC (), 14.0 }, // Study UPC collisions only
10331037 }};
10341038
10351039 for (const auto & check : checks) {
@@ -1041,6 +1045,12 @@ struct Derivedupcanalysis {
10411045 }
10421046 }
10431047
1048+ if (studyUPConly && !collision.isUPC ()) {
1049+ return false ;
1050+ } else if (collision.isUPC ()) {
1051+ histos.fill (HIST (" eventQA/hEventSelection" ), 14.0 ); // is UPC compatible
1052+ }
1053+
10441054 // Additional check for UPC collision flag
10451055 if (useUPCflag && collision.flags () < 1 ) {
10461056 return false ;
@@ -1511,27 +1521,27 @@ struct Derivedupcanalysis {
15111521 const int pdgV0 = v0.pdgCode ();
15121522 const bool isPhysPrim = v0.isPhysicalPrimary ();
15131523
1514- const bool isPositiveProton = (pdgPos == 2212 );
1515- const bool isPositivePion = (pdgPos == 211 ) || (doTreatPiToMuon && pdgPos == - 13 );
1516- const bool isNegativeProton = (pdgNeg == - 2212 );
1517- const bool isNegativePion = (pdgNeg == - 211 ) || (doTreatPiToMuon && pdgNeg == 13 );
1524+ const bool isPositiveProton = (pdgPos == PDG_t:: kProton );
1525+ const bool isPositivePion = (pdgPos == PDG_t:: kPiPlus ) || (doTreatPiToMuon && pdgPos == PDG_t:: kMuonPlus );
1526+ const bool isNegativeProton = (pdgNeg == kProtonBar );
1527+ const bool isNegativePion = (pdgNeg == PDG_t:: kPiMinus ) || (doTreatPiToMuon && pdgNeg == PDG_t:: kMuonMinus );
15181528
15191529 switch (pdgV0) {
1520- case 310 : // K0Short
1530+ case PDG_t:: kK0Short : // K0Short
15211531 if (isPositivePion && isNegativePion) {
15221532 bitMap.set (selConsiderK0Short);
15231533 if (isPhysPrim)
15241534 bitMap.set (selPhysPrimK0Short);
15251535 }
15261536 break ;
1527- case 3122 : // Lambda
1537+ case PDG_t:: kLambda0 : // Lambda
15281538 if (isPositiveProton && isNegativePion) {
15291539 bitMap.set (selConsiderLambda);
15301540 if (isPhysPrim)
15311541 bitMap.set (selPhysPrimLambda);
15321542 }
15331543 break ;
1534- case - 3122 : // AntiLambda
1544+ case PDG_t:: kLambda0Bar : // AntiLambda
15351545 if (isPositivePion && isNegativeProton) {
15361546 bitMap.set (selConsiderAntiLambda);
15371547 if (isPhysPrim)
@@ -1541,6 +1551,60 @@ struct Derivedupcanalysis {
15411551 }
15421552 }
15431553
1554+ template <typename TCasc>
1555+ void computeCascadeMCAssociation (const TCasc& casc, std::bitset<kSelNum >& bitMap)
1556+ {
1557+ const int pdgPos = casc.pdgCodePositive ();
1558+ const int pdgNeg = casc.pdgCodeNegative ();
1559+ const int pdgBach = casc.pdgCodeBachelor ();
1560+ const int pdgCasc = casc.pdgCode ();
1561+ const bool isPhysPrim = casc.isPhysicalPrimary ();
1562+
1563+ const bool isPositiveProton = (pdgPos == PDG_t::kProton );
1564+
1565+ const bool isPositivePion = (pdgPos == PDG_t::kPiPlus );
1566+ const bool isBachelorPositivePion = (pdgBach == PDG_t::kPiPlus );
1567+
1568+ const bool isNegativeProton = (pdgNeg == kProtonBar );
1569+
1570+ const bool isNegativePion = (pdgNeg == PDG_t::kPiMinus );
1571+ const bool isBachelorNegativePion = (pdgBach == PDG_t::kPiMinus );
1572+
1573+ const bool isBachelorPositiveKaon = (pdgBach == PDG_t::kKPlus );
1574+ const bool isBachelorNegativeKaon = (pdgBach == PDG_t::kKMinus );
1575+
1576+ switch (pdgCasc) {
1577+ case PDG_t::kXiMinus : // Xi
1578+ if (isPositiveProton && isNegativePion && isBachelorNegativePion) {
1579+ bitMap.set (selConsiderXi);
1580+ if (isPhysPrim)
1581+ bitMap.set (selPhysPrimXi);
1582+ }
1583+ break ;
1584+ case PDG_t::kXiPlusBar : // Anti-Xi
1585+ if (isNegativeProton && isPositivePion && isBachelorPositivePion) {
1586+ bitMap.set (selConsiderAntiXi);
1587+ if (isPhysPrim)
1588+ bitMap.set (selPhysPrimAntiXi);
1589+ }
1590+ break ;
1591+ case PDG_t::kOmegaMinus : // Omega
1592+ if (isPositiveProton && isNegativePion && isBachelorNegativeKaon) {
1593+ bitMap.set (selConsiderOmega);
1594+ if (isPhysPrim)
1595+ bitMap.set (selPhysPrimOmega);
1596+ }
1597+ break ;
1598+ case PDG_t::kOmegaPlusBar : // Anti-Omega
1599+ if (isNegativeProton && isPositivePion && isBachelorPositiveKaon) {
1600+ bitMap.set (selConsiderAntiOmega);
1601+ if (isPhysPrim)
1602+ bitMap.set (selPhysPrimAntiOmega);
1603+ }
1604+ break ;
1605+ }
1606+ }
1607+
15441608 template <typename TV0, typename TCollision>
15451609 void analyseV0Candidate (TV0 const & v0, TCollision const & coll, int const & gap, std::bitset<kSelNum > const & selMap)
15461610 {
@@ -1585,6 +1649,9 @@ struct Derivedupcanalysis {
15851649 std::vector<int > listBestCollisionIds (mcCollisions.size (), -1 );
15861650
15871651 for (auto const & mcCollision : mcCollisions) {
1652+ if (std::find (generatorIds->begin (), generatorIds->end (), mcCollision.generatorsID ()) == generatorIds->end ()) {
1653+ continue ;
1654+ }
15881655 auto groupedCollisions = collisions.sliceBy (perMcCollision, mcCollision.globalIndex ());
15891656 // Find the collision with the biggest nbr of PV contributors
15901657 // Follows what was done here: https://github.com/AliceO2Group/O2Physics/blob/master/Common/TableProducer/mcCollsExtra.cxx#L93
@@ -1595,7 +1662,7 @@ struct Derivedupcanalysis {
15951662 continue ;
15961663 }
15971664
1598- int selGapSide = collision.isUPC () ? getGapSide (collision) : -1 ;
1665+ int selGapSide = collision.isUPC () ? getGapSide (collision, false ) : -1 ;
15991666 if (studyUPConly && (selGapSide < -0.5 ))
16001667 continue ;
16011668
@@ -1613,6 +1680,11 @@ struct Derivedupcanalysis {
16131680 void fillGenMCHistogramsQA (StraMCCollisionsFull const & mcCollisions, StraCollisonsFullMC const & collisions)
16141681 {
16151682 for (auto const & mcCollision : mcCollisions) {
1683+ // LOGF(info, "Generator ID is %i", mcCollision.generatorsID());
1684+ if (std::find (generatorIds->begin (), generatorIds->end (), mcCollision.generatorsID ()) == generatorIds->end ()) {
1685+ continue ;
1686+ }
1687+
16161688 histos.fill (HIST (" eventQA/mc/hEventSelectionMC" ), 0.0 );
16171689 histos.fill (HIST (" eventQA/mc/hMCNParticlesEta10" ), mcCollision.multMCNParticlesEta10 (), 0 /* all gen. events*/ );
16181690
@@ -1636,7 +1708,7 @@ struct Derivedupcanalysis {
16361708 continue ;
16371709 }
16381710
1639- int selGapSide = collision.isUPC () ? getGapSide (collision) : -1 ;
1711+ int selGapSide = collision.isUPC () ? getGapSide (collision, false ) : -1 ;
16401712 if (studyUPConly && (selGapSide < -0.5 ))
16411713 continue ;
16421714
@@ -1670,7 +1742,9 @@ struct Derivedupcanalysis {
16701742 return ;
16711743 } // event is accepted
16721744
1673- int selGapSide = collision.isUPC () ? getGapSide (collision) : -1 ;
1745+ histos.fill (HIST (" eventQA/hRawGapSide" ), collision.gapSide ());
1746+
1747+ int selGapSide = collision.isUPC () ? getGapSide (collision, true ) : -1 ;
16741748 if (studyUPConly && (selGapSide < -0.5 ))
16751749 return ;
16761750
@@ -1697,15 +1771,22 @@ struct Derivedupcanalysis {
16971771 StraMCCollisionsFull const &,
16981772 V0MCCoresFull const &)
16991773 {
1774+ if (!collision.has_straMCCollision ()) {
1775+ histos.fill (HIST (" eventQA/mc/hFakeEvents" ), 0 ); // no assoc. MC collisions
1776+ } else {
1777+ const auto & mcCollision = collision.straMCCollision_as <StraMCCollisionsFull>();
1778+ if (std::find (generatorIds->begin (), generatorIds->end (), mcCollision.generatorsID ()) == generatorIds->end ()) {
1779+ return ;
1780+ }
1781+ }
1782+
17001783 if (!acceptEvent (collision, true )) {
17011784 return ;
17021785 } // event is accepted
17031786
1704- if (!collision.has_straMCCollision ()) {
1705- histos.fill (HIST (" eventQA/mc/hFakeEvents" ), 0 ); // no assoc. MC collisions
1706- }
1787+ histos.fill (HIST (" eventQA/hRawGapSide" ), collision.gapSide ());
17071788
1708- int selGapSide = collision.isUPC () ? getGapSide (collision) : -1 ;
1789+ int selGapSide = collision.isUPC () ? getGapSide (collision, true ) : -1 ;
17091790 if (studyUPConly && (selGapSide < -0.5 ))
17101791 return ;
17111792
@@ -1748,7 +1829,9 @@ struct Derivedupcanalysis {
17481829 return ;
17491830 } // event is accepted
17501831
1751- int selGapSide = collision.isUPC () ? getGapSide (collision) : -1 ;
1832+ histos.fill (HIST (" eventQA/hRawGapSide" ), collision.gapSide ());
1833+
1834+ int selGapSide = collision.isUPC () ? getGapSide (collision, true ) : -1 ;
17521835 if (studyUPConly && (selGapSide < -0.5 ))
17531836 return ;
17541837
@@ -1764,6 +1847,60 @@ struct Derivedupcanalysis {
17641847 } // end casc loop
17651848 }
17661849
1850+ void processCascadesMC (StraCollisonFullMC const & collision,
1851+ CascadeCandidatesMC const & fullCascades,
1852+ DauTracks const &,
1853+ aod::MotherMCParts const &,
1854+ StraMCCollisionsFull const &,
1855+ CascMCCoresFull const &)
1856+ {
1857+ if (!collision.has_straMCCollision ()) {
1858+ histos.fill (HIST (" eventQA/mc/hFakeEvents" ), 0 ); // no assoc. MC collisions
1859+ } else {
1860+ const auto & mcCollision = collision.straMCCollision_as <StraMCCollisionsFull>();
1861+ if (std::find (generatorIds->begin (), generatorIds->end (), mcCollision.generatorsID ()) == generatorIds->end ()) {
1862+ return ;
1863+ }
1864+ }
1865+
1866+ if (!acceptEvent (collision, true )) {
1867+ return ;
1868+ } // event is accepted
1869+
1870+ histos.fill (HIST (" eventQA/hRawGapSide" ), collision.gapSide ());
1871+
1872+ int selGapSide = collision.isUPC () ? getGapSide (collision, true ) : -1 ;
1873+ if (studyUPConly && (selGapSide < -0.5 ))
1874+ return ;
1875+
1876+ fillHistogramsQA (collision, selGapSide);
1877+
1878+ if (collision.has_straMCCollision ()) {
1879+ const auto & mcCollision = collision.straMCCollision_as <StraMCCollisionsFull>();
1880+ histos.fill (HIST (" eventQA/mc/hNTracksGlobalvsMCNParticlesEta10rec" ), collision.multNTracksGlobal (), mcCollision.multMCNParticlesEta10 ());
1881+ histos.fill (HIST (" eventQA/mc/hNTracksPVeta1vsMCNParticlesEta10rec" ), collision.multNTracksPVeta1 (), mcCollision.multMCNParticlesEta10 ());
1882+ histos.fill (HIST (" eventQA/mc/hNTracksGlobalvstotalMultMCParticles" ), collision.multNTracksGlobal (), mcCollision.totalMultMCParticles ());
1883+ histos.fill (HIST (" eventQA/mc/hNTracksPVeta1vstotalMultMCParticles" ), collision.multNTracksPVeta1 (), mcCollision.totalMultMCParticles ());
1884+ }
1885+
1886+ for (const auto & casc : fullCascades) {
1887+ std::bitset<kSelNum > selMap = computeBitmapCascade (casc, collision);
1888+
1889+ if (doMCAssociation) {
1890+ if (casc.has_cascMCCore ()) {
1891+ const auto & cascMC = casc.cascMCCore_as <CascMCCoresFull>();
1892+ computeCascadeMCAssociation (cascMC, selMap);
1893+ }
1894+ } else {
1895+ // the candidate may belong to any particle species
1896+ setBits (selMap, {selConsiderXi, selConsiderAntiXi, selConsiderOmega, selConsiderAntiOmega,
1897+ selPhysPrimXi, selPhysPrimAntiXi, selPhysPrimOmega, selPhysPrimAntiOmega});
1898+ }
1899+
1900+ analyseCascCandidate (casc, collision, selGapSide, selMap);
1901+ } // end casc loop
1902+ }
1903+
17671904 void processGenerated (StraMCCollisionsFull const & mcCollisions,
17681905 V0MCCoresFull const & V0MCCores,
17691906 CascMCCoresFull const & CascMCCores,
@@ -1780,9 +1917,9 @@ struct Derivedupcanalysis {
17801917 // Kinematics (|y| < rapidityCut)
17811918 float pTmc = v0MC.ptMC ();
17821919 float ymc = 1e3 ;
1783- if (v0MC.pdgCode () == 310 )
1920+ if (v0MC.pdgCode () == PDG_t:: kK0Short )
17841921 ymc = v0MC.rapidityMC (0 );
1785- else if (std::abs (v0MC.pdgCode ()) == 3122 )
1922+ else if (std::abs (v0MC.pdgCode ()) == PDG_t:: kLambda0 )
17861923 ymc = v0MC.rapidityMC (1 );
17871924 if (std::abs (ymc) > rapidityCut)
17881925 continue ;
@@ -1792,6 +1929,11 @@ struct Derivedupcanalysis {
17921929 if (std::abs (mcCollision.posZ ()) > maxZVtxPosition)
17931930 continue ;
17941931
1932+ // Collision is of the proccess of interest
1933+ if (std::find (generatorIds->begin (), generatorIds->end (), mcCollision.generatorsID ()) == generatorIds->end ()) {
1934+ continue ;
1935+ }
1936+
17951937 float centrality = -1 .f ;
17961938 int nTracksGlobal = -1 ;
17971939
@@ -1802,13 +1944,13 @@ struct Derivedupcanalysis {
18021944 }
18031945
18041946 // Fill histograms
1805- if (v0MC.pdgCode () == 310 ) {
1947+ if (v0MC.pdgCode () == PDG_t:: kK0Short ) {
18061948 histos.fill (HIST (kParticlenames [0 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18071949 }
1808- if (v0MC.pdgCode () == 3122 ) {
1950+ if (v0MC.pdgCode () == PDG_t:: kLambda0 ) {
18091951 histos.fill (HIST (kParticlenames [1 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18101952 }
1811- if (v0MC.pdgCode () == - 3122 ) {
1953+ if (v0MC.pdgCode () == PDG_t:: kLambda0Bar ) {
18121954 histos.fill (HIST (kParticlenames [2 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18131955 }
18141956 } // V0 end
@@ -1821,9 +1963,9 @@ struct Derivedupcanalysis {
18211963 // Kinematics (|y| < rapidityCut)
18221964 float pTmc = cascMC.ptMC ();
18231965 float ymc = 1e3 ;
1824- if (std::abs (cascMC.pdgCode ()) == 3312 )
1966+ if (std::abs (cascMC.pdgCode ()) == PDG_t:: kXiMinus )
18251967 ymc = RecoDecay::y (std::array{cascMC.pxMC (), cascMC.pyMC (), cascMC.pzMC ()}, o2::constants::physics::MassXiMinus);
1826- else if (std::abs (cascMC.pdgCode ()) == 3334 )
1968+ else if (std::abs (cascMC.pdgCode ()) == PDG_t:: kOmegaMinus )
18271969 ymc = RecoDecay::y (std::array{cascMC.pxMC (), cascMC.pyMC (), cascMC.pzMC ()}, o2::constants::physics::MassOmegaMinus);
18281970 if (std::abs (ymc) > rapidityCut)
18291971 continue ;
@@ -1842,16 +1984,16 @@ struct Derivedupcanalysis {
18421984 }
18431985
18441986 // Fill histograms
1845- if (cascMC.pdgCode () == 3312 ) {
1987+ if (cascMC.pdgCode () == PDG_t:: kXiMinus ) {
18461988 histos.fill (HIST (kParticlenames [3 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18471989 }
1848- if (cascMC.pdgCode () == - 3312 ) {
1990+ if (cascMC.pdgCode () == PDG_t:: kXiPlusBar ) {
18491991 histos.fill (HIST (kParticlenames [4 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18501992 }
1851- if (cascMC.pdgCode () == 3334 ) {
1993+ if (cascMC.pdgCode () == PDG_t:: kOmegaMinus ) {
18521994 histos.fill (HIST (kParticlenames [5 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18531995 }
1854- if (cascMC.pdgCode () == - 3334 ) {
1996+ if (cascMC.pdgCode () == PDG_t:: kOmegaPlusBar ) {
18551997 histos.fill (HIST (kParticlenames [6 ]) + HIST (" /mc/h6dGen" ), centrality, nTracksGlobal, mcCollision.multMCNParticlesEta10 (), pTmc, static_cast <int >(upcCuts.genGapSide ), ymc);
18561998 }
18571999 } // Cascade end
0 commit comments