Skip to content

[PWGJE,EMCAL] Add new clusterizer definitions, add new dispersion tables and QA p… #17948

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

Open
mhemmer-cern wants to merge 3 commits into
base: master
Choose a base branch
from
+74 −35
Open

[PWGJE,EMCAL] Add new clusterizer definitions, add new dispersion tables and QA p… #17948

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
ba45bb1
[PWGJE,EMCAL] Add new clusterizer definitions, add new table and QA p…
mhemmer-cern Sep 12, 2026
9a1e78c
Remove local debug output
mhemmer-cern Sep 17, 2026
da7656d
Use table and iterator concepts where applicable in templates
mhemmer-cern Sep 17, 2026
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
19 changes: 19 additions & 0 deletions PWGJE/DataModel/EMCALClusters.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,9 @@ inline const EMCALClusterDefinition kV3SmallestTimeDiff(ClusterAlgorithm::kV3, 4
inline const EMCALClusterDefinition kV3MostSplitSmallTimeDiff(ClusterAlgorithm::kV3, 43, 1, "kV3MostSplitSmallTimeDiff", 0.5, 0.1, -10000, 10000, 500, true, 0., false);
inline const EMCALClusterDefinition kV3MostSplitSmallerTimeDiff(ClusterAlgorithm::kV3, 44, 1, "kV3MostSplitSmallerTimeDiff", 0.5, 0.1, -10000, 10000, 100, true, 0., false);
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiff(ClusterAlgorithm::kV3, 45, 1, "kV3MostSplitSmallestTimeDiff", 0.5, 0.1, -10000, 10000, 50, true, 0., false);
inline const EMCALClusterDefinition kV3SmallTimeDiffLowestSeed(ClusterAlgorithm::kV3, 46, 1, "kV3SmallTimeDiffLowestSeed", 0.1, 0.1, -10000, 10000, 500, true, 0.03, false);
inline const EMCALClusterDefinition kV3SmallerTimeDiffLowSeed(ClusterAlgorithm::kV3, 47, 1, "kV3SmallerTimeDiffLowSeed", 0.3, 0.1, -10000, 10000, 100, true, 0.03, false);
inline const EMCALClusterDefinition kV3SmallestTimeDiffLowerSeed(ClusterAlgorithm::kV3, 48, 1, "kV3SmallestTimeDiffLowerSeed", 0.2, 0.1, -10000, 10000, 50, true, 0.03, false);
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowestSeed(ClusterAlgorithm::kV3, 50, 1, "kV3MostSplitSmallestTimeDiffLowestSeed", 0.1, 0.1, -10000, 10000, 50, true, 0., false);
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowSeed(ClusterAlgorithm::kV3, 51, 1, "kV3MostSplitSmallestTimeDiffLowSeed", 0.3, 0.1, -10000, 10000, 50, true, 0., false);
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowerSeed(ClusterAlgorithm::kV3, 52, 1, "kV3MostSplitSmallestTimeDiffLowerSeed", 0.2, 0.1, -10000, 10000, 50, true, 0., false);
Expand Down Expand Up @@ -94,6 +97,12 @@ inline const EMCALClusterDefinition& getClusterDefinitionFromID(int storageID)
return kV3MostSplitSmallerTimeDiff;
case 45:
return kV3MostSplitSmallestTimeDiff;
case 46:
return kV3SmallTimeDiffLowestSeed;
case 47:
return kV3SmallerTimeDiffLowSeed;
case 48:
return kV3SmallestTimeDiffLowerSeed;
case 50:
return kV3MostSplitSmallestTimeDiffLowestSeed;
case 51:
Expand Down Expand Up @@ -128,6 +137,9 @@ inline const EMCALClusterDefinition& getClusterDefinitionFromString(const std::s
{"kV3MostSplitSmallTimeDiff", 43},
{"kV3MostSplitSmallerTimeDiff", 44},
{"kV3MostSplitSmallestTimeDiff", 45},
{"kV3SmallTimeDiffLowestSeed", 46},
{"kV3SmallerTimeDiffLowSeed", 47},
{"kV3SmallestTimeDiffLowerSeed", 48},
{"kV3MostSplitSmallestTimeDiffLowestSeed", 50},
{"kV3MostSplitSmallestTimeDiffLowSeed", 51},
{"kV3MostSplitSmallestTimeDiffLowerSeed", 52},
Expand Down Expand Up @@ -156,6 +168,7 @@ DECLARE_SOA_COLUMN(IsExotic, isExotic, bool); //! flag
DECLARE_SOA_COLUMN(DistanceToBadChannel, distanceToBadChannel, float); //! distance to bad channel
DECLARE_SOA_COLUMN(NLM, nlm, int); //! number of local maxima
DECLARE_SOA_COLUMN(Definition, definition, int); //! cluster definition, see EMCALClusterDefinition.h
DECLARE_SOA_COLUMN(Dispersion, dispersion, float); //! cluster energy dispersion (cells)

} // namespace emcalcluster
// table of clusters that could be matched to a collision
Expand All @@ -174,6 +187,12 @@ DECLARE_SOA_TABLE(EMCALAmbiguousClusters, "AOD", "EMCALAMBCLUS", //!
using EMCALCluster = EMCALClusters::iterator;
using EMCALAmbiguousCluster = EMCALAmbiguousClusters::iterator;

DECLARE_SOA_TABLE(Dispersions, "AOD", "DISPERSIONS", //!
o2::soa::Index<>, emcalcluster::Dispersion);

DECLARE_SOA_TABLE(AmbigousDispersions, "AOD", "AMBDISPERSIONS", //!
o2::soa::Index<>, emcalcluster::Dispersion);

namespace emcalclustermc
{
DECLARE_SOA_ARRAY_INDEX_COLUMN(McParticle, mcParticle); //! Array of MC particles that deposited energy in this calo cell
Expand Down
31 changes: 22 additions & 9 deletions PWGJE/TableProducer/emcalCorrectionTask.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,7 @@
#include <Framework/AnalysisDataModel.h>
#include <Framework/AnalysisHelpers.h>
#include <Framework/AnalysisTask.h>
#include <Framework/Concepts.h>
#include <Framework/Configurable.h>
#include <Framework/Expressions.h>
#include <Framework/HistogramRegistry.h>
Expand Down Expand Up @@ -96,8 +97,10 @@ enum CellScaleMode {
struct EmcalCorrectionTask {
Produces<o2::aod::EMCALClusters> clusters;
Produces<o2::aod::EMCALMCClusters> mcclusters;
Produces<o2::aod::Dispersions> dispersions;
Produces<o2::aod::EMCALAmbiguousClusters> clustersAmbiguous;
Produces<o2::aod::EMCALAmbiguousMCClusters> mcclustersAmbiguous;
Produces<o2::aod::AmbigousDispersions> ambigousDispersions;
Produces<o2::aod::EMCALClusterCells> clustercells; // cells belonging to given cluster
Produces<o2::aod::EMCALAmbiguousClusterCells> clustercellsambiguous;
Produces<o2::aod::EMCALMatchedTracks> matchedTracks;
Expand Down Expand Up @@ -404,7 +407,7 @@ struct EmcalCorrectionTask {
}
}

template <typename BCType>
template <o2::soa::is_iterator BCType>
void initZorroCCDB(const BCType& bc)
{
if (applySoftwareTriggerSelection) {
Expand All @@ -422,7 +425,9 @@ struct EmcalCorrectionTask {
{
LOG(debug) << "Starting process full.";
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());

Expand Down Expand Up @@ -529,7 +534,7 @@ struct EmcalCorrectionTask {

// Store the clusters in the table where a matching collision could
// be identified.
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex);
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex, nullptr, nullptr);
} else {
mHistManager.fill(HIST("hBCMatchErrors"), 2);
}
Expand Down Expand Up @@ -582,7 +587,9 @@ struct EmcalCorrectionTask {
LOG(debug) << "Starting process full.";

clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());

Expand Down Expand Up @@ -748,8 +755,10 @@ struct EmcalCorrectionTask {

clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
mcclusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
mcclustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());

Expand Down Expand Up @@ -943,8 +952,10 @@ struct EmcalCorrectionTask {

clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
mcclusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
mcclustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());

Expand Down Expand Up @@ -1139,7 +1150,9 @@ struct EmcalCorrectionTask {
LOG(debug) << "Starting process standalone.";

clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());

Expand Down Expand Up @@ -1303,12 +1316,11 @@ struct EmcalCorrectionTask {
LOG(debug) << "Converted to analysis clusters.";
}

template <typename Collision>
template <o2::soa::is_iterator Collision>
void fillClusterTable(Collision const& col, math_utils::Point3D<float> const& vertexPos, size_t iClusterizer, const gsl::span<int64_t> cellIndicesBC, MatchResult* indexMapPair = nullptr, const std::vector<int64_t>* trackGlobalIndex = nullptr, MatchResult* indexMapPairSecondaries = nullptr, const std::vector<int64_t>* secondariesGlobalIndex = nullptr)
{
// get the clusterType once
const auto clusterType = static_cast<int>(mClusterDefinitions[iClusterizer]);

int cellindex = -1;
unsigned int iCluster = 0;
float energy = 0.f;
Expand Down Expand Up @@ -1341,11 +1353,11 @@ struct EmcalCorrectionTask {
cluster.getClusterTime(), cluster.getIsExotic(),
cluster.getDistanceToBadChannel(), cluster.getNExMax(),
clusterType);
dispersions(cluster.getDispersion());
++nCluster;
if (!mClusterLabels.empty()) {
mcclusters(mClusterLabels[iCluster].getLabels(), mClusterLabels[iCluster].getEnergyFractions());
}
// loop over cells in cluster and save to table
for (int ncell = 0; ncell < cluster.getNCells(); ncell++) {
cellindex = cluster.getCellIndex(ncell);
LOG(debug) << "trying to find cell index " << cellindex << " in map";
Expand Down Expand Up @@ -1391,7 +1403,7 @@ struct EmcalCorrectionTask {
} // end of cluster loop
}

template <typename BC>
template <o2::soa::is_iterator BC>
void fillAmbigousClusterTable(BC const& bc, size_t iClusterizer, const gsl::span<int64_t> cellIndicesBC, bool hasCollision)
{
int cellindex = -1;
Expand Down Expand Up @@ -1426,6 +1438,7 @@ struct EmcalCorrectionTask {
cluster.getM20(), cluster.getNCells(), cluster.getClusterTime(),
cluster.getIsExotic(), cluster.getDistanceToBadChannel(),
cluster.getNExMax(), static_cast<int>(mClusterDefinitions.at(iClusterizer)));
ambigousDispersions(cluster.getDispersion());
++nClusterAmb;
if (!mClusterLabels.empty()) {
mcclustersAmbiguous(mClusterLabels[iCluster].getLabels(), mClusterLabels[iCluster].getEnergyFractions());
Expand All @@ -1439,7 +1452,7 @@ struct EmcalCorrectionTask {
} // end of cluster loop
}

template <typename Collision>
template <o2::soa::is_iterator Collision>
void doTrackMatching(Collision const& col, MyGlobTracks const& tracks, MatchResult& indexMapPair, std::vector<int64_t>& trackGlobalIndex)
{
auto groupedTracks = tracks.sliceBy(perCollision, col.globalIndex());
Expand All @@ -1455,7 +1468,7 @@ struct EmcalCorrectionTask {
indexMapPair = matchTracksToCluster(mClusterPhi, mClusterEta, trackPhi, trackEta, maxMatchingDistance, MaxMatchesPerCluster);
}

template <typename Collision>
template <o2::soa::is_iterator Collision>
void doSecondaryTrackMatching(Collision const& col, EMV0Legs const& v0legs, MatchResult& indexMapPair, std::vector<int64_t>& trackGlobalIndex, MyGlobTracks const& tracks)
{
auto groupedV0Legs = v0legs.sliceBy(perCollisionEMV0Legs, col.globalIndex());
Expand Down Expand Up @@ -1491,7 +1504,7 @@ struct EmcalCorrectionTask {
indexMapPair = matchTracksToCluster(mClusterPhi, mClusterEta, trackPhi, trackEta, maxMatchingDistance, MaxMatchesPerCluster);
}

template <typename Tracks>
template <o2::soa::is_table Tracks>
void fillTrackInfo(Tracks const& tracks, std::vector<float>& trackPhi, std::vector<float>& trackEta, std::vector<int64_t>& trackGlobalIndex)
{
for (const auto& track : tracks) {
Expand Down
Loading
Loading