1414// / \since Apr/2/2026
1515// / \brief flow efficiency analysis on UPC MC
1616
17- #include " PWGUD/Core/SGSelector.h"
1817#include " PWGUD/DataModel/UDTables.h"
1918
2019#include " Common/Core/RecoDecay.h"
@@ -81,6 +80,8 @@ struct FlowMcUpc {
8180 histos.add <TH1>(" mcEventCounter" , " Monte Carlo Truth EventCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
8281 histos.add <TH1>(" RecoProcessEventCounter" , " Reconstruction EventCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
8382 histos.add <TH1>(" hImpactParameter" , " hImpactParameter" , HistType::kTH1D , {axisB});
83+ histos.add <TH1>(" RecoProcessTrackCounter" , " Reconstruction TrackCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
84+ histos.add <TH1>(" numberOfRecoCollisions" , " numberOfRecoCollisions" , kTH1F , {{100 , -0 .5f , 99 .5f }});
8485
8586 histos.add <TH1>(" hPtMCGen" , " Monte Carlo Truth; pT (GeV/c);" , {HistType::kTH1D , {axisPt}});
8687 histos.add <TH3>(" hEtaPtVtxzMCGen" , " Monte Carlo Truth; #eta; p_{T} (GeV/c); V_{z} (cm);" , {HistType::kTH3D , {axisEta, axisPt, axisVertex}});
@@ -97,6 +98,10 @@ struct FlowMcUpc {
9798 template <typename TTrack>
9899 bool trackSelected (TTrack const & track)
99100 {
101+ if (!track.hasTPC ())
102+ return false ;
103+ if (!track.isPVContributor ())
104+ return false ;
100105 // auto momentum = std::array<double, 3>{track.px(), track.py(), track.pz()};
101106 auto pt = track.pt ();
102107 if (pt < cfgPtCutMin || pt > cfgPtCutMax) {
@@ -151,51 +156,71 @@ struct FlowMcUpc {
151156 PROCESS_SWITCH (FlowMcUpc, processMCTrue, " process pure simulation information" , true );
152157
153158 using MCRecoTracks = soa::Join<aod::UDTracks, aod::UDTracksPID, aod::UDTracksExtra, aod::UDTracksFlags, aod::UDTracksDCA, aod::UDMcTrackLabels>;
154- using MCRecoCollisions = soa::Join<aod::UDCollisions, aod::UDCollisionSelExtras, aod:: UDCollisionsSels, aod::UDZdcsReduced , aod::UDMcCollsLabels>;
159+ using MCRecoCollisions = soa::Join<aod::UDCollisions, aod::UDCollisionsSels, aod::UDMcCollsLabels>;
155160
156161 // PresliceUnsorted<MCRecoTracks> trackPerMcParticle = aod::udmctracklabel::udMcParticleId;
157162 Preslice<MCRecoTracks> trackPerCollision = aod::udtrack::udCollisionId; // sorted preslice used because the pair track-collision is already sorted in processDataSG function
158163
159- void processReco (MCRecoCollisions const & collisions, MCRecoTracks const & tracks)
164+ void processReco (MCRecoCollisions const & collisions, MCRecoTracks const & tracks, aod::UDMcParticles const & mcParticles )
160165 {
166+ histos.fill (HIST (" numberOfRecoCollisions" ), mcParticles.size ()); // number of times coll was reco-ed
167+ // std::cout << "process reco" << std::endl;
161168 for (const auto & collision : collisions) {
169+ Partition<MCRecoTracks> pvContributors = aod::udtrack::isPVContributor == true ;
170+ pvContributors.bindTable (tracks);
171+ // std::cout << "collision loop" << std::endl;
162172 histos.fill (HIST (" RecoProcessEventCounter" ), 0.5 );
163173 // if (!eventSelected(collision))
164174 // return;
165175 histos.fill (HIST (" RecoProcessEventCounter" ), 1.5 );
166- if (!collision.has_udMcCollision ())
167- return ;
176+ // if (!collision.has_udMcCollision())
177+ // return;
168178 histos.fill (HIST (" RecoProcessEventCounter" ), 2.5 );
169- if (tracks.size () < 1 )
170- return ;
171179 histos.fill (HIST (" RecoProcessEventCounter" ), 3.5 );
172180
173181 float vtxz = collision.posZ ();
174182
175- auto const & tempTracks = tracks.sliceBy (trackPerCollision, static_cast <int64_t >(collision.globalIndex ()));
183+ // auto const& tempTracks = tracks.sliceBy(trackPerCollision, static_cast<int64_t>(collision.globalIndex()));
184+ // std::cout << "sliced" << std::endl;
176185
177- for (const auto & track : tempTracks) {
186+ for (const auto & track : tracks) {
187+ histos.fill (HIST (" RecoProcessTrackCounter" ), 0.5 );
188+ // std::cout << "track loop" << std::endl;
178189 // focus on bulk: e, mu, pi, k, p
179190 auto momentum = std::array<double , 3 >{track.px (), track.py (), track.pz ()};
180191 double pt = RecoDecay::pt (momentum);
181192 double eta = RecoDecay::eta (momentum);
182193 // double phi = RecoDecay::phi(momentum);
183194 if (!trackSelected (track) || (!track.has_udMcParticle ()))
184195 continue ;
196+ histos.fill (HIST (" RecoProcessTrackCounter" ), 1.5 );
197+ // std::cout << "track selected" << std::endl;
185198 auto mcParticle = track.udMcParticle ();
199+ // std::cout << "mc particle" << std::endl;
186200 int pdgCode = std::abs (mcParticle.pdgCode ());
201+ // std::cout << "pdg code" << std::endl;
187202
188203 // double pt = recoMC.Pt();
189204 // double eta = recoMC.Eta();
190- if (pdgCode != PDG_t::kElectron && pdgCode != PDG_t::kMuonMinus && pdgCode != PDG_t::kPiPlus && pdgCode != PDG_t::kKPlus && pdgCode != PDG_t::kProton )
205+ if (pdgCode != PDG_t::kElectron && pdgCode != PDG_t::kMuonMinus && pdgCode != PDG_t::kPiPlus && pdgCode != PDG_t::kKPlus && pdgCode != PDG_t::kProton ) {
206+ // std::cout << "pdg code not in list" << std::endl;
191207 continue ;
192- if (std::fabs (eta) > cfgCutEta) // main acceptance
208+ }
209+ histos.fill (HIST (" RecoProcessTrackCounter" ), 2.5 );
210+ if (std::fabs (eta) > cfgCutEta) {
211+ // std::cout << "cfgcuteta" << std::endl;
193212 continue ;
194- if (!mcParticle.isPhysicalPrimary ())
213+ } // main acceptance
214+ histos.fill (HIST (" RecoProcessTrackCounter" ), 3.5 );
215+ if (!mcParticle.isPhysicalPrimary ()) {
216+ // std::cout << "not physical primary" << std::endl;
195217 continue ;
218+ }
219+ histos.fill (HIST (" RecoProcessTrackCounter" ), 4.5 );
196220
197221 histos.fill (HIST (" hPtReco" ), pt);
198222 histos.fill (HIST (" hEtaPtVtxzMCReco" ), eta, pt, vtxz);
223+ // std::cout << "first loop end" << std::endl;
199224 }
200225 }
201226 }
0 commit comments