Skip to content

Commit f729ed0

Browse files
fjonasALICEalibuild
fjonasALICE
and
alibuild
authored
[PWGJE] additional features for GammaJetTrees (#8565)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent 8fc23d1 commit f729ed0

File tree

3 files changed

+84
-40
lines changed

3 files changed

+84
-40
lines changed

‎PWGJE/DataModel/GammaJetAnalysisTree.h‎

Lines changed: 6 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -31,9 +31,10 @@ DECLARE_SOA_COLUMN(Multiplicity, multiplicity, float);
3131
DECLARE_SOA_COLUMN(Centrality, centrality, float);
3232
DECLARE_SOA_COLUMN(Rho, rho, float);
3333
DECLARE_SOA_COLUMN(EventSel, eventSel, uint8_t);
34+
DECLARE_SOA_COLUMN(Occupancy, occupancy, int);
3435
DECLARE_SOA_BITMAP_COLUMN(Alias, alias, 32);
3536
} // namespace gjevent
36-
DECLARE_SOA_TABLE(GjEvents, "AOD", "GJEVENT", o2::soa::Index<>, gjevent::Multiplicity, gjevent::Centrality, gjevent::Rho, gjevent::EventSel, gjevent::Alias)
37+
DECLARE_SOA_TABLE(GjEvents, "AOD", "GJEVENT", o2::soa::Index<>, gjevent::Multiplicity, gjevent::Centrality, gjevent::Rho, gjevent::EventSel, gjevent::Occupancy, gjevent::Alias)
3738

3839
using GjEvent = GjEvents::iterator;
3940

@@ -64,12 +65,15 @@ DECLARE_SOA_INDEX_COLUMN(GjEvent, gjevent);
6465
DECLARE_SOA_COLUMN(Pt, pt, float);
6566
DECLARE_SOA_COLUMN(Eta, eta, float);
6667
DECLARE_SOA_COLUMN(Phi, phi, float);
68+
DECLARE_SOA_COLUMN(Radius, radius, float);
6769
DECLARE_SOA_COLUMN(Energy, energy, float);
6870
DECLARE_SOA_COLUMN(Mass, mass, float);
6971
DECLARE_SOA_COLUMN(Area, area, float);
72+
DECLARE_SOA_COLUMN(LeadingTrackPt, leadingtrackpt, float);
73+
DECLARE_SOA_COLUMN(PerpConeRho, perpconerho, float);
7074
DECLARE_SOA_COLUMN(NConstituents, nConstituents, ushort);
7175
} // namespace gjchjet
72-
DECLARE_SOA_TABLE(GjChargedJets, "AOD", "GJCHJET", gjchjet::GjEventId, gjchjet::Pt, gjchjet::Eta, gjchjet::Phi, gjchjet::Energy, gjchjet::Mass, gjchjet::Area, gjchjet::NConstituents)
76+
DECLARE_SOA_TABLE(GjChargedJets, "AOD", "GJCHJET", gjchjet::GjEventId, gjchjet::Pt, gjchjet::Eta, gjchjet::Phi, gjchjet::Radius, gjchjet::Energy, gjchjet::Mass, gjchjet::Area, gjchjet::LeadingTrackPt, gjchjet::PerpConeRho, gjchjet::NConstituents)
7377
} // namespace o2::aod
7478

7579
#endif // PWGJE_DATAMODEL_GAMMAJETANALYSISTREE_H_

‎PWGJE/TableProducer/emcalCorrectionTask.cxx‎

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -184,13 +184,15 @@ struct EmcalCorrectionTask {
184184
o2Axis energyAxis{200, 0., 100., "E (GeV)"},
185185
timeAxis{300, -100, 200., "t (ns)"},
186186
etaAxis{160, -0.8, 0.8, "#eta"},
187-
phiAxis{72, 0, 2 * 3.14159, "phi"};
187+
phiAxis{72, 0, 2 * 3.14159, "phi"},
188+
nlmAxis{50, -0.5, 49.5, "NLM"};
188189
mHistManager.add("hCellE", "hCellE", o2HistType::kTH1F, {energyAxis});
189190
mHistManager.add("hCellTowerID", "hCellTowerID", o2HistType::kTH1D, {{20000, 0, 20000}});
190191
mHistManager.add("hCellEtaPhi", "hCellEtaPhi", o2HistType::kTH2F, {etaAxis, phiAxis});
191192
// NOTE: Reversed column and row because it's more natural for presentation.
192193
mHistManager.add("hCellRowCol", "hCellRowCol;Column;Row", o2HistType::kTH2D, {{97, 0, 97}, {600, 0, 600}});
193194
mHistManager.add("hClusterE", "hClusterE", o2HistType::kTH1F, {energyAxis});
195+
mHistManager.add("hClusterNLM", "hClusterNLM", o2HistType::kTH1F, {nlmAxis});
194196
mHistManager.add("hClusterEtaPhi", "hClusterEtaPhi", o2HistType::kTH2F, {etaAxis, phiAxis});
195197
mHistManager.add("hClusterTime", "hClusterTime", o2HistType::kTH1F, {timeAxis});
196198
mHistManager.add("hGlobalTrackEtaPhi", "hGlobalTrackEtaPhi", o2HistType::kTH2F, {etaAxis, phiAxis});
@@ -627,6 +629,7 @@ struct EmcalCorrectionTask {
627629
} // end of cells of cluser loop
628630
// fill histograms
629631
mHistManager.fill(HIST("hClusterE"), cluster.E());
632+
mHistManager.fill(HIST("hClusterNLM"), cluster.getNExMax());
630633
mHistManager.fill(HIST("hClusterTime"), cluster.getClusterTime());
631634
mHistManager.fill(HIST("hClusterEtaPhi"), pos.Eta(), TVector2::Phi_0_2pi(pos.Phi()));
632635
if (IndexMapPair && trackGlobalIndex) {

‎PWGJE/Tasks/gammajettreeproducer.cxx‎

Lines changed: 74 additions & 37 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,12 @@
99
// granted to it by virtue of its status as an Intergovernmental Organization
1010
// or submit itself to any jurisdiction.
1111

12+
// C++ system headers first
13+
#include <string>
14+
#include <unordered_map>
15+
#include <vector>
16+
17+
// Framework and other headers after
1218
#include "Framework/ASoA.h"
1319
#include "Framework/AnalysisDataModel.h"
1420
#include "Framework/AnalysisTask.h"
@@ -44,9 +50,10 @@
4450
/// \since 02.08.2024
4551
///
4652
using namespace o2;
53+
using namespace o2::aod;
4754
using namespace o2::framework;
4855
using namespace o2::framework::expressions;
49-
using selectedClusters = o2::soa::Filtered<o2::aod::JClusters>;
56+
using selectedClusters = o2::soa::Filtered<o2::soa::Join<o2::aod::JClusters, o2::aod::JClusterTracks>>;
5057

5158
#include "Framework/runDataProcessing.h"
5259

@@ -72,12 +79,11 @@ struct GammaJetTreeProducer {
7279
trackSelections{"trackSelections", "globalTracks", "set track selections"};
7380
Configurable<float> trackMinPt{"trackMinPt", 0.15, "minimum track pT cut"};
7481
Configurable<float> jetPtMin{"jetPtMin", 5.0, "minimum jet pT cut"};
75-
Configurable<float> jetR{"jetR", 0.4, "jet resolution parameter"};
7682
Configurable<float> isoR{"isoR", 0.4, "isolation cone radius"};
77-
83+
Configurable<float> perpConeJetR{"perpConeJetR", 0.4, "perpendicular cone radius used to calculate perp cone rho for jet"};
84+
Configurable<float> trackMatchingEoverP{"trackMatchingEoverP", 2.0, "closest track is required to have E/p < value"};
7885
// cluster cuts
7986
Configurable<int> mClusterDefinition{"clusterDefinition", 10, "cluster definition to be selected, e.g. 10=kV3Default"};
80-
// Preslice<o2::aod::JClusterTracks> perClusterMatchedTracks = o2::aod::jcluster::clusterId;
8187

8288
int mRunNumber = 0;
8389
int eventSelection = -1;
@@ -98,14 +104,21 @@ struct GammaJetTreeProducer {
98104
// create histograms
99105
LOG(info) << "Creating histograms";
100106

101-
const o2Axis ptAxis{100, 0, 100, "p_{T} (GeV/c)"};
107+
const o2Axis ptAxis{100, 0, 200, "p_{T} (GeV/c)"};
102108
const o2Axis energyAxis{100, 0, 100, "E (GeV)"};
103109
const o2Axis m02Axis{100, 0, 3, "m02"};
104-
110+
const o2Axis etaAxis{100, -1, 1, "#eta"};
111+
const o2Axis phiAxis{100, 0, 2 * TMath::Pi(), "#phi"};
112+
const o2Axis occupancyAxis{300, 0, 30000, "occupancy"};
105113
mHistograms.add("clusterE", "Energy of cluster", o2HistType::kTH1F, {energyAxis});
106114
mHistograms.add("trackPt", "pT of track", o2HistType::kTH1F, {ptAxis});
107115
mHistograms.add("chjetPt", "pT of charged jet", o2HistType::kTH1F, {ptAxis});
116+
mHistograms.add("chjetPtEtaPhi", "pT of charged jet", o2HistType::kTHnSparseF, {ptAxis, etaAxis, phiAxis});
108117
mHistograms.add("chjetpt_vs_constpt", "pT of charged jet vs pT of constituents", o2HistType::kTH2F, {ptAxis, ptAxis});
118+
119+
// track QA THnSparse
120+
mHistograms.add("trackPtEtaPhi", "Track QA", o2HistType::kTHnSparseF, {ptAxis, etaAxis, phiAxis});
121+
mHistograms.add("trackPtEtaOccupancy", "Track QA vs occupancy", o2HistType::kTHnSparseF, {ptAxis, etaAxis, occupancyAxis});
109122
}
110123

111124
// ---------------------
@@ -156,12 +169,25 @@ struct GammaJetTreeProducer {
156169
}
157170
return iso;
158171
}
159-
double ch_perp_cone_rho(const auto& cluster, aod::JetTracks const& tracks, float radius = 0.4)
172+
173+
void runTrackQA(const auto& collision, aod::JetTracks const& tracks)
174+
{
175+
for (auto track : tracks) {
176+
if (!isTrackSelected(track)) {
177+
continue;
178+
}
179+
mHistograms.fill(HIST("trackPt"), track.pt());
180+
mHistograms.fill(HIST("trackPtEtaPhi"), track.pt(), track.eta(), track.phi());
181+
mHistograms.fill(HIST("trackPtEtaOccupancy"), track.pt(), track.eta(), collision.trackOccupancyInTimeRange());
182+
}
183+
}
184+
185+
double ch_perp_cone_rho(const auto& object, aod::JetTracks const& tracks, float radius = 0.4)
160186
{
161187
double ptSumLeft = 0;
162188
double ptSumRight = 0;
163189

164-
double cPhi = TVector2::Phi_0_2pi(cluster.phi());
190+
double cPhi = TVector2::Phi_0_2pi(object.phi());
165191

166192
// rotate cone left by 90 degrees
167193
float cPhiLeft = cPhi - TMath::Pi() / 2;
@@ -173,8 +199,8 @@ struct GammaJetTreeProducer {
173199
if (!isTrackSelected(track)) {
174200
continue;
175201
}
176-
dRLeft = jetutilities::deltaR(cluster.eta(), cPhiLeft, track.eta(), track.phi());
177-
dRRight = jetutilities::deltaR(cluster.eta(), cPhiRight, track.eta(), track.phi());
202+
dRLeft = jetutilities::deltaR(object.eta(), cPhiLeft, track.eta(), track.phi());
203+
dRRight = jetutilities::deltaR(object.eta(), cPhiRight, track.eta(), track.phi());
178204

179205
if (dRLeft < radius) {
180206
ptSumLeft += track.pt();
@@ -201,16 +227,21 @@ struct GammaJetTreeProducer {
201227
// sadly passing of the string at runtime is not possible for technical region so cluster definition is
202228
// an integer instead
203229
Filter clusterDefinitionSelection = (o2::aod::jcluster::definition == mClusterDefinition);
230+
PresliceUnsorted<aod::JEMCTracks> EMCTrackPerTrack = aod::jemctrack::trackId;
231+
204232
// Process clusters
205-
void processClusters(soa::Join<aod::JetCollisions, aod::BkgChargedRhos, aod::JCollisionBCs>::iterator const& collision, selectedClusters const& clusters, aod::JetTracks const& tracks)
233+
void processClusters(soa::Join<aod::JetCollisions, aod::BkgChargedRhos, aod::JCollisionBCs>::iterator const& collision, selectedClusters const& clusters, aod::JetTracks const& tracks, aod::JEMCTracks const& emctracks)
206234
{
207235
if (!isEventAccepted(collision)) {
208236
return;
209237
}
210238

211-
eventsTable(collision.multiplicity(), collision.centrality(), collision.rho(), collision.eventSel(), collision.alias_raw());
239+
eventsTable(collision.multiplicity(), collision.centrality(), collision.rho(), collision.eventSel(), collision.trackOccupancyInTimeRange(), collision.alias_raw());
212240
collisionMapping[collision.globalIndex()] = eventsTable.lastIndex();
213241

242+
// loop over tracks one time for QA
243+
runTrackQA(collision, tracks);
244+
214245
// loop over clusters
215246
for (auto cluster : clusters) {
216247

@@ -226,26 +257,24 @@ struct GammaJetTreeProducer {
226257
// double dRMin = 100;
227258
double p = -1;
228259

229-
// auto tracksofcluster = matchedtracks.sliceBy(perClusterMatchedTracks, cluster.globalIndex());
230-
// for (const auto& match : tracksofcluster) {
231-
// // ask the jtracks table for track with ID trackID
232-
// double dR = deltaR(cluster.eta(), cluster.phi(), match.tracks_as<o2::aod::JTracks>().Eta(), match.tracks_as<o2::aod::JTracks>().Phi());
233-
// if (dR < dRMin) {
234-
// dRMin = dR;
235-
// dEta = cluster.eta() - match.tracks_as<o2::aod::JTracks>().eta();
236-
// dPhi = TVector2::Phi_0_2pi(cluster.phi()) - TVector2::Phi_0_2pi(match.tracks_as<o2::aod::JTracks>().phi());
237-
// if (abs(dPhi) > M_PI) {
238-
// dPhi = 2 * M_PI - abs(dPhi);
239-
// }
240-
// p = match.tracks_as<o2::aod::JTracks>().p();
241-
// }
242-
// }
243-
244-
// // for compression reasons make dPhi and dEta 0 if no match is found
245-
// if (p == -1) {
246-
// dPhi = 0;
247-
// dEta = 0;
248-
// }
260+
// do track matching
261+
auto tracksofcluster = cluster.matchedTracks_as<aod::JetTracks>();
262+
for (auto track : tracksofcluster) {
263+
if (!isTrackSelected(track)) {
264+
continue;
265+
}
266+
auto emcTracksPerTrack = emctracks.sliceBy(EMCTrackPerTrack, track.globalIndex());
267+
auto emcTrack = emcTracksPerTrack.iteratorAt(0);
268+
// find closest track that still has E/p < trackMatchingEoverP
269+
if (cluster.energy() / track.p() > trackMatchingEoverP) {
270+
continue;
271+
} else {
272+
dEta = cluster.eta() - emcTrack.etaEmcal();
273+
dPhi = RecoDecay::constrainAngle(RecoDecay::constrainAngle(emcTrack.phiEmcal(), -M_PI) - RecoDecay::constrainAngle(cluster.phi(), -M_PI), -M_PI);
274+
p = track.p();
275+
break;
276+
}
277+
}
249278

250279
gammasTable(eventsTable.lastIndex(), cluster.energy(), cluster.eta(), cluster.phi(), cluster.m02(), cluster.m20(), cluster.nCells(), cluster.time(), cluster.isExotic(), cluster.distanceToBadChannel(), cluster.nlm(), isoraw, perpconerho, dPhi, dEta, p);
251280
}
@@ -257,30 +286,38 @@ struct GammaJetTreeProducer {
257286
}
258287
PROCESS_SWITCH(GammaJetTreeProducer, processClusters, "Process EMCal clusters", true);
259288

260-
Filter jetCuts = aod::jet::pt > jetPtMin&& aod::jet::r == nround(jetR.node() * 100.0f);
289+
Filter jetCuts = aod::jet::pt > jetPtMin;
261290
// Process charged jets
262-
void processChargedJets(soa::Join<aod::JetCollisions, aod::BkgChargedRhos, aod::JCollisionBCs>::iterator const& collision, soa::Filtered<soa::Join<aod::ChargedJets, aod::ChargedJetConstituents>> const& chargedJets, aod::JetTracks const&)
291+
void processChargedJets(soa::Join<aod::JetCollisions, aod::BkgChargedRhos, aod::JCollisionBCs>::iterator const& collision, soa::Filtered<soa::Join<aod::ChargedJets, aod::ChargedJetConstituents>> const& chargedJets, aod::JetTracks const& tracks)
263292
{
264293
// event selection
265294
if (!isEventAccepted(collision)) {
266295
return;
267296
}
268-
297+
float leadingTrackPt = 0;
298+
ushort nconst = 0;
269299
// loop over charged jets
270300
for (auto jet : chargedJets) {
271301
if (jet.pt() < jetPtMin)
272302
continue;
273-
ushort nconst = 0;
303+
nconst = 0;
304+
leadingTrackPt = 0;
274305
// loop over constituents
275306
for (auto& constituent : jet.template tracks_as<aod::JetTracks>()) {
276307
mHistograms.fill(HIST("chjetpt_vs_constpt"), jet.pt(), constituent.pt());
277308
nconst++;
309+
if (constituent.pt() > leadingTrackPt) {
310+
leadingTrackPt = constituent.pt();
311+
}
278312
}
279313
int32_t storedColIndex = -1;
280314
if (auto foundCol = collisionMapping.find(collision.globalIndex()); foundCol != collisionMapping.end()) {
281315
storedColIndex = foundCol->second;
282316
}
283-
chargedJetsTable(storedColIndex, jet.pt(), jet.eta(), jet.phi(), jet.energy(), jet.mass(), jet.area(), nconst);
317+
// calculate perp cone rho
318+
double perpconerho = ch_perp_cone_rho(jet, tracks, perpConeJetR);
319+
mHistograms.fill(HIST("chjetPtEtaPhi"), jet.pt(), jet.eta(), jet.phi());
320+
chargedJetsTable(storedColIndex, jet.pt(), jet.eta(), jet.phi(), jet.r(), jet.energy(), jet.mass(), jet.area(), leadingTrackPt, perpconerho, nconst);
284321
// fill histograms
285322
mHistograms.fill(HIST("chjetPt"), jet.pt());
286323
}

0 commit comments

Comments
 (0)