Skip to content

Commit ec60e16

Browse files
creetz16Barbara Chytla
authored andcommitted
DPG: Lambda daughters momentum resolution study (AliceO2Group#4596)
* Lambda daughter track momentum resolution study * Add daughter track covariance matrices to V0Datas
1 parent c1f133c commit ec60e16

File tree

5 files changed

+439
-2
lines changed

5 files changed

+439
-2
lines changed

DPG/Tasks/AOTTrack/V0Cascades/CMakeLists.txt

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -29,3 +29,8 @@ o2physics_add_dpl_workflow(qa-cascades
2929
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
3030
COMPONENT_NAME Analysis)
3131

32+
o2physics_add_dpl_workflow(qa-lam-momentum-resolution
33+
SOURCES qaLamMomResolution.cxx
34+
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
35+
COMPONENT_NAME Analysis)
36+
Lines changed: 284 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,284 @@
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+
/// \file qaLamMomResolution.cxx
14+
/// \author Carolina Reetz c.reetz@cern.ch
15+
/// \brief QA task to study momentum resolution of Lambda daughter tracks
16+
17+
#include "Framework/runDataProcessing.h"
18+
#include "Framework/AnalysisTask.h"
19+
#include "Framework/AnalysisDataModel.h"
20+
#include "Framework/ASoAHelpers.h"
21+
#include "ReconstructionDataFormats/Track.h"
22+
#include "ReconstructionDataFormats/PID.h"
23+
#include "Common/Core/trackUtilities.h"
24+
#include "PWGLF/DataModel/LFStrangenessTables.h"
25+
#include "DPG/Tasks/AOTTrack/V0Cascades/qaLamMomResolution.h"
26+
#include "Common/Core/TrackSelection.h"
27+
#include "Common/DataModel/TrackSelectionTables.h"
28+
#include "Common/DataModel/EventSelection.h"
29+
#include "Common/DataModel/PIDResponse.h"
30+
31+
using namespace o2;
32+
using namespace o2::framework;
33+
using namespace o2::framework::expressions;
34+
35+
using V0DatasLabeled = soa::Join<aod::V0Datas, aod::V0Covs, aod::McV0Labels>;
36+
using MyTracks = soa::Join<aod::TracksIU, aod::TracksCovIU, aod::TracksExtra, aod::TracksDCA, aod::TracksDCACov, aod::McTrackLabels>;
37+
38+
struct qaLamMomResolution {
39+
40+
Produces<aod::LamDaughters> lamdaughters;
41+
42+
Configurable<bool> collSelection{"collSelection", true, "collision selection"};
43+
44+
HistogramRegistry hist{"Histograms"};
45+
46+
int LambdaPDG = 3122;
47+
int AntiLambdaPDG = -3122;
48+
int ProtonPDG = 2212;
49+
int AntiProtonPDG = -2212;
50+
int PosPionPDG = 211;
51+
int NegPionPDG = -211;
52+
53+
float massLambda = -1.0f;
54+
float radiusLambda = -1.0f;
55+
float ptLambda = -1.0f;
56+
float etaProton = -1.0f, etaPion = -1.0f;
57+
int tpcNClsProton = 0, tpcNClsPion = 0;
58+
int chargeProton = 0, chargePion = 0;
59+
// daughter momenta
60+
std::array<float, 3> momProtonRecIU;
61+
std::array<float, 3> momPionRecIU;
62+
std::array<float, 3> momProtonRecIUErr;
63+
std::array<float, 3> momPionRecIUErr;
64+
std::array<float, 3> momProtonRec;
65+
std::array<float, 3> momPionRec;
66+
std::array<float, 3> momProtonRecErr;
67+
std::array<float, 3> momPionRecErr;
68+
float sigma1PtProtonIU = -1.0f, sigma1PtPionIU = -1.0f;
69+
// daughter DCA
70+
std::array<float, 2> DCAProtonRec; // 0: xy, 1: z
71+
std::array<float, 2> DCAPionRec; // 0: xy, 1: z
72+
std::array<float, 2> DCAProtonRecErr; // 0: xy, 1: z
73+
std::array<float, 2> DCAPionRecErr; // 0: xy, 1: z
74+
// MC info
75+
std::array<float, 3> momProtonGen;
76+
std::array<float, 3> momPionGen;
77+
std::array<float, 2> DCAProtonGen; // 0: xy, 1: z
78+
std::array<float, 2> DCAPionGen; // 0: xy, 1: z
79+
80+
void init(InitContext const&)
81+
{
82+
auto hEventSelectionFlow = hist.add<TH1>("hEventSelectionFlow", "Event selection flow", kTH1F, {{2, 0.5f, 2.5f}});
83+
hEventSelectionFlow->GetXaxis()->SetBinLabel(hEventSelectionFlow->FindBin(1), "Sel8");
84+
hEventSelectionFlow->GetXaxis()->SetBinLabel(hEventSelectionFlow->FindBin(2), "|Vtx_{z}|<10cm");
85+
}
86+
87+
template <typename TCollision>
88+
void fillTable(TCollision const& collision)
89+
{
90+
lamdaughters(collision.globalIndex(),
91+
massLambda, radiusLambda, ptLambda,
92+
chargeProton, chargePion,
93+
etaProton, etaPion,
94+
tpcNClsProton, tpcNClsPion,
95+
momProtonRec[0], momProtonRec[1], momProtonRec[2],
96+
momProtonRecErr[0], momProtonRecErr[1], momProtonRecErr[2],
97+
momPionRec[0], momPionRec[1], momPionRec[2],
98+
momPionRecErr[0], momPionRecErr[1], momPionRecErr[2],
99+
momProtonRecIU[0], momProtonRecIU[1], momProtonRecIU[2],
100+
momProtonRecIUErr[0], momProtonRecIUErr[1], momProtonRecIUErr[2],
101+
momPionRecIU[0], momPionRecIU[1], momPionRecIU[2],
102+
momPionRecIUErr[0], momPionRecIUErr[1], momPionRecIUErr[2],
103+
momProtonGen[0], momProtonGen[1], momProtonGen[2],
104+
momPionGen[0], momPionGen[1], momPionGen[2],
105+
sigma1PtProtonIU, sigma1PtPionIU,
106+
DCAProtonRec[0], DCAProtonRec[1],
107+
DCAProtonRecErr[0], DCAProtonRecErr[1],
108+
DCAPionRec[0], DCAPionRec[1],
109+
DCAPionRecErr[0], DCAPionRecErr[1]);
110+
}
111+
112+
template <typename TProton, typename TPion>
113+
void getTrackInfo(TProton const& protonTrackIU, TPion const& pionTrackIU)
114+
{
115+
// daughter momenta at IU
116+
o2::track::TrackParCov protonTrackParCov, pionTrackParCov;
117+
std::array<float, 21> protoncv, pioncv;
118+
protonTrackParCov = getTrackParCov(protonTrackIU);
119+
pionTrackParCov = getTrackParCov(pionTrackIU);
120+
protonTrackParCov.getCovXYZPxPyPzGlo(protoncv);
121+
pionTrackParCov.getCovXYZPxPyPzGlo(pioncv);
122+
// proton
123+
momProtonRecIU[0] = protonTrackIU.px();
124+
momProtonRecIU[1] = protonTrackIU.py();
125+
momProtonRecIU[2] = protonTrackIU.pz();
126+
momProtonRecIUErr[0] = sqrt(protoncv[9]);
127+
momProtonRecIUErr[1] = sqrt(protoncv[14]);
128+
momProtonRecIUErr[2] = sqrt(protoncv[20]);
129+
sigma1PtProtonIU = protonTrackIU.sigma1Pt();
130+
// pion
131+
momPionRecIU[0] = pionTrackIU.px();
132+
momPionRecIU[1] = pionTrackIU.py();
133+
momPionRecIU[2] = pionTrackIU.pz();
134+
momPionRecIUErr[0] = sqrt(pioncv[9]);
135+
momPionRecIUErr[1] = sqrt(pioncv[14]);
136+
momPionRecIUErr[2] = sqrt(pioncv[20]);
137+
sigma1PtPionIU = pionTrackIU.sigma1Pt();
138+
139+
// daughter DCA
140+
DCAProtonRec[0] = protonTrackIU.dcaXY();
141+
DCAProtonRec[1] = protonTrackIU.dcaZ();
142+
DCAProtonRecErr[0] = sqrt(protonTrackIU.sigmaDcaXY2());
143+
DCAProtonRecErr[1] = sqrt(protonTrackIU.sigmaDcaZ2());
144+
DCAPionRec[0] = pionTrackIU.dcaXY();
145+
DCAPionRec[1] = pionTrackIU.dcaZ();
146+
DCAPionRecErr[0] = sqrt(pionTrackIU.sigmaDcaXY2());
147+
DCAPionRecErr[1] = sqrt(pionTrackIU.sigmaDcaZ2());
148+
149+
// daughter charges, eta, nTPCclusters
150+
chargeProton = protonTrackIU.sign();
151+
chargePion = pionTrackIU.sign();
152+
etaProton = protonTrackIU.eta();
153+
etaPion = pionTrackIU.eta();
154+
tpcNClsProton = protonTrackIU.tpcNClsFound();
155+
tpcNClsPion = pionTrackIU.tpcNClsFound();
156+
}
157+
158+
void processMC(soa::Join<aod::Collisions, aod::EvSels>::iterator const& collision,
159+
V0DatasLabeled const& V0Datas,
160+
aod::McParticles const& mcparticles,
161+
MyTracks const&)
162+
{
163+
164+
// event selection
165+
if (collSelection && !collision.sel8()) {
166+
return;
167+
}
168+
hist.fill(HIST("hEventSelectionFlow"), 1.f);
169+
if (collSelection && (abs(collision.posZ()) >= 10.)) {
170+
return;
171+
}
172+
hist.fill(HIST("hEventSelectionFlow"), 2.f);
173+
174+
for (auto& v0data : V0Datas) {
175+
176+
if (v0data.has_mcParticle() && v0data.mcParticleId() > -1 && v0data.mcParticleId() <= mcparticles.size()) {
177+
auto MCv0 = v0data.mcParticle_as<aod::McParticles>();
178+
179+
// Lambda
180+
if (MCv0.has_daughters() && MCv0.pdgCode() == LambdaPDG) {
181+
LOG(debug) << "V0 is a Lambda.";
182+
const auto& protonTrackIU = v0data.posTrack_as<MyTracks>();
183+
const auto& pionTrackIU = v0data.negTrack_as<MyTracks>();
184+
185+
if (protonTrackIU.has_mcParticle() && pionTrackIU.has_mcParticle()) {
186+
187+
const auto& MCproton = protonTrackIU.mcParticle_as<aod::McParticles>();
188+
const auto& MCpion = pionTrackIU.mcParticle_as<aod::McParticles>();
189+
190+
if (MCproton.pdgCode() == ProtonPDG && MCpion.pdgCode() == NegPionPDG) {
191+
192+
// lambda
193+
massLambda = v0data.mLambda();
194+
radiusLambda = v0data.v0radius();
195+
ptLambda = v0data.pt();
196+
/// daughter momenta at Lambda vertex
197+
// proton
198+
momProtonRec[0] = v0data.pxpos();
199+
momProtonRec[1] = v0data.pypos();
200+
momProtonRec[2] = v0data.pzpos();
201+
momProtonRecErr[0] = sqrt(v0data.covmatposdau()[9]);
202+
momProtonRecErr[1] = sqrt(v0data.covmatposdau()[14]);
203+
momProtonRecErr[2] = sqrt(v0data.covmatposdau()[20]);
204+
momProtonGen[0] = MCproton.px();
205+
momProtonGen[1] = MCproton.py();
206+
momProtonGen[2] = MCproton.pz();
207+
// pion
208+
momPionRec[0] = v0data.pxneg();
209+
momPionRec[1] = v0data.pyneg();
210+
momPionRec[2] = v0data.pzneg();
211+
momPionRecErr[0] = sqrt(v0data.covmatnegdau()[9]);
212+
momPionRecErr[1] = sqrt(v0data.covmatnegdau()[14]);
213+
momPionRecErr[2] = sqrt(v0data.covmatnegdau()[20]);
214+
momPionGen[0] = MCpion.px();
215+
momPionGen[1] = MCpion.py();
216+
momPionGen[2] = MCpion.pz();
217+
218+
// get daughter momenta at IU, charge, eta, nTPCclusters
219+
getTrackInfo(protonTrackIU, pionTrackIU);
220+
221+
// fill table
222+
fillTable(collision);
223+
}
224+
}
225+
} // end Lambda
226+
227+
// Anti-Lambda
228+
if (MCv0.pdgCode() == AntiLambdaPDG) {
229+
LOG(debug) << "V0 is an Anti-Lambda.";
230+
const auto& protonTrackIU = v0data.negTrack_as<MyTracks>();
231+
const auto& pionTrackIU = v0data.posTrack_as<MyTracks>();
232+
233+
if (protonTrackIU.has_mcParticle() && pionTrackIU.has_mcParticle()) {
234+
235+
const auto& MCproton = protonTrackIU.mcParticle_as<aod::McParticles>();
236+
const auto& MCpion = pionTrackIU.mcParticle_as<aod::McParticles>();
237+
238+
if (MCproton.pdgCode() == AntiProtonPDG && MCpion.pdgCode() == PosPionPDG) {
239+
240+
// lambda mass and radius
241+
massLambda = v0data.mAntiLambda();
242+
radiusLambda = v0data.v0radius();
243+
ptLambda = v0data.pt();
244+
/// daughter momenta at Lambda vertex
245+
// proton
246+
momProtonRec[0] = v0data.pxneg();
247+
momProtonRec[1] = v0data.pyneg();
248+
momProtonRec[2] = v0data.pzneg();
249+
momProtonRecErr[0] = sqrt(v0data.covmatnegdau()[9]);
250+
momProtonRecErr[1] = sqrt(v0data.covmatnegdau()[14]);
251+
momProtonRecErr[2] = sqrt(v0data.covmatnegdau()[20]);
252+
momProtonGen[0] = MCproton.px();
253+
momProtonGen[1] = MCproton.py();
254+
momProtonGen[2] = MCproton.pz();
255+
// pion
256+
momPionRec[0] = v0data.pxpos();
257+
momPionRec[1] = v0data.pypos();
258+
momPionRec[2] = v0data.pzpos();
259+
momPionRecErr[0] = sqrt(v0data.covmatposdau()[9]);
260+
momPionRecErr[1] = sqrt(v0data.covmatposdau()[14]);
261+
momPionRecErr[2] = sqrt(v0data.covmatposdau()[20]);
262+
momPionGen[0] = MCpion.px();
263+
momPionGen[1] = MCpion.py();
264+
momPionGen[2] = MCpion.pz();
265+
266+
// get daughter momenta at IU, charge, eta, nTPCclusters
267+
getTrackInfo(protonTrackIU, pionTrackIU);
268+
269+
// fill table
270+
fillTable(collision);
271+
}
272+
}
273+
} // end Anti-Lambda
274+
} // end MC
275+
} // end V0 loop
276+
}
277+
PROCESS_SWITCH(qaLamMomResolution, processMC, "Process MC", true);
278+
};
279+
280+
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
281+
{
282+
return WorkflowSpec{
283+
adaptAnalysisTask<qaLamMomResolution>(cfgc)};
284+
}

0 commit comments

Comments
 (0)