@@ -82,6 +82,8 @@ struct FlowMcUpc {
8282 histos.add <TH1>(" mcEventCounter" , " Monte Carlo Truth EventCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
8383 histos.add <TH1>(" RecoProcessEventCounter" , " Reconstruction EventCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
8484 histos.add <TH1>(" hImpactParameter" , " hImpactParameter" , HistType::kTH1D , {axisB});
85+ histos.add <TH1>(" RecoProcessTrackCounter" , " Reconstruction TrackCounter" , HistType::kTH1F , {{5 , 0 , 5 }});
86+ histos.add <TH1>(" numberOfRecoCollisions" , " numberOfRecoCollisions" , kTH1F , {{100 , -0 .5f , 99 .5f }});
8587
8688 histos.add <TH1>(" hPtMCGen" , " Monte Carlo Truth; pT (GeV/c);" , {HistType::kTH1D , {axisPt}});
8789 histos.add <TH3>(" hEtaPtVtxzMCGen" , " Monte Carlo Truth; #eta; p_{T} (GeV/c); V_{z} (cm);" , {HistType::kTH3D , {axisEta, axisPt, axisVertex}});
@@ -98,6 +100,10 @@ struct FlowMcUpc {
98100 template <typename TTrack>
99101 bool trackSelected (TTrack const & track)
100102 {
103+ if (!track.hasTPC ())
104+ return false ;
105+ if (!track.isPVContributor ())
106+ return false ;
101107 // auto momentum = std::array<double, 3>{track.px(), track.py(), track.pz()};
102108 auto pt = track.pt ();
103109 if (pt < cfgPtCutMin || pt > cfgPtCutMax) {
@@ -152,51 +158,71 @@ struct FlowMcUpc {
152158 PROCESS_SWITCH (FlowMcUpc, processMCTrue, " process pure simulation information" , true );
153159
154160 using MCRecoTracks = soa::Join<aod::UDTracks, aod::UDTracksPID, aod::UDTracksExtra, aod::UDTracksFlags, aod::UDTracksDCA, aod::UDMcTrackLabels>;
155- using MCRecoCollisions = soa::Join<aod::UDCollisions, aod::UDCollisionSelExtras, aod:: UDCollisionsSels, aod::UDZdcsReduced , aod::UDMcCollsLabels>;
161+ using MCRecoCollisions = soa::Join<aod::UDCollisions, aod::UDCollisionsSels, aod::UDMcCollsLabels>;
156162
157163 // PresliceUnsorted<MCRecoTracks> trackPerMcParticle = aod::udmctracklabel::udMcParticleId;
158164 Preslice<MCRecoTracks> trackPerCollision = aod::udtrack::udCollisionId; // sorted preslice used because the pair track-collision is already sorted in processDataSG function
159165
160- void processReco (MCRecoCollisions const & collisions, MCRecoTracks const & tracks)
166+ void processReco (MCRecoCollisions const & collisions, MCRecoTracks const & tracks, aod::UDMcParticles const & mcParticles )
161167 {
168+ histos.fill (HIST (" numberOfRecoCollisions" ), mcParticles.size ()); // number of times coll was reco-ed
169+ // std::cout << "process reco" << std::endl;
162170 for (const auto & collision : collisions) {
171+ Partition<MCRecoTracks> pvContributors = aod::udtrack::isPVContributor == true ;
172+ pvContributors.bindTable (tracks);
173+ // std::cout << "collision loop" << std::endl;
163174 histos.fill (HIST (" RecoProcessEventCounter" ), 0.5 );
164175 // if (!eventSelected(collision))
165176 // return;
166177 histos.fill (HIST (" RecoProcessEventCounter" ), 1.5 );
167- if (!collision.has_udMcCollision ())
168- return ;
178+ // if (!collision.has_udMcCollision())
179+ // return;
169180 histos.fill (HIST (" RecoProcessEventCounter" ), 2.5 );
170- if (tracks.size () < 1 )
171- return ;
172181 histos.fill (HIST (" RecoProcessEventCounter" ), 3.5 );
173182
174183 float vtxz = collision.posZ ();
175184
176- auto const & tempTracks = tracks.sliceBy (trackPerCollision, static_cast <int64_t >(collision.globalIndex ()));
185+ // auto const& tempTracks = tracks.sliceBy(trackPerCollision, static_cast<int64_t>(collision.globalIndex()));
186+ // std::cout << "sliced" << std::endl;
177187
178- for (const auto & track : tempTracks) {
188+ for (const auto & track : tracks) {
189+ histos.fill (HIST (" RecoProcessTrackCounter" ), 0.5 );
190+ // std::cout << "track loop" << std::endl;
179191 // focus on bulk: e, mu, pi, k, p
180192 auto momentum = std::array<double , 3 >{track.px (), track.py (), track.pz ()};
181193 double pt = RecoDecay::pt (momentum);
182194 double eta = RecoDecay::eta (momentum);
183195 // double phi = RecoDecay::phi(momentum);
184196 if (!trackSelected (track) || (!track.has_udMcParticle ()))
185197 continue ;
198+ histos.fill (HIST (" RecoProcessTrackCounter" ), 1.5 );
199+ // std::cout << "track selected" << std::endl;
186200 auto mcParticle = track.udMcParticle ();
201+ // std::cout << "mc particle" << std::endl;
187202 int pdgCode = std::abs (mcParticle.pdgCode ());
203+ // std::cout << "pdg code" << std::endl;
188204
189205 // double pt = recoMC.Pt();
190206 // double eta = recoMC.Eta();
191- if (pdgCode != PDG_t::kElectron && pdgCode != PDG_t::kMuonMinus && pdgCode != PDG_t::kPiPlus && pdgCode != PDG_t::kKPlus && pdgCode != PDG_t::kProton )
207+ if (pdgCode != PDG_t::kElectron && pdgCode != PDG_t::kMuonMinus && pdgCode != PDG_t::kPiPlus && pdgCode != PDG_t::kKPlus && pdgCode != PDG_t::kProton ) {
208+ // std::cout << "pdg code not in list" << std::endl;
192209 continue ;
193- if (std::fabs (eta) > cfgCutEta) // main acceptance
210+ }
211+ histos.fill (HIST (" RecoProcessTrackCounter" ), 2.5 );
212+ if (std::fabs (eta) > cfgCutEta) {
213+ // std::cout << "cfgcuteta" << std::endl;
194214 continue ;
195- if (!mcParticle.isPhysicalPrimary ())
215+ } // main acceptance
216+ histos.fill (HIST (" RecoProcessTrackCounter" ), 3.5 );
217+ if (!mcParticle.isPhysicalPrimary ()) {
218+ // std::cout << "not physical primary" << std::endl;
196219 continue ;
220+ }
221+ histos.fill (HIST (" RecoProcessTrackCounter" ), 4.5 );
197222
198223 histos.fill (HIST (" hPtReco" ), pt);
199224 histos.fill (HIST (" hEtaPtVtxzMCReco" ), eta, pt, vtxz);
225+ // std::cout << "first loop end" << std::endl;
200226 }
201227 }
202228 }
0 commit comments