Skip to content

Commit 343514d

Browse files
wefeng1110Wenhui Feng
wefeng1110
and
Wenhui Feng
authored
[PWGJE] Add workflows for charged jet spectra analysis (#8884)
Co-authored-by: Wenhui Feng <whfeng@localhost.localdomain>
1 parent 9316a78 commit 343514d

File tree

3 files changed

+858
-0
lines changed

3 files changed

+858
-0
lines changed

‎PWGJE/Tasks/CMakeLists.txt‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,10 @@ o2physics_add_dpl_workflow(photon-isolation-qa
4848
COMPONENT_NAME Analysis)
4949

5050
if(FastJet_FOUND)
51+
o2physics_add_dpl_workflow(jet-background-analysis
52+
SOURCES jetBackgroundAnalysis.cxx
53+
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
54+
COMPONENT_NAME Analysis)
5155
o2physics_add_dpl_workflow(jet-substructure
5256
SOURCES jetSubstructure.cxx
5357
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
@@ -132,6 +136,10 @@ o2physics_add_dpl_workflow(jet-substructure-bplus
132136
SOURCES jetFinderV0QA.cxx
133137
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
134138
COMPONENT_NAME Analysis)
139+
o2physics_add_dpl_workflow(jet-charged-spectra
140+
SOURCES jetSpectraCharged.cxx
141+
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore
142+
COMPONENT_NAME Analysis)
135143
o2physics_add_dpl_workflow(trigger-correlations
136144
SOURCES triggerCorrelations.cxx
137145
PUBLIC_LINK_LIBRARIES O2::Framework O2::EMCALBase O2::EMCALCalib O2Physics::PWGJECore O2Physics::AnalysisCore

‎PWGJE/Tasks/jetBackgroundAnalysis.cxx‎

Lines changed: 235 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,235 @@
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 finder QA task
13+
//
14+
/// \author Aimeric Landou <aimeric.landou@cern.ch>
15+
/// \author Nima Zardoshti <nima.zardoshti@cern.ch>
16+
17+
#include <cmath>
18+
#include <TRandom3.h>
19+
#include <string>
20+
#include "TLorentzVector.h"
21+
22+
#include "Framework/ASoA.h"
23+
#include "Framework/AnalysisDataModel.h"
24+
#include "Framework/AnalysisTask.h"
25+
#include "Framework/O2DatabasePDGPlugin.h"
26+
#include "Framework/HistogramRegistry.h"
27+
#include "Framework/runDataProcessing.h"
28+
29+
#include "Common/Core/TrackSelection.h"
30+
#include "Common/Core/TrackSelectionDefaults.h"
31+
32+
#include "Common/DataModel/EventSelection.h"
33+
#include "Common/DataModel/TrackSelectionTables.h"
34+
35+
#include "PWGJE/Core/FastJetUtilities.h"
36+
#include "PWGJE/Core/JetFinder.h"
37+
#include "PWGJE/Core/JetFindingUtilities.h"
38+
#include "PWGJE/DataModel/Jet.h"
39+
40+
#include "PWGJE/Core/JetDerivedDataUtilities.h"
41+
42+
#include "EventFiltering/filterTables.h"
43+
44+
using namespace o2;
45+
using namespace o2::framework;
46+
using namespace o2::framework::expressions;
47+
48+
struct JetBackgroundAnalysisTask {
49+
HistogramRegistry registry;
50+
51+
Configurable<std::string> eventSelections{"eventSelections", "sel8", "choose event selection"};
52+
Configurable<float> vertexZCut{"vertexZCut", 10.0f, "Accepted z-vertex range for collisions"};
53+
Configurable<float> centralityMin{"centralityMin", -999.0, "minimum centrality for collisions"};
54+
Configurable<float> centralityMax{"centralityMax", 999.0, "maximum centrality for collisions"};
55+
Configurable<int> trackOccupancyInTimeRangeMax{"trackOccupancyInTimeRangeMax", 999999, "maximum track occupancy of collisions in neighbouring collisions in a given time range; only applied to reconstructed collisions (data and mcd jets), not mc collisions (mcp jets)"};
56+
Configurable<int> trackOccupancyInTimeRangeMin{"trackOccupancyInTimeRangeMin", -999999, "minimum track occupancy of collisions in neighbouring collisions in a given time range; only applied to reconstructed collisions (data and mcd jets), not mc collisions (mcp jets)"};
57+
58+
Configurable<float> trackEtaMin{"trackEtaMin", -0.9, "minimum eta acceptance for tracks"};
59+
Configurable<float> trackEtaMax{"trackEtaMax", 0.9, "maximum eta acceptance for tracks"};
60+
Configurable<float> trackPtMin{"trackPtMin", 0.15, "minimum pT acceptance for tracks"};
61+
Configurable<float> trackPtMax{"trackPtMax", 100.0, "maximum pT acceptance for tracks"};
62+
Configurable<std::string> trackSelections{"trackSelections", "globalTracks", "set track selections"};
63+
64+
// // If weights are implemented for MCD bkg checks, the cut on pTHatMax of jets should be introduced
65+
// Configurable<float> pTHatMaxMCD{"pTHatMaxMCD", 999.0, "maximum fraction of hard scattering for jet acceptance in detector MC"};
66+
// Configurable<float> pTHatExponent{"pTHatExponent", 6.0, "exponent of the event weight for the calculation of pTHat"};
67+
68+
Configurable<float> randomConeR{"randomConeR", 0.4, "size of random Cone for estimating background fluctuations"};
69+
Configurable<float> randomConeLeadJetDeltaR{"randomConeLeadJetDeltaR", -99.0, "min distance between leading jet axis and random cone (RC) axis; if negative, min distance is set to automatic value of R_leadJet+R_RC "};
70+
71+
int eventSelection = -1;
72+
int trackSelection = -1;
73+
74+
void init(o2::framework::InitContext&)
75+
{
76+
// selection settings initialisation
77+
eventSelection = jetderiveddatautilities::initialiseEventSelection(static_cast<std::string>(eventSelections));
78+
trackSelection = jetderiveddatautilities::initialiseTrackSelection(static_cast<std::string>(trackSelections));
79+
80+
// histogram definitions
81+
82+
if (doprocessRho) {
83+
registry.add("h2_centrality_ntracks", "; centrality; N_{tracks};", {HistType::kTH2F, {{1100, 0., 110.0}, {10000, 0.0, 10000.0}}});
84+
registry.add("h2_ntracks_rho", "; N_{tracks}; #it{rho} (GeV/area);", {HistType::kTH2F, {{10000, 0.0, 10000.0}, {400, 0.0, 400.0}}});
85+
registry.add("h2_ntracks_rhom", "; N_{tracks}; #it{rho}_{m} (GeV/area);", {HistType::kTH2F, {{10000, 0.0, 10000.0}, {100, 0.0, 100.0}}});
86+
registry.add("h2_centrality_rho", "; centrality; #it{rho} (GeV/area);", {HistType::kTH2F, {{1100, 0., 110.}, {400, 0., 400.0}}});
87+
registry.add("h2_centrality_rhom", ";centrality; #it{rho}_{m} (GeV/area)", {HistType::kTH2F, {{1100, 0., 110.}, {100, 0., 100.0}}});
88+
}
89+
90+
if (doprocessBkgFluctuationsData || doprocessBkgFluctuationsMCD) {
91+
registry.add("h2_centrality_rhorandomcone", "; centrality; #it{p}_{T,random cone} - #it{area, random cone} * #it{rho} (GeV/c);", {HistType::kTH2F, {{1100, 0., 110.}, {800, -400.0, 400.0}}});
92+
registry.add("h2_centrality_rhorandomconerandomtrackdirection", "; centrality; #it{p}_{T,random cone} - #it{area, random cone} * #it{rho} (GeV/c);", {HistType::kTH2F, {{1100, 0., 110.}, {800, -400.0, 400.0}}});
93+
registry.add("h2_centrality_rhorandomconewithoutleadingjet", "; centrality; #it{p}_{T,random cone} - #it{area, random cone} * #it{rho} (GeV/c);", {HistType::kTH2F, {{1100, 0., 110.}, {800, -400.0, 400.0}}});
94+
registry.add("h2_centrality_rhorandomconerandomtrackdirectionwithoutoneleadingjets", "; centrality; #it{p}_{T,random cone} - #it{area, random cone} * #it{rho} (GeV/c);", {HistType::kTH2F, {{1100, 0., 110.}, {800, -400.0, 400.0}}});
95+
registry.add("h2_centrality_rhorandomconerandomtrackdirectionwithouttwoleadingjets", "; centrality; #it{p}_{T,random cone} - #it{area, random cone} * #it{rho} (GeV/c);", {HistType::kTH2F, {{1100, 0., 110.}, {800, -400.0, 400.0}}});
96+
}
97+
}
98+
99+
Filter trackCuts = (aod::jtrack::pt >= trackPtMin && aod::jtrack::pt < trackPtMax && aod::jtrack::eta > trackEtaMin && aod::jtrack::eta < trackEtaMax);
100+
Filter eventCuts = (nabs(aod::jcollision::posZ) < vertexZCut && aod::jcollision::centrality >= centralityMin && aod::jcollision::centrality < centralityMax);
101+
102+
template <typename TTracks, typename TJets>
103+
bool trackIsInJet(TTracks const& track, TJets const& jet)
104+
{
105+
for (auto const& constituentId : jet.tracksIds()) {
106+
if (constituentId == track.globalIndex()) {
107+
return true;
108+
}
109+
}
110+
return false;
111+
}
112+
113+
template <typename TCollisions, typename TJets, typename TTracks>
114+
void bkgFluctuationsRandomCone(TCollisions const& collision, TJets const& jets, TTracks const& tracks)
115+
{
116+
if (!jetderiveddatautilities::selectCollision(collision, eventSelection)) {
117+
return;
118+
}
119+
if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) {
120+
return;
121+
}
122+
TRandom3 randomNumber(0);
123+
float randomConeEta = randomNumber.Uniform(trackEtaMin + randomConeR, trackEtaMax - randomConeR);
124+
float randomConePhi = randomNumber.Uniform(0.0, 2 * M_PI);
125+
float randomConePt = 0;
126+
for (auto const& track : tracks) {
127+
if (jetderiveddatautilities::selectTrack(track, trackSelection)) {
128+
float dPhi = RecoDecay::constrainAngle(track.phi() - randomConePhi, static_cast<float>(-M_PI));
129+
float dEta = track.eta() - randomConeEta;
130+
if (TMath::Sqrt(dEta * dEta + dPhi * dPhi) < randomConeR) {
131+
randomConePt += track.pt();
132+
}
133+
}
134+
}
135+
registry.fill(HIST("h2_centrality_rhorandomcone"), collision.centrality(), randomConePt - M_PI * randomConeR * randomConeR * collision.rho());
136+
137+
// randomised eta,phi for tracks, to assess part of fluctuations coming from statistically independently emitted particles
138+
randomConePt = 0;
139+
for (auto const& track : tracks) {
140+
if (jetderiveddatautilities::selectTrack(track, trackSelection)) {
141+
float dPhi = RecoDecay::constrainAngle(randomNumber.Uniform(0.0, 2 * M_PI) - randomConePhi, static_cast<float>(-M_PI)); // ignores actual phi of track
142+
float dEta = randomNumber.Uniform(trackEtaMin, trackEtaMax) - randomConeEta; // ignores actual eta of track
143+
if (TMath::Sqrt(dEta * dEta + dPhi * dPhi) < randomConeR) {
144+
randomConePt += track.pt();
145+
}
146+
}
147+
}
148+
registry.fill(HIST("h2_centrality_rhorandomconerandomtrackdirection"), collision.centrality(), randomConePt - M_PI * randomConeR * randomConeR * collision.rho());
149+
150+
// removing the leading jet from the random cone
151+
if (jets.size() > 0) { // if there are no jets in the acceptance (from the jetfinder cuts) then there can be no leading jet
152+
float dPhiLeadingJet = RecoDecay::constrainAngle(jets.iteratorAt(0).phi() - randomConePhi, static_cast<float>(-M_PI));
153+
float dEtaLeadingJet = jets.iteratorAt(0).eta() - randomConeEta;
154+
155+
bool jetWasInCone = false;
156+
while ((randomConeLeadJetDeltaR <= 0 && (TMath::Sqrt(dEtaLeadingJet * dEtaLeadingJet + dPhiLeadingJet * dPhiLeadingJet) < jets.iteratorAt(0).r() / 100.0 + randomConeR)) || (randomConeLeadJetDeltaR > 0 && (TMath::Sqrt(dEtaLeadingJet * dEtaLeadingJet + dPhiLeadingJet * dPhiLeadingJet) < randomConeLeadJetDeltaR))) {
157+
jetWasInCone = true;
158+
randomConeEta = randomNumber.Uniform(trackEtaMin + randomConeR, trackEtaMax - randomConeR);
159+
randomConePhi = randomNumber.Uniform(0.0, 2 * M_PI);
160+
dPhiLeadingJet = RecoDecay::constrainAngle(jets.iteratorAt(0).phi() - randomConePhi, static_cast<float>(-M_PI));
161+
dEtaLeadingJet = jets.iteratorAt(0).eta() - randomConeEta;
162+
}
163+
if (jetWasInCone) {
164+
randomConePt = 0.0;
165+
for (auto const& track : tracks) {
166+
if (jetderiveddatautilities::selectTrack(track, trackSelection)) { // if track selection is uniformTrack, dcaXY and dcaZ cuts need to be added as they aren't in the selection so that they can be studied here
167+
float dPhi = RecoDecay::constrainAngle(track.phi() - randomConePhi, static_cast<float>(-M_PI));
168+
float dEta = track.eta() - randomConeEta;
169+
if (TMath::Sqrt(dEta * dEta + dPhi * dPhi) < randomConeR) {
170+
randomConePt += track.pt();
171+
}
172+
}
173+
}
174+
}
175+
}
176+
177+
registry.fill(HIST("h2_centrality_rhorandomconewithoutleadingjet"), collision.centrality(), randomConePt - M_PI * randomConeR * randomConeR * collision.rho());
178+
179+
// randomised eta,phi for tracks, to assess part of fluctuations coming from statistically independently emitted particles, removing tracks from 2 leading jets
180+
double randomConePtWithoutOneLeadJet = 0;
181+
double randomConePtWithoutTwoLeadJet = 0;
182+
for (auto const& track : tracks) {
183+
if (jetderiveddatautilities::selectTrack(track, trackSelection)) {
184+
float dPhi = RecoDecay::constrainAngle(randomNumber.Uniform(0.0, 2 * M_PI) - randomConePhi, static_cast<float>(-M_PI)); // ignores actual phi of track
185+
float dEta = randomNumber.Uniform(trackEtaMin, trackEtaMax) - randomConeEta; // ignores actual eta of track
186+
if (TMath::Sqrt(dEta * dEta + dPhi * dPhi) < randomConeR) {
187+
if (!trackIsInJet(track, jets.iteratorAt(0))) {
188+
randomConePtWithoutOneLeadJet += track.pt();
189+
if (!trackIsInJet(track, jets.iteratorAt(1))) {
190+
randomConePtWithoutTwoLeadJet += track.pt();
191+
}
192+
}
193+
}
194+
}
195+
}
196+
registry.fill(HIST("h2_centrality_rhorandomconerandomtrackdirectionwithoutoneleadingjets"), collision.centrality(), randomConePtWithoutOneLeadJet - M_PI * randomConeR * randomConeR * collision.rho());
197+
registry.fill(HIST("h2_centrality_rhorandomconerandomtrackdirectionwithouttwoleadingjets"), collision.centrality(), randomConePtWithoutTwoLeadJet - M_PI * randomConeR * randomConeR * collision.rho());
198+
}
199+
200+
void processRho(soa::Filtered<soa::Join<aod::JetCollisions, aod::BkgChargedRhos>>::iterator const& collision, soa::Filtered<aod::JetTracks> const& tracks)
201+
{
202+
if (!jetderiveddatautilities::selectCollision(collision, eventSelection)) {
203+
return;
204+
}
205+
if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) {
206+
return;
207+
}
208+
int nTracks = 0;
209+
for (auto const& track : tracks) {
210+
if (jetderiveddatautilities::selectTrack(track, trackSelection)) {
211+
nTracks++;
212+
}
213+
}
214+
registry.fill(HIST("h2_centrality_ntracks"), collision.centrality(), nTracks);
215+
registry.fill(HIST("h2_ntracks_rho"), nTracks, collision.rho());
216+
registry.fill(HIST("h2_ntracks_rhom"), nTracks, collision.rhoM());
217+
registry.fill(HIST("h2_centrality_rho"), collision.centrality(), collision.rho());
218+
registry.fill(HIST("h2_centrality_rhom"), collision.centrality(), collision.rhoM());
219+
}
220+
PROCESS_SWITCH(JetBackgroundAnalysisTask, processRho, "QA for rho-area subtracted jets", false);
221+
222+
void processBkgFluctuationsData(soa::Filtered<soa::Join<aod::JetCollisions, aod::BkgChargedRhos>>::iterator const& collision, soa::Join<aod::ChargedJets, aod::ChargedJetConstituents> const& jets, soa::Filtered<aod::JetTracks> const& tracks)
223+
{
224+
bkgFluctuationsRandomCone(collision, jets, tracks);
225+
}
226+
PROCESS_SWITCH(JetBackgroundAnalysisTask, processBkgFluctuationsData, "QA for random cone estimation of background fluctuations in data", false);
227+
228+
void processBkgFluctuationsMCD(soa::Filtered<soa::Join<aod::JetCollisions, aod::BkgChargedRhos>>::iterator const& collision, soa::Join<aod::ChargedMCDetectorLevelJets, aod::ChargedMCDetectorLevelJetConstituents> const& jets, soa::Filtered<aod::JetTracks> const& tracks)
229+
{
230+
bkgFluctuationsRandomCone(collision, jets, tracks);
231+
}
232+
PROCESS_SWITCH(JetBackgroundAnalysisTask, processBkgFluctuationsMCD, "QA for random cone estimation of background fluctuations in mcd", false);
233+
};
234+
235+
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{adaptAnalysisTask<JetBackgroundAnalysisTask>(cfgc, TaskName{"jet-background-analysis"})}; }

0 commit comments

Comments
 (0)