Skip to content

Commit a762576

Browse files
blacwblacw
blacw
authored and
blacw
committed
pi-hypertriton ME process
1 parent 6aa00f1 commit a762576

File tree

1 file changed

+227
-40
lines changed

1 file changed

+227
-40
lines changed

‎PWGCF/Femto/FemtoNuclei/TableProducer/PiNucleiFemto.cxx‎

Lines changed: 227 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -204,12 +204,13 @@ struct PiNucleiFemto {
204204
PresliceUnsorted<o2::aod::DataHypCandsWColl> hypPerCol = o2::aod::hyperrec::collisionId;
205205

206206
// binning for EM background
207-
ConfigurableAxis axisVertex{"axisVertex", {30, -10, 10}, "Binning for multiplicity"};
207+
ConfigurableAxis axisVertex{"axisVertex", {30, -10, 10}, "Binning for vtxz"};
208208
ConfigurableAxis axisCentrality{"axisCentrality", {40, 0, 100}, "Binning for centrality"};
209209
using BinningType = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0C>;
210210
BinningType binningPolicy{{axisVertex, axisCentrality}, true};
211211
SliceCache cache;
212212
SameKindPair<CollisionsFull, TrackCandidates, BinningType> mPair{binningPolicy, settingNoMixedEvents, -1, &cache};
213+
//Pair<CollisionsFull, TrackCandidates, o2::aod::DataHypCandsWColl, BinningType> hyperPair{binningPolicy, settingNoMixedEvents, -1, &cache};
213214

214215
std::array<float, 6> mBBparamsDe;
215216

@@ -229,6 +230,7 @@ struct PiNucleiFemto {
229230
"QA",
230231
{{"hVtxZ", "Vertex distribution in Z;Z (cm)", {HistType::kTH1F, {{400, -20.0, 20.0}}}},
231232
{"hNcontributor", "Number of primary vertex contributor", {HistType::kTH1F, {{2000, 0.0f, 2000.0f}}}},
233+
{"hCentrality", "Centrality", {HistType::kTH1F, {{200, 0.0f, 100.0f}}}},
232234
{"hTrackSel", "Accepted tracks", {HistType::kTH1F, {{Selections::kAll, -0.5, static_cast<double>(Selections::kAll) - 0.5}}}},
233235
{"hEvents", "; Events;", {HistType::kTH1F, {{3, -0.5, 2.5}}}},
234236
{"hEmptyPool", "svPoolCreator did not find track pairs false/true", {HistType::kTH1F, {{2, -0.5, 1.5}}}},
@@ -240,6 +242,9 @@ struct PiNucleiFemto {
240242
{"hNuPitInvMass", "; M(Nu + p) (GeV/#it{c}^{2})", {HistType::kTH1F, {{300, 3.74f, 4.34f}}}},
241243
{"hNuPt", "#it{p}_{T} distribution; #it{p}_{T} (GeV/#it{c})", {HistType::kTH1F, {{240, -6.0f, 6.0f}}}},
242244
{"hPiPt", "Pt distribution; #it{p}_{T} (GeV/#it{c})", {HistType::kTH1F, {{120, -3.0f, 3.0f}}}},
245+
{"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}}}},
246+
{"hHe3P", "Pin distribution; p (GeV/#it{c})", {HistType::kTH1F, {{120, -3.0f, 3.0f}}}},
247+
{"hHe3P_preselected", "Pin distribution_preselected; p (GeV/#it{c})", {HistType::kTH1F, {{120, -3.0f, 3.0f}}}},
243248
{"hNuEta", "eta distribution; #eta(Nu)", {HistType::kTH1F, {{200, -1.0f, 1.0f}}}},
244249
{"hPiEta", "eta distribution; #eta(#pi)", {HistType::kTH1F, {{200, -1.0f, 1.0f}}}},
245250
{"hNuPhi", "phi distribution; phi(Nu)", {HistType::kTH1F, {{600, -4.0f, 4.0f}}}},
@@ -268,6 +273,49 @@ struct PiNucleiFemto {
268273
false,
269274
true};
270275

276+
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) iy = 0;
309+
if (iy >= numOfCentBins) iy = numOfCentBins - 1;
310+
311+
int ix = static_cast<int>(std::floor((vz - Vz_low) / Vz_step));
312+
if (ix < 0) ix = 0;
313+
if (ix >= numOfVertexZBins) ix = numOfVertexZBins - 1;
314+
315+
int bin = ix + numOfVertexZBins * iy;
316+
return bin;
317+
}
318+
271319
void init(o2::framework::InitContext&)
272320
{
273321
mZorroSummary.setObject(mZorro.getZorroSummary());
@@ -488,6 +536,41 @@ struct PiNucleiFemto {
488536
return false;
489537
}
490538

539+
float averageClusterSizeCosl(uint32_t itsClusterSizes, float eta)
540+
{
541+
float average = 0;
542+
int nclusters = 0;
543+
const float cosl = 1. / std::cosh(eta);
544+
const int nlayerITS = 7;
545+
546+
for (int layer = 0; layer < nlayerITS; layer++) {
547+
if ((itsClusterSizes >> (layer * 4)) & 0xf) {
548+
nclusters++;
549+
average += (itsClusterSizes >> (layer * 4)) & 0xf;
550+
}
551+
}
552+
if (nclusters == 0) {
553+
return 0;
554+
}
555+
return average * cosl / nclusters;
556+
};
557+
558+
bool selectionPIDHyper(const aod::DataHypCandsWColl::iterator& V0Hyper)
559+
{
560+
mQaRegistry.fill(HIST("hHe3P_preselected"), V0Hyper.tpcMomHe());
561+
float averClusSizeHe = averageClusterSizeCosl(V0Hyper.itsClusterSizesHe(), V0Hyper.etaHe3());
562+
if (averClusSizeHe <= 4) {
563+
return false;
564+
}
565+
if (V0Hyper.tpcChi2He() <= 0.5) {
566+
return false;
567+
}
568+
mQaRegistry.fill(HIST("hHe3P"), V0Hyper.tpcMomHe());
569+
mQaRegistry.fill(HIST("hHe3TPCnsigma"), V0Hyper.nSigmaHe());
570+
571+
return true;
572+
}
573+
491574
// ==================================================================================================================
492575

493576
template <typename Ttrack, typename Tcollisions, typename Ttracks>
@@ -606,22 +689,22 @@ struct PiNucleiFemto {
606689
bool fillCandidateInfoHyper(const aod::DataHypCandsWColl::iterator& V0Hyper, const Ttrack& trackPi, PiNucandidate& piHypercand, bool isMixedEvent)
607690
{
608691
piHypercand.collisionID = V0Hyper.collisionId();
609-
// get hypertriton information
610-
// constexpr double mHe3 = o2::constants::physics::MassHelium3;
611-
// constexpr double mPi = o2::constants::physics::MassPiPlus;
612-
// --- He3
692+
//get hypertriton information
693+
//constexpr double mHe3 = o2::constants::physics::MassHelium3;
694+
//constexpr double mPi = o2::constants::physics::MassPiPlus;
695+
// --- He3
613696
float pxHe3 = V0Hyper.ptHe3() * std::cos(V0Hyper.phiHe3());
614697
float pyHe3 = V0Hyper.ptHe3() * std::sin(V0Hyper.phiHe3());
615698
float pzHe3 = V0Hyper.ptHe3() * std::sinh(V0Hyper.etaHe3());
616-
// float pHe3 = V0Hyper.ptHe3() * std::cosh(V0Hyper.etaHe3());
617-
// float enHe3 = std::sqrt(pHe3 * pHe3 + mHe3 * mHe3);
618-
// --- pi
699+
//float pHe3 = V0Hyper.ptHe3() * std::cosh(V0Hyper.etaHe3());
700+
//float enHe3 = std::sqrt(pHe3 * pHe3 + mHe3 * mHe3);
701+
// --- pi
619702
float pxPi = V0Hyper.ptPi() * std::cos(V0Hyper.phiPi());
620703
float pyPi = V0Hyper.ptPi() * std::sin(V0Hyper.phiPi());
621704
float pzPi = V0Hyper.ptPi() * std::sinh(V0Hyper.etaPi());
622-
// float pPi = V0Hyper.ptPi() * std::cosh(V0Hyper.etaPi());
623-
// float enPi = std::sqrt(pPi * pPi + mPi * mPi);
624-
// --- hypertriton
705+
//float pPi = V0Hyper.ptPi() * std::cosh(V0Hyper.etaPi());
706+
//float enPi = std::sqrt(pPi * pPi + mPi * mPi);
707+
// --- hypertriton
625708
float px = pxHe3 + pxPi;
626709
float py = pyHe3 + pyPi;
627710
float pz = pzHe3 + pzPi;
@@ -639,12 +722,12 @@ struct PiNucleiFemto {
639722
}
640723

641724
piHypercand.signPi = trackPi.sign();
642-
if (V0Hyper.isMatter()) {
725+
if(V0Hyper.isMatter()){
643726
piHypercand.signNu = 1;
644-
} else {
727+
}else{
645728
piHypercand.signNu = -1;
646729
}
647-
730+
648731
piHypercand.dcaxyPi = trackPi.dcaXY();
649732
piHypercand.dcazPi = trackPi.dcaZ();
650733
piHypercand.tpcSignalPi = trackPi.tpcSignal();
@@ -720,29 +803,31 @@ struct PiNucleiFemto {
720803
template <typename Ttrack, typename Thypers>
721804
void pairTracksSameEventHyper(const Ttrack& piTracks, const Thypers& V0Hypers)
722805
{
723-
for (const auto& piTrack : piTracks) {
724-
725-
mQaRegistry.fill(HIST("hTrackSel"), Selections::kNoCuts);
726-
727-
if (!selectTrack(piTrack)) {
728-
continue;
729-
}
730-
mQaRegistry.fill(HIST("hTrackSel"), Selections::kTrackCuts);
731-
732-
if (!selectionPIDPion(piTrack)) {
733-
continue;
734-
}
735-
mQaRegistry.fill(HIST("hTrackSel"), Selections::kPID);
736-
737-
for (const auto& V0Hyper : V0Hypers) {
738-
739-
SVCand pair;
740-
pair.tr0Idx = piTrack.globalIndex();
741-
pair.tr1Idx = V0Hyper.globalIndex();
742-
const int collIdx = V0Hyper.collisionId();
743-
CollBracket collBracket{collIdx, collIdx};
744-
pair.collBracket = collBracket;
745-
mTrackHypPairs.push_back(pair);
806+
for (const auto& V0Hyper : V0Hypers) {
807+
if (!selectionPIDHyper(V0Hyper)) {
808+
continue;
809+
}
810+
for (const auto& piTrack : piTracks) {
811+
812+
mQaRegistry.fill(HIST("hTrackSel"), Selections::kNoCuts);
813+
814+
if (!selectTrack(piTrack)) {
815+
continue;
816+
}
817+
mQaRegistry.fill(HIST("hTrackSel"), Selections::kTrackCuts);
818+
819+
if (!selectionPIDPion(piTrack)) {
820+
continue;
821+
}
822+
mQaRegistry.fill(HIST("hTrackSel"), Selections::kPID);
823+
824+
SVCand pair;
825+
pair.tr0Idx = piTrack.globalIndex();
826+
pair.tr1Idx = V0Hyper.globalIndex();
827+
const int collIdx = V0Hyper.collisionId();
828+
CollBracket collBracket{collIdx, collIdx};
829+
pair.collBracket = collBracket;
830+
mTrackHypPairs.push_back(pair);
746831
}
747832
}
748833
}
@@ -770,6 +855,29 @@ struct PiNucleiFemto {
770855
}
771856
}
772857

858+
template <typename T1, typename T2>
859+
void pairHyperEventMixing(T1& pionCands, T2& hypCands)
860+
{
861+
for (const auto& hypCand : hypCands) {
862+
if (!selectionPIDHyper(hypCand)) {
863+
continue;
864+
}
865+
for (const auto& pionCand : pionCands) {
866+
if (!selectTrack(pionCand) || !selectionPIDPion(pionCand)) {
867+
continue;
868+
}
869+
870+
SVCand pair;
871+
pair.tr0Idx = hypCand.globalIndex();
872+
pair.tr1Idx = pionCand.globalIndex();
873+
const int collIdx = hypCand.collisionId();
874+
CollBracket collBracket{collIdx, collIdx};
875+
pair.collBracket = collBracket;
876+
mTrackHypPairs.push_back(pair);
877+
}
878+
}
879+
}
880+
773881
template <typename Tcoll>
774882
void fillTable(const PiNucandidate& piNucand, const Tcoll& collision)
775883
{
@@ -971,12 +1079,21 @@ struct PiNucleiFemto {
9711079

9721080
auto v0hyper = V0Hypers.rawIteratorAt(trackPair.tr1Idx);
9731081
auto piTrack = piTracks.rawIteratorAt(trackPair.tr0Idx);
974-
// auto collBracket = trackPair.collBracket;
1082+
//auto collBracket = trackPair.collBracket;
9751083

9761084
PiNucandidate piNucand;
9771085
if (!fillCandidateInfoHyper(v0hyper, piTrack, piNucand, isMixedEvent)) {
9781086
continue;
9791087
}
1088+
mQaRegistry.fill(HIST("hNuPt"), piNucand.recoPtNu());
1089+
mQaRegistry.fill(HIST("hPiPt"), piNucand.recoPtPi());
1090+
mQaRegistry.fill(HIST("hNuEta"), piNucand.recoEtaNu());
1091+
mQaRegistry.fill(HIST("hPiEta"), piNucand.recoEtaPi());
1092+
mQaRegistry.fill(HIST("hNuPhi"), piNucand.recoPhiNu());
1093+
mQaRegistry.fill(HIST("hPiPhi"), piNucand.recoPhiPi());
1094+
mQaRegistry.fill(HIST("hNuPitInvMass"), piNucand.invMass);
1095+
mQaRegistry.fill(HIST("hNClsPiITS"), piNucand.nClsItsPi);
1096+
mQaRegistry.fill(HIST("hisBkgEM"), piNucand.isBkgEM);
9801097

9811098
auto collision = collisions.rawIteratorAt(piNucand.collisionID);
9821099

@@ -1021,6 +1138,7 @@ struct PiNucleiFemto {
10211138
{
10221139
mGoodCollisions.clear();
10231140
mGoodCollisions.resize(collisions.size(), false);
1141+
//LOG(info) << "Number of hyperCandidates read = " << V0Hypers.size();
10241142

10251143
for (const auto& collision : collisions) {
10261144

@@ -1033,7 +1151,7 @@ struct PiNucleiFemto {
10331151
mGoodCollisions[collision.globalIndex()] = true;
10341152
const uint64_t collIdx = collision.globalIndex();
10351153
auto trackTableThisCollision = pitracks.sliceBy(mPerCol, collIdx);
1036-
auto hypdTableThisCollision = V0Hypers.sliceBy(hypPerCol, collIdx);
1154+
auto hypdTableThisCollision = V0Hypers.sliceBy(hypPerCol, collIdx);
10371155
trackTableThisCollision.bindExternalIndices(&pitracks);
10381156
hypdTableThisCollision.bindExternalIndices(&V0Hypers);
10391157

@@ -1043,7 +1161,7 @@ struct PiNucleiFemto {
10431161
continue;
10441162
}
10451163

1046-
fillPairsHyper(collisions, pitracks, V0Hypers, /*isMixedEvent*/ false);
1164+
fillPairsHyper(collisions, pitracks, V0Hypers,/*isMixedEvent*/ false);
10471165
}
10481166
}
10491167
PROCESS_SWITCH(PiNucleiFemto, processSameEventHyper, "Process Same event", false);
@@ -1068,6 +1186,75 @@ struct PiNucleiFemto {
10681186
fillPairs(collisions, tracks, /*isMixedEvent*/ true);
10691187
}
10701188
PROCESS_SWITCH(PiNucleiFemto, processMixedEvent, "Process Mixed event", false);
1189+
1190+
/*void processMixedEventHyper(const CollisionsFull& collisions, o2::aod::DataHypCandsWColl const& V0Hypers, const TrackCandidates& pitracks)
1191+
{
1192+
LOG(debug) << "Processing mixed event for hypertriton";
1193+
mTrackHypPairs.clear();
1194+
1195+
for (const auto& [c1, tracks1, c2, V0Hypers2] : hyperPair) {
1196+
if (!c1.sel8() || !c2.sel8()) {
1197+
continue;
1198+
}
1199+
1200+
mQaRegistry.fill(HIST("hNcontributor"), c2.numContrib());
1201+
//mQaRegistry.fill(HIST("hCentrality"), c2.centFT0C());
1202+
mQaRegistry.fill(HIST("hVtxZ"), c2.posZ());
1203+
1204+
pairHyperEventMixing(tracks1, V0Hypers2);
1205+
}
1206+
1207+
fillPairsHyper(collisions, pitracks, V0Hypers,/*isMixedEvent*/ /*false);
1208+
}
1209+
PROCESS_SWITCH(PiNucleiFemto, processMixedEventHyper, "Process Mixed event", false);*/
1210+
1211+
void processMixedEventHyperPool(const CollisionsFull& collisions, o2::aod::DataHypCandsWColl const& V0Hypers, const TrackCandidates& pitracks)
1212+
{
1213+
LOG(info) << "Processing mixed event for hypertriton";
1214+
1215+
mTrackHypPairs.clear();
1216+
if (!isInitialized) {
1217+
initializePools();
1218+
LOG(info) << "Initialized event pool with size = " << All_Event_pool.size();
1219+
}
1220+
for (auto const& collision : collisions) {
1221+
int poolIndexPi = where_pool(collision.posZ(), collision.centFT0C());
1222+
auto& pool = All_Event_pool[poolIndexPi];
1223+
1224+
if (poolIndexPi < 0 || poolIndexPi >= All_Event_pool.size()) {
1225+
continue;
1226+
}
1227+
1228+
for (auto const& storedEvent : pool.events) {
1229+
auto c1 = collisions.iteratorAt(storedEvent.collisionId);
1230+
const auto& c2 = collision;
1231+
if (!c1.sel8() || !c2.sel8()) continue;
1232+
1233+
std::vector<TrackCandidates::iterator> tracks1;
1234+
for (auto const& t : pitracks) {
1235+
if (t.collisionId() == c1.globalIndex()) {
1236+
tracks1.push_back(t);
1237+
}
1238+
}
1239+
1240+
std::vector<o2::aod::DataHypCandsWColl::iterator> hypers2;
1241+
for (auto const& h : V0Hypers) {
1242+
if (h.collisionId() != c2.globalIndex()) continue;
1243+
int poolIndexHyp = where_pool(h.zPrimVtx(), h.centralityFT0C());
1244+
if (poolIndexHyp != poolIndexPi) continue;
1245+
hypers2.push_back(h);
1246+
}
1247+
pairHyperEventMixing(tracks1, hypers2);
1248+
}
1249+
fillPairsHyper(collisions, pitracks, V0Hypers,/*isMixedEvent*/ true);
1250+
1251+
if (static_cast<int>(pool.events.size()) >= settingNoMixedEvents) {
1252+
pool.events.pop_front();
1253+
}
1254+
pool.events.push_back({collision.globalIndex()});
1255+
}
1256+
}
1257+
PROCESS_SWITCH(PiNucleiFemto, processMixedEventHyperPool, "Process Mixed event", false);
10711258
};
10721259

10731260
WorkflowSpec defineDataProcessing(const ConfigContext& cfgc)

0 commit comments

Comments
 (0)