88// In applying this license CERN does not waive the privileges and immunities
99// granted to it by virtue of its status as an Intergovernmental Organization
1010// 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>
1811//
12+ // / \file jetDsSpecSubs.cxx
13+ // / \brief Ds-tagged jet analysis with substructure histogram outputs
14+ // / \author Monalisa Melo <monalisa.melo@cern.ch>, Universidade de São Paulo
1915
2016#include " PWGJE/Core/JetDerivedDataUtilities.h"
2117#include " PWGJE/Core/JetUtilities.h"
2218#include " PWGJE/DataModel/Jet.h"
23- #include " PWGJE/DataModel/JetReducedData.h"
2419
2520#include < Framework/ASoA.h>
2621#include < Framework/AnalysisDataModel.h>
27- #include < Framework/AnalysisHelpers.h>
2822#include < Framework/AnalysisTask.h>
2923#include < Framework/ConfigContext.h>
3024#include < Framework/Configurable.h>
25+ #include < Framework/DataProcessorSpec.h>
3126#include < Framework/HistogramRegistry.h>
3227#include < Framework/HistogramSpec.h>
3328#include < Framework/InitContext.h>
@@ -43,51 +38,6 @@ using namespace o2;
4338using namespace o2 ::framework;
4439using namespace o2 ::framework::expressions;
4540
46- namespace o2 ::aod
47- {
48- namespace jet_distance
49- {
50- DECLARE_SOA_COLUMN (JetHfDist, jetHfDist, float );
51- DECLARE_SOA_COLUMN (JetPt, jetPt, float );
52- DECLARE_SOA_COLUMN (JetEta, jetEta, float );
53- DECLARE_SOA_COLUMN (JetPhi, jetPhi, float );
54- DECLARE_SOA_COLUMN (JetNConst, jetNConst, int );
55- DECLARE_SOA_COLUMN (HfPt, hfPt, float );
56- DECLARE_SOA_COLUMN (HfEta, hfEta, float );
57- DECLARE_SOA_COLUMN (HfPhi, hfPhi, float );
58- DECLARE_SOA_COLUMN (HfMass, hfMass, float );
59- DECLARE_SOA_COLUMN (HfY, hfY, float );
60- DECLARE_SOA_COLUMN (HfMlScore0, hfMlScore0, float );
61- DECLARE_SOA_COLUMN (HfMlScore1, hfMlScore1, float );
62- DECLARE_SOA_COLUMN (HfMlScore2, hfMlScore2, float );
63-
64- // extra
65- DECLARE_SOA_COLUMN (JetMass, jetMass, float );
66- DECLARE_SOA_COLUMN (JetGirth, jetGirth, float );
67- DECLARE_SOA_COLUMN (JetThrust, jetThrust, float ); // lambda_2^1
68- DECLARE_SOA_COLUMN (JetLambda11, jetLambda11, float ); // lambda_1^1
69- } // namespace jet_distance
70-
71- DECLARE_SOA_TABLE (JetDistanceTable, " AOD" , " JETDISTTABLE" ,
72- jet_distance::JetHfDist,
73- jet_distance::JetPt,
74- jet_distance::JetEta,
75- jet_distance::JetPhi,
76- jet_distance::JetNConst,
77- jet_distance::HfPt,
78- jet_distance::HfEta,
79- jet_distance::HfPhi,
80- jet_distance::HfMass,
81- jet_distance::HfY,
82- jet_distance::HfMlScore0,
83- jet_distance::HfMlScore1,
84- jet_distance::HfMlScore2,
85- jet_distance::JetMass,
86- jet_distance::JetGirth,
87- jet_distance::JetThrust,
88- jet_distance::JetLambda11);
89- } // namespace o2::aod
90-
9141struct JetDsSpecSubs {
9242 HistogramRegistry registry{
9343 " registry" ,
@@ -100,13 +50,14 @@ struct JetDsSpecSubs {
10050 {" h_jet_eta" , " jet #eta;#eta_{jet};entries" , {HistType::kTH1F , {{100 , -1.0 , 1.0 }}}},
10151 {" h_jet_phi" , " jet #phi;#phi_{jet};entries" , {HistType::kTH1F , {{80 , -1.0 , 7 .}}}},
10252 {" h_collision_counter" , " # of collisions;" , {HistType::kTH1F , {{200 , 0 ., 200 .}}}},
103- {" h_jet_counter" , " ;# of D_{S} jets; " , {HistType::kTH1F , {{6 , 0 ., 3.0 }}}},
53+ {" h_jet_counter" , " ;type;counts " , {HistType::kTH1F , {{2 , 0 ., 2 . }}}},
10454 {" h_ds_jet_projection" , " ;z^{D_{S},jet}_{||};dN/dz^{D_{S},jet}_{||}" , {HistType::kTH1F , {{1000 , 0 ., 2 .}}}},
10555 {" h_ds_jet_distance_vs_projection" , " ;#DeltaR_{D_{S},jet};z^{D_{S},jet}_{||}" , {HistType::kTH2F , {{1000 , 0 ., 1 .}, {1000 , 0 ., 2 .}}}},
10656 {" h_ds_jet_distance" , " ;#DeltaR_{D_{S},jet};dN/d(#DeltaR)" , {HistType::kTH1F , {{1000 , 0 ., 1 .}}}},
10757 {" h_ds_jet_pt" , " ;p_{T,D_{S} jet};dN/dp_{T,D_{S} jet}" , {HistType::kTH1F , {{1000 , 0 ., 100 .}}}},
10858 {" h_ds_jet_eta" , " ;#eta_{D_{S} jet};entries" , {HistType::kTH1F , {{250 , -1 ., 1 .}}}},
10959 {" h_ds_jet_phi" , " ;#phi_{D_{S} jet};entries" , {HistType::kTH1F , {{250 , -1 ., 7 .}}}},
60+ {" hSparse_ds" , " ;m_{D_{S}};p_{T,D_{S}};z^{D_{S},jet}_{||};#DeltaR_{D_{S},jet}" , {HistType::kTHnSparseF , {{60 , 1.7 , 2.1 }, {60 , 0 ., 100 .}, {60 , 0 ., 2 .}, {60 , 0 ., 1.0 }}}},
11061 {" h_ds_mass" , " ;m_{D_{S}} (GeV/c^{2});entries" , {HistType::kTH1F , {{1000 , 0 ., 6 .}}}},
11162 {" h_ds_eta" , " ;#eta_{D_{S}};entries" , {HistType::kTH1F , {{250 , -1 ., 1 .}}}},
11263 {" h_ds_phi" , " ;#phi_{D_{S}};entries" , {HistType::kTH1F , {{250 , -1 ., 7 .}}}},
@@ -118,8 +69,6 @@ struct JetDsSpecSubs {
11869 {" h2_dsjet_pt_lambda12" , " ;#it{p}_{T,jet} (GeV/#it{c});#lambda_{2}^{1}" , {HistType::kTH2F , {{100 , 0 ., 100 .}, {200 , 0 ., 1.0 }}}},
11970 {" h2_dsjet_pt_mass" , " ;#it{p}_{T,jet} (GeV/#it{c});m_{jet}^{ch} (GeV/#it{c}^{2})" , {HistType::kTH2F , {{100 , 0 ., 100 .}, {200 , 0 ., 50.0 }}}},
12071 {" h2_dsjet_pt_girth" , " ;#it{p}_{T,jet} (GeV/#it{c});g" , {HistType::kTH2F , {{100 , 0 ., 100 .}, {200 , 0 ., 0.5 }}}},
121- {" h_ds_jet_lambda_extra" , " ;#lambda_{#alpha}^{#kappa};entries" , {HistType::kTH1F , {{200 , 0 ., 1.0 }}}},
122- {" h2_dsjet_pt_lambda_extra" , " ;#it{p}_{T,jet} (GeV/#it{c});#lambda_{#alpha}^{#kappa}" , {HistType::kTH2F , {{100 , 0 ., 100 .}, {200 , 0 ., 1.0 }}}},
12372 }};
12473
12574 Configurable<float > vertexZCut{" vertexZCut" , 10 .0f , " Accepted z-vertex range" };
@@ -129,17 +78,9 @@ struct JetDsSpecSubs {
12978 Configurable<std::string> eventSelections{" eventSelections" , " sel8" , " choose event selection" };
13079 Configurable<std::string> trackSelections{" trackSelections" , " globalTracks" , " set track selections" };
13180
132- // extra angularity knob
133- Configurable<float > kappa{" kappa" , 1 .0f , " angularity kappa" };
134- Configurable<float > alpha{" alpha" , 1 .0f , " angularity alpha" };
135-
136- bool doExtraAngularity = false ;
137-
13881 std::vector<int > eventSelectionBits;
13982 int trackSelection = -1 ;
14083
141- Produces<aod::JetDistanceTable> distJetTable;
142-
14384 template <typename JET, typename TRACKS>
14485 float computeLambda (JET const & jet, TRACKS const & tracks, float a, float k)
14586 {
@@ -151,8 +92,8 @@ struct JetDsSpecSubs {
15192 const float dr = jetutilities::deltaR (jet, trk);
15293 sum += std::pow (trk.pt (), k) * std::pow (dr, a);
15394 }
154- const float R = jet.r () / 100 .f ;
155- const float denom = std::pow (jet.pt (), k) * std::pow (R , a);
95+ const float jetR = jet.r () / 100 .f ;
96+ const float denom = std::pow (jet.pt (), k) * std::pow (jetR , a);
15697 if (denom <= 0 .f ) {
15798 return -1 .f ;
15899 }
@@ -189,9 +130,9 @@ struct JetDsSpecSubs {
189130 eventSelectionBits = jetderiveddatautilities::initialiseEventSelectionBits (static_cast <std::string>(eventSelections));
190131 trackSelection = jetderiveddatautilities::initialiseTrackSelection (static_cast <std::string>(trackSelections));
191132
192- const bool is11 = ( std::abs (kappa. value - 1 . f ) < 1e- 6f ) && ( std::abs (alpha. value - 1 . f ) < 1e- 6f );
193- const bool is12 = ( std::abs (kappa. value - 1 . f ) < 1e- 6f ) && ( std::abs (alpha. value - 2 . f ) < 1e- 6f );
194- doExtraAngularity = !(is11 || is12 );
133+ auto h = registry. get <TH1>( HIST ( " h_jet_counter " ) );
134+ h-> GetXaxis ()-> SetBinLabel ( 1 , " All jets " );
135+ h-> GetXaxis ()-> SetBinLabel ( 2 , " Ds-tagged jets " );
195136 }
196137
197138 Filter jetCuts = aod::jet::pt > jetPtMin&& aod::jet::r == nround(jetR.node() * 100 .0f );
@@ -203,6 +144,7 @@ struct JetDsSpecSubs {
203144 if (!jetderiveddatautilities::selectCollision (collision, eventSelectionBits)) {
204145 return ;
205146 }
147+
206148 registry.fill (HIST (" h_collisions" ), 1.5 );
207149
208150 for (auto const & track : tracks) {
@@ -222,7 +164,7 @@ struct JetDsSpecSubs {
222164 if (!jetderiveddatautilities::selectCollision (collision, eventSelectionBits)) {
223165 return ;
224166 }
225- for (auto & jet : jets) {
167+ for (const auto & jet : jets) {
226168 registry.fill (HIST (" h_jet_pt" ), jet.pt ());
227169 registry.fill (HIST (" h_jet_eta" ), jet.eta ());
228170 registry.fill (HIST (" h_jet_phi" ), jet.phi ());
@@ -243,77 +185,85 @@ struct JetDsSpecSubs {
243185 registry.fill (HIST (" h_collision_counter" ), 3.0 );
244186
245187 for (const auto & jet : jets) {
188+
246189 registry.fill (HIST (" h_jet_counter" ), 0.5 );
247190
191+ bool hasDsCandidate = false ;
192+
248193 TVector3 jetVector (jet.px (), jet.py (), jet.pz ());
249194
195+ // Compute jet-level quantities once (independent of Ds candidates)
196+ auto jetTracks = jet.tracks_as <aod::JetTracks>();
197+
198+ const float lambda11 = computeLambda (jet, jetTracks, 1 .f , 1 .f );
199+ const float lambda12 = computeLambda (jet, jetTracks, 2 .f , 1 .f );
200+ const float mjet = computeJetMassFromTracksMass (jetTracks);
201+
202+ const float jetR = jet.r () / 100 .f ;
203+ const float girth = (lambda11 >= 0 .f ) ? (lambda11 * jetR) : -1 .f ;
204+
205+ // Loop over Ds candidates (particle level)
250206 for (const auto & dsCandidate : jet.candidates_as <aod::CandidatesDsData>()) {
207+
208+ hasDsCandidate = true ;
209+
251210 TVector3 dsVector (dsCandidate.px (), dsCandidate.py (), dsCandidate.pz ());
252211
212+ // zParallel defined as longitudinal momentum fraction along the jet axis
253213 const double zParallel = (jetVector * dsVector) / (jetVector * jetVector);
254- const double axisDistance = jetutilities::deltaR (jet, dsCandidate);
214+ const float axisDistance = jetutilities::deltaR (jet, dsCandidate);
255215
216+ // --- Ds-level observables ---
256217 registry.fill (HIST (" h_ds_jet_projection" ), zParallel);
257218 registry.fill (HIST (" h_ds_jet_distance_vs_projection" ), axisDistance, zParallel);
258219 registry.fill (HIST (" h_ds_jet_distance" ), axisDistance);
259- registry.fill (HIST (" h_ds_jet_pt" ), jet.pt ());
260- registry.fill (HIST (" h_ds_jet_eta" ), jet.eta ());
261- registry.fill (HIST (" h_ds_jet_phi" ), jet.phi ());
220+
262221 registry.fill (HIST (" h_ds_mass" ), dsCandidate.m ());
263222 registry.fill (HIST (" h_ds_eta" ), dsCandidate.eta ());
264223 registry.fill (HIST (" h_ds_phi" ), dsCandidate.phi ());
265224
266- auto jetTracks = jet.tracks_as <aod::JetTracks>();
225+ // Main THnSparse: invariant mass, pT, z, and DeltaR
226+ registry.fill (HIST (" hSparse_ds" ),
227+ dsCandidate.m (),
228+ dsCandidate.pt (),
229+ zParallel,
230+ axisDistance);
231+ }
267232
268- const float lambda11 = computeLambda (jet, jetTracks, 1 .f , 1 .f );
269- const float lambda12 = computeLambda (jet, jetTracks, 2 .f , 1 .f ); // thrust = λ_2^1
270- const float mjet = computeJetMassFromTracksMass (jetTracks);
233+ // Jet-level quantities (filled once per jet containing at least one Ds)
234+ if (hasDsCandidate) {
271235
272- const float R = jet.r () / 100 .f ;
273- const float girth = (lambda11 >= 0 .f ) ? (lambda11 * R) : -1 .f ;
236+ registry.fill (HIST (" h_jet_counter" ), 1.5 );
237+
238+ // Jet properties
239+ registry.fill (HIST (" h_ds_jet_pt" ), jet.pt ());
240+ registry.fill (HIST (" h_ds_jet_eta" ), jet.eta ());
241+ registry.fill (HIST (" h_ds_jet_phi" ), jet.phi ());
274242
243+ // Jet substructure observables
275244 if (lambda11 >= 0 .f ) {
276245 registry.fill (HIST (" h_ds_jet_lambda11" ), lambda11);
277246 registry.fill (HIST (" h2_dsjet_pt_lambda11" ), jet.pt (), lambda11);
278247 }
248+
279249 if (lambda12 >= 0 .f ) {
280250 registry.fill (HIST (" h_ds_jet_lambda12" ), lambda12);
281251 registry.fill (HIST (" h2_dsjet_pt_lambda12" ), jet.pt (), lambda12);
282252 }
253+
283254 registry.fill (HIST (" h_ds_jet_mass" ), mjet);
284255 registry.fill (HIST (" h2_dsjet_pt_mass" ), jet.pt (), mjet);
285256
286257 if (girth >= 0 .f ) {
287258 registry.fill (HIST (" h_ds_jet_girth" ), girth);
288259 registry.fill (HIST (" h2_dsjet_pt_girth" ), jet.pt (), girth);
289260 }
290-
291- if (doExtraAngularity) {
292- const float lambdaExtra = computeLambda (jet, jetTracks, alpha.value , kappa.value );
293- if (lambdaExtra >= 0 .f ) {
294- registry.fill (HIST (" h_ds_jet_lambda_extra" ), lambdaExtra);
295- registry.fill (HIST (" h2_dsjet_pt_lambda_extra" ), jet.pt (), lambdaExtra);
296- }
297- }
298-
299- auto scores = dsCandidate.mlScores ();
300- const float s0 = (scores.size () > 0 ) ? scores[0 ] : -999 .f ;
301- const float s1 = (scores.size () > 1 ) ? scores[1 ] : -999 .f ;
302- const float s2 = (scores.size () > 2 ) ? scores[2 ] : -999 .f ;
303-
304- distJetTable (static_cast <float >(axisDistance),
305- jet.pt (), jet.eta (), jet.phi (),
306- static_cast <int >(jetTracks.size ()),
307- dsCandidate.pt (), dsCandidate.eta (), dsCandidate.phi (),
308- dsCandidate.m (), dsCandidate.y (),
309- s0, s1, s2,
310- mjet, girth, lambda12, lambda11);
311-
312- break ; // only first Ds per jet
313261 }
314262 }
315263 }
316264 PROCESS_SWITCH (JetDsSpecSubs, processDataChargedSubstructure, " charged HF jet substructure" , false );
317265};
318-
319- WorkflowSpec defineDataProcessing (ConfigContext const & cfgc) { return WorkflowSpec{adaptAnalysisTask<JetDsSpecSubs>(cfgc, TaskName{" jet-ds-spectrum-subs" })}; }
266+ WorkflowSpec defineDataProcessing (ConfigContext const & cfgc)
267+ {
268+ return WorkflowSpec{adaptAnalysisTask<JetDsSpecSubs>(cfgc)};
269+ }
0 commit comments