1111
1212// /
1313// / \file pidTPC.cxx
14- // / \author Nicolo' Jacazio
14+ // / \author Nicolò Jacazio
1515// / \brief Task to produce PID tables for TPC split for each particle with only the Nsigma information.
1616// / Only the tables for the mass hypotheses requested are filled, the others are sent empty.
1717// /
2525#include " Common/Core/PID/PIDResponse.h"
2626#include " Common/Core/PID/PIDTPC.h"
2727#include " Common/DataModel/TrackSelectionTables.h"
28+ #include " Common/DataModel/EventSelection.h"
29+ #include " TableHelper.h"
30+ #include " Framework/StaticFor.h"
2831
2932using namespace o2 ;
3033using namespace o2 ::framework;
@@ -40,6 +43,7 @@ void customize(std::vector<o2::framework::ConfigParamSpec>& workflowOptions)
4043
4144#include " Framework/runDataProcessing.h"
4245
46+ // / Task to produce the TPC response table
4347struct tpcPid {
4448 using Trks = soa::Join<aod::Tracks, aod::TracksExtra>;
4549 using Coll = aod::Collisions;
@@ -76,34 +80,31 @@ struct tpcPid {
7680 void init (o2::framework::InitContext& initContext)
7781 {
7882 // Checking the tables are requested in the workflow and enabling them
79- auto & workflows = initContext.services ().get <RunningWorkflowInfo const >();
80- for (DeviceSpec device : workflows.devices ) {
81- for (auto input : device.inputs ) {
82- auto enableFlag = [&input](const std::string particle, Configurable<int >& flag) {
83- const std::string table = " pidTPC" + particle;
84- if (input.matcher .binding == table) {
85- if (flag < 0 ) {
86- flag.value = 1 ;
87- LOG (info) << " Auto-enabling table: " + table;
88- } else if (flag > 0 ) {
89- flag.value = 1 ;
90- LOG (info) << " Table enabled: " + table;
91- } else {
92- LOG (info) << " Table disabled: " + table;
93- }
94- }
95- };
96- enableFlag (" El" , pidEl);
97- enableFlag (" Mu" , pidMu);
98- enableFlag (" Pi" , pidPi);
99- enableFlag (" Ka" , pidKa);
100- enableFlag (" Pr" , pidPr);
101- enableFlag (" De" , pidDe);
102- enableFlag (" Tr" , pidTr);
103- enableFlag (" He" , pidHe);
104- enableFlag (" Al" , pidAl);
83+ auto enableFlag = [&](const std::string particle, Configurable<int >& flag) {
84+ const std::string table = " pidTPC" + particle;
85+ if (isTableRequiredInWorkflow (initContext, table)) {
86+ if (flag < 0 ) {
87+ flag.value = 1 ;
88+ LOG (info) << " Auto-enabling table: " + table;
89+ } else if (flag > 0 ) {
90+ flag.value = 1 ;
91+ LOG (info) << " Table enabled: " + table;
92+ } else {
93+ LOG (info) << " Table disabled: " + table;
94+ }
10595 }
106- }
96+ };
97+
98+ enableFlag (" El" , pidEl);
99+ enableFlag (" Mu" , pidMu);
100+ enableFlag (" Pi" , pidPi);
101+ enableFlag (" Ka" , pidKa);
102+ enableFlag (" Pr" , pidPr);
103+ enableFlag (" De" , pidDe);
104+ enableFlag (" Tr" , pidTr);
105+ enableFlag (" He" , pidHe);
106+ enableFlag (" Al" , pidAl);
107+
107108 // Getting the parametrization parameters
108109 ccdb->setURL (url.value );
109110 ccdb->setTimestamp (timestamp.value );
@@ -172,12 +173,23 @@ struct tpcPid {
172173 }
173174};
174175
176+ // / Task to produce the TPC QA plots
175177struct tpcPidQa {
176178 static constexpr int Np = 9 ;
177179 static constexpr const char * pT[Np] = {" e" , " #mu" , " #pi" , " K" , " p" , " d" , " t" , " ^{3}He" , " #alpha" };
178180 static constexpr std::string_view hnsigma[Np] = {" nsigma/El" , " nsigma/Mu" , " nsigma/Pi" ,
179181 " nsigma/Ka" , " nsigma/Pr" , " nsigma/De" ,
180182 " nsigma/Tr" , " nsigma/He" , " nsigma/Al" };
183+ static constexpr std::string_view hnsigmapt[Np] = {" nsigmapt/El" , " nsigmapt/Mu" , " nsigmapt/Pi" ,
184+ " nsigmapt/Ka" , " nsigmapt/Pr" , " nsigmapt/De" ,
185+ " nsigmapt/Tr" , " nsigmapt/He" , " nsigmapt/Al" };
186+ static constexpr std::string_view hnsigmapospt[Np] = {" nsigmapospt/El" , " nsigmapospt/Mu" , " nsigmapospt/Pi" ,
187+ " nsigmapospt/Ka" , " nsigmapospt/Pr" , " nsigmapospt/De" ,
188+ " nsigmapospt/Tr" , " nsigmapospt/He" , " nsigmapospt/Al" };
189+ static constexpr std::string_view hnsigmanegpt[Np] = {" nsigmanegpt/El" , " nsigmanegpt/Mu" , " nsigmanegpt/Pi" ,
190+ " nsigmanegpt/Ka" , " nsigmanegpt/Pr" , " nsigmanegpt/De" ,
191+ " nsigmanegpt/Tr" , " nsigmanegpt/He" , " nsigmanegpt/Al" };
192+
181193 HistogramRegistry histos{" Histos" , {}, OutputObjHandlingPolicy::QAObject};
182194
183195 Configurable<int > logAxis{" logAxis" , 0 , " Flag to use a log momentum axis" };
@@ -187,73 +199,125 @@ struct tpcPidQa {
187199 Configurable<int > nBinsNSigma{" nBinsNSigma" , 200 , " Number of bins for the NSigma" };
188200 Configurable<float > minNSigma{" minNSigma" , -10 .f , " Minimum NSigma in range" };
189201 Configurable<float > maxNSigma{" maxNSigma" , 10 .f , " Maximum NSigma in range" };
202+ Configurable<int > applyEvSel{" applyEvSel" , 0 , " Flag to apply rapidity cut: 0 -> no event selection, 1 -> Run 2 event selection, 2 -> Run 3 event selection" };
203+ Configurable<bool > applyTrackCut{" applyTrackCut" , false , " Flag to apply standard track cuts" };
204+ Configurable<bool > applyRapidityCut{" applyRapidityCut" , false , " Flag to apply rapidity cut" };
190205
191206 template <uint8_t i>
192- void addParticleHistos ()
207+ void addParticleHistos (const AxisSpec& pAxis, const AxisSpec& ptAxis )
193208 {
194- AxisSpec pAxis{nBinsP, minP, maxP, " #it{p} (GeV/#it{c})" };
195- if (logAxis) {
196- pAxis.makeLogaritmic ();
197- }
198- const AxisSpec nSigmaAxis{nBinsNSigma, minNSigma, maxNSigma, Form (" N_{#sigma}^{TPC}(%s)" , pT[i])};
199-
200209 // NSigma
210+ const AxisSpec nSigmaAxis{nBinsNSigma, minNSigma, maxNSigma, Form (" N_{#sigma}^{TPC}(%s)" , pT[i])};
201211 histos.add (hnsigma[i].data (), Form (" N_{#sigma}^{TPC}(%s)" , pT[i]), kTH2F , {pAxis, nSigmaAxis});
212+ histos.add (hnsigmapt[i].data (), Form (" N_{#sigma}^{TPC}(%s)" , pT[i]), kTH2F , {ptAxis, nSigmaAxis});
213+ histos.add (hnsigmapospt[i].data (), Form (" N_{#sigma}^{TPC}(%s)" , pT[i]), kTH2F , {ptAxis, nSigmaAxis});
214+ histos.add (hnsigmanegpt[i].data (), Form (" N_{#sigma}^{TPC}(%s)" , pT[i]), kTH2F , {ptAxis, nSigmaAxis});
202215 }
203216
204217 void init (o2::framework::InitContext&)
205218 {
206-
219+ const AxisSpec multAxis{1000 , 0 .f , 1000 .f , " Track multiplicity" };
220+ const AxisSpec vtxZAxis{100 , -20 , 20 , " Vtx_{z} (cm)" };
221+ const AxisSpec pAxisPosNeg{nBinsP, -maxP, maxP, " Signed #it{p} (GeV/#it{c})" };
207222 AxisSpec pAxis{nBinsP, minP, maxP, " #it{p} (GeV/#it{c})" };
223+ AxisSpec ptAxis{nBinsP, minP, maxP, " #it{p}_{T} (GeV/#it{c})" };
208224 if (logAxis) {
225+ ptAxis.makeLogaritmic ();
209226 pAxis.makeLogaritmic ();
210227 }
211- const AxisSpec vtxZAxis{100 , -20 , 20 , " Vtx_{z} (cm)" };
212228 const AxisSpec dedxAxis{1000 , 0 , 1000 , " d#it{E}/d#it{x} A.U." };
213229
214230 // Event properties
231+ auto h = histos.add <TH1 >(" event/evsel" , " " , kTH1F , {{10 , 0.5 , 10.5 , " Ev. Sel." }});
232+ h->GetXaxis ()->SetBinLabel (1 , " Events read" );
233+ h->GetXaxis ()->SetBinLabel (2 , " Passed ev. sel." );
234+ h->GetXaxis ()->SetBinLabel (3 , " Passed mult." );
235+ h->GetXaxis ()->SetBinLabel (4 , " Passed vtx Z" );
236+
215237 histos.add (" event/vertexz" , " " , kTH1F , {vtxZAxis});
238+ histos.add (" event/multiplicity" , " " , kTH1F , {multAxis});
216239 histos.add (" event/tpcsignal" , " " , kTH2F , {pAxis, dedxAxis});
240+ histos.add (" event/signedtpcsignal" , " " , kTH2F , {pAxisPosNeg, dedxAxis});
217241
218- addParticleHistos<0 >();
219- addParticleHistos<1 >();
220- addParticleHistos<2 >();
221- addParticleHistos<3 >();
222- addParticleHistos<4 >();
223- addParticleHistos<5 >();
224- addParticleHistos<6 >();
225- addParticleHistos<7 >();
226- addParticleHistos<8 >();
242+ static_for<0 , 8 >([&](auto i) {
243+ addParticleHistos<i>(pAxis, ptAxis);
244+ });
227245 }
228246
229- template <uint8_t i , typename T>
230- void fillParticleHistos (const T& t, const float nsigma)
247+ template <o2::track:: PID :: ID id , typename T>
248+ void fillParticleHistos (const T& t, const float & nsigma)
231249 {
232- histos.fill (HIST (hnsigma[i]), t.p (), nsigma);
250+ if (applyRapidityCut) {
251+ const float y = TMath::ASinH (t.pt () / TMath::Sqrt (PID::getMass2 (id) + t.pt () * t.pt ()) * TMath::SinH (t.eta ()));
252+ if (abs (y) > 0.5 ) {
253+ return ;
254+ }
255+ }
256+ // Fill histograms
257+ histos.fill (HIST (hnsigma[id]), t.p (), nsigma);
258+ histos.fill (HIST (hnsigmapt[id]), t.pt (), nsigma);
259+ if (t.sign () > 0 ) {
260+ histos.fill (HIST (hnsigmapospt[id]), t.pt (), nsigma);
261+ } else {
262+ histos.fill (HIST (hnsigmanegpt[id]), t.pt (), nsigma);
263+ }
233264 }
234265
235- void process (aod::Collision const & collision, soa::Join<aod::Tracks, aod::TracksExtra,
236- aod::pidTPCEl, aod::pidTPCMu, aod::pidTPCPi,
237- aod::pidTPCKa, aod::pidTPCPr, aod::pidTPCDe,
238- aod::pidTPCTr, aod::pidTPCHe, aod::pidTPCAl,
239- aod::TrackSelection> const & tracks)
266+ void process (soa::Join<aod::Collisions, aod::EvSels>::iterator const & collision,
267+ soa::Join<aod::Tracks, aod::TracksExtra,
268+ aod::pidTPCEl, aod::pidTPCMu, aod::pidTPCPi,
269+ aod::pidTPCKa, aod::pidTPCPr, aod::pidTPCDe,
270+ aod::pidTPCTr, aod::pidTPCHe, aod::pidTPCAl,
271+ aod::TrackSelection> const & tracks)
240272 {
273+ histos.fill (HIST (" event/evsel" ), 1 );
274+ if (applyEvSel == 1 ) {
275+ if (!collision.sel7 ()) {
276+ return ;
277+ }
278+ } else if (applyEvSel == 2 ) {
279+ if (!collision.sel8 ()) {
280+ return ;
281+ }
282+ }
283+ histos.fill (HIST (" event/evsel" ), 2 );
284+
285+ float ntracks = 0 ;
286+ for (auto t : tracks) {
287+ if (applyTrackCut && !t.isGlobalTrack ()) {
288+ continue ;
289+ }
290+ ntracks += 1 ;
291+ }
292+ // if (0 && ntracks < 1) {
293+ // return;
294+ // }
295+ histos.fill (HIST (" event/evsel" ), 3 );
296+ if (abs (collision.posZ ()) > 10 .f ) {
297+ return ;
298+ }
299+ histos.fill (HIST (" event/evsel" ), 4 );
241300 histos.fill (HIST (" event/vertexz" ), collision.posZ ());
301+ histos.fill (HIST (" event/multiplicity" ), ntracks);
242302
243303 for (auto t : tracks) {
304+ if (applyTrackCut && !t.isGlobalTrack ()) {
305+ continue ;
306+ }
244307 // const float mom = t.p();
245308 const float mom = t.tpcInnerParam ();
246309 histos.fill (HIST (" event/tpcsignal" ), mom, t.tpcSignal ());
310+ histos.fill (HIST (" event/signedtpcsignal" ), mom * t.sign (), t.tpcSignal ());
247311 //
248- fillParticleHistos<0 >(t, t.tpcNSigmaEl ());
249- fillParticleHistos<1 >(t, t.tpcNSigmaMu ());
250- fillParticleHistos<2 >(t, t.tpcNSigmaPi ());
251- fillParticleHistos<3 >(t, t.tpcNSigmaKa ());
252- fillParticleHistos<4 >(t, t.tpcNSigmaPr ());
253- fillParticleHistos<5 >(t, t.tpcNSigmaDe ());
254- fillParticleHistos<6 >(t, t.tpcNSigmaTr ());
255- fillParticleHistos<7 >(t, t.tpcNSigmaHe ());
256- fillParticleHistos<8 >(t, t.tpcNSigmaAl ());
312+ fillParticleHistos<PID ::Electron >(t, t.tpcNSigmaEl ());
313+ fillParticleHistos<PID ::Muon >(t, t.tpcNSigmaMu ());
314+ fillParticleHistos<PID ::Pion >(t, t.tpcNSigmaPi ());
315+ fillParticleHistos<PID ::Kaon >(t, t.tpcNSigmaKa ());
316+ fillParticleHistos<PID ::Proton >(t, t.tpcNSigmaPr ());
317+ fillParticleHistos<PID ::Deuteron >(t, t.tpcNSigmaDe ());
318+ fillParticleHistos<PID ::Triton >(t, t.tpcNSigmaTr ());
319+ fillParticleHistos<PID ::Helium3 >(t, t.tpcNSigmaHe ());
320+ fillParticleHistos<PID ::Alpha >(t, t.tpcNSigmaAl ());
257321 }
258322 }
259323};
0 commit comments