Skip to content

Commit 598928a

Browse files
mhemmer-cernmhemmer-cern
committed
[PWGEM] PhotonMeson: change more cutlibs to use concepts
- Use is_table and is_iterator concepts in DalitzEE and V0Photon Cut libs. - Add PCM to flow task
1 parent 8218812 commit 598928a

File tree

3 files changed

+166
-12
lines changed

3 files changed

+166
-12
lines changed

‎PWGEM/PhotonMeson/Core/DalitzEECut.h‎

Lines changed: 8 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -20,6 +20,7 @@
2020
#include "PWGEM/Dilepton/Utils/PairUtilities.h"
2121

2222
#include <CommonConstants/PhysicsConstants.h>
23+
#include <Framework/ASoA.h>
2324

2425
#include <Math/Vector4D.h> // IWYU pragma: keep
2526
#include <Math/Vector4Dfwd.h>
@@ -72,7 +73,7 @@ class DalitzEECut : public TNamed
7273
kTPConly = 1,
7374
};
7475

75-
template <typename TTrack1, typename TTrack2>
76+
template <o2::soa::is_iterator TTrack1, o2::soa::is_iterator TTrack2>
7677
bool IsSelected(TTrack1 const& t1, TTrack2 const& t2, float bz) const
7778
{
7879
if (!IsSelectedTrack(t1) || !IsSelectedTrack(t2)) {
@@ -86,7 +87,7 @@ class DalitzEECut : public TNamed
8687
return true;
8788
}
8889

89-
template <typename TTrack1, typename TTrack2>
90+
template <o2::soa::is_iterator TTrack1, o2::soa::is_iterator TTrack2>
9091
bool IsSelectedPair(TTrack1 const& t1, TTrack2 const& t2, const float bz) const
9192
{
9293
ROOT::Math::PtEtaPhiMVector v1(t1.pt(), t1.eta(), t1.phi(), o2::constants::physics::MassElectron);
@@ -108,7 +109,7 @@ class DalitzEECut : public TNamed
108109
return true;
109110
}
110111

111-
template <bool isML = false, typename TTrack, typename TCollision = int>
112+
template <bool isML = false, o2::soa::is_iterator TTrack, typename TCollision = int>
112113
bool IsSelectedTrack(TTrack const& track, TCollision const& = 0) const
113114
{
114115
if (!track.hasITS()) {
@@ -194,7 +195,7 @@ class DalitzEECut : public TNamed
194195
return true;
195196
}
196197

197-
template <typename T>
198+
template <o2::soa::is_iterator T>
198199
bool PassPID(T const& track) const
199200
{
200201
switch (mPIDScheme) {
@@ -212,15 +213,15 @@ class DalitzEECut : public TNamed
212213
}
213214
}
214215

215-
template <typename T>
216+
template <o2::soa::is_iterator T>
216217
bool PassTPConly(T const& track) const
217218
{
218219
bool is_el_included_TPC = mMinTPCNsigmaEl < track.tpcNSigmaEl() && track.tpcNSigmaEl() < mMaxTPCNsigmaEl;
219220
bool is_pi_excluded_TPC = track.tpcNSigmaPi() < mMinTPCNsigmaPi || mMaxTPCNsigmaPi < track.tpcNSigmaPi();
220221
return is_el_included_TPC && is_pi_excluded_TPC;
221222
}
222223

223-
template <typename T>
224+
template <o2::soa::is_iterator T>
224225
bool PassTOFif(T const& track) const
225226
{
226227
bool is_el_included_TPC = mMinTPCNsigmaEl < track.tpcNSigmaEl() && track.tpcNSigmaEl() < mMaxTPCNsigmaEl;
@@ -229,7 +230,7 @@ class DalitzEECut : public TNamed
229230
return is_el_included_TPC && is_pi_excluded_TPC && is_el_included_TOF;
230231
}
231232

232-
template <typename T>
233+
template <o2::soa::is_iterator T>
233234
bool IsSelectedTrack(T const& track, const DalitzEECuts& cut) const
234235
{
235236
switch (cut) {

‎PWGEM/PhotonMeson/Core/V0PhotonCut.h‎

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -16,6 +16,7 @@
1616
#ifndef PWGEM_PHOTONMESON_CORE_V0PHOTONCUT_H_
1717
#define PWGEM_PHOTONMESON_CORE_V0PHOTONCUT_H_
1818

19+
#include "PWGEM/PhotonMeson/Core/EMBitFlags.h"
1920
#include "PWGEM/PhotonMeson/Core/EmMlResponsePCM.h"
2021
#include "PWGEM/PhotonMeson/Core/V0PhotonCandidate.h"
2122
#include "PWGEM/PhotonMeson/Utils/TrackSelection.h"
@@ -193,6 +194,21 @@ class V0PhotonCut : public TNamed
193194
kNCuts
194195
};
195196

197+
/// \brief check if given v0 photon survives all cuts
198+
/// \param flags EMBitFlags where results will be stored
199+
/// \param v0s v0 photon table to check
200+
template <o2::soa::is_table TV0, typename TLeg>
201+
void AreSelectedRunning(EMBitFlags& flags, TV0 const& v0s) const
202+
{
203+
size_t iV0 = 0;
204+
for (const auto& v0 : v0s) {
205+
if (!IsSelected<decltype(v0), TLeg>(v0)) {
206+
flags.set(iV0);
207+
}
208+
++iV0;
209+
}
210+
}
211+
196212
template <o2::soa::is_iterator TV0, typename TLeg>
197213
bool IsSelected(TV0 const& v0) const
198214
{

‎PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx‎

Lines changed: 142 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1240,6 +1240,9 @@ struct TaskPi0FlowEMC {
12401240
EMBitFlags emcFlags(clusters.size());
12411241
fEMCCut.AreSelectedRunning(emcFlags, clusters, matchedPrims, matchedSeconds);
12421242

1243+
EMBitFlags v0flags(photons.size());
1244+
fV0PhotonCut.AreSelectedRunning<decltype(photons), aod::V0Legs>(v0flags, photons);
1245+
12431246
if (cfgDoReverseScaling.value) {
12441247
energyCorrectionFactor = 1.0505f;
12451248
}
@@ -1267,7 +1270,7 @@ struct TaskPi0FlowEMC {
12671270
}
12681271
}
12691272
for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photonsEMCPerCollision, photonsPCMPerCollision))) {
1270-
if (!(emcFlags.test(g1.globalIndex())) || !(fV0PhotonCut.IsSelected<decltype(g2), aod::V0Legs>(g2))) {
1273+
if (!(emcFlags.test(g1.globalIndex())) || !(v0flags.test(g2.globalIndex()))) {
12711274
continue;
12721275
}
12731276

@@ -1304,7 +1307,7 @@ struct TaskPi0FlowEMC {
13041307
registry.fill(HIST("hMesonCuts"), 4);
13051308
continue;
13061309
}
1307-
if (mesonConfig.cfgEnableQA) {
1310+
if (mesonConfig.cfgEnableQA.value) {
13081311
registry.fill(HIST("mesonQA/hInvMassPt"), vMeson.M(), vMeson.Pt());
13091312
registry.fill(HIST("mesonQA/hTanThetaPhi"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi)));
13101313
registry.fill(HIST("mesonQA/hAlphaPt"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt());
@@ -1332,6 +1335,10 @@ struct TaskPi0FlowEMC {
13321335
float energyCorrectionFactor = 1.f;
13331336
EMBitFlags emcFlags(clusters.size());
13341337
fEMCCut.AreSelectedRunning(emcFlags, clusters, matchedPrims, matchedSeconds);
1338+
1339+
EMBitFlags v0flags(pcmPhotons.size());
1340+
fV0PhotonCut.AreSelectedRunning<decltype(pcmPhotons), aod::V0Legs>(v0flags, pcmPhotons);
1341+
13351342
for (const auto& [c1, photonEMC, c2, photonPCM] : pairPCMEMC) {
13361343
if (!(fEMEventCut.IsSelected(c1)) || !(fEMEventCut.IsSelected(c2))) {
13371344
// general event selection
@@ -1352,7 +1359,7 @@ struct TaskPi0FlowEMC {
13521359
}
13531360
registry.fill(HIST("h3DMixingCount"), c1.posZ(), getCentrality(c1), c1.ep2ft0m());
13541361
for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photonEMC, photonPCM))) {
1355-
if (!(emcFlags.test(g1.globalIndex())) || !(fV0PhotonCut.IsSelected<decltype(g2), aod::V0Legs>(g2))) {
1362+
if (!(emcFlags.test(g1.globalIndex())) || !(v0flags.test(g2.globalIndex()))) {
13561363
continue;
13571364
}
13581365
// Cut edge clusters away, similar to rotation method to ensure same acceptance is used
@@ -1383,7 +1390,7 @@ struct TaskPi0FlowEMC {
13831390
registry.fill(HIST("hMesonCutsMixed"), 4);
13841391
continue;
13851392
}
1386-
if (mesonConfig.cfgEnableQA) {
1393+
if (mesonConfig.cfgEnableQA.value) {
13871394
registry.fill(HIST("mesonQA/hInvMassPtMixed"), vMeson.M(), vMeson.Pt());
13881395
registry.fill(HIST("mesonQA/hTanThetaPhiMixed"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi)));
13891396
registry.fill(HIST("mesonQA/hAlphaPtMixed"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt());
@@ -1402,6 +1409,9 @@ struct TaskPi0FlowEMC {
14021409
// Pi0 from EMCal
14031410
void processM02(CollsWithQvecs const& collisions, EMCalPhotons const& clusters, MinMTracks const& matchedPrims, MinMSTracks const& matchedSeconds)
14041411
{
1412+
EMBitFlags emcFlags(clusters.size());
1413+
fEMCCut.AreSelectedRunning(emcFlags, clusters, matchedPrims, matchedSeconds);
1414+
14051415
for (const auto& collision : collisions) {
14061416
o2::aod::pwgem::photonmeson::utils::eventhistogram::fillEventInfo<0>(&registry, collision);
14071417
if (!(fEMEventCut.IsSelected(collision))) {
@@ -1439,7 +1449,7 @@ struct TaskPi0FlowEMC {
14391449
}
14401450
auto matchedPrimsPerCluster = matchedPrims.sliceBy(perEMCClusterMT, photon.globalIndex());
14411451
auto matchedSecondsPerCluster = matchedSeconds.sliceBy(perEMCClusterMS, photon.globalIndex());
1442-
if (!(fEMCCut.IsSelected(photon, matchedPrimsPerCluster, matchedSecondsPerCluster))) {
1452+
if (!(emcFlags.test(photon.globalIndex()))) {
14431453
continue;
14441454
}
14451455
if (cfgDistanceToEdge.value && (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelSameEvent.value)) {
@@ -1472,6 +1482,133 @@ struct TaskPi0FlowEMC {
14721482
} // processM02
14731483
PROCESS_SWITCH(TaskPi0FlowEMC, processM02, "Process single EMCal clusters as function of M02", false);
14741484

1485+
// Pi0 from EMCal
1486+
void processPCM(CollsWithQvecs const& collisions, PCMPhotons const& photons, aod::V0Legs const&)
1487+
{
1488+
EMBitFlags v0flags(photons.size());
1489+
fV0PhotonCut.AreSelectedRunning<decltype(photons), aod::V0Legs>(v0flags, photons);
1490+
for (const auto& collision : collisions) {
1491+
1492+
if (!isFullEventSelected(collision)) {
1493+
continue;
1494+
}
1495+
runNow = collision.runNumber();
1496+
if (runNow != runBefore) {
1497+
initCCDB(collision);
1498+
runBefore = runNow;
1499+
}
1500+
1501+
auto photonsPerCollision = photons.sliceBy(perCollisionPCM, collision.globalIndex());
1502+
for (const auto& [g1, g2] : combinations(CombinationsStrictlyUpperIndexPolicy(photonsPerCollision, photonsPerCollision))) {
1503+
if (!(v0flags.test(g1.globalIndex())) || !(v0flags.test(g2.globalIndex()))) {
1504+
continue;
1505+
}
1506+
if (correctionConfig.cfgEnableNonLin.value) {
1507+
energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy);
1508+
}
1509+
if (correctionConfig.cfgEnableNonLin.value) {
1510+
energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g2.e() > MinEnergy ? g2.e() : MinEnergy);
1511+
}
1512+
1513+
ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.);
1514+
ROOT::Math::PtEtaPhiMVector v2(energyCorrectionFactor * g2.pt(), g2.eta(), g2.phi(), 0.);
1515+
ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2;
1516+
1517+
float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P()));
1518+
1519+
registry.fill(HIST("hMesonCuts"), 1);
1520+
if (openingAngle <= mesonConfig.minOpenAngle) {
1521+
registry.fill(HIST("hMesonCuts"), 2);
1522+
continue;
1523+
}
1524+
if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) {
1525+
registry.fill(HIST("hMesonCuts"), 3);
1526+
continue;
1527+
}
1528+
if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) {
1529+
registry.fill(HIST("hMesonCuts"), 4);
1530+
continue;
1531+
}
1532+
if (mesonConfig.cfgEnableQA.value) {
1533+
registry.fill(HIST("mesonQA/hInvMassPt"), vMeson.M(), vMeson.Pt());
1534+
registry.fill(HIST("mesonQA/hAlphaPt"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt());
1535+
}
1536+
registry.fill(HIST("hMesonCuts"), 6);
1537+
runFlowAnalysis<0>(collision, vMeson);
1538+
}
1539+
} // end of loop over collisions
1540+
}
1541+
PROCESS_SWITCH(TaskPi0FlowEMC, processPCM, "Process PCM Pi0 candidates", false);
1542+
1543+
// PCM-EMCal mixed event
1544+
void processPCMMixed(FilteredCollsWithQvecs const& collisions, PCMPhotons const& pcmPhotons, aod::V0Legs const&)
1545+
{
1546+
1547+
using BinningType = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0C, emevent::EP2FT0C<emevent::Q2xFT0C, emevent::Q2yFT0C>>;
1548+
BinningType binningMixedEvent{{mixingConfig.cfgVtxBins, mixingConfig.cfgCentBins, mixingConfig.cfgEPBins}, true};
1549+
1550+
auto pcmPhotonTuple = std::make_tuple(pcmPhotons);
1551+
SameKindPair<FilteredCollsWithQvecs, PCMPhotons, BinningType> pair{binningMixedEvent, mixingConfig.cfgMixingDepth, -1, collisions, pcmPhotonTuple, &cache}; // indicates that 5 events should be mixed and under/overflow (-1) to be ignored
1552+
1553+
float energyCorrectionFactor = 1.f;
1554+
1555+
EMBitFlags v0flags(pcmPhotons.size());
1556+
fV0PhotonCut.AreSelectedRunning<decltype(pcmPhotons), aod::V0Legs>(v0flags, pcmPhotons);
1557+
1558+
for (const auto& [c1, photon1, c2, photon2] : pair) {
1559+
if (!(fEMEventCut.IsSelected(c1)) || !(fEMEventCut.IsSelected(c2))) {
1560+
// general event selection
1561+
continue;
1562+
}
1563+
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
1564+
// occupancy selection
1565+
continue;
1566+
}
1567+
if (getCentrality(c1) < eventcuts.cfgMinCent || getCentrality(c1) > eventcuts.cfgMaxCent || getCentrality(c2) < eventcuts.cfgMinCent || getCentrality(c2) > eventcuts.cfgMaxCent) {
1568+
// event selection
1569+
continue;
1570+
}
1571+
runNow = c1.runNumber();
1572+
if (runNow != runBefore) {
1573+
initCCDB(c1);
1574+
runBefore = runNow;
1575+
}
1576+
registry.fill(HIST("h3DMixingCount"), c1.posZ(), getCentrality(c1), c1.ep2ft0m());
1577+
for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photon1, photon2))) {
1578+
if (!(v0flags.test(g1.globalIndex())) || !(v0flags.test(g2.globalIndex()))) {
1579+
continue;
1580+
}
1581+
energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy);
1582+
ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.);
1583+
ROOT::Math::PtEtaPhiMVector v2(g2.pt(), g2.eta(), g2.phi(), 0.);
1584+
ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2;
1585+
1586+
float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P()));
1587+
1588+
registry.fill(HIST("hMesonCutsMixed"), 1);
1589+
if (openingAngle <= mesonConfig.minOpenAngle) {
1590+
registry.fill(HIST("hMesonCutsMixed"), 2);
1591+
continue;
1592+
}
1593+
if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) {
1594+
registry.fill(HIST("hMesonCutsMixed"), 3);
1595+
continue;
1596+
}
1597+
if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) {
1598+
registry.fill(HIST("hMesonCutsMixed"), 4);
1599+
continue;
1600+
}
1601+
if (mesonConfig.cfgEnableQA.value) {
1602+
registry.fill(HIST("mesonQA/hInvMassPtMixed"), vMeson.M(), vMeson.Pt());
1603+
registry.fill(HIST("mesonQA/hAlphaPtMixed"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt());
1604+
}
1605+
registry.fill(HIST("hMesonCutsMixed"), 6);
1606+
runFlowAnalysis<2>(c1, vMeson);
1607+
}
1608+
}
1609+
}
1610+
PROCESS_SWITCH(TaskPi0FlowEMC, processPCMMixed, "Process neutral meson flow using PCM-EMC mixed event", false);
1611+
14751612
}; // End struct TaskPi0FlowEMC
14761613

14771614
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)

0 commit comments

Comments
 (0)