Skip to content

Commit ba45bb1

Browse files
committed
[PWGJE,EMCAL] Add new clusterizer definitions, add new table and QA plots for cluster dispersion and update M02 and M20 QA plots to 2D plots
- Add new clusterizer definitions - Add new tables for the clusters and the ambigous clusters to store the cluster dispersion - Add new QA plot for dispersion as function of energy and add energy as new dimension to the M02 and M20 plots
1 parent 939e5c7 commit ba45bb1

3 files changed

Lines changed: 101 additions & 34 deletions

File tree

PWGJE/DataModel/EMCALClusters.h

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,6 +50,9 @@ inline const EMCALClusterDefinition kV3SmallestTimeDiff(ClusterAlgorithm::kV3, 4
5050
inline const EMCALClusterDefinition kV3MostSplitSmallTimeDiff(ClusterAlgorithm::kV3, 43, 1, "kV3MostSplitSmallTimeDiff", 0.5, 0.1, -10000, 10000, 500, true, 0., false);
5151
inline const EMCALClusterDefinition kV3MostSplitSmallerTimeDiff(ClusterAlgorithm::kV3, 44, 1, "kV3MostSplitSmallerTimeDiff", 0.5, 0.1, -10000, 10000, 100, true, 0., false);
5252
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiff(ClusterAlgorithm::kV3, 45, 1, "kV3MostSplitSmallestTimeDiff", 0.5, 0.1, -10000, 10000, 50, true, 0., false);
53+
inline const EMCALClusterDefinition kV3SmallTimeDiffLowestSeed(ClusterAlgorithm::kV3, 46, 1, "kV3SmallTimeDiffLowestSeed", 0.1, 0.1, -10000, 10000, 500, true, 0.03, false);
54+
inline const EMCALClusterDefinition kV3SmallerTimeDiffLowSeed(ClusterAlgorithm::kV3, 47, 1, "kV3SmallerTimeDiffLowSeed", 0.3, 0.1, -10000, 10000, 100, true, 0.03, false);
55+
inline const EMCALClusterDefinition kV3SmallestTimeDiffLowerSeed(ClusterAlgorithm::kV3, 48, 1, "kV3SmallestTimeDiffLowerSeed", 0.2, 0.1, -10000, 10000, 50, true, 0.03, false);
5356
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowestSeed(ClusterAlgorithm::kV3, 50, 1, "kV3MostSplitSmallestTimeDiffLowestSeed", 0.1, 0.1, -10000, 10000, 50, true, 0., false);
5457
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowSeed(ClusterAlgorithm::kV3, 51, 1, "kV3MostSplitSmallestTimeDiffLowSeed", 0.3, 0.1, -10000, 10000, 50, true, 0., false);
5558
inline const EMCALClusterDefinition kV3MostSplitSmallestTimeDiffLowerSeed(ClusterAlgorithm::kV3, 52, 1, "kV3MostSplitSmallestTimeDiffLowerSeed", 0.2, 0.1, -10000, 10000, 50, true, 0., false);
@@ -94,6 +97,12 @@ inline const EMCALClusterDefinition& getClusterDefinitionFromID(int storageID)
9497
return kV3MostSplitSmallerTimeDiff;
9598
case 45:
9699
return kV3MostSplitSmallestTimeDiff;
100+
case 46:
101+
return kV3SmallTimeDiffLowestSeed;
102+
case 47:
103+
return kV3SmallerTimeDiffLowSeed;
104+
case 48:
105+
return kV3SmallestTimeDiffLowerSeed;
97106
case 50:
98107
return kV3MostSplitSmallestTimeDiffLowestSeed;
99108
case 51:
@@ -128,6 +137,9 @@ inline const EMCALClusterDefinition& getClusterDefinitionFromString(const std::s
128137
{"kV3MostSplitSmallTimeDiff", 43},
129138
{"kV3MostSplitSmallerTimeDiff", 44},
130139
{"kV3MostSplitSmallestTimeDiff", 45},
140+
{"kV3SmallTimeDiffLowestSeed", 46},
141+
{"kV3SmallerTimeDiffLowSeed", 47},
142+
{"kV3SmallestTimeDiffLowerSeed", 48},
131143
{"kV3MostSplitSmallestTimeDiffLowestSeed", 50},
132144
{"kV3MostSplitSmallestTimeDiffLowSeed", 51},
133145
{"kV3MostSplitSmallestTimeDiffLowerSeed", 52},
@@ -156,6 +168,7 @@ DECLARE_SOA_COLUMN(IsExotic, isExotic, bool); //! flag
156168
DECLARE_SOA_COLUMN(DistanceToBadChannel, distanceToBadChannel, float); //! distance to bad channel
157169
DECLARE_SOA_COLUMN(NLM, nlm, int); //! number of local maxima
158170
DECLARE_SOA_COLUMN(Definition, definition, int); //! cluster definition, see EMCALClusterDefinition.h
171+
DECLARE_SOA_COLUMN(Dispersion, dispersion, float); //! cluster energy dispersion (cells)
159172

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

190+
DECLARE_SOA_TABLE(Dispersions, "AOD", "DISPERSIONS", //!
191+
o2::soa::Index<>, emcalcluster::Dispersion);
192+
193+
DECLARE_SOA_TABLE(AmbigousDispersions, "AOD", "AMBDISPERSIONS", //!
194+
o2::soa::Index<>, emcalcluster::Dispersion);
195+
177196
namespace emcalclustermc
178197
{
179198
DECLARE_SOA_ARRAY_INDEX_COLUMN(McParticle, mcParticle); //! Array of MC particles that deposited energy in this calo cell

PWGJE/TableProducer/emcalCorrectionTask.cxx

Lines changed: 49 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -45,6 +45,7 @@
4545
#include <Framework/AnalysisDataModel.h>
4646
#include <Framework/AnalysisHelpers.h>
4747
#include <Framework/AnalysisTask.h>
48+
#include <Framework/Concepts.h>
4849
#include <Framework/Configurable.h>
4950
#include <Framework/Expressions.h>
5051
#include <Framework/HistogramRegistry.h>
@@ -96,8 +97,10 @@ enum CellScaleMode {
9697
struct EmcalCorrectionTask {
9798
Produces<o2::aod::EMCALClusters> clusters;
9899
Produces<o2::aod::EMCALMCClusters> mcclusters;
100+
Produces<o2::aod::Dispersions> dispersions;
99101
Produces<o2::aod::EMCALAmbiguousClusters> clustersAmbiguous;
100102
Produces<o2::aod::EMCALAmbiguousMCClusters> mcclustersAmbiguous;
103+
Produces<o2::aod::AmbigousDispersions> ambigousDispersions;
101104
Produces<o2::aod::EMCALClusterCells> clustercells; // cells belonging to given cluster
102105
Produces<o2::aod::EMCALAmbiguousClusterCells> clustercellsambiguous;
103106
Produces<o2::aod::EMCALMatchedTracks> matchedTracks;
@@ -387,6 +390,10 @@ struct EmcalCorrectionTask {
387390
}
388391
}
389392

393+
auto hClusters = mHistManager.add<TH1>("hClusters", "hClusters", O2HistType::kTH1D, {{2, -0.5, 1.5}});
394+
hClusters->GetXaxis()->SetBinLabel(1, "all clusters");
395+
hClusters->GetXaxis()->SetBinLabel(2, "cross boundary cluster");
396+
390397
// For some runs, LG cells require an extra time shift of 2 * 8.8ns due to problems in the time calibration
391398
// Affected run ranges (inclusive) are initialised here (min,max)
392399
mExtraTimeShiftRunRanges.emplace_back(535365, 535645); // LHC23g-LHC23h
@@ -422,7 +429,9 @@ struct EmcalCorrectionTask {
422429
{
423430
LOG(debug) << "Starting process full.";
424431
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
432+
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
425433
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
434+
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
426435
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
427436
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());
428437

@@ -529,7 +538,7 @@ struct EmcalCorrectionTask {
529538

530539
// Store the clusters in the table where a matching collision could
531540
// be identified.
532-
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex);
541+
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, cells, &indexMapPair, &trackGlobalIndex, nullptr, nullptr);
533542
} else {
534543
mHistManager.fill(HIST("hBCMatchErrors"), 2);
535544
}
@@ -582,7 +591,9 @@ struct EmcalCorrectionTask {
582591
LOG(debug) << "Starting process full.";
583592

584593
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
594+
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
585595
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
596+
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
586597
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
587598
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());
588599

@@ -694,7 +705,7 @@ struct EmcalCorrectionTask {
694705

695706
// Store the clusters in the table where a matching collision could
696707
// be identified.
697-
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex, &indexMapPairSecondary, &secondaryGlobalIndex);
708+
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, cells, &indexMapPair, &trackGlobalIndex, &indexMapPairSecondary, &secondaryGlobalIndex);
698709
} else {
699710
mHistManager.fill(HIST("hBCMatchErrors"), 2);
700711
}
@@ -748,8 +759,10 @@ struct EmcalCorrectionTask {
748759

749760
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
750761
mcclusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
762+
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
751763
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
752764
mcclustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
765+
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
753766
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
754767
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());
755768

@@ -891,7 +904,7 @@ struct EmcalCorrectionTask {
891904

892905
// Store the clusters in the table where a matching collision could
893906
// be identified.
894-
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex);
907+
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, cells, &indexMapPair, &trackGlobalIndex);
895908
} else {
896909
mHistManager.fill(HIST("hBCMatchErrors"), 2);
897910
}
@@ -943,8 +956,10 @@ struct EmcalCorrectionTask {
943956

944957
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
945958
mcclusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
959+
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
946960
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
947961
mcclustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
962+
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
948963
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
949964
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());
950965

@@ -1088,7 +1103,7 @@ struct EmcalCorrectionTask {
10881103

10891104
// Store the clusters in the table where a matching collision could
10901105
// be identified.
1091-
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, &indexMapPair, &trackGlobalIndex, &indexMapPairSecondary, &secondaryGlobalIndex);
1106+
fillClusterTable<CollEventSels::filtered_iterator>(col, vertexPos, iClusterizer, cellIndicesBC, cells, &indexMapPair, &trackGlobalIndex, &indexMapPairSecondary, &secondaryGlobalIndex);
10921107
} else {
10931108
mHistManager.fill(HIST("hBCMatchErrors"), 2);
10941109
}
@@ -1139,7 +1154,9 @@ struct EmcalCorrectionTask {
11391154
LOG(debug) << "Starting process standalone.";
11401155

11411156
clusters.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
1157+
dispersions.reserve(MaxClusterPerDFPerClusterizer * mClusterizers.size());
11421158
clustersAmbiguous.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
1159+
ambigousDispersions.reserve(MaxAmbClusterPerDFPerClusterizer * mClusterizers.size());
11431160
clustercells.reserve(MaxCellsPerClusterPerDFPerClusterizer * mClusterizers.size());
11441161
clustercellsambiguous.reserve(MaxCellsPerAmbClusterPerDFPerClusterizer * mClusterizers.size());
11451162

@@ -1234,7 +1251,7 @@ struct EmcalCorrectionTask {
12341251

12351252
// Store the clusters in the table where a matching collision could
12361253
// be identified.
1237-
fillClusterTable<aod::Collision>(col, vertexPos, iClusterizer, cellIndicesBC);
1254+
fillClusterTable<aod::Collision>(col, vertexPos, iClusterizer, cellIndicesBC, cells);
12381255
}
12391256
} else { // ambiguous
12401257
// LOG(warning) << "No vertex found for event. Assuming (0,0,0).";
@@ -1303,12 +1320,11 @@ struct EmcalCorrectionTask {
13031320
LOG(debug) << "Converted to analysis clusters.";
13041321
}
13051322

1306-
template <typename Collision>
1307-
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)
1323+
template <o2::soa::is_iterator Collision, o2::soa::is_table Cells>
1324+
void fillClusterTable(Collision const& col, math_utils::Point3D<float> const& vertexPos, size_t iClusterizer, const gsl::span<int64_t> cellIndicesBC, Cells const& cells, MatchResult* indexMapPair = nullptr, const std::vector<int64_t>* trackGlobalIndex = nullptr, MatchResult* indexMapPairSecondaries = nullptr, const std::vector<int64_t>* secondariesGlobalIndex = nullptr)
13081325
{
13091326
// get the clusterType once
13101327
const auto clusterType = static_cast<int>(mClusterDefinitions[iClusterizer]);
1311-
13121328
int cellindex = -1;
13131329
unsigned int iCluster = 0;
13141330
float energy = 0.f;
@@ -1341,16 +1357,40 @@ struct EmcalCorrectionTask {
13411357
cluster.getClusterTime(), cluster.getIsExotic(),
13421358
cluster.getDistanceToBadChannel(), cluster.getNExMax(),
13431359
clusterType);
1360+
dispersions(cluster.getDispersion());
13441361
++nCluster;
1362+
mHistManager.fill(HIST("hClusters"), 0);
13451363
if (!mClusterLabels.empty()) {
13461364
mcclusters(mClusterLabels[iCluster].getLabels(), mClusterLabels[iCluster].getEnergyFractions());
13471365
}
13481366
// loop over cells in cluster and save to table
1367+
bool hasLargeDispersion = cluster.getDispersion() > 8;
1368+
if (hasLargeDispersion) {
1369+
mHistManager.fill(HIST("hClusters"), 1);
1370+
LOG(info) << "Found cluster with large dispersion = " << cluster.getDispersion() << "\t M02 = " << cluster.getM02() << "\t NCells = " << cluster.getNCells();
1371+
}
13491372
for (int ncell = 0; ncell < cluster.getNCells(); ncell++) {
13501373
cellindex = cluster.getCellIndex(ncell);
13511374
LOG(debug) << "trying to find cell index " << cellindex << " in map";
13521375
if (cellIndicesBC[cellindex] >= 0) {
13531376
clustercells(clusters.lastIndex(), cellIndicesBC[cellindex]);
1377+
auto cellGlobalIndex = cellIndicesBC[cellindex];
1378+
1379+
if (hasLargeDispersion) {
1380+
auto theCell = cells.rawIteratorAt(cellGlobalIndex);
1381+
auto towerId = theCell.cellNumber();
1382+
1383+
auto [nSupMod, nModule, nIphi, nIeta] = geometry->GetCellIndex(towerId);
1384+
auto [iphiLocal, ietaLocal] = geometry->GetCellPhiEtaIndexInSModule(nSupMod, nModule, nIphi, nIeta);
1385+
auto [rowGlobal, colGlobal] = geometry->GlobalRowColFromIndex(towerId);
1386+
1387+
LOG(info) << " Cell globalIndex = " << cellGlobalIndex
1388+
<< " towerId = " << towerId
1389+
<< " SM = " << nSupMod
1390+
<< " local(eta, phi) = (" << ietaLocal << ", " << iphiLocal << ")"
1391+
<< " global(eta, phi) = (" << colGlobal << ", " << rowGlobal << ")"
1392+
<< " E = " << theCell.amplitude();
1393+
}
13541394
++nCells;
13551395
}
13561396
} // end of cells of cluser loop
@@ -1426,6 +1466,7 @@ struct EmcalCorrectionTask {
14261466
cluster.getM20(), cluster.getNCells(), cluster.getClusterTime(),
14271467
cluster.getIsExotic(), cluster.getDistanceToBadChannel(),
14281468
cluster.getNExMax(), static_cast<int>(mClusterDefinitions.at(iClusterizer)));
1469+
ambigousDispersions(cluster.getDispersion());
14291470
++nClusterAmb;
14301471
if (!mClusterLabels.empty()) {
14311472
mcclustersAmbiguous(mClusterLabels[iCluster].getLabels(), mClusterLabels[iCluster].getEnergyFractions());

0 commit comments

Comments
 (0)