Skip to content

Commit 786854d

Browse files
monamelopmonamelop
authored
[PWGJE] Add jet-ds-spectrum-subs task (#14723)
1 parent efc3474 commit 786854d

File tree

2 files changed

+223
-1
lines changed

2 files changed

+223
-1
lines changed

‎PWGJE/Tasks/CMakeLists.txt‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -358,5 +358,9 @@ if(FastJet_FOUND)
358358
o2physics_add_dpl_workflow(substructure-debug
359359
SOURCES substructureDebug.cxx
360360
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
361+
COMPONENT_NAME Analysis)
362+
o2physics_add_dpl_workflow(jet-ds-spectrum-subs
363+
SOURCES jetDsSpectrumAndSubstructure.cxx
364+
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
361365
COMPONENT_NAME Analysis)
362-
endif()
366+
endif()

‎PWGJE/Tasks/jetDsSpectrumAndSubstructure.cxx‎

Lines changed: 218 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,218 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
12+
// Jet substructure and spectrum task for D_s mesons
13+
//
14+
// This task is used to reconstruct and analyse jets containing charged D_s
15+
// mesons
16+
//
17+
/// \author Monalisa Melo <monalisa.melo@cern.ch>
18+
//
19+
20+
#include "PWGHF/Core/DecayChannels.h"
21+
#include "PWGJE/Core/JetDerivedDataUtilities.h"
22+
#include "PWGJE/Core/JetUtilities.h"
23+
#include "PWGJE/DataModel/Jet.h"
24+
#include "PWGJE/DataModel/JetReducedData.h"
25+
#include "PWGJE/DataModel/JetSubtraction.h"
26+
27+
#include "Common/Core/RecoDecay.h"
28+
#include "Common/DataModel/TrackSelectionTables.h"
29+
30+
#include "Framework/ASoA.h"
31+
#include "Framework/AnalysisDataModel.h"
32+
#include "Framework/HistogramRegistry.h"
33+
#include <CommonConstants/MathConstants.h>
34+
#include <Framework/AnalysisTask.h>
35+
#include <Framework/ConfigContext.h>
36+
#include <Framework/Configurable.h>
37+
#include <Framework/DataProcessorSpec.h>
38+
#include <Framework/HistogramSpec.h>
39+
#include <Framework/InitContext.h>
40+
#include <Framework/runDataProcessing.h>
41+
42+
#include "TVector3.h"
43+
44+
#include <cmath>
45+
#include <cstdlib>
46+
#include <string>
47+
#include <vector>
48+
49+
using namespace o2;
50+
using namespace o2::framework;
51+
using namespace o2::framework::expressions;
52+
53+
namespace o2::aod
54+
{
55+
namespace jet_distance
56+
{
57+
DECLARE_SOA_COLUMN(JetHfDist, jetHfDist, float);
58+
DECLARE_SOA_COLUMN(JetPt, jetPt, float);
59+
DECLARE_SOA_COLUMN(JetEta, jetEta, float);
60+
DECLARE_SOA_COLUMN(JetPhi, jetPhi, float);
61+
DECLARE_SOA_COLUMN(JetNConst, jetNConst, int);
62+
DECLARE_SOA_COLUMN(HfPt, hfPt, float);
63+
DECLARE_SOA_COLUMN(HfEta, hfEta, float);
64+
DECLARE_SOA_COLUMN(HfPhi, hfPhi, float);
65+
DECLARE_SOA_COLUMN(HfMass, hfMass, float);
66+
DECLARE_SOA_COLUMN(HfY, hfY, float);
67+
DECLARE_SOA_COLUMN(HfMlScore0, hfMlScore0, float);
68+
DECLARE_SOA_COLUMN(HfMlScore1, hfMlScore1, float);
69+
DECLARE_SOA_COLUMN(HfMlScore2, hfMlScore2, float);
70+
} // namespace jet_distance
71+
72+
DECLARE_SOA_TABLE(JetDistanceTable, "AOD", "JETDISTTABLE",
73+
jet_distance::JetHfDist,
74+
jet_distance::JetPt,
75+
jet_distance::JetEta,
76+
jet_distance::JetPhi,
77+
jet_distance::JetNConst,
78+
jet_distance::HfPt,
79+
jet_distance::HfEta,
80+
jet_distance::HfPhi,
81+
jet_distance::HfMass,
82+
jet_distance::HfY,
83+
jet_distance::HfMlScore0,
84+
jet_distance::HfMlScore1,
85+
jet_distance::HfMlScore2);
86+
} // namespace o2::aod
87+
88+
struct JetDsSpecSubs {
89+
HistogramRegistry registry{"registry",
90+
{{"h_collisions", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}}},
91+
{"h_track_pt", "track pT;#it{p}_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}}},
92+
{"h_track_eta", "track #eta;#eta_{track};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}}},
93+
{"h_track_phi", "track #varphi;#varphi_{track};entries", {HistType::kTH1F, {{80, -1.0, 7.}}}},
94+
{"h_jet_pt", "jet pT;#it{p}_{T,jet} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0., 200.}}}},
95+
{"h_jet_eta", "jet #eta;#eta_{jet};entries", {HistType::kTH1F, {{100, -1.0, 1.0}}}},
96+
{"h_jet_phi", "jet #phi;#phi_{jet};entries", {HistType::kTH1F, {{80, -1.0, 7.}}}},
97+
{"h_collision_counter", "# of collisions;", {HistType::kTH1F, {{200, 0., 200.}}}},
98+
{"h_jet_counter", ";# of D_{S} jets;", {HistType::kTH1F, {{6, 0., 3.0}}}},
99+
{"h_ds_jet_projection", ";z^{D_{S},jet}_{||};dN/dz^{D_{S},jet}_{||}", {HistType::kTH1F, {{1000, 0., 10.}}}},
100+
{"h_ds_jet_distance_vs_projection", ";#DeltaR_{D_{S},jet};z^{D_{S},jet}_{||}", {HistType::kTH2F, {{1000, 0., 10.}, {1000, 0., 10.}}}},
101+
{"h_ds_jet_distance", ";#DeltaR_{D_{S},jet};dN/d(#DeltaR)", {HistType::kTH1F, {{1000, 0., 10.}}}},
102+
{"h_ds_jet_pt", ";p_{T,D_{S} jet};dN/dp_{T,D_{S} jet}", {HistType::kTH1F, {{200, 0., 10.}}}},
103+
{"h_ds_jet_eta", ";#eta_{T,D_{S} jet};dN/d#eta_{D_{S} jet}", {HistType::kTH1F, {{250, -5., 5.}}}},
104+
{"h_ds_jet_phi", ";#phi_{T,D_{S} jet};dN/d#phi_{D_{S} jet}", {HistType::kTH1F, {{250, -10., 10.}}}},
105+
{"h_ds_mass", ";m_{D_{S}} (GeV/c^{2});dN/dm_{D_{S}}", {HistType::kTH1F, {{1000, 0., 10.}}}},
106+
{"h_ds_eta", ";#eta_{D_{S}} (GeV/c^{2});dN/d#eta_{D_{S}}", {HistType::kTH1F, {{250, -5., 5.}}}},
107+
{"h_ds_phi", ";#phi_{D_{S}} (GeV/c^{2});dN/d#phi_{D_{S}}", {HistType::kTH1F, {{250, -10., 10.}}}}}};
108+
109+
Configurable<float> vertexZCut{"vertexZCut", 10.0f, "Accepted z-vertex range"};
110+
111+
Configurable<float> jetPtMin{"jetPtMin", 5.0, "minimum jet pT cut"};
112+
Configurable<float> jetR{"jetR", 0.4, "jet resolution parameter"};
113+
114+
Configurable<std::string> eventSelections{"eventSelections", "sel8", "choose event selection"};
115+
Configurable<std::string> trackSelections{"trackSelections", "globalTracks", "set track selections"};
116+
117+
std::vector<int> eventSelectionBits;
118+
int trackSelection = -1;
119+
120+
Produces<aod::JetDistanceTable> distJetTable;
121+
122+
void init(o2::framework::InitContext&)
123+
{
124+
eventSelectionBits = jetderiveddatautilities::initialiseEventSelectionBits(static_cast<std::string>(eventSelections));
125+
trackSelection = jetderiveddatautilities::initialiseTrackSelection(static_cast<std::string>(trackSelections));
126+
}
127+
128+
Filter jetCuts = aod::jet::pt > jetPtMin&& aod::jet::r == nround(jetR.node() * 100.0f);
129+
Filter collisionFilter = nabs(aod::jcollision::posZ) < vertexZCut;
130+
131+
void processCollisions(aod::JetCollision const& collision, aod::JetTracks const& tracks)
132+
{
133+
134+
registry.fill(HIST("h_collisions"), 0.5);
135+
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits)) {
136+
return;
137+
}
138+
registry.fill(HIST("h_collisions"), 1.5);
139+
for (auto const& track : tracks) {
140+
if (!jetderiveddatautilities::selectTrack(track, trackSelection)) {
141+
continue;
142+
}
143+
registry.fill(HIST("h_track_pt"), track.pt());
144+
registry.fill(HIST("h_track_eta"), track.eta());
145+
registry.fill(HIST("h_track_phi"), track.phi());
146+
}
147+
}
148+
PROCESS_SWITCH(JetDsSpecSubs, processCollisions, "process JE collisions", false);
149+
150+
void processDataCharged(soa::Filtered<aod::JetCollisions>::iterator const& collision, soa::Filtered<aod::ChargedJets> const& jets)
151+
{
152+
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits)) {
153+
return;
154+
}
155+
// jets -> charged jet
156+
for (auto& jet : jets) {
157+
registry.fill(HIST("h_jet_pt"), jet.pt());
158+
registry.fill(HIST("h_jet_eta"), jet.eta());
159+
registry.fill(HIST("h_jet_phi"), jet.phi());
160+
}
161+
}
162+
PROCESS_SWITCH(JetDsSpecSubs, processDataCharged, "charged jets in data", false);
163+
164+
void processDataChargedSubstructure(aod::JetCollision const& collision,
165+
soa::Join<aod::DsChargedJets, aod::DsChargedJetConstituents> const& jets,
166+
aod::CandidatesDsData const&,
167+
aod::JetTracks const&)
168+
{
169+
// apply event selection and fill histograms for sanity check
170+
registry.fill(HIST("h_collision_counter"), 2.0);
171+
if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits) || !(std::abs(collision.posZ()) < vertexZCut)) {
172+
return;
173+
}
174+
registry.fill(HIST("h_collision_counter"), 3.0);
175+
// jets -> charged jet with Ds
176+
for (const auto& jet : jets) {
177+
// number of charged jets with Ds
178+
registry.fill(HIST("h_jet_counter"), 0.5);
179+
// obtaining jet 3-vector
180+
TVector3 jetVector(jet.px(), jet.py(), jet.pz());
181+
182+
for (const auto& dsCandidate : jet.candidates_as<aod::CandidatesDsData>()) {
183+
184+
// obtaining jet 3-vector
185+
TVector3 dsVector(dsCandidate.px(), dsCandidate.py(), dsCandidate.pz());
186+
187+
// calculating fraction of the jet momentum carried by the Ds along the direction of the jet axis
188+
double zParallel = (jetVector * dsVector) / (jetVector * jetVector);
189+
190+
// calculating angular distance in eta-phi plane
191+
double axisDistance = jetutilities::deltaR(jet, dsCandidate);
192+
193+
// filling histograms
194+
registry.fill(HIST("h_ds_jet_projection"), zParallel);
195+
registry.fill(HIST("h_ds_jet_distance_vs_projection"), axisDistance, zParallel);
196+
registry.fill(HIST("h_ds_jet_distance"), axisDistance);
197+
registry.fill(HIST("h_ds_jet_pt"), jet.pt());
198+
registry.fill(HIST("h_ds_jet_eta"), jet.eta());
199+
registry.fill(HIST("h_ds_jet_phi"), jet.phi());
200+
registry.fill(HIST("h_ds_mass"), dsCandidate.m());
201+
registry.fill(HIST("h_ds_eta"), dsCandidate.eta());
202+
registry.fill(HIST("h_ds_phi"), dsCandidate.phi());
203+
204+
// filling table
205+
distJetTable(axisDistance,
206+
jet.pt(), jet.eta(), jet.phi(), jet.tracks_as<aod::JetTracks>().size(),
207+
dsCandidate.pt(), dsCandidate.eta(), dsCandidate.phi(), dsCandidate.m(), dsCandidate.y(), dsCandidate.mlScores()[0], dsCandidate.mlScores()[1], dsCandidate.mlScores()[2]);
208+
209+
break; // get out of candidates' loop after first HF particle is found in jet
210+
} // end of DS candidates loop
211+
212+
} // end of jets loop
213+
214+
} // end of process function
215+
PROCESS_SWITCH(JetDsSpecSubs, processDataChargedSubstructure, "charged HF jet substructure", false);
216+
};
217+
218+
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{adaptAnalysisTask<JetDsSpecSubs>(cfgc, TaskName{"jet-ds-spectrum-subs"})}; }

0 commit comments

Comments
 (0)