@@ -222,7 +222,6 @@ struct statPromptPhoton {
222222 histos.add (" REC_cluster_both_phiQA" , " REC_cluster_both_phiQA" , kTH1F , {{640 * 2 , 0 , 2 * TMath::Pi ()}});
223223 histos.add (" REC_cluster_both_etaQA" , " REC_cluster_both_etaQA" , kTH1F , {{100 , -1 , 1 }});
224224 histos.add (" REC_cluster_both_energyQA" , " REC_cluster_both_energyQA" , kTH1F , {{82 , -1.0 , 40.0 }});
225-
226225 }
227226 if (cfgGenHistograms) {
228227 histos.add (" GEN_prompt_phiQA" , " GEN_prompt_phiQA" , kTH1F , {{640 * 2 , 0 , 2 * TMath::Pi ()}});
@@ -465,7 +464,7 @@ struct statPromptPhoton {
465464 }; // end of track selection
466465 // ///////////////////////////////////////////////////////////////////////////
467466 // ///////////////////////////////////////////////////////////////////////////
468- // Below is shamelessly stolen from Florian's gammetreeproducer code.
467+ // Below is shamelessly stolen from Florian's gammetreeproducer code.
469468 template <typename T>
470469 T iTopCopy (const T& particle) const
471470 {
@@ -1474,7 +1473,7 @@ struct statPromptPhoton {
14741473
14751474 PROCESS_SWITCH (statPromptPhoton, processData, " processJE data" , false );
14761475
1477- int nEventsGenMC_Simple = 0 ;
1476+ int nEventsGenMC_Simple = 0 ;
14781477 void processMCGen_simple (filteredMCCollisions::iterator const & collision, soa::SmallGroups<soa::Join<aod::JMcCollisionLbs, jfilteredCollisions>> const & recocolls, aod::JMcParticles const & mcParticles, jfilteredMCClusters const &)
14791478 {
14801479 nEventsGenMC_Simple++;
@@ -1487,28 +1486,29 @@ struct statPromptPhoton {
14871486 if (fabs (collision.posZ ()) > cfgVtxCut)
14881487 return ;
14891488
1490- if (cfgGenReqRec){
1489+ if (cfgGenReqRec) {
14911490 if (recocolls.size () <= 0 ) // not reconstructed
1492- return ;
1491+ return ;
14931492 for (auto & recocoll : recocolls) { // poorly reconstructed
1494- if (!recocoll.sel8 ())
1495- return ;
1496- if (fabs (recocoll.posZ ()) > cfgVtxCut)
1497-
1498- return ;
1499- histos.fill (HIST (" GEN_nEvents_simple" ), 1.5 );
1500-
1501- if (cfgEmcTrigger) {
1502- if (!recocoll.isEmcalReadout ())
1503- return ;
1504- }
1505- histos.fill (HIST (" GEN_nEvents_simple" ), 2.5 );
1493+ if (!recocoll.sel8 ())
1494+ return ;
1495+ if (fabs (recocoll.posZ ()) > cfgVtxCut)
1496+
1497+ return ;
1498+ histos.fill (HIST (" GEN_nEvents_simple" ), 1.5 );
1499+
1500+ if (cfgEmcTrigger) {
1501+ if (!recocoll.isEmcalReadout ())
1502+ return ;
1503+ }
1504+ histos.fill (HIST (" GEN_nEvents_simple" ), 2.5 );
15061505 }
15071506 }
15081507
15091508 // First pass: find all status -23 particles
15101509 for (auto & hardParticle : mcParticles) {
1511- if (hardParticle.getGenStatusCode () != -23 ) continue ;
1510+ if (hardParticle.getGenStatusCode () != -23 )
1511+ continue ;
15121512
15131513 bool isPhoton23 = (hardParticle.pdgCode () == 22 );
15141514
@@ -1517,17 +1517,21 @@ struct statPromptPhoton {
15171517 // We search all final-state photons and check if they trace back here//
15181518
15191519 for (auto & mcParticle : mcParticles) {
1520- if (mcParticle.pdgCode () != 22 ) continue ;
1521- if (mcParticle.getGenStatusCode () < 0 ) continue ;
1522- if (std::fabs (mcParticle.getGenStatusCode ()) >= 81 || !mcParticle.isPhysicalPrimary ()) continue ;
1520+ if (mcParticle.pdgCode () != 22 )
1521+ continue ;
1522+ if (mcParticle.getGenStatusCode () < 0 )
1523+ continue ;
1524+ if (std::fabs (mcParticle.getGenStatusCode ()) >= 81 || !mcParticle.isPhysicalPrimary ())
1525+ continue ;
15231526
15241527 // Chase this final-state photon upward
15251528 int chaseindex = -1 ;
15261529 for (auto & mom : mcParticle.mothers_as <aod::JMcParticles>()) {
15271530 chaseindex = mom.globalIndex ();
15281531 break ;
15291532 }
1530- if (chaseindex < 0 ) continue ;
1533+ if (chaseindex < 0 )
1534+ continue ;
15311535
15321536 std::set<int > visited;
15331537 bool chase = true ;
@@ -1536,11 +1540,15 @@ struct statPromptPhoton {
15361540 bool cleanPhotonChain = true ; // all intermediates are photons
15371541
15381542 while (chase) {
1539- if (visited.count (chaseindex)) { chase = false ; break ; }
1543+ if (visited.count (chaseindex)) {
1544+ chase = false ;
1545+ break ;
1546+ }
15401547 visited.insert (chaseindex);
15411548
15421549 for (auto & particle : mcParticles) {
1543- if (particle.globalIndex () != chaseindex) continue ;
1550+ if (particle.globalIndex () != chaseindex)
1551+ continue ;
15441552
15451553 if (particle.globalIndex () == hardParticle.globalIndex ()) {
15461554 reachedThisHard = true ;
@@ -1549,21 +1557,27 @@ struct statPromptPhoton {
15491557 }
15501558
15511559 int abspdg = std::abs (particle.pdgCode ());
1552- if (abspdg > 100 ) hadronInChain = true ;
1553- if (abspdg != 22 ) cleanPhotonChain = false ;
1560+ if (abspdg > 100 )
1561+ hadronInChain = true ;
1562+ if (abspdg != 22 )
1563+ cleanPhotonChain = false ;
15541564
15551565 int nextindex = -1 ;
15561566 for (auto & mom : particle.mothers_as <aod::JMcParticles>()) {
15571567 nextindex = mom.globalIndex ();
15581568 break ;
15591569 }
1560- if (nextindex < 0 ) { chase = false ; }
1561- else { chaseindex = nextindex; }
1570+ if (nextindex < 0 ) {
1571+ chase = false ;
1572+ } else {
1573+ chaseindex = nextindex;
1574+ }
15621575 break ;
15631576 }
15641577 }
15651578
1566- if (!reachedThisHard) continue ;
1579+ if (!reachedThisHard)
1580+ continue ;
15671581
15681582 if (isPhoton23 && cleanPhotonChain) {
15691583 // Case 1: -23 photon, clean photon chain — direct prompt
@@ -1586,7 +1600,7 @@ struct statPromptPhoton {
15861600 nEventsRecMC_simple++;
15871601 if (cfgDebug) {
15881602 if ((nEventsRecMC_simple + 1 ) % 10000 == 0 ) {
1589- std::cout << " Processed JE Rec MC Events: " << nEventsRecMC_simple << std::endl;
1603+ std::cout << " Processed JE Rec MC Events: " << nEventsRecMC_simple << std::endl;
15901604 }
15911605 }
15921606 histos.fill (HIST (" REC_nEvents" ), 0.5 );
@@ -1597,7 +1611,7 @@ struct statPromptPhoton {
15971611 histos.fill (HIST (" REC_nEvents" ), 1.5 );
15981612 if (cfgEmcTrigger) {
15991613 if (!collision.isEmcalReadout ())
1600- return ;
1614+ return ;
16011615 }
16021616 histos.fill (HIST (" REC_nEvents" ), 2.5 );
16031617 if (!jetderiveddatautilities::selectTrigger (collision, triggerMaskBits))
@@ -1606,103 +1620,116 @@ struct statPromptPhoton {
16061620 for (auto & mccluster : mcclusters) {
16071621 histos.fill (HIST (" REC_M02_BC" ), mccluster.m02 ());
16081622 if (mccluster.m02 () < cfgLowM02)
1609- continue ;
1623+ continue ;
16101624 if (mccluster.m02 () > cfgHighM02)
1611- continue ;
1625+ continue ;
16121626 if (mccluster.energy () < cfgLowClusterE)
1613- continue ;
1627+ continue ;
16141628 if (mccluster.energy () > cfgHighClusterE)
1615- continue ;
1629+ continue ;
16161630 if (fabs (mccluster.eta ()) > cfgtrkMaxEta)
1617- continue ;
1618- int ClusterHasDirectPhoton =0 ;
1631+ continue ;
1632+ int ClusterHasDirectPhoton = 0 ;
16191633 int ClusterHasFragPhoton = 0 ;
1620- auto ClusterParticles = mccluster.mcParticles_as <aod::JMcParticles>();
1621- for (auto & clusterparticle : ClusterParticles) {
1622- if (clusterparticle.pdgCode () != 22 && std::fabs (clusterparticle.pdgCode ()) != 13 ) continue ;
1623- if (clusterparticle.getGenStatusCode () < 0 ) continue ;
1624- if (std::fabs (clusterparticle.getGenStatusCode ()) >= 81 ) continue ;
1625-
1626- int chaseindex = -1 ;
1627- for (auto & mom : clusterparticle.mothers_as <aod::JMcParticles>()) {
1628- chaseindex = mom.globalIndex ();
1629- break ;
1630- }
1631- if (chaseindex < 0 ) continue ;
1632-
1633- std::set<int > visited;
1634- bool chase = true ;
1635- bool hadronInChain = false ;
1636- bool cleanPhotonChain = true ;
1637- bool adrianprompt = false ;
1638- bool adrianfrag = false ;
1639-
1640- while (chase) {
1641- if (visited.count (chaseindex)) { chase = false ; break ; }
1642- visited.insert (chaseindex);
1643-
1644- for (auto & particle : mcparticles) {
1645- if (particle.globalIndex () != chaseindex) continue ;
1646-
1647- if (particle.getGenStatusCode () == -23 ) {
1648- if (particle.pdgCode () == 22 && cleanPhotonChain) {
1649- adrianprompt = true ;
1650- } else if (particle.pdgCode () != 22 && !hadronInChain) {
1651- adrianfrag = true ;
1652- }
1634+ auto ClusterParticles = mccluster.mcParticles_as <aod::JMcParticles>();
1635+ for (auto & clusterparticle : ClusterParticles) {
1636+ if (clusterparticle.pdgCode () != 22 && std::fabs (clusterparticle.pdgCode ()) != 13 )
1637+ continue ;
1638+ if (clusterparticle.getGenStatusCode () < 0 )
1639+ continue ;
1640+ if (std::fabs (clusterparticle.getGenStatusCode ()) >= 81 )
1641+ continue ;
1642+
1643+ int chaseindex = -1 ;
1644+ for (auto & mom : clusterparticle.mothers_as <aod::JMcParticles>()) {
1645+ chaseindex = mom.globalIndex ();
1646+ break ;
1647+ }
1648+ if (chaseindex < 0 )
1649+ continue ;
1650+
1651+ std::set<int > visited;
1652+ bool chase = true ;
1653+ bool hadronInChain = false ;
1654+ bool cleanPhotonChain = true ;
1655+ bool adrianprompt = false ;
1656+ bool adrianfrag = false ;
1657+
1658+ while (chase) {
1659+ if (visited.count (chaseindex)) {
16531660 chase = false ;
16541661 break ;
16551662 }
1663+ visited.insert (chaseindex);
16561664
1657- int abspdg = std::abs (particle. pdgCode ());
1658- if (abspdg > 100 ) hadronInChain = true ;
1659- if (abspdg != 22 ) cleanPhotonChain = false ;
1665+ for ( auto & particle : mcparticles) {
1666+ if (particle. globalIndex () != chaseindex)
1667+ continue ;
16601668
1661- int nextindex = -1 ;
1662- for (auto & mom : particle.mothers_as <aod::JMcParticles>()) {
1663- nextindex = mom.globalIndex ();
1669+ if (particle.getGenStatusCode () == -23 ) {
1670+ if (particle.pdgCode () == 22 && cleanPhotonChain) {
1671+ adrianprompt = true ;
1672+ } else if (particle.pdgCode () != 22 && !hadronInChain) {
1673+ adrianfrag = true ;
1674+ }
1675+ chase = false ;
1676+ break ;
1677+ }
1678+
1679+ int abspdg = std::abs (particle.pdgCode ());
1680+ if (abspdg > 100 )
1681+ hadronInChain = true ;
1682+ if (abspdg != 22 )
1683+ cleanPhotonChain = false ;
1684+
1685+ int nextindex = -1 ;
1686+ for (auto & mom : particle.mothers_as <aod::JMcParticles>()) {
1687+ nextindex = mom.globalIndex ();
1688+ break ;
1689+ }
1690+ if (nextindex < 0 ) {
1691+ chase = false ;
1692+ } else {
1693+ chaseindex = nextindex;
1694+ }
16641695 break ;
16651696 }
1666- if (nextindex < 0 ) { chase = false ; }
1667- else { chaseindex = nextindex; }
1668- break ;
1697+ } // chase
1698+
1699+ if (adrianprompt) {
1700+ ClusterHasDirectPhoton++;
1701+ histos.fill (HIST (" REC_direct_phiQA" ), clusterparticle.phi ());
1702+ histos.fill (HIST (" REC_direct_etaQA" ), clusterparticle.eta ());
1703+ histos.fill (HIST (" REC_direct_ptQA" ), clusterparticle.pt ());
16691704 }
1670- } // chase
1705+ if (adrianfrag) {
1706+ ClusterHasFragPhoton++;
1707+ histos.fill (HIST (" REC_frag_phiQA" ), clusterparticle.phi ());
1708+ histos.fill (HIST (" REC_frag_etaQA" ), clusterparticle.eta ());
1709+ histos.fill (HIST (" REC_frag_ptQA" ), clusterparticle.pt ());
1710+ }
1711+ } // clusterparticle loop
16711712
1672- if (adrianprompt) {
1673- ClusterHasDirectPhoton++;
1674- histos.fill (HIST (" REC_direct_phiQA" ), clusterparticle.phi ());
1675- histos.fill (HIST (" REC_direct_etaQA" ), clusterparticle.eta ());
1676- histos.fill (HIST (" REC_direct_ptQA" ), clusterparticle.pt ());
1713+ if (ClusterHasFragPhoton > 0 ) {
1714+ histos.fill (HIST (" REC_cluster_frag_phiQA" ), mccluster.phi ());
1715+ histos.fill (HIST (" REC_cluster_frag_etaQA" ), mccluster.eta ());
1716+ histos.fill (HIST (" REC_cluster_frag_energyQA" ), mccluster.energy ());
16771717 }
1678- if (adrianfrag) {
1679- ClusterHasFragPhoton++;
1680- histos.fill (HIST (" REC_frag_phiQA" ), clusterparticle.phi ());
1681- histos.fill (HIST (" REC_frag_etaQA" ), clusterparticle.eta ());
1682- histos.fill (HIST (" REC_frag_ptQA" ), clusterparticle.pt ());
1718+ if (ClusterHasDirectPhoton > 0 ) {
1719+ histos.fill (HIST (" REC_cluster_direct_phiQA" ), mccluster.phi ());
1720+ histos.fill (HIST (" REC_cluster_direct_etaQA" ), mccluster.eta ());
1721+ histos.fill (HIST (" REC_cluster_direct_energyQA" ), mccluster.energy ());
1722+ }
1723+ if (ClusterHasDirectPhoton > 0 && ClusterHasFragPhoton > 0 ) {
1724+ histos.fill (HIST (" REC_cluster_both_phiQA" ), mccluster.phi ());
1725+ histos.fill (HIST (" REC_cluster_both_etaQA" ), mccluster.eta ());
1726+ histos.fill (HIST (" REC_cluster_both_energyQA" ), mccluster.energy ());
16831727 }
1684- } // clusterparticle loop
1685-
1686- if (ClusterHasFragPhoton>0 ){
1687- histos.fill (HIST (" REC_cluster_frag_phiQA" ), mccluster.phi ());
1688- histos.fill (HIST (" REC_cluster_frag_etaQA" ), mccluster.eta ());
1689- histos.fill (HIST (" REC_cluster_frag_energyQA" ), mccluster.energy ());
1690- }
1691- if (ClusterHasDirectPhoton>0 ){
1692- histos.fill (HIST (" REC_cluster_direct_phiQA" ), mccluster.phi ());
1693- histos.fill (HIST (" REC_cluster_direct_etaQA" ), mccluster.eta ());
1694- histos.fill (HIST (" REC_cluster_direct_energyQA" ), mccluster.energy ());
1695- }
1696- if (ClusterHasDirectPhoton>0 && ClusterHasFragPhoton>0 ){
1697- histos.fill (HIST (" REC_cluster_both_phiQA" ), mccluster.phi ());
1698- histos.fill (HIST (" REC_cluster_both_etaQA" ), mccluster.eta ());
1699- histos.fill (HIST (" REC_cluster_both_energyQA" ), mccluster.energy ());
1700- }
17011728
1702- // now we do cluster tracks
1703- bool photontrigger = false ; // is a neutral cluster
1704- bool chargetrigger = false ; // is definitely not a neutral cluster
1705- auto tracksofcluster = mccluster.matchedTracks_as <soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs>>();
1729+ // now we do cluster tracks
1730+ bool photontrigger = false ; // is a neutral cluster
1731+ bool chargetrigger = false ; // is definitely not a neutral cluster
1732+ auto tracksofcluster = mccluster.matchedTracks_as <soa::Join<aod::JTracks, aod::JTrackExtras, aod::JTrackPIs>>();
17061733 // first, we check if veto is required
17071734 double sumptT = 0 ;
17081735 bool clusterqa = false ;
@@ -1767,7 +1794,7 @@ struct statPromptPhoton {
17671794 histos.fill (HIST (" REC_Track_v_Cluster_Phi" ), phidiff);
17681795 histos.fill (HIST (" REC_Track_v_Cluster_Eta" ), etadiff);
17691796 histos.fill (HIST (" REC_Track_v_Cluster_Phi_Eta" ), phidiff, etadiff);
1770- histos.fill (HIST (" REC_track_phiQA" ), ctrack.phi ());
1797+ histos.fill (HIST (" REC_track_phiQA" ), ctrack.phi ());
17711798 histos.fill (HIST (" REC_track_etaQA" ), ctrack.eta ());
17721799 histos.fill (HIST (" REC_track_ptQA" ), ctrack.pt ());
17731800 } // track of cluster loop
@@ -1798,8 +1825,8 @@ struct statPromptPhoton {
17981825 histos.fill (HIST (" REC_cluster_etaQA" ), mccluster.eta ());
17991826 histos.fill (HIST (" REC_cluster_energyQA" ), mccluster.energy ());
18001827 }
1801- }// clusters
1802- }// main function
1828+ } // clusters
1829+ } // main function
18031830 PROCESS_SWITCH (statPromptPhoton, processMCRec_simple, " processMC_QA_Rce" , false );
18041831
18051832}; // end of main struct
@@ -1808,4 +1835,3 @@ WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
18081835{
18091836 return WorkflowSpec{adaptAnalysisTask<statPromptPhoton>(cfgc)};
18101837};
1811-
0 commit comments