Project
Loading...
Searching...
No Matches
StrangenessTracker.cxx
Go to the documentation of this file.
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.
13
14#include <numeric>
16#include "ITStracking/IOUtils.h"
18
19#ifdef ENABLE_UPGRADES
21#endif
22
23namespace o2
24{
25namespace strangeness_tracking
26{
27
29{
30 clear();
31
32 mInputV0tracks = recoData.getV0s();
33 mInputV0Indices = recoData.getV0sIdx();
38 if (mInputV0Indices.size() != mInputV0tracks.size() || mInputCascadeIndices.size() != mInputCascadeTracks.size() || mInput3BodyIndices.size() != mInput3BodyTracks.size()) {
39 LOGP(fatal, "Mismatch between input SVertices indices and kinematics (not requested?): V0: {}/{} Cascades: {}/{} Decay3Bodys: {}/{}",
41 }
42 mInputITStracks = recoData.getITSTracks();
44
45 auto compClus = recoData.getITSClusters();
46 auto clusPatt = recoData.getITSClustersPatterns();
47 auto pattIt = clusPatt.begin();
48 auto pattIt2 = clusPatt.begin();
49 mInputITSclusters.reserve(compClus.size());
50 mInputClusterSizes.resize(compClus.size());
51#ifdef ENABLE_UPGRADES
52 if (o2::GlobalParams::Instance().withITS3) {
54 getClusterSizesIT3(mInputClusterSizes, compClus, pattIt2, mIT3Dict);
55 } else {
58 }
59#else
62#endif
63
64 mITSvtxBrackets.resize(mInputITStracks.size());
65 for (int i = 0; i < mInputITStracks.size(); i++) {
66 mITSvtxBrackets[i] = {-1, -1};
67 }
68
69 // build time bracket for each ITS track
70 auto trackIndex = recoData.getPrimaryVertexMatchedTracks(); // Global ID's for associated tracks
71 auto vtxRefs = recoData.getPrimaryVertexMatchedTrackRefs(); // references from vertex to these track IDs
72
74 int nv = vtxRefs.size();
75 for (int iv = 0; iv < nv; iv++) {
76 const auto& vtref = vtxRefs[iv];
77 int it = vtref.getFirstEntry(), itLim = it + vtref.getEntries();
78 for (; it < itLim; it++) {
79 auto tvid = trackIndex[it];
80 if (!recoData.isTrackSourceLoaded(tvid.getSource()) || tvid.getSource() != GIndex::ITS) {
81 continue;
82 }
83 if (mITSvtxBrackets[tvid.getIndex()].getMin() == -1) {
84 mITSvtxBrackets[tvid.getIndex()].setMin(iv);
85 mITSvtxBrackets[tvid.getIndex()].setMax(iv);
86 } else {
87 mITSvtxBrackets[tvid.getIndex()].setMax(iv);
88 }
89 }
90 }
91 }
92
93 if (mMCTruthON) {
94 mITSClsLabels = recoData.mcITSClusters.get();
96 }
97
98 LOG(debug) << "V0 tracks size: " << mInputV0tracks.size();
99 LOG(debug) << "Cascade tracks size: " << mInputCascadeTracks.size();
100 LOG(debug) << "Decay3Body tracks size: " << mInput3BodyTracks.size();
101 LOG(debug) << "ITS tracks size: " << mInputITStracks.size();
102 LOG(debug) << "ITS idxs size: " << mInputITSidxs.size();
103 LOG(debug) << "ITS clusters size: " << mInputITSclusters.size();
104 LOG(debug) << "VtxRefs size: " << vtxRefs.size();
105
106 return true;
107}
108
109void StrangenessTracker::prepareITStracks() // sort tracks by eta and phi and select only tracks with vertex matching
110{
111
112 for (int iTrack{0}; iTrack < mInputITStracks.size(); iTrack++) {
113 if (mStrParams->mVertexMatching && mITSvtxBrackets[iTrack].getMin() == -1) {
114 continue;
115 }
116 mSortedITStracks.push_back(mInputITStracks[iTrack]);
117 mSortedITSindexes.push_back(iTrack);
118 }
119
121 std::sort(mSortedITStracks.begin(), mSortedITStracks.end(), [&](o2::its::TrackITS& a, o2::its::TrackITS& b) { return mUtils.getBinIndex(a.getEta(), a.getPhi()) < mUtils.getBinIndex(b.getEta(), b.getPhi()); });
122 std::sort(mSortedITSindexes.begin(), mSortedITSindexes.end(), [&](int i, int j) { return mUtils.getBinIndex(mInputITStracks[i].getEta(), mInputITStracks[i].getPhi()) < mUtils.getBinIndex(mInputITStracks[j].getEta(), mInputITStracks[j].getPhi()); });
123
124 for (auto& track : mSortedITStracks) {
125 mTracksIdxTable[mUtils.getBinIndex(track.getEta(), track.getPhi())]++;
126 }
127 std::exclusive_scan(mTracksIdxTable.begin(), mTracksIdxTable.begin() + mUtils.mPhiBins * mUtils.mEtaBins, mTracksIdxTable.begin(), 0);
129}
130
131void StrangenessTracker::processV0(int iv0, const V0& v0, const V0Index& v0Idx, int iThread)
132{
133 if (mStrParams->mSkipTPC && ((v0Idx.getProngID(kV0DauNeg).getSource() == GIndex::TPC) || (v0Idx.getProngID(kV0DauPos).getSource() == GIndex::TPC))) {
134 return;
135 }
136 ClusAttachments structClus;
137 auto& daughterTracks = mDaughterTracks[iThread];
138 daughterTracks.resize(2); // resize to 2 prongs: first positive second negative
139 auto posTrack = v0.getProng(kV0DauPos);
140 auto negTrack = v0.getProng(kV0DauNeg);
141 auto alphaV0 = calcV0alpha(v0);
142 if (alphaV0 > 0) {
143 if (posTrack.getPID() != PID::Alpha) {
144 posTrack.setPID(PID::Helium3, true);
145 }
146 } else {
147 if (negTrack.getPID() != PID::Alpha) {
148 negTrack.setPID(PID::Helium3, true);
149 }
150 }
151 V0 correctedV0; // recompute V0 for Hypertriton
152 if (!recreateV0(posTrack, negTrack, correctedV0, iThread)) {
153 return;
154 }
155 StrangeTrack strangeTrack;
156 strangeTrack.mPartType = dataformats::kStrkV0;
157 auto v0R = std::sqrt(v0.calcR2());
158 auto iBinsV0 = mUtils.getBinRect(correctedV0.getEta(), correctedV0.getPhi(), mStrParams->mEtaBinSize, mStrParams->mPhiBinSize);
159 for (int& iBinV0 : iBinsV0) {
160 for (int iTrack{mTracksIdxTable[iBinV0]}; iTrack < TMath::Min(mTracksIdxTable[iBinV0 + 1], int(mSortedITStracks.size())); iTrack++) {
161 strangeTrack.mMother = (o2::track::TrackParCovF)correctedV0;
162 daughterTracks[kV0DauPos] = correctedV0.getProng(kV0DauPos);
163 daughterTracks[kV0DauNeg] = correctedV0.getProng(kV0DauNeg);
164 const auto& itsTrack = mSortedITStracks[iTrack];
165 const auto& ITSindexRef = mSortedITSindexes[iTrack];
166 if (mStrParams->mVertexMatching && (mITSvtxBrackets[ITSindexRef].getMin() > v0Idx.getVertexID() || mITSvtxBrackets[ITSindexRef].getMax() < v0Idx.getVertexID())) {
167 continue;
168 }
169 if (matchDecayToITStrack(v0R, strangeTrack, structClus, itsTrack, daughterTracks, iThread)) {
170 auto propInstance = o2::base::Propagator::Instance();
171 o2::track::TrackParCov decayVtxTrackClone = strangeTrack.mMother; // clone track and propagate to decay vertex
172 if (!propInstance->propagateToX(decayVtxTrackClone, strangeTrack.mDecayVtx[0], getBz(), o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, mCorrType)) {
173 LOG(debug) << "Mother propagation to decay vertex failed";
174 continue;
175 }
176 decayVtxTrackClone.getPxPyPzGlo(strangeTrack.mDecayMom);
177 std::array<float, 3> momPos, momNeg;
178 mFitter3Body[iThread].getTrack(kV0DauPos).getPxPyPzGlo(momPos);
179 mFitter3Body[iThread].getTrack(kV0DauNeg).getPxPyPzGlo(momNeg);
180 if (alphaV0 > 0) {
181 strangeTrack.mMasses[0] = calcMotherMass(momPos, momNeg, PID::Helium3, PID::Pion); // Hypertriton invariant mass at decay vertex
182 strangeTrack.mMasses[1] = calcMotherMass(momPos, momNeg, PID::Alpha, PID::Pion); // Hyperhydrogen4Lam invariant mass at decay vertex
183 } else {
184 strangeTrack.mMasses[0] = calcMotherMass(momPos, momNeg, PID::Pion, PID::Helium3); // Anti-Hypertriton invariant mass at decay vertex
185 strangeTrack.mMasses[1] = calcMotherMass(momPos, momNeg, PID::Pion, PID::Alpha); // Anti-Hyperhydrogen4Lam invariant mass at decay vertex
186 }
187
188 LOG(debug) << "ITS Track matched with a V0 decay topology ....";
189 LOG(debug) << "Number of ITS track clusters attached: " << itsTrack.getNumberOfClusters();
190 strangeTrack.mDecayRef = iv0;
191 strangeTrack.mITSRef = mSortedITSindexes[iTrack];
192 mStrangeTrackVec[iThread].push_back(strangeTrack);
193 mClusAttachments[iThread].push_back(structClus);
194 if (mMCTruthON) {
195 auto lab = getStrangeTrackLabel(itsTrack, strangeTrack, structClus);
196 mStrangeTrackLabels[iThread].push_back(lab);
197 }
198 }
199 }
200 }
201}
202
203void StrangenessTracker::processCascade(int iCasc, const Cascade& casc, const CascadeIndex& cascIdx, const V0& cascV0, int iThread)
204{
205 ClusAttachments structClus;
206 auto& daughterTracks = mDaughterTracks[iThread];
207 daughterTracks.resize(3); // resize to 3 prongs: first bachelor, second V0 pos, third V0 neg
208
209 StrangeTrack strangeTrack;
211 // first: bachelor, second: V0 pos, third: V0 neg
212 auto cascR = std::sqrt(casc.calcR2());
213 auto iBinsCasc = mUtils.getBinRect(casc.getEta(), casc.getPhi(), mStrParams->mEtaBinSize, mStrParams->mPhiBinSize);
214 for (int& iBinCasc : iBinsCasc) {
215 for (int iTrack{mTracksIdxTable[iBinCasc]}; iTrack < TMath::Min(mTracksIdxTable[iBinCasc + 1], int(mSortedITStracks.size())); iTrack++) {
216 strangeTrack.mMother = (o2::track::TrackParCovF)casc;
217 daughterTracks[kV0DauPos] = cascV0.getProng(kV0DauPos);
218 daughterTracks[kV0DauNeg] = cascV0.getProng(kV0DauNeg);
219 daughterTracks[kBach] = casc.getBachelorTrack();
220 const auto& itsTrack = mSortedITStracks[iTrack];
221 const auto& ITSindexRef = mSortedITSindexes[iTrack];
222 LOG(debug) << "----------------------";
223 LOG(debug) << "CascV0: " << cascIdx.getV0ID() << ", Bach ID: " << cascIdx.getBachelorID() << ", ITS track ref: " << mSortedITSindexes[iTrack];
224 if (mStrParams->mVertexMatching && (mITSvtxBrackets[ITSindexRef].getMin() > cascIdx.getVertexID() || mITSvtxBrackets[ITSindexRef].getMax() < cascIdx.getVertexID())) {
225 LOG(debug) << "Vertex ID mismatch: " << mITSvtxBrackets[ITSindexRef].getMin() << " < " << cascIdx.getVertexID() << " < " << mITSvtxBrackets[ITSindexRef].getMax();
226 continue;
227 }
228 if (matchDecayToITStrack(cascR, strangeTrack, structClus, itsTrack, daughterTracks, iThread)) {
229 auto propInstance = o2::base::Propagator::Instance();
230 o2::track::TrackParCov decayVtxTrackClone = strangeTrack.mMother; // clone track and propagate to decay vertex
231 if (!propInstance->propagateToX(decayVtxTrackClone, strangeTrack.mDecayVtx[0], getBz(), o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, mCorrType)) {
232 LOG(debug) << "Mother propagation to decay vertex failed";
233 continue;
234 }
235 decayVtxTrackClone.getPxPyPzGlo(strangeTrack.mDecayMom);
236 std::array<float, 3> momV0, momBach;
237 mFitter3Body[iThread].getTrack(0).getPxPyPzGlo(momV0); // V0 momentum at decay vertex
238 mFitter3Body[iThread].getTrack(1).getPxPyPzGlo(momBach); // bachelor momentum at decay vertex
239 strangeTrack.mMasses[0] = calcMotherMass(momV0, momBach, PID::Lambda, PID::Pion); // Xi invariant mass at decay vertex
240 strangeTrack.mMasses[1] = calcMotherMass(momV0, momBach, PID::Lambda, PID::Kaon); // Omega invariant mass at decay vertex
241
242 LOG(debug) << "ITS Track matched with a Cascade decay topology ....";
243 LOG(debug) << "Number of ITS track clusters attached: " << itsTrack.getNumberOfClusters();
244
245 strangeTrack.mDecayRef = iCasc;
246 strangeTrack.mITSRef = mSortedITSindexes[iTrack];
247 mStrangeTrackVec[iThread].push_back(strangeTrack);
248 mClusAttachments[iThread].push_back(mStructClus);
249 if (mMCTruthON) {
250 auto lab = getStrangeTrackLabel(itsTrack, strangeTrack, structClus);
251 mStrangeTrackLabels[iThread].push_back(lab);
252 }
253 }
254 }
255 }
256}
257
258void StrangenessTracker::process3Body(int i3Body, const Decay3Body& dec3body, const Decay3BodyIndex& dec3bodyIdx, int iThread)
259{
260 if (!mStrParams->mSkip3Body) {
261 ClusAttachments structClus;
262 auto& daughterTracks = mDaughterTracks[iThread];
263 daughterTracks.resize(3); // resize to 3 prongs
264
265 StrangeTrack strangeTrack;
267 auto dec3bodyR = std::sqrt(dec3body.calcR2());
268 auto iBins3Body = mUtils.getBinRect(dec3body.getEta(), dec3body.getPhi(), mStrParams->mEtaBinSize, mStrParams->mPhiBinSize);
269 for (int& iBin3Body : iBins3Body) {
270 for (int iTrack{mTracksIdxTable[iBin3Body]}; iTrack < TMath::Min(mTracksIdxTable[iBin3Body + 1], int(mSortedITStracks.size())); iTrack++) {
271 strangeTrack.mMother = (o2::track::TrackParCovF)dec3body;
273 daughterTracks[kV0DauPos] = dec3body.getProng(kV0DauPos); // proton
274 daughterTracks[kV0DauNeg] = dec3body.getProng(kV0DauNeg); // pion
275 daughterTracks[kBach] = dec3body.getProng(kBach); // deuteron
276 const auto& itsTrack = mSortedITStracks[iTrack];
277 const auto& ITSindexRef = mSortedITSindexes[iTrack];
278 if (mStrParams->mVertexMatching && (mITSvtxBrackets[ITSindexRef].getMin() > dec3bodyIdx.getVertexID() || mITSvtxBrackets[ITSindexRef].getMax() < dec3bodyIdx.getVertexID())) {
279 continue;
280 }
281 if (matchDecayToITStrack(dec3bodyR, strangeTrack, structClus, itsTrack, daughterTracks, iThread)) {
282 auto propInstance = o2::base::Propagator::Instance();
283 o2::track::TrackParCov decayVtxTrackClone = strangeTrack.mMother; // clone track and propagate to decay vertex
284 if (!propInstance->propagateToX(decayVtxTrackClone, strangeTrack.mDecayVtx[0], getBz(), o2::base::PropagatorImpl<float>::MAX_SIN_PHI, o2::base::PropagatorImpl<float>::MAX_STEP, mCorrType)) {
285 LOG(debug) << "Mother propagation to decay vertex failed";
286 continue;
287 }
288 decayVtxTrackClone.getPxPyPzGlo(strangeTrack.mDecayMom);
289 std::array<float, 3> momPos, momNeg, momBach;
290 mFitter4Body[iThread].propagateTracksToVertex();
291 mFitter4Body[iThread].getTrack(kV0DauPos).getPxPyPzGlo(momPos);
292 mFitter4Body[iThread].getTrack(kV0DauNeg).getPxPyPzGlo(momNeg);
293 mFitter4Body[iThread].getTrack(kBach).getPxPyPzGlo(momBach);
295 if (daughterTracks[kBach].getCharge() > 0) {
296 strangeTrack.mMasses[0] = calcMotherMass3body(momPos, momNeg, momBach, PID::Proton, PID::Pion, PID::Deuteron);
297 } else {
298 strangeTrack.mMasses[0] = calcMotherMass3body(momPos, momNeg, momBach, PID::Pion, PID::Proton, PID::Deuteron);
299 }
300
301 LOG(debug) << "ITS Track matched with a dec3body decay topology ....";
302 LOG(debug) << "Number of ITS track clusters attached: " << itsTrack.getNumberOfClusters();
303 strangeTrack.mDecayRef = i3Body;
304 strangeTrack.mITSRef = mSortedITSindexes[iTrack];
305 mStrangeTrackVec[iThread].push_back(strangeTrack);
306 mClusAttachments[iThread].push_back(structClus);
307 if (mMCTruthON) {
308 auto lab = getStrangeTrackLabel(itsTrack, strangeTrack, structClus);
309 mStrangeTrackLabels[iThread].push_back(lab);
310 }
311 }
312 }
313 }
314 } else {
315 return;
316 }
317}
318
320{
321 // Loop over V0s
322 for (int iV0{0}; iV0 < mInputV0tracks.size(); iV0++) {
323 LOG(debug) << "Analysing V0: " << iV0 + 1 << "/" << mInputV0tracks.size();
325 }
326
327 // Loop over Cascades
328 for (int iCasc{0}; iCasc < mInputCascadeTracks.size(); iCasc++) {
329 LOG(debug) << "Analysing Cascade: " << iCasc + 1 << "/" << mInputCascadeTracks.size();
331 }
332
333 // Loop over 3bodys
334 if (!mStrParams->mSkip3Body) {
335 for (int i3Body{0}; i3Body < mInput3BodyTracks.size(); i3Body++) {
336 LOG(debug) << "Analysing 3-Body: " << i3Body + 1 << "/" << mInput3BodyTracks.size();
337 process3Body(i3Body, mInput3BodyTracks[i3Body], mInput3BodyIndices[i3Body]);
338 }
339 }
340}
341
342bool StrangenessTracker::matchDecayToITStrack(float decayR, StrangeTrack& strangeTrack, ClusAttachments& structClus, const TrackITS& itsTrack, std::vector<o2::track::TrackParCovF>& daughterTracks, int iThread)
343{
345 auto trackClusters = getTrackClusters(itsTrack);
346 auto trackClusSizes = getTrackClusterSizes(itsTrack);
347 auto& lastClus = trackClusters[0];
348
349 auto radTol = decayR < 4 ? mStrParams->mRadiusTolIB : mStrParams->mRadiusTolOB;
350 auto nMinClusMother = trackClusters.size() < 4 ? 2 : mStrParams->mMinMotherClus;
351
352 std::vector<ITSCluster> motherClusters;
353 std::array<unsigned int, 7> nAttachments;
354 nAttachments.fill(-1); // fill arr with -1
355
356 int nUpdates = 0;
357 bool isMotherUpdated = false;
358
359 for (int iClus{0}; iClus < trackClusters.size(); iClus++) {
360 auto& clus = trackClusters[iClus];
361 auto& compClus = trackClusSizes[iClus];
362 int nUpdOld = nUpdates;
363 double clusRad = sqrt(clus.getX() * clus.getX() - clus.getY() * clus.getY());
364 auto diffR = decayR - clusRad;
365 auto relDiffR = diffR / decayR;
366 auto lay = geom->getLayer(clus.getSensorID());
367 // Look for the Mother if the Decay radius allows for it, within a tolerance
368 LOG(debug) << "decayR: " << decayR << ", diffR: " << diffR << ", clus rad: " << clusRad << ", radTol: " << radTol;
369 if (relDiffR > -radTol) {
370 LOG(debug) << "Try to attach cluster to Mother, layer: " << lay;
371 if (updateTrack(clus, strangeTrack.mMother)) {
372 motherClusters.push_back(clus);
373 strangeTrack.setClusterSize(lay, compClus);
374 nAttachments[lay] = 0;
375 isMotherUpdated = true;
376 nUpdates++;
377 LOG(debug) << "Cluster attached to Mother";
378 continue; // if the cluster is attached to the mother, skip the rest of the loop
379 }
380 }
381
382 // if Mother is not found, check for V0 daughters compatibility
383 if (relDiffR < radTol && !isMotherUpdated) {
384 bool isDauUpdated = false;
385 LOG(debug) << "Try to attach cluster to Daughters, layer: " << lay;
386 for (int iDau{0}; iDau < daughterTracks.size(); iDau++) {
387 auto& dauTrack = daughterTracks[iDau];
388 if (updateTrack(clus, dauTrack)) {
389 nAttachments[lay] = iDau + 1;
390 isDauUpdated = true;
391 break;
392 }
393 }
394 if (!isDauUpdated) {
395 break; // no daughter track updated, stop the loop
396 }
397 nUpdates++;
398 }
399 if (nUpdates == nUpdOld) {
400 break; // no track updated, stop the loop
401 }
402 }
403
404 if (nUpdates < trackClusters.size() || motherClusters.size() < nMinClusMother) {
405 return false;
406 }
407
408 o2::track::TrackParCov motherTrackClone = strangeTrack.mMother; // clone and reset covariance for final topology refit
409 motherTrackClone.resetCovariance();
410
411 LOG(debug) << "Clusters attached, starting inward-outward refit";
412
413 std::reverse(motherClusters.begin(), motherClusters.end());
414
415 mGlobalChi2 = -1;
416 for (auto& clus : motherClusters) {
417 if (!updateTrack(clus, motherTrackClone)) {
418 break;
419 }
420 }
421 strangeTrack.mMatchChi2 = mGlobalChi2;
422
423 LOG(debug) << "Inward-outward refit finished, starting final topology refit";
424 // final Topology refit
425
426 int cand = 0; // best V0 candidate
427 int nCand;
428
429 // refit cascade
430 if (strangeTrack.mPartType == dataformats::kStrkCascade) {
431 V0 cascV0Upd;
432 if (!recreateV0(daughterTracks[kV0DauPos], daughterTracks[kV0DauNeg], cascV0Upd, iThread)) {
433 LOG(debug) << "Cascade V0 refit failed";
434 return false;
435 }
436 try {
437 nCand = mFitter3Body[iThread].process(cascV0Upd, daughterTracks[kBach], motherTrackClone);
438 } catch (std::runtime_error& e) {
439 LOG(debug) << "Fitter3Body failed: " << e.what();
440 return false;
441 }
442 if (!nCand || !mFitter3Body[iThread].propagateTracksToVertex()) {
443 LOG(debug) << "Fitter3Body failed: propagation to vertex failed";
444 return false;
445 }
446 }
447
448 // refit V0
449 else if (strangeTrack.mPartType == dataformats::kStrkV0) {
450 try {
451 nCand = mFitter3Body[iThread].process(daughterTracks[kV0DauPos], daughterTracks[kV0DauNeg], motherTrackClone);
452 } catch (std::runtime_error& e) {
453 LOG(debug) << "Fitter3Body failed: " << e.what();
454 return false;
455 }
456 if (!nCand || !mFitter3Body[iThread].propagateTracksToVertex()) {
457 LOG(debug) << "Fitter3Body failed: propagation to vertex failed";
458 return false;
459 }
460 }
461
462 // refit 3body
463 else if (strangeTrack.mPartType == dataformats::kStrkThreeBody) {
464 try {
465 nCand = mFitter4Body[iThread].process(daughterTracks[kV0DauPos], daughterTracks[kV0DauNeg], daughterTracks[kBach], motherTrackClone);
466 } catch (std::runtime_error& e) {
467 LOG(debug) << "Fitter4Body failed: " << e.what();
468 return false;
469 }
470 if (!nCand || !mFitter4Body[iThread].propagateTracksToVertex()) {
471 LOG(debug) << "Fitter4Body failed: propagation to vertex failed";
472 return false;
473 }
474 }
475
476 if (strangeTrack.mPartType == dataformats::kStrkThreeBody) {
477 strangeTrack.mDecayVtx = mFitter4Body[iThread].getPCACandidatePos();
478 strangeTrack.mTopoChi2 = mFitter4Body[iThread].getChi2AtPCACandidate();
479 } else {
480 strangeTrack.mDecayVtx = mFitter3Body[iThread].getPCACandidatePos();
481 strangeTrack.mTopoChi2 = mFitter3Body[iThread].getChi2AtPCACandidate();
482 }
483 structClus.arr = nAttachments;
484
485 return true;
486}
487
489{
491 auto propInstance = o2::base::Propagator::Instance();
492 float alpha = geom->getSensorRefAlpha(clus.getSensorID()), x = clus.getX();
493 int layer{geom->getLayer(clus.getSensorID())};
494
495 if (!track.rotate(alpha)) {
496 return false;
497 }
498
500 return false;
501 }
502
503 if (mCorrType == o2::base::PropagatorF::MatCorrType::USEMatCorrNONE) {
504 float thick = layer < 3 ? 0.005 : 0.01;
505 constexpr float radl = 9.36f; // Radiation length of Si [cm]
506 constexpr float rho = 2.33f; // Density of Si [g/cm^3]
507 if (!track.correctForMaterial(thick, thick * rho * radl)) {
508 return false;
509 }
510 }
511 auto chi2 = std::abs(track.getPredictedChi2Quiet(clus)); // abs to be understood
512 LOG(debug) << "Chi2: " << chi2;
513 if (chi2 > mStrParams->mMaxChi2 || chi2 < 0) {
514 return false;
515 }
516
517 if (!track.update(clus)) {
518 return false;
519 }
520
521 mGlobalChi2 += chi2;
522 return true;
523}
524
525} // namespace strangeness_tracking
526} // namespace o2
int32_t i
uint32_t j
Definition RawData.h:0
std::ostringstream debug
T getX() const
Definition BaseCluster.h:62
std::int16_t getSensorID() const
Definition BaseCluster.h:81
GPUd() value_type estimateLTFast(o2 static GPUd() float estimateLTIncrement(const o2 PropagatorImpl * Instance(bool uninitialized=false)
Definition Propagator.h:143
const Track & getBachelorTrack() const
Definition Cascade.h:40
TO BE DONE: extend to generic N body vertex.
Definition Decay3Body.h:26
const Track & getProng(int i) const
Definition Decay3Body.h:34
const Track & getProng(int i) const
Definition V0.h:39
float calcR2() const
Definition V0.h:63
static GeometryTGeo * Instance()
std::vector< std::vector< o2::track::TrackParCovF > > mDaughterTracks
o2::MCCompLabel getStrangeTrackLabel(const TrackITS &itsTrack, const StrangeTrack &strangeTrack, const ClusAttachments &structClus)
const StrangenessTrackingParamConfig * mStrParams
double calcMotherMass(const std::array< float, 3 > &pDauFirst, const std::array< float, 3 > &pDauSecond, PID pidDauFirst, PID pidDauSecond)
gsl::span< const TrackITS > mInputITStracks
global topology matching chi2
std::vector< std::vector< StrangeTrack > > mStrangeTrackVec
void processV0(int iv0, const V0 &v0, const V0Index &v0Idx, int iThread=0)
void processCascade(int icasc, const Cascade &casc, const CascadeIndex &cascIdx, const V0 &cascV0, int iThread=0)
std::vector< ITSCluster > getTrackClusters(const TrackITS &itsTrack)
o2::base::PropagatorImpl< float >::MatCorrType mCorrType
std::vector< int > getTrackClusterSizes(const TrackITS &itsTrack)
void getClusterSizesITS(std::vector< int > &clusSizeVec, const gsl::span< const o2::itsmft::CompClusterExt > ITSclus, gsl::span< const unsigned char >::iterator &pattIt, const o2::itsmft::TopologyDictionary *mdict)
std::vector< std::vector< o2::MCCompLabel > > mStrangeTrackLabels
bool matchDecayToITStrack(float decayR, StrangeTrack &strangeTrack, ClusAttachments &structClus, const TrackITS &itsTrack, std::vector< o2::track::TrackParCovF > &daughterTracks, int iThread=0)
std::vector< o2::its::TrackITS > mSortedITStracks
input ITS Track MC labels
gsl::span< const Decay3Body > mInput3BodyTracks
bool recreateV0(const o2::track::TrackParCov &posTrack, const o2::track::TrackParCov &negTrack, V0 &newV0, int iThread=0)
std::vector< std::vector< ClusAttachments > > mClusAttachments
double calcMotherMass3body(const std::array< float, 3 > &pDauFirst, const std::array< float, 3 > &pDauSecond, const std::array< float, 3 > &pDauThird, PID pidDauFirst, PID pidDauSecond, PID pidDauThird)
float mGlobalChi2
number of threads (externally driven)
bool updateTrack(const ITSCluster &clus, o2::track::TrackParCov &track)
void process3Body(int i3body, const Decay3Body &dec3body, const Decay3BodyIndex &dec3bodyIdx, int iThread=0)
gsl::span< const Decay3BodyIndex > mInput3BodyIndices
const o2::itsmft::TopologyDictionary * mITSDict
gsl::span< const CascadeIndex > mInputCascadeIndices
bool loadData(const o2::globaltracking::RecoContainer &recoData)
MCLabSpan mITSTrkLabels
input ITS Cluster MC labels
static constexpr ID Lambda
Definition PID.h:112
static constexpr ID Helium3
Definition PID.h:101
static constexpr ID Kaon
Definition PID.h:97
static constexpr ID Pion
Definition PID.h:96
static constexpr ID Deuteron
Definition PID.h:99
static constexpr ID Proton
Definition PID.h:98
static constexpr ID Alpha
Definition PID.h:102
GLfloat GLfloat GLfloat alpha
Definition glcorearb.h:279
GLint GLenum GLint x
Definition glcorearb.h:403
GLboolean GLboolean GLboolean b
Definition glcorearb.h:1233
GLfloat v0
Definition glcorearb.h:811
GLenum GLuint GLint GLint layer
Definition glcorearb.h:1310
GLboolean GLboolean GLboolean GLboolean a
Definition glcorearb.h:1233
void convertCompactClusters(gsl::span< const itsmft::CompClusterExt > clusters, gsl::span< const unsigned char >::iterator &pattIt, std::vector< o2::BaseCluster< float > > &output, const its3::TopologyDictionary *dict)
convert compact clusters to 3D spacepoints
Definition IOUtils.cxx:31
void convertCompactClusters(gsl::span< const itsmft::CompClusterExt > clusters, gsl::span< const unsigned char >::iterator &pattIt, std::vector< o2::BaseCluster< float > > &output, const itsmft::TopologyDictionary *dict)
convert compact clusters to 3D spacepoints
Definition IOUtils.cxx:49
TrackParCovF TrackParCov
Definition Track.h:33
TrackParametrizationWithError< float > TrackParCovF
Definition Track.h:31
a couple of static helper functions to create timestamp values for CCDB queries or override obsolete ...
o2::track::TrackParCovF mMother
std::array< float, 3 > mDecayVtx
std::array< float, 3 > mDecayMom
std::array< float, 2 > mMasses
void setClusterSize(int l, int size)
std::unique_ptr< const o2::dataformats::MCTruthContainer< o2::MCCompLabel > > mcITSClusters
std::vector< int > getBinRect(float eta, float phi, float deltaEta, float deltaPhi)
LOG(info)<< "Compressed in "<< sw.CpuTime()<< " s"