Skip to content

Commit 2c300f5

Browse files
PushingMyRocheLimitsPushingMyRocheLimits
authored
[PWGUD] Created task for studying jspi using CB
1 parent 1b0014f commit 2c300f5

File tree

2 files changed

+332
-0
lines changed

2 files changed

+332
-0
lines changed

‎PWGUD/Tasks/CMakeLists.txt‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -263,3 +263,8 @@ o2physics_add_dpl_workflow(upc-test-rct-tables
263263
SOURCES upcTestRctTables.cxx
264264
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
265265
COMPONENT_NAME Analysis)
266+
267+
o2physics_add_dpl_workflow(upc-exclusive-jpsi-cb
268+
SOURCES upcExclusiveJpsiCB.cxx
269+
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
270+
COMPONENT_NAME Analysis)

‎PWGUD/Tasks/upcExclusiveJpsiCB.cxx‎

Lines changed: 327 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,327 @@
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+
13+
#include <iostream>
14+
#include <cmath>
15+
#include <cstdlib>
16+
#include <string>
17+
#include <vector>
18+
#include "Math/Vector4D.h"
19+
#include "TDatabasePDG.h"
20+
#include "TMath.h"
21+
22+
#include "Framework/runDataProcessing.h"
23+
#include "Framework/AnalysisTask.h"
24+
#include "Framework/AnalysisDataModel.h"
25+
#include "PWGUD/DataModel/UDTables.h"
26+
#include "PWGUD/Core/SGSelector.h"
27+
#include "PWGUD/Core/SGTrackSelector.h"
28+
29+
using namespace o2;
30+
using namespace o2::aod;
31+
using namespace o2::framework;
32+
using namespace o2::framework::expressions;
33+
34+
using LorentzVector = ROOT::Math::PxPyPzMVector;
35+
36+
/// \brief Task to study two pions candidates using SG derive data
37+
/// \author Levi Van Ryder (based on Anisa Khatun's UD Tutorial 5 example)
38+
/// \date 02/12/2025
39+
40+
// Struct to define the analysis task
41+
struct upcExclusiveJpsiCB {
42+
// SGSelector object to manage track and collision selections
43+
SGSelector sgSelector;
44+
45+
Configurable<float> FV0_cut{"FV0", 100., "FV0A threshold"};
46+
Configurable<float> FT0A_cut{"FT0A", 100., "FT0A threshold"};
47+
Configurable<float> FT0C_cut{"FT0C", 50., "FT0C threshold"};
48+
Configurable<float> FDDA_cut{"FDDA", 10000., "FDDA threshold"};
49+
Configurable<float> FDDC_cut{"FDDC", 10000., "FDDC threshold"};
50+
Configurable<float> ZDC_cut{"ZDC", 10., "ZDC threshold"};
51+
Configurable<float> gap_Side{"gap", 2, "gap selection"};
52+
53+
// Track Selections
54+
Configurable<float> PV_cut{"PV_cut", 1.0, "Use Only PV tracks"};
55+
Configurable<float> dcaZ_cut{"dcaZ_cut", 2.0, "dcaZ cut"};
56+
Configurable<float> dcaXY_cut{"dcaXY_cut", 0.0, "dcaXY cut (0 for Pt-function)"};
57+
Configurable<float> tpcChi2_cut{"tpcChi2_cut", 4, "Max tpcChi2NCl"};
58+
Configurable<float> tpcNClsFindable_cut{"tpcNClsFindable_cut", 70, "Min tpcNClsFindable"};
59+
Configurable<float> itsChi2_cut{"itsChi2_cut", 36, "Max itsChi2NCl"};
60+
Configurable<float> eta_cut{"eta_cut", 0.9, "Track Pseudorapidity"};
61+
Configurable<float> pt_cut{"pt_cut", 0.01, "Track Pt"};
62+
Configurable<float> TPC_cluster{"TPC_cluster", 50, "No.of TPC cluster"};
63+
64+
// Kinematic cuts
65+
Configurable<float> PID_cut{"PID_cut", 5, "TPC PID"};
66+
Configurable<float> Rap_cut{"Rap_cut", 0.9, "Track rapidity"};
67+
Configurable<float> Mass_Max{"Mass_Max", 10, "Invariant Mass range high"};
68+
Configurable<float> Mass_Min{"Mass_Min", 0, "Invariant Mass range low"};
69+
Configurable<float> Pt_coherent{"Pt_coherent", 0.15, "Coherent selection"};
70+
71+
// defining histograms using histogram registry
72+
HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject};
73+
74+
//-----------------------------------------------------------------------------------------------------------------------
75+
void init(o2::framework::InitContext&)
76+
{
77+
78+
registry.add("GapSide", "Gap Side; Entries", kTH1F, {{4, -1.5, 2.5}});
79+
registry.add("TrueGapSide", "Gap Side; Entries", kTH1F, {{4, -1.5, 2.5}});
80+
81+
// Fill counter to see effect of each selection criteria
82+
auto hSelectionCounter = registry.add<TH1>("hSelectionCounter", "hSelectionCounter;;NEvents", HistType::kTH1I, {{20, 0., 20.}});
83+
84+
TString SelectionCuts[16] = {"NoSelection", "gapside", "goodtracks", "truegap", "2collcontrib", "2goodtrk", "TPCPID", "Rap_cut", "unlikesign", "mass_cut", "coherent", "incoherent", "likesign", "mass_cut", "coherent", "incoherent"};
85+
86+
for (int i = 0; i < 16; i++) {
87+
hSelectionCounter->GetXaxis()->SetBinLabel(i + 1, SelectionCuts[i].Data());
88+
}
89+
// tracks
90+
registry.add("hTracks", "N_{tracks}", kTH1F, {{100, -0.5, 99.5}});
91+
registry.add("hTracksPions", "N_{tracks}", kTH1F, {{100, -0.5, 99.5}});
92+
registry.add("TwoPion/h2TracksPions", "N_{tracks}", kTH1F, {{100, -0.5, 99.5}});
93+
94+
registry.add("hdEdx", "p vs dE/dx Signal", kTH2F, {{100, 0.0, 3.0}, {100, 0.0, 200.0}});
95+
registry.add("hdEdxPion", "p_{#pi} vs dE/dx Signal", kTH2F, {{100, 0.0, 3.0}, {100, 0.0, 200.0}});
96+
97+
registry.add("TwoPion/hNsigPi1vsPi2", "NSigmaPi(t1) vs NSigmapi (t2);n#sigma_{1};n#sigma_{2}", kTH2F, {{100, -15., 15.}, {100, -15., 15}});
98+
registry.add("TwoPion/hNsigEl1vsEl2", "NSigmaEl(t1) vs NSigmaEl (t2);n#sigma_{1};n#sigma_{2}", kTH2F, {{100, -15., 15.}, {100, -15., 15}});
99+
registry.add("TwoPion/hNsigPivsPt1", "Pt vs NSigmaPi (t1);#it{p_{t}}, GeV/c;n#sigma_{#pi}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
100+
registry.add("TwoPion/hNsigPivsPt2", "Pt vs NSigmaPi (t2);#it{p_{t}}, GeV/c;n#sigma_{#pi}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
101+
registry.add("TwoPion/hNsigElvsPt1", "Pt vs NSigmaEl (t1);#it{p_{t}}, GeV/c;n#sigma_{#e}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
102+
registry.add("TwoPion/hNsigElvsPt2", "Pt vs NSigmaEl (t2);#it{p_{t}}, GeV/c;n#sigma_{#e}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
103+
registry.add("TwoPion/hNsigMuvsPt1", "Pt vs NSigmaMu (t1);#it{p_{t}}, GeV/c;n#sigma_{#pi}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
104+
registry.add("TwoPion/hNsigMuvsPt2", "Pt vs NSigmaMu (t2);#it{p_{t}}, GeV/c;n#sigma_{#pi}", kTH2F, {{100, 0., 2.5}, {100, -15., 15}});
105+
106+
registry.add("TwoPion/hPtsingle_track1", "Pt t1;#it{p_{t}}, GeV/c;", kTH1F, {{600, 0., 3.}});
107+
registry.add("TwoPion/hPtsingle_track2", "Pt t2;#it{p_{t}}, GeV/c;", kTH1F, {{600, 0., 3.}});
108+
registry.add("TwoPion/hEta_t1", "Eta of t1;#it{#eta};", kTH1F, {{100, -5., 5.}});
109+
registry.add("TwoPion/hEta_t2", "Eta of t2;#it{#eta};", kTH1F, {{100, -5., 5.}});
110+
registry.add("TwoPion/hP1", "P vs TPC signal;#it{P_{track}}, GeV/c; signal_{TPC} t1", kTH2F, {{100, 0., 2.}, {300, 0, 150}});
111+
registry.add("TwoPion/hTPCsig", "TPC signal;signal_{TPC} t2; signal_{TPC} t2", kTH2F, {{300, 0., 150.}, {300, 0, 150}});
112+
registry.add("TwoPion/hP2", "P vs TPC signal;#it{P_{track}}, GeV/c; signal_{TPC} t1", kTH2F, {{100, 0., 2.}, {300, 0, 150}});
113+
registry.add("TwoPion/hTPCsig1", "TPC signal;signal_{TPC} t2; signal_{TPC} t2", kTH2F, {{300, 0, 150.}, {300, 0, 150}});
114+
115+
registry.add("TwoPion/hMassLike", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
116+
registry.add("TwoPion/hMassUnlike", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
117+
registry.add("TwoPion/Coherent/hMassUnlikeCoherent", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
118+
registry.add("TwoPion/Coherent/hMassLikeCoherent", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
119+
registry.add("TwoPion/Incoherent/hMassUnlikeInCoherent", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
120+
registry.add("TwoPion/Incoherent/hMassLikeInCoherent", "m_{#pi#pi} [GeV/#it{c}^{2}]", kTH1F, {{20000, 0., 20.}});
121+
122+
registry.add("TwoPion/hPt", "Pt;#it{p_{t}}, GeV/c;", kTH1D, {{1000, 0., 10.}});
123+
registry.add("TwoPion/hPtLike", "Pt;#it{p_{t}}, GeV/c;", kTH1D, {{1000, 0., 10.}});
124+
registry.add("TwoPion/hEta", "Eta;#it{#eta};", kTH1F, {{500, -10., 10.}});
125+
registry.add("TwoPion/hRap", "Rapidity;#it{y};", kTH1F, {{500, -10., 10.}});
126+
registry.add("TwoPion/hPhiSystem", "Phi;#it{#Phi};", kTH1F, {{250, 0. * TMath::Pi(), 2. * TMath::Pi()}});
127+
registry.add("TwoPion/hMPt", "Inv.M vs Pt;M, GeV/c^{2};#it{P_{t}}, GeV/c;", kTH2F, {{100, 0., 10.}, {100, 0., 10.}});
128+
}
129+
130+
using udtracks = soa::Join<aod::UDTracks, aod::UDTracksExtra, aod::UDTracksPID>;
131+
using udtracksfull = soa::Join<aod::UDTracks, aod::UDTracksPID, aod::UDTracksExtra, aod::UDTracksFlags, aod::UDTracksDCA>;
132+
133+
using UDCollisionsFull = soa::Join<aod::UDCollisions, aod::SGCollisions, aod::UDCollisionsSels, aod::UDZdcsReduced>;
134+
135+
//__________________________________________________________________________
136+
// Main process
137+
void process(UDCollisionsFull::iterator const& collision, udtracksfull const& tracks)
138+
{
139+
// No selection criteria
140+
registry.fill(HIST("hSelectionCounter"), 0);
141+
142+
// Accessing gap sides
143+
int gapSide = collision.gapSide();
144+
if (gapSide < 0 || gapSide > 2)
145+
return;
146+
147+
registry.fill(HIST("hSelectionCounter"), 1);
148+
149+
// Accessing FIT information for further exclusivity and/or inclusivity
150+
float FIT_cut[5] = {FV0_cut, FT0A_cut, FT0C_cut, FDDA_cut, FDDC_cut};
151+
int truegapSide = sgSelector.trueGap(collision, FIT_cut[0], FIT_cut[1], FIT_cut[2], ZDC_cut);
152+
153+
// Initiating track parameters to select good tracks, values to be optimized in the configurables, parameters will be taken from SGtrackselector.h task included in the header
154+
std::vector<float> parameters = {PV_cut, dcaZ_cut, dcaXY_cut, tpcChi2_cut, tpcNClsFindable_cut, itsChi2_cut, eta_cut, pt_cut};
155+
156+
registry.fill(HIST("GapSide"), gapSide);
157+
registry.fill(HIST("TrueGapSide"), truegapSide);
158+
159+
// Gap side to be selected in the configurables
160+
gapSide = truegapSide;
161+
162+
if (gapSide == gap_Side) {
163+
164+
registry.fill(HIST("hSelectionCounter"), 2);
165+
166+
if(collision.flags()!=1) return; //UPC setting vs std setting
167+
//____________________________________________________________________________________
168+
169+
// Create LorentzVector to store all tracks, Pion tracks and TPC Pion PID
170+
std::vector<LorentzVector> allTracks;
171+
std::vector<LorentzVector> onlyPionTracks;
172+
std::vector<float> onlyPionSigma;
173+
std::vector<decltype(tracks.begin())> rawPionTracks;
174+
175+
// initialize pair 4-vector to zero before accumulation
176+
LorentzVector p(0.0, 0.0, 0.0, 0.0);
177+
178+
registry.fill(HIST("hTracks"), tracks.size());
179+
180+
for (auto t : tracks) {
181+
// Apply good track selection criteria
182+
if (!trackselector(t, parameters))
183+
continue;
184+
185+
double dEdx = t.tpcSignal();
186+
187+
registry.fill(HIST("hdEdx"), t.tpcInnerParam() / t.sign(), dEdx);
188+
189+
// Create Lorentz vector for this track (use constructor, portable)
190+
LorentzVector a(t.px(), t.py(), t.pz(), o2::constants::physics::MassPionCharged);
191+
allTracks.push_back(a);
192+
193+
// Apply TPC pion sigma
194+
auto nSigmaPi = t.tpcNSigmaPi();
195+
if (fabs(nSigmaPi) < PID_cut) {
196+
onlyPionTracks.push_back(a);
197+
onlyPionSigma.push_back(nSigmaPi);
198+
rawPionTracks.push_back(t);
199+
registry.fill(HIST("hdEdxPion"), t.tpcInnerParam() / t.sign(), dEdx);
200+
}
201+
}
202+
203+
registry.fill(HIST("hTracksPions"), onlyPionTracks.size());
204+
205+
//_____________________________________
206+
// Add all onlyPionTracks into p
207+
for (const auto &pion : onlyPionTracks) {
208+
p += pion;
209+
}
210+
211+
//_____________________________________
212+
// Selecting collisions with Two PV contributors
213+
if (collision.numContrib() == 2) {
214+
215+
registry.fill(HIST("hSelectionCounter"), 3);
216+
217+
// Selecting only Two good tracks
218+
if ((rawPionTracks.size() == 2) && (onlyPionTracks.size() == 2)) {
219+
220+
registry.fill(HIST("hSelectionCounter"), 4);
221+
222+
registry.fill(HIST("TwoPion/h2TracksPions"), onlyPionTracks.size());
223+
224+
registry.fill(HIST("TwoPion/hNsigPivsPt1"), onlyPionTracks[0].Pt(), rawPionTracks[0].tpcNSigmaPi());
225+
registry.fill(HIST("TwoPion/hNsigPivsPt2"), onlyPionTracks[1].Pt(), rawPionTracks[1].tpcNSigmaPi());
226+
registry.fill(HIST("TwoPion/hTPCsig"), rawPionTracks[0].tpcSignal(), rawPionTracks[1].tpcSignal());
227+
registry.fill(HIST("TwoPion/hNsigPi1vsPi2"), rawPionTracks[0].tpcNSigmaPi(), rawPionTracks[1].tpcNSigmaPi());
228+
229+
// Make sure two good tracks are within TPC pion sigma limit
230+
if ((onlyPionSigma[0] * onlyPionSigma[0] + onlyPionSigma[1] * onlyPionSigma[1]) > (PID_cut * PID_cut)) {
231+
return;
232+
}
233+
234+
registry.fill(HIST("hSelectionCounter"), 5);
235+
236+
// Rapidity of midrapidity acceptance
237+
if ((p.Rapidity() < -Rap_cut) || (p.Rapidity() > Rap_cut)) {
238+
return;
239+
}
240+
241+
registry.fill(HIST("hSelectionCounter"), 6);
242+
243+
// Two opposite sign tracks
244+
if (rawPionTracks[0].sign() != rawPionTracks[1].sign()) {
245+
246+
registry.fill(HIST("hSelectionCounter"), 7);
247+
registry.fill(HIST("TwoPion/hMassUnlike"), p.M());
248+
249+
// Flexible mass limits, can be selected in the configurable
250+
if ((p.M() > Mass_Min) && (p.M() < Mass_Max)) {
251+
252+
registry.fill(HIST("hSelectionCounter"), 8);
253+
254+
registry.fill(HIST("TwoPion/hPt"), p.Pt());
255+
registry.fill(HIST("TwoPion/hEta"), p.Eta());
256+
registry.fill(HIST("TwoPion/hRap"), p.Rapidity());
257+
registry.fill(HIST("TwoPion/hPhiSystem"), p.Phi());
258+
registry.fill(HIST("TwoPion/hMPt"), p.M(), p.Pt());
259+
260+
// flexible pt limit for selecting coherent Rho(0)
261+
if (p.Pt() < Pt_coherent) {
262+
263+
registry.fill(HIST("hSelectionCounter"), 9);
264+
265+
// Quality Control plots after coherent Rho(0) selection
266+
registry.fill(HIST("TwoPion/hEta_t1"), onlyPionTracks[0].Eta());
267+
registry.fill(HIST("TwoPion/hEta_t2"), onlyPionTracks[1].Eta());
268+
registry.fill(HIST("TwoPion/hPtsingle_track1"), onlyPionTracks[0].Pt());
269+
registry.fill(HIST("TwoPion/hPtsingle_track2"), onlyPionTracks[1].Pt());
270+
271+
registry.fill(HIST("TwoPion/hNsigMuvsPt1"), onlyPionTracks[0].Pt(), rawPionTracks[0].tpcNSigmaPi());
272+
registry.fill(HIST("TwoPion/hNsigMuvsPt2"), onlyPionTracks[1].Pt(), rawPionTracks[1].tpcNSigmaPi());
273+
registry.fill(HIST("TwoPion/hNsigElvsPt1"), onlyPionTracks[0].Pt(), rawPionTracks[0].tpcNSigmaEl());
274+
registry.fill(HIST("TwoPion/hNsigElvsPt2"), onlyPionTracks[1].Pt(), rawPionTracks[1].tpcNSigmaEl());
275+
registry.fill(HIST("TwoPion/hNsigEl1vsEl2"), rawPionTracks[0].tpcNSigmaPi(), rawPionTracks[1].tpcNSigmaPi());
276+
277+
registry.fill(HIST("TwoPion/hP1"), onlyPionTracks[0].P(), rawPionTracks[0].tpcSignal());
278+
registry.fill(HIST("TwoPion/hP2"), onlyPionTracks[1].P(), rawPionTracks[1].tpcSignal());
279+
registry.fill(HIST("TwoPion/hTPCsig1"), rawPionTracks[0].tpcSignal(), rawPionTracks[1].tpcSignal());
280+
281+
registry.fill(HIST("TwoPion/Coherent/hMassUnlikeCoherent"), p.M());
282+
}
283+
// Incoherent Rho(0) selection
284+
if (p.Pt() > Pt_coherent) {
285+
registry.fill(HIST("hSelectionCounter"), 10);
286+
registry.fill(HIST("TwoPion/Incoherent/hMassUnlikeInCoherent"), p.M());
287+
}
288+
}
289+
}
290+
291+
// Same charge particles
292+
if (rawPionTracks[0].sign() == rawPionTracks[1].sign()) {
293+
294+
registry.fill(HIST("hSelectionCounter"), 11);
295+
registry.fill(HIST("TwoPion/hMassLike"), p.M());
296+
297+
// Mass limit
298+
if ((p.M() > Mass_Min) && (p.M() < Mass_Max)) {
299+
300+
registry.fill(HIST("hSelectionCounter"), 12);
301+
registry.fill(HIST("TwoPion/hPtLike"), p.Pt());
302+
303+
// Coherent Rho(0) selection
304+
if (p.Pt() < Pt_coherent) {
305+
306+
registry.fill(HIST("hSelectionCounter"), 13);
307+
registry.fill(HIST("TwoPion/Coherent/hMassLikeCoherent"), p.M());
308+
}
309+
// Incoherent Rho(0) selection
310+
if (p.Pt() > Pt_coherent) {
311+
312+
registry.fill(HIST("hSelectionCounter"), 14);
313+
registry.fill(HIST("TwoPion/Incoherent/hMassLikeInCoherent"), p.M());
314+
}
315+
}
316+
}
317+
}
318+
}
319+
}
320+
}
321+
};
322+
323+
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
324+
{
325+
return WorkflowSpec{
326+
adaptAnalysisTask<upcExclusiveJpsiCB>(cfgc)};
327+
}

0 commit comments

Comments
 (0)