6363#include < array>
6464#include < cmath>
6565#include < cstdlib>
66+ #include < deque>
6667#include < iterator> // std::prev
6768#include < string>
6869#include < vector>
@@ -204,12 +205,13 @@ struct PiNucleiFemto {
204205 PresliceUnsorted<o2::aod::DataHypCandsWColl> hypPerCol = o2::aod::hyperrec::collisionId;
205206
206207 // binning for EM background
207- ConfigurableAxis axisVertex{" axisVertex" , {30 , -10 , 10 }, " Binning for multiplicity " };
208+ ConfigurableAxis axisVertex{" axisVertex" , {30 , -10 , 10 }, " Binning for vtxz " };
208209 ConfigurableAxis axisCentrality{" axisCentrality" , {40 , 0 , 100 }, " Binning for centrality" };
209210 using BinningType = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0C>;
210211 BinningType binningPolicy{{axisVertex, axisCentrality}, true };
211212 SliceCache cache;
212213 SameKindPair<CollisionsFull, TrackCandidates, BinningType> mPair {binningPolicy, settingNoMixedEvents, -1 , &cache};
214+ // Pair<CollisionsFull, TrackCandidates, o2::aod::DataHypCandsWColl, BinningType> hyperPair{binningPolicy, settingNoMixedEvents, -1, &cache};
213215
214216 std::array<float , 6 > mBBparamsDe ;
215217
@@ -229,6 +231,7 @@ struct PiNucleiFemto {
229231 " QA" ,
230232 {{" hVtxZ" , " Vertex distribution in Z;Z (cm)" , {HistType::kTH1F , {{400 , -20.0 , 20.0 }}}},
231233 {" hNcontributor" , " Number of primary vertex contributor" , {HistType::kTH1F , {{2000 , 0 .0f , 2000 .0f }}}},
234+ {" hCentrality" , " Centrality" , {HistType::kTH1F , {{200 , 0 .0f , 100 .0f }}}},
232235 {" hTrackSel" , " Accepted tracks" , {HistType::kTH1F , {{Selections::kAll , -0.5 , static_cast <double >(Selections::kAll ) - 0.5 }}}},
233236 {" hEvents" , " ; Events;" , {HistType::kTH1F , {{3 , -0.5 , 2.5 }}}},
234237 {" hEmptyPool" , " svPoolCreator did not find track pairs false/true" , {HistType::kTH1F , {{2 , -0.5 , 1.5 }}}},
@@ -240,6 +243,9 @@ struct PiNucleiFemto {
240243 {" hNuPitInvMass" , " ; M(Nu + p) (GeV/#it{c}^{2})" , {HistType::kTH1F , {{300 , 3 .74f , 4 .34f }}}},
241244 {" hNuPt" , " #it{p}_{T} distribution; #it{p}_{T} (GeV/#it{c})" , {HistType::kTH1F , {{240 , -6 .0f , 6 .0f }}}},
242245 {" hPiPt" , " Pt distribution; #it{p}_{T} (GeV/#it{c})" , {HistType::kTH1F , {{120 , -3 .0f , 3 .0f }}}},
246+ {" hHe3TPCnsigma" , " NsigmaHe3 TPC distribution; #it{p}_{T} (GeV/#it{c}); n#sigma_{TPC}(He3)" , {HistType::kTH2F , {{100 , -2 .0f , 2 .0f }, {200 , -5 .0f , 5 .0f }}}},
247+ {" hHe3P" , " Pin distribution; p (GeV/#it{c})" , {HistType::kTH1F , {{120 , -3 .0f , 3 .0f }}}},
248+ {" hHe3P_preselected" , " Pin distribution_preselected; p (GeV/#it{c})" , {HistType::kTH1F , {{120 , -3 .0f , 3 .0f }}}},
243249 {" hNuEta" , " eta distribution; #eta(Nu)" , {HistType::kTH1F , {{200 , -1 .0f , 1 .0f }}}},
244250 {" hPiEta" , " eta distribution; #eta(#pi)" , {HistType::kTH1F , {{200 , -1 .0f , 1 .0f }}}},
245251 {" hNuPhi" , " phi distribution; phi(Nu)" , {HistType::kTH1F , {{600 , -4 .0f , 4 .0f }}}},
@@ -268,6 +274,52 @@ struct PiNucleiFemto {
268274 false ,
269275 true };
270276
277+ int numOfCentBins = 40 ;
278+ int numOfVertexZBins = 30 ;
279+ float Vz_low = -10 .0f ;
280+ float Vz_high = 10 .0f ;
281+ float Vz_step = (Vz_high - Vz_low) / numOfVertexZBins;
282+
283+ struct EventRef {
284+ int64_t collisionId;
285+ };
286+
287+ struct PoolBin {
288+ std::deque<EventRef> events;
289+ };
290+
291+ std::vector<PoolBin> All_Event_pool;
292+ bool isInitialized = false ;
293+
294+ int nPoolBins () const { return numOfVertexZBins * numOfCentBins; }
295+
296+ void initializePools ()
297+ {
298+ All_Event_pool.clear ();
299+ All_Event_pool.resize (nPoolBins ());
300+ isInitialized = true ;
301+ }
302+
303+ int where_pool (float vz, float v0Centr) const
304+ {
305+ float CentBinWidth = 100.0 / numOfCentBins; // = 2.5
306+
307+ int iy = static_cast <int >(std::floor (v0Centr / CentBinWidth));
308+ if (iy < 0 )
309+ iy = 0 ;
310+ if (iy >= numOfCentBins)
311+ iy = numOfCentBins - 1 ;
312+
313+ int ix = static_cast <int >(std::floor ((vz - Vz_low) / Vz_step));
314+ if (ix < 0 )
315+ ix = 0 ;
316+ if (ix >= numOfVertexZBins)
317+ ix = numOfVertexZBins - 1 ;
318+
319+ int bin = ix + numOfVertexZBins * iy;
320+ return bin;
321+ }
322+
271323 void init (o2::framework::InitContext&)
272324 {
273325 mZorroSummary .setObject (mZorro .getZorroSummary ());
@@ -488,6 +540,41 @@ struct PiNucleiFemto {
488540 return false ;
489541 }
490542
543+ float averageClusterSizeCosl (uint32_t itsClusterSizes, float eta)
544+ {
545+ float average = 0 ;
546+ int nclusters = 0 ;
547+ const float cosl = 1 . / std::cosh (eta);
548+ const int nlayerITS = 7 ;
549+
550+ for (int layer = 0 ; layer < nlayerITS; layer++) {
551+ if ((itsClusterSizes >> (layer * 4 )) & 0xf ) {
552+ nclusters++;
553+ average += (itsClusterSizes >> (layer * 4 )) & 0xf ;
554+ }
555+ }
556+ if (nclusters == 0 ) {
557+ return 0 ;
558+ }
559+ return average * cosl / nclusters;
560+ };
561+
562+ bool selectionPIDHyper (const aod::DataHypCandsWColl::iterator& V0Hyper)
563+ {
564+ mQaRegistry .fill (HIST (" hHe3P_preselected" ), V0Hyper.tpcMomHe ());
565+ float averClusSizeHe = averageClusterSizeCosl (V0Hyper.itsClusterSizesHe (), V0Hyper.etaHe3 ());
566+ if (averClusSizeHe <= 4 ) {
567+ return false ;
568+ }
569+ if (V0Hyper.tpcChi2He () <= 0.5 ) {
570+ return false ;
571+ }
572+ mQaRegistry .fill (HIST (" hHe3P" ), V0Hyper.tpcMomHe ());
573+ mQaRegistry .fill (HIST (" hHe3TPCnsigma" ), V0Hyper.nSigmaHe ());
574+
575+ return true ;
576+ }
577+
491578 // ==================================================================================================================
492579
493580 template <typename Ttrack, typename Tcollisions, typename Ttracks>
@@ -720,21 +807,23 @@ struct PiNucleiFemto {
720807 template <typename Ttrack, typename Thypers>
721808 void pairTracksSameEventHyper (const Ttrack& piTracks, const Thypers& V0Hypers)
722809 {
723- for (const auto & piTrack : piTracks) {
724-
725- mQaRegistry .fill (HIST (" hTrackSel" ), Selections::kNoCuts );
726-
727- if (!selectTrack (piTrack)) {
810+ for (const auto & V0Hyper : V0Hypers) {
811+ if (!selectionPIDHyper (V0Hyper)) {
728812 continue ;
729813 }
730- mQaRegistry . fill ( HIST ( " hTrackSel " ), Selections:: kTrackCuts );
814+ for ( const auto & piTrack : piTracks) {
731815
732- if (!selectionPIDPion (piTrack)) {
733- continue ;
734- }
735- mQaRegistry .fill (HIST (" hTrackSel" ), Selections::kPID );
816+ mQaRegistry .fill (HIST (" hTrackSel" ), Selections::kNoCuts );
736817
737- for (const auto & V0Hyper : V0Hypers) {
818+ if (!selectTrack (piTrack)) {
819+ continue ;
820+ }
821+ mQaRegistry .fill (HIST (" hTrackSel" ), Selections::kTrackCuts );
822+
823+ if (!selectionPIDPion (piTrack)) {
824+ continue ;
825+ }
826+ mQaRegistry .fill (HIST (" hTrackSel" ), Selections::kPID );
738827
739828 SVCand pair;
740829 pair.tr0Idx = piTrack.globalIndex ();
@@ -770,6 +859,29 @@ struct PiNucleiFemto {
770859 }
771860 }
772861
862+ template <typename T1, typename T2>
863+ void pairHyperEventMixing (T1& pionCands, T2& hypCands)
864+ {
865+ for (const auto & hypCand : hypCands) {
866+ if (!selectionPIDHyper (hypCand)) {
867+ continue ;
868+ }
869+ for (const auto & pionCand : pionCands) {
870+ if (!selectTrack (pionCand) || !selectionPIDPion (pionCand)) {
871+ continue ;
872+ }
873+
874+ SVCand pair;
875+ pair.tr0Idx = hypCand.globalIndex ();
876+ pair.tr1Idx = pionCand.globalIndex ();
877+ const int collIdx = hypCand.collisionId ();
878+ CollBracket collBracket{collIdx, collIdx};
879+ pair.collBracket = collBracket;
880+ mTrackHypPairs .push_back (pair);
881+ }
882+ }
883+ }
884+
773885 template <typename Tcoll>
774886 void fillTable (const PiNucandidate& piNucand, const Tcoll& collision)
775887 {
@@ -977,6 +1089,15 @@ struct PiNucleiFemto {
9771089 if (!fillCandidateInfoHyper (v0hyper, piTrack, piNucand, isMixedEvent)) {
9781090 continue ;
9791091 }
1092+ mQaRegistry .fill (HIST (" hNuPt" ), piNucand.recoPtNu ());
1093+ mQaRegistry .fill (HIST (" hPiPt" ), piNucand.recoPtPi ());
1094+ mQaRegistry .fill (HIST (" hNuEta" ), piNucand.recoEtaNu ());
1095+ mQaRegistry .fill (HIST (" hPiEta" ), piNucand.recoEtaPi ());
1096+ mQaRegistry .fill (HIST (" hNuPhi" ), piNucand.recoPhiNu ());
1097+ mQaRegistry .fill (HIST (" hPiPhi" ), piNucand.recoPhiPi ());
1098+ mQaRegistry .fill (HIST (" hNuPitInvMass" ), piNucand.invMass );
1099+ mQaRegistry .fill (HIST (" hNClsPiITS" ), piNucand.nClsItsPi );
1100+ mQaRegistry .fill (HIST (" hisBkgEM" ), piNucand.isBkgEM );
9801101
9811102 auto collision = collisions.rawIteratorAt (piNucand.collisionID );
9821103
@@ -1021,6 +1142,7 @@ struct PiNucleiFemto {
10211142 {
10221143 mGoodCollisions .clear ();
10231144 mGoodCollisions .resize (collisions.size (), false );
1145+ // LOG(info) << "Number of hyperCandidates read = " << V0Hypers.size();
10241146
10251147 for (const auto & collision : collisions) {
10261148
@@ -1068,6 +1190,76 @@ struct PiNucleiFemto {
10681190 fillPairs (collisions, tracks, /* isMixedEvent*/ true );
10691191 }
10701192 PROCESS_SWITCH (PiNucleiFemto, processMixedEvent, " Process Mixed event" , false );
1193+
1194+ /* void processMixedEventHyper(const CollisionsFull& collisions, o2::aod::DataHypCandsWColl const& V0Hypers, const TrackCandidates& pitracks)
1195+ {
1196+ LOG(debug) << "Processing mixed event for hypertriton";
1197+ mTrackHypPairs.clear();
1198+
1199+ for (const auto& [c1, tracks1, c2, V0Hypers2] : hyperPair) {
1200+ if (!c1.sel8() || !c2.sel8()) {
1201+ continue;
1202+ }
1203+
1204+ mQaRegistry.fill(HIST("hNcontributor"), c2.numContrib());
1205+ //mQaRegistry.fill(HIST("hCentrality"), c2.centFT0C());
1206+ mQaRegistry.fill(HIST("hVtxZ"), c2.posZ());
1207+
1208+ pairHyperEventMixing(tracks1, V0Hypers2);
1209+ }
1210+ }
1211+ PROCESS_SWITCH(PiNucleiFemto, processMixedEventHyper, "Process Mixed event", false);*/
1212+
1213+ void processMixedEventHyperPool (const CollisionsFull& collisions, o2::aod::DataHypCandsWColl const & V0Hypers, const TrackCandidates& pitracks)
1214+ {
1215+ LOG (info) << " Processing mixed event for hypertriton" ;
1216+
1217+ mTrackHypPairs .clear ();
1218+ if (!isInitialized) {
1219+ initializePools ();
1220+ LOG (info) << " Initialized event pool with size = " << All_Event_pool.size ();
1221+ }
1222+ for (auto const & collision : collisions) {
1223+ int poolIndexPi = where_pool (collision.posZ (), collision.centFT0C ());
1224+ auto & pool = All_Event_pool[poolIndexPi];
1225+
1226+ if (poolIndexPi < 0 || static_cast <size_t >(poolIndexPi) >= All_Event_pool.size ()) {
1227+ continue ;
1228+ }
1229+
1230+ for (auto const & storedEvent : pool.events ) {
1231+ auto c1 = collisions.iteratorAt (storedEvent.collisionId );
1232+ const auto & c2 = collision;
1233+ if (!c1.sel8 () || !c2.sel8 ())
1234+ continue ;
1235+
1236+ std::vector<TrackCandidates::iterator> tracks1;
1237+ for (auto const & t : pitracks) {
1238+ if (t.collisionId () == c1.globalIndex ()) {
1239+ tracks1.push_back (t);
1240+ }
1241+ }
1242+
1243+ std::vector<o2::aod::DataHypCandsWColl::iterator> hypers2;
1244+ for (auto const & h : V0Hypers) {
1245+ if (h.collisionId () != c2.globalIndex ())
1246+ continue ;
1247+ int poolIndexHyp = where_pool (h.zPrimVtx (), h.centralityFT0C ());
1248+ if (poolIndexHyp != poolIndexPi)
1249+ continue ;
1250+ hypers2.push_back (h);
1251+ }
1252+ pairHyperEventMixing (tracks1, hypers2);
1253+ }
1254+ fillPairsHyper (collisions, pitracks, V0Hypers, /* isMixedEvent*/ true );
1255+
1256+ if (static_cast <int >(pool.events .size ()) >= settingNoMixedEvents) {
1257+ pool.events .pop_front ();
1258+ }
1259+ pool.events .push_back ({collision.globalIndex ()});
1260+ }
1261+ }
1262+ PROCESS_SWITCH (PiNucleiFemto, processMixedEventHyperPool, " Process Mixed event" , false );
10711263};
10721264
10731265WorkflowSpec defineDataProcessing (const ConfigContext& cfgc)
0 commit comments