Skip to content

[PWGEM,PWGEM-36] PhotonMeson: Fix issues found in emcalQC and Pi0 Flo… #16905

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
mhemmer-cern merged 1 commit into from
Jul 2, 2026
+67 −45
Merged

[PWGEM,PWGEM-36] PhotonMeson: Fix issues found in emcalQC and Pi0 Flo… #16905

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 39 additions & 13 deletions PWGEM/PhotonMeson/Tasks/emcalQC.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -97,6 +97,11 @@ struct EmcalQC {
Configurable<bool> emcUseSecondaryTM{"emcUseSecondaryTM", false, "flag to use EMCal secondary track matching cut or not"};
} emccuts;

struct : ConfigurableGroup {
std::string prefix = "axis_group";
ConfigurableAxis thnConfigAxisPt{"thnConfigAxisPt", {100, 0., 20.}, "pT axis for photon candidates"};
} axisGroup;

void defineEMEventCut()
{
fEMEventCut = EMPhotonEventCut("fEMEventCut", "fEMEventCut");
Expand Down Expand Up @@ -146,6 +151,10 @@ struct EmcalQC {
defineEMCCut();
defineEMEventCut();

const AxisSpec thnAxisPt{axisGroup.thnConfigAxisPt, "#it{p}_{T} (GeV/#it{c})"};
const AxisSpec thAxisdEta{200, -0.06, 0.06, "#Delta#eta"};
const AxisSpec thAxisdPhi{200, -0.06, 0.06, "#Delta#varphi (rad)"};

o2::aod::pwgem::photonmeson::utils::eventhistogram::addEventHistograms(&fRegistry);
auto hEMCCollisionCounter = fRegistry.add<TH1>("Event/hEMCCollisionCounter", "Number of collisions after event cuts", HistType::kTH1D, {{7, 0.5, 7.5}}, false);
hEMCCollisionCounter->GetXaxis()->SetBinLabel(1, "all");
Expand All @@ -156,6 +165,13 @@ struct EmcalQC {
hEMCCollisionCounter->GetXaxis()->SetBinLabel(6, "+unique"); // TVX with z < 10cm and Sel8 and good z xertex and unique (only collision in the BC)
hEMCCollisionCounter->GetXaxis()->SetBinLabel(7, "+EMC readout"); // TVX with z < 10cm and Sel8 and good z xertex and unique (only collision in the BC) and kTVXinEMC
o2::aod::pwgem::photonmeson::utils::clusterhistogram::addClusterHistograms(&fRegistry, cfgDo2DQA);

if (doprocessQCTM) {
fRegistry.add("Cluster/hDeltaEtaPhiAllPrimTracks", "#Delta#eta #Delta#varphi distribution of all matched prim. tracks", HistType::kTH3D, {thAxisdEta, thAxisdPhi, thnAxisPt}, false);
fRegistry.add("Cluster/hDeltaEtaPhiClosestPrimTracks", "#Delta#eta #Delta#varphi distribution of the Closest matched prim. track", HistType::kTH3D, {thAxisdEta, thAxisdPhi, thnAxisPt}, false);
fRegistry.add("Cluster/hDeltaEtaPhiAllSecTracks", "#Delta#eta #Delta#varphi distribution of all matched sec. tracks", HistType::kTH3D, {thAxisdEta, thAxisdPhi, thnAxisPt}, false);
fRegistry.add("Cluster/hDeltaEtaPhiClosestSecTracks", "#Delta#eta #Delta#varphi distribution of the Closest matched sec. track", HistType::kTH3D, {thAxisdEta, thAxisdPhi, thnAxisPt}, false);
}
}

bool isEventGood(MyCollision const& collision)
Expand Down Expand Up @@ -189,7 +205,7 @@ struct EmcalQC {
if (!fEMEventCut.IsSelected(collision)) {
return false;
}
if (!(eventcuts.cfgOccupancyMin <= collision.trackOccupancyInTimeRange() && collision.trackOccupancyInTimeRange() < eventcuts.cfgOccupancyMax)) {
if (!(eventcuts.cfgOccupancyMin <= collision.trackOccupancyInTimeRange()) || !(collision.trackOccupancyInTimeRange() < eventcuts.cfgOccupancyMax)) {
return false;
}

Expand All @@ -216,11 +232,11 @@ struct EmcalQC {

// Define two boleans to see, whether the cluster "survives" the EMC cluster cuts to later check, whether the cuts in this task align with the ones in EMCPhotonCut.h:
bool survivesIsSelectedEMCalCuts = true; // Survives "manual" cuts listed in this task
bool survivesIsSelectedCuts = false;
bool survivesIsSelectedCuts = true;

survivesIsSelectedCuts = fEMCCut.IsSelected(cluster);
for (int icut = 0; icut < static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts); icut++) { // Loop through different cut observables
EMCPhotonCut::EMCPhotonCuts specificcut = static_cast<EMCPhotonCut::EMCPhotonCuts>(icut);
auto specificcut = static_cast<EMCPhotonCut::EMCPhotonCuts>(icut);
if (!fEMCCut.IsSelectedEMCal(specificcut, cluster)) { // Check whether cluster passes this cluster requirement, if not, fill why in the next row
fRegistry.fill(HIST("Cluster/hClusterQualityCuts"), icut + 1, cluster.e(), collision.weight());
survivesIsSelectedEMCalCuts = false;
Expand All @@ -240,12 +256,10 @@ struct EmcalQC {
}
}
fRegistry.fill(HIST("Cluster/after/hNgamma"), ngAfter, collision.weight());

return;
}

template <o2::soa::is_table TClusters, o2::soa::is_table TMatchedTracks, o2::soa::is_table TMatchedSecondaries>
void doClusterQA(MyCollision const& collision, TClusters const& clusters, TMatchedTracks const& primTracks, TMatchedSecondaries const& secTracks)
void doClusterQAWithTM(MyCollision const& collision, TClusters const& clusters, TMatchedTracks const& primTracks, TMatchedSecondaries const& secTracks)
{
fRegistry.fill(HIST("Cluster/before/hNgamma"), clusters.size(), collision.weight());
int ngAfter = 0;
Expand All @@ -263,11 +277,11 @@ struct EmcalQC {

// Define two boleans to see, whether the cluster "survives" the EMC cluster cuts to later check, whether the cuts in this task align with the ones in EMCPhotonCut.h:
bool survivesIsSelectedEMCalCuts = true; // Survives "manual" cuts listed in this task
bool survivesIsSelectedCuts = false;
bool survivesIsSelectedCuts = true;

survivesIsSelectedCuts = fEMCCut.IsSelected(cluster, primTracksPerCluster, secTracksPerCluster);
for (int icut = 0; icut < static_cast<int>(EMCPhotonCut::EMCPhotonCuts::kNCuts); icut++) { // Loop through different cut observables
EMCPhotonCut::EMCPhotonCuts specificcut = static_cast<EMCPhotonCut::EMCPhotonCuts>(icut);
auto specificcut = static_cast<EMCPhotonCut::EMCPhotonCuts>(icut);
if (specificcut == EMCPhotonCut::EMCPhotonCuts::kTM || specificcut == EMCPhotonCut::EMCPhotonCuts::kSecondaryTM) {
// will do track matching cuts extra later or never depending on chosen process function
continue;
Expand All @@ -277,6 +291,20 @@ struct EmcalQC {
survivesIsSelectedEMCalCuts = false;
}
}
if (primTracksPerCluster.size() > 0) {
const auto closestPrimTrack = primTracksPerCluster.begin();
fRegistry.fill(HIST("Cluster/hDeltaEtaPhiClosestPrimTracks"), closestPrimTrack.deltaEta(), closestPrimTrack.deltaPhi(), cluster.pt());
}
for (const auto& matchedPrimTrack : primTracksPerCluster) {
fRegistry.fill(HIST("Cluster/hDeltaEtaPhiAllPrimTracks"), matchedPrimTrack.deltaEta(), matchedPrimTrack.deltaPhi(), cluster.pt());
}
if (secTracksPerCluster.size() > 0) {
const auto closestSecTrack = secTracksPerCluster.begin();
fRegistry.fill(HIST("Cluster/hDeltaEtaPhiClosestSecTracks"), closestSecTrack.deltaEta(), closestSecTrack.deltaPhi(), cluster.pt());
}
for (const auto& matchedSecTrack : secTracksPerCluster) {
fRegistry.fill(HIST("Cluster/hDeltaEtaPhiAllSecTracks"), matchedSecTrack.deltaEta(), matchedSecTrack.deltaPhi(), cluster.pt());
}

if (!fEMCCut.IsSelectedEMCal(EMCPhotonCut::EMCPhotonCuts::kTM, cluster, primTracksPerCluster)) { // Check whether cluster passes this cluster requirement, if not, fill why in the next row
fRegistry.fill(HIST("Cluster/hClusterQualityCuts"), static_cast<double>(EMCPhotonCut::EMCPhotonCuts::kTM) + 1, cluster.e(), collision.weight());
Expand All @@ -300,8 +328,6 @@ struct EmcalQC {
}
}
fRegistry.fill(HIST("Cluster/after/hNgamma"), ngAfter, collision.weight());

return;
}

void processQC(MyCollisions const& collisions, EMCalPhotons const& clusters)
Expand All @@ -327,7 +353,7 @@ struct EmcalQC {
}

auto clustersPerColl = clusters.sliceBy(perCollisionEMC, collision.collisionId());
doClusterQA(collision, clustersPerColl, matchedPrims, matchedSeconds);
doClusterQAWithTM(collision, clustersPerColl, matchedPrims, matchedSeconds);

} // end of collision loop
} // end of process
Expand All @@ -339,7 +365,7 @@ struct EmcalQC {
PROCESS_SWITCH(EmcalQC, processDummy, "Dummy function", true);
};

WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
WorkflowSpec defineDataProcessing(ConfigContext const& context)
{
return WorkflowSpec{adaptAnalysisTask<EmcalQC>(cfgc)};
return WorkflowSpec{adaptAnalysisTask<EmcalQC>(context)};
}
60 changes: 28 additions & 32 deletions PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,19 +76,19 @@ using namespace o2::aod::pwgem::photon;
enum QvecEstimator {
FT0M = 0,
FT0A = 1,
FT0C,
TPCPos,
TPCNeg,
TPCTot,
FV0A
FT0C = 2,
TPCPos = 3,
TPCNeg = 4,
TPCTot = 5,
FV0A = 6
};

enum CentralityEstimator {
None = 0,
CFT0A = 1,
CFT0C,
CFT0M,
NCentralityEstimators
CFT0C = 2,
CFT0M = 3,
NCentralityEstimators = 4
};

enum Harmonics {
Expand Down Expand Up @@ -233,21 +233,21 @@ struct TaskPi0FlowEMC {
std::string prefix = "rotationConfig";
Configurable<bool> cfgDoRotation{"cfgDoRotation", false, "Flag to enable rotation background method."};
Configurable<int> cfgDownsampling{"cfgDownsampling", 1, "Calculate rotation background only for every <value> collision."};
Configurable<float> cfgRotAngle{"cfgRotAngle", std::move(const_cast<float&>(o2::constants::math::PIHalf)), "Angle used for the rotation method."};
Configurable<float> cfgRotAngle{"cfgRotAngle", static_cast<float>(o2::constants::math::PIHalf), "Angle used for the rotation method."};
Configurable<bool> cfgUseWeights{"cfgUseWeights", false, "Flag to enable weights for rotation background method."};
} rotationConfig;

struct : ConfigurableGroup {
std::string prefix = "correctionConfig";
Configurable<std::string> cfgSpresoPath{"cfgSpresoPath", "Users/m/mhemmer/EM/Flow/Resolution", "Path to SP resolution file"};
Configurable<int> cfgApplySPresolution{"cfgApplySPresolution", 0, "Apply resolution correction"};
Configurable<bool> cfgApplySPresolution{"cfgApplySPresolution", 0, "Apply resolution correction"};
Comment thread
mhemmer-cern marked this conversation as resolved.
Configurable<bool> doEMCalCalib{"doEMCalCalib", 0, "Produce output for EMCal calibration"};
Configurable<bool> cfgEnableNonLin{"cfgEnableNonLin", false, "flag to turn extra non linear energy calibration on/off"};
} correctionConfig;

SliceCache cache;
EventPlaneHelper epHelper;
o2::framework::Service<o2::ccdb::BasicCCDBManager> ccdb;
o2::framework::Service<o2::ccdb::BasicCCDBManager> ccdb{};
int runNow = 0;
int runBefore = -1;

Expand All @@ -269,8 +269,8 @@ struct TaskPi0FlowEMC {

HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject, false, false};

o2::emcal::Geometry* emcalGeom;
o2::emcal::BadChannelMap* mBadChannels;
o2::emcal::Geometry* emcalGeom = nullptr;
o2::emcal::BadChannelMap* mBadChannels = nullptr;
TH1D* h1SPResolution = nullptr;
// Constants for eta and phi ranges for the look up table
static constexpr double EtaMin = -0.75, etaMax = 0.75;
Expand All @@ -279,7 +279,7 @@ struct TaskPi0FlowEMC {
static constexpr double PhiMin = 1.35, phiMax = 5.75;
static constexpr int NBinsPhi = 440; // (440 bins = 0.01 step size covering most regions)

std::array<int8_t, NBinsEta * NBinsPhi> lookupTable1D;
std::array<int8_t, NBinsEta * NBinsPhi> lookupTable1D{};
float epsilon = 1.e-8;

// static constexpr
Expand Down Expand Up @@ -506,8 +506,8 @@ struct TaskPi0FlowEMC {
template <const int histType>
void fillThn(const float mass, const float pt, const float cent, const float sp)
{
static constexpr std::string_view FlowHistTypes[3] = {"hSparsePi0Flow", "hSparseBkgRotFlow", "hSparseBkgMixFlow"};
static constexpr std::string_view HistTypes[3] = {"hSparsePi0", "hSparseBkgRot", "hSparseBkgMix"};
static constexpr std::array<std::string_view, 3> FlowHistTypes = {"hSparsePi0Flow", "hSparseBkgRotFlow", "hSparseBkgMixFlow"};
static constexpr std::array<std::string_view, 3> HistTypes = {"hSparsePi0", "hSparseBkgRot", "hSparseBkgMix"};
registry.fill(HIST(FlowHistTypes[histType]), mass, pt, cent, sp);
registry.fill(HIST(HistTypes[histType]), mass, pt, cent);
}
Expand Down Expand Up @@ -668,9 +668,8 @@ struct TaskPi0FlowEMC {
int iRowLast = 24;
if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::EMCAL_HALF) {
iRowLast /= 2; // 2/3 sm case
} else if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::EMCAL_THIRD) {
iRowLast /= 3; // 1/3 sm case
} else if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::DCAL_EXT) {
} else if ((emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::EMCAL_THIRD) ||
(emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::DCAL_EXT)) {
iRowLast /= 3; // 1/3 sm case
}

Expand Down Expand Up @@ -836,7 +835,6 @@ struct TaskPi0FlowEMC {
} // end of loop over third photon
}
}
return;
}

/// \brief Calculate background using rotation background method
Expand Down Expand Up @@ -952,7 +950,6 @@ struct TaskPi0FlowEMC {
}
} // end of loop over PCM photons
} // if(iCellIDphotonEMC > -1)
return;
}

/// Compute the scalar product
Expand Down Expand Up @@ -983,7 +980,6 @@ struct TaskPi0FlowEMC {
}

fillThn<histType>(massCand, ptCand, cent, scalprodCand);
return;
}

/// \brief check if standard event cuts + FT0 occupancy + centrality + QVec good is
Expand Down Expand Up @@ -1105,7 +1101,7 @@ struct TaskPi0FlowEMC {
if (!(fEMCCut.IsSelected(photon))) {
continue;
}
if (cfgDistanceToEdge.value && (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelSameEvent.value)) {
if (cfgDistanceToEdge.value > 0 && (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelSameEvent.value)) {
continue;
}
registry.fill(HIST("clusterQA/hEClusterAfter"), photon.corrE()); // accepted after cuts
Expand Down Expand Up @@ -1145,7 +1141,7 @@ struct TaskPi0FlowEMC {
// general event selection
continue;
}
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange()) || !(c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange()) || !(c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
// occupancy selection
continue;
}
Expand All @@ -1168,7 +1164,7 @@ struct TaskPi0FlowEMC {
continue;
}
// Cut edge clusters away, similar to rotation method to ensure same acceptance is used
if (cfgDistanceToEdge.value) {
if (cfgDistanceToEdge.value > 0) {
if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelBackground.value) {
continue;
}
Expand Down Expand Up @@ -1259,7 +1255,7 @@ struct TaskPi0FlowEMC {
}

// Cut edge clusters away, similar to rotation method to ensure same acceptance is used
if (cfgDistanceToEdge.value) {
if (cfgDistanceToEdge.value > 0) {
if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelSameEvent.value) {
continue;
}
Expand Down Expand Up @@ -1335,7 +1331,7 @@ struct TaskPi0FlowEMC {
// general event selection
continue;
}
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange()) || !(c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange()) || !(c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
// occupancy selection
continue;
}
Expand All @@ -1358,7 +1354,7 @@ struct TaskPi0FlowEMC {
continue;
}
// Cut edge clusters away, similar to rotation method to ensure same acceptance is used
if (cfgDistanceToEdge.value) {
if (cfgDistanceToEdge.value > 0) {
if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelBackground.value) {
continue;
}
Expand Down Expand Up @@ -1431,7 +1427,7 @@ struct TaskPi0FlowEMC {
if (!(emcFlags.test(photon.globalIndex()))) {
continue;
}
if (cfgDistanceToEdge.value && (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelSameEvent.value)) {
if (cfgDistanceToEdge.value > 0 && (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelSameEvent.value)) {
continue;
}
if (emccuts.cfgEnableQA.value) {
Expand Down Expand Up @@ -1536,7 +1532,7 @@ struct TaskPi0FlowEMC {
// general event selection
continue;
}
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange()) || !(c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange()) || !(c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) {
// occupancy selection
continue;
}
Expand Down Expand Up @@ -1586,7 +1582,7 @@ struct TaskPi0FlowEMC {

}; // End struct TaskPi0FlowEMC

WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
WorkflowSpec defineDataProcessing(ConfigContext const& context)
{
return WorkflowSpec{adaptAnalysisTask<TaskPi0FlowEMC>(cfgc)};
return WorkflowSpec{adaptAnalysisTask<TaskPi0FlowEMC>(context)};
}
Loading