Project
Loading...
Searching...
No Matches
GPUTPCGMO2Output.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.
11
14
15#include "GPUTPCDef.h"
16#include "GPUTPCGMO2Output.h"
17#include "GPUCommonAlgorithm.h"
21#include "TPCFastTransform.h"
23#include "GPUGetConstexpr.h"
24
25#ifndef GPUCA_GPUCODE
28#include "GPUQAHelper.h"
29#endif
30
31using namespace o2::gpu;
32using namespace o2::tpc;
33using namespace o2::tpc::constants;
34
35GPUdi() static constexpr uint8_t getFlagsReject() { return GPUTPCGMMergedTrackHit::flagReject | GPUTPCGMMergedTrackHit::flagHighIncl; }
36GPUdi() static uint32_t getFlagsRequired(const GPUSettingsRec& rec) { return gputpcgmmergertypes::attachGoodLeg; }
37
38namespace o2::gpu::internal
39{
40
42 GPUd() bool operator()(const GPUTPCGMMerger::tmpSort& a, const GPUTPCGMMerger::tmpSort& b)
43 {
44 return GPUCA_DETERMINISTIC_CODE(a.y != b.y ? a.y > b.y : a.x > b.x, a.y > b.y);
45 }
46};
47
48} // namespace o2::gpu::internal
49
50template <>
51GPUdii() void GPUTPCGMO2Output::Thread<GPUTPCGMO2Output::prepare>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUsharedref() GPUSharedMemory& smem, processorType& GPUrestrict() merger)
52{
53 const GPUTPCGMMergedTrack* tracks = merger.MergedTracks();
54 const uint32_t nTracks = merger.NMergedTracks();
55 const GPUTPCGMMergedTrackHit* trackClusters = merger.Clusters();
56 const GPUdEdxInfo* tracksdEdx = merger.MergedTracksdEdx();
57
58 constexpr uint8_t flagsReject = getFlagsReject();
59 const uint32_t flagsRequired = getFlagsRequired(merger.Param().rec);
60 bool cutOnTrackdEdx = merger.Param().par.dodEdx && merger.Param().dodEdxEnabled && merger.Param().rec.tpc.minTrackdEdxMax2Tot > 0.f;
61
62 GPUTPCGMMerger::tmpSort* GPUrestrict() trackSort = merger.TrackSortO2();
63 uint2* GPUrestrict() tmpData = merger.ClusRefTmp();
64 for (uint32_t i = get_global_id(0); i < nTracks; i += get_global_size(0)) {
65 if (!tracks[i].OK()) {
66 continue;
67 }
68 if (tracks[i].MergedLooper()) {
69 continue;
70 }
71
72 uint32_t nCl = 0;
73 for (uint32_t j = 0; j < tracks[i].NClusters(); j++) {
74 if ((trackClusters[tracks[i].FirstClusterRef() + j].state & flagsReject) || (merger.ClusterAttachment()[trackClusters[tracks[i].FirstClusterRef() + j].num] & flagsRequired) != flagsRequired) {
75 continue;
76 }
77 nCl++;
78 }
79 if (nCl == 0) {
80 continue;
81 }
82 if (nCl + 2 < GPUCA_TRACKLET_SELECTOR_MIN_HITS_B5(tracks[i].GetParam().GetQPt() * merger.Param().qptB5Scaler)) { // Give 2 hits tolerance in the primary leg, compared to the full fit of the looper
83 continue;
84 }
85 if (merger.Param().rec.tpc.minNClustersFinalTrack != -1 && nCl < (uint32_t)merger.Param().rec.tpc.minNClustersFinalTrack) {
86 continue;
87 }
88 if (cutOnTrackdEdx && (tracksdEdx[i].dEdxMaxTPC < merger.Param().rec.tpc.minTrackdEdxMax || tracksdEdx[i].dEdxMaxTPC < tracksdEdx[i].dEdxTotTPC * merger.Param().rec.tpc.minTrackdEdxMax2Tot) && !(tracksdEdx[i].dEdxMaxTPC == 0 && CAMath::Abs(tracks[i].GetParam().GetDzDs()) > 0.03f)) {
89 continue;
90 }
91 uint32_t myId = CAMath::AtomicAdd(&merger.Memory()->nO2Tracks, 1u);
92 tmpData[i] = {nCl, CAMath::AtomicAdd(&merger.Memory()->nO2ClusRefs, nCl + (nCl + 1) / 2)};
93 trackSort[myId] = {i, tracks[i].CSide() ? tracks[i].GetParam().GetTOffset() : -tracks[i].GetParam().GetTOffset()};
94 }
95}
96
97template <>
98GPUdii() void GPUTPCGMO2Output::Thread<GPUTPCGMO2Output::sort>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUsharedref() GPUSharedMemory& smem, processorType& GPUrestrict() merger)
99{
100#ifndef GPUCA_SPECIALIZE_THRUST_SORTS
101 if (iThread == 0 && iBlock == 0) {
102 GPUTPCGMMerger::tmpSort* GPUrestrict() trackSort = merger.TrackSortO2();
103 GPUCommonAlgorithm::sortDeviceDynamic(trackSort, trackSort + merger.Memory()->nO2Tracks, internal::GPUTPCGMO2OutputSort_comp());
104 }
105#endif
106}
107
108template <>
109GPUdii() void GPUTPCGMO2Output::Thread<GPUTPCGMO2Output::output>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUsharedref() GPUSharedMemory& smem, processorType& GPUrestrict() merger)
110{
111 constexpr float MinDelta = 0.1f;
112 const GPUTPCGMMergedTrack* tracks = merger.MergedTracks();
113 GPUdEdxInfo* tracksdEdx = merger.MergedTracksdEdx();
114 GPUdEdxInfo* tracksdEdxAlt = merger.MergedTracksdEdxAlt();
115 const int32_t nTracks = merger.NOutputTracksTPCO2();
116 const GPUTPCGMMergedTrackHit* trackClusters = merger.Clusters();
117 constexpr uint8_t flagsReject = getFlagsReject();
118 const uint32_t flagsRequired = getFlagsRequired(merger.Param().rec);
119 TrackTPC* outputTracks = merger.OutputTracksTPCO2();
120 uint32_t* clusRefs = merger.OutputClusRefsTPCO2();
121 const auto& param = merger.Param();
122
123 GPUTPCGMMerger::tmpSort* GPUrestrict() trackSort = merger.TrackSortO2();
124 uint2* GPUrestrict() tmpData = merger.ClusRefTmp();
125 float const SNPThresh = 0.999990f;
126
127 for (int32_t iTmp = get_global_id(0); iTmp < nTracks; iTmp += get_global_size(0)) {
128 TrackTPC oTrack;
129 const int32_t i = trackSort[iTmp].x;
130 const auto& track = tracks[i];
131 auto snpIn = track.GetParam().GetSinPhi();
132 if (snpIn > SNPThresh) {
133 snpIn = SNPThresh;
134 } else if (snpIn < -SNPThresh) {
135 snpIn = -SNPThresh;
136 }
137 oTrack.set(track.GetParam().GetX(), track.GetAlpha(),
138 {track.GetParam().GetY(), track.GetParam().GetZ(), snpIn, track.GetParam().GetDzDs(), track.GetParam().GetQPt()},
139 {track.GetParam().GetCov(0),
140 track.GetParam().GetCov(1), track.GetParam().GetCov(2),
141 track.GetParam().GetCov(3), track.GetParam().GetCov(4), track.GetParam().GetCov(5),
142 track.GetParam().GetCov(6), track.GetParam().GetCov(7), track.GetParam().GetCov(8), track.GetParam().GetCov(9),
143 track.GetParam().GetCov(10), track.GetParam().GetCov(11), track.GetParam().GetCov(12), track.GetParam().GetCov(13), track.GetParam().GetCov(14)});
144
145 oTrack.setChi2(track.GetParam().GetChi2());
146 auto& outerPar = track.OuterParam();
148 if (param.dodEdxEnabled) {
149 oTrack.setdEdx(tracksdEdx[i]);
150 if GPUCA_RTC_CONSTEXPR (GPUCA_GET_CONSTEXPR(param.rec.tpc, dEdxClusterRejectionFlagMask) != GPUCA_GET_CONSTEXPR(param.rec.tpc, dEdxClusterRejectionFlagMaskAlt)) {
151 oTrack.setdEdxAlt(tracksdEdxAlt[i]);
152 } else {
153 oTrack.setdEdxAlt(tracksdEdx[i]);
154 }
155 }
156 }
157
158 auto snpOut = outerPar.P[2];
159 if (snpOut > SNPThresh) {
160 snpOut = SNPThresh;
161 } else if (snpOut < -SNPThresh) {
162 snpOut = -SNPThresh;
163 }
164 oTrack.setOuterParam(o2::track::TrackParCov(
165 outerPar.X, outerPar.alpha,
166 {outerPar.P[0], outerPar.P[1], snpOut, outerPar.P[3], outerPar.P[4]},
167 {outerPar.C[0], outerPar.C[1], outerPar.C[2], outerPar.C[3], outerPar.C[4], outerPar.C[5],
168 outerPar.C[6], outerPar.C[7], outerPar.C[8], outerPar.C[9], outerPar.C[10], outerPar.C[11],
169 outerPar.C[12], outerPar.C[13], outerPar.C[14]}));
170
171 if (param.par.dodEdx && param.dodEdxEnabled && param.rec.tpc.enablePID) {
172 PIDResponse pidResponse{};
173 auto pid = pidResponse.getMostProbablePID(oTrack, param.rec.tpc.PID_EKrangeMin, param.rec.tpc.PID_EKrangeMax, param.rec.tpc.PID_EPrangeMin, param.rec.tpc.PID_EPrangeMax, param.rec.tpc.PID_EDrangeMin, param.rec.tpc.PID_EDrangeMax, param.rec.tpc.PID_ETrangeMin, merger.Param().rec.tpc.PID_ETrangeMax, merger.Param().rec.tpc.PID_useNsigma, merger.Param().rec.tpc.PID_sigma);
174 auto pidRemap = merger.Param().rec.tpc.PID_remap[pid];
175 if (pidRemap >= 0) {
176 pid = pidRemap;
177 }
178 oTrack.setPID(pid, true);
179 oTrack.getParamOut().setPID(pid, true);
180 }
181
182 uint32_t nOutCl = tmpData[i].x;
183 uint32_t clBuff = tmpData[i].y;
184 oTrack.setClusterRef(clBuff, nOutCl);
185 uint32_t* clIndArr = &clusRefs[clBuff];
186 uint8_t* sectorIndexArr = reinterpret_cast<uint8_t*>(clIndArr + nOutCl);
187 uint8_t* rowIndexArr = sectorIndexArr + nOutCl;
188
189 uint32_t nOutCl2 = 0;
190 float t1 = 0, t2 = 0;
191 int32_t sector1 = 0, sector2 = 0;
192 const o2::tpc::ClusterNativeAccess* GPUrestrict() clusters = merger.GetConstantMem()->ioPtrs.clustersNative;
193 for (uint32_t j = 0; j < track.NClusters(); j++) {
194 if ((trackClusters[track.FirstClusterRef() + j].state & flagsReject) || (merger.ClusterAttachment()[trackClusters[track.FirstClusterRef() + j].num] & flagsRequired) != flagsRequired) {
195 continue;
196 }
197 int32_t clusterIdGlobal = trackClusters[track.FirstClusterRef() + j].num;
198 int32_t sector = trackClusters[track.FirstClusterRef() + j].sector;
199 int32_t globalRow = trackClusters[track.FirstClusterRef() + j].row;
200 int32_t clusterIdInRow = clusterIdGlobal - clusters->clusterOffset[sector][globalRow];
201 clIndArr[nOutCl2] = clusterIdInRow;
202 sectorIndexArr[nOutCl2] = sector;
203 rowIndexArr[nOutCl2] = globalRow;
204 if (nOutCl2 == 0) {
205 t1 = clusters->clustersLinear[clusterIdGlobal].getTime();
206 sector1 = sector;
207 }
208 if (++nOutCl2 == nOutCl) {
209 t2 = clusters->clustersLinear[clusterIdGlobal].getTime();
210 sector2 = sector;
211 }
212 }
213
214 if (track.PrevSegment() >= 0) {
215 const GPUTPCGMMergedTrack* chkTrk = track.GetFirstSegment(tracks, merger.Param().rec.enableCyclicGraphWorkarounds);
216 const auto& firstPrevCluster = trackClusters[chkTrk->FirstClusterRef()];
217 t1 = clusters->clustersLinear[firstPrevCluster.num].getTime();
218 sector1 = firstPrevCluster.sector;
219 }
220
221 bool cce = track.CCE() && ((sector1 < MAXSECTOR / 2) ^ (sector2 < MAXSECTOR / 2));
222 float time0 = 0.f, tFwd = 0.f, tBwd = 0.f;
223 if (merger.Param().par.continuousTracking) {
224 time0 = track.GetParam().GetTOffset();
225 if (cce) {
226 bool lastSide = trackClusters[track.FirstClusterRef()].sector < MAXSECTOR / 2;
227 float delta = 0.f;
228 for (uint32_t iCl = 1; iCl < track.NClusters(); iCl++) {
229 auto& cacl1 = trackClusters[track.FirstClusterRef() + iCl];
230 if (lastSide ^ (cacl1.sector < MAXSECTOR / 2)) {
231 auto& cl1 = clusters->clustersLinear[cacl1.num];
232 auto& cl2 = clusters->clustersLinear[trackClusters[track.FirstClusterRef() + iCl - 1].num];
233 delta = CAMath::Abs(cl1.getTime() - cl2.getTime()) * 0.5f;
234 break;
235 }
236 }
237 if (delta < MinDelta) {
238 delta = MinDelta;
239 }
240 tFwd = tBwd = delta;
241 } else {
242 // estimate max/min time increments which still keep track in the physical limits of the TPC
243 const float tmin = CAMath::Min(t1, t2);
244 const float maxDriftTime = merger.GetConstantMem()->calibObjects.fastTransformHelper->getCorrMap()->getMaxDriftTime(t1 > t2 ? sector1 : sector2);
245 const float clusterT0 = merger.GetConstantMem()->calibObjects.fastTransformHelper->getCorrMap()->getT0();
246 const float tmax = CAMath::Min(tmin + maxDriftTime, CAMath::Max(t1, t2));
247 float delta = 0.f;
248 if (time0 + maxDriftTime < tmax) {
249 delta = tmax - time0 - maxDriftTime;
250 }
251 if (tmin - clusterT0 < time0 + delta) {
252 delta = tmin - clusterT0 - time0;
253 }
254 if (delta != 0.f) {
255 time0 += delta;
256 const float deltaZ = merger.GetConstantMem()->calibObjects.fastTransformHelper->getCorrMap()->convDeltaTimeToDeltaZinTimeFrame(sector2, delta);
257 oTrack.setZ(oTrack.getZ() + deltaZ);
258 }
259 tFwd = tmin - clusterT0 - time0;
260 tBwd = time0 - tmax + maxDriftTime;
261 }
262 }
263 if (tBwd < 0.f) {
264 tBwd = 0.f;
265 }
266 oTrack.setTime0(time0);
267 oTrack.setDeltaTBwd(tBwd);
268 oTrack.setDeltaTFwd(tFwd);
269 if (cce) {
270 oTrack.setHasCSideClusters();
271 oTrack.setHasASideClusters();
272 } else if (track.CSide()) {
273 oTrack.setHasCSideClusters();
274 } else {
275 oTrack.setHasASideClusters();
276 }
277 outputTracks[iTmp] = oTrack;
278 }
279}
280
281template <>
282GPUdii() void GPUTPCGMO2Output::Thread<GPUTPCGMO2Output::mc>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUsharedref() GPUSharedMemory& smem, processorType& GPUrestrict() merger)
283{
284#ifndef GPUCA_GPUCODE
285 const o2::tpc::ClusterNativeAccess* GPUrestrict() clusters = merger.GetConstantMem()->ioPtrs.clustersNative;
286 if (clusters == nullptr || clusters->clustersMCTruth == nullptr) {
287 return;
288 }
289 if (merger.OutputTracksTPCO2MC() == nullptr) {
290 return;
291 }
292
293 auto labelAssigner = GPUTPCTrkLbl(clusters->clustersMCTruth, 0.1f);
294 uint32_t* clusRefs = merger.OutputClusRefsTPCO2();
295 for (uint32_t i = get_global_id(0); i < merger.NOutputTracksTPCO2(); i += get_global_size(0)) {
296 labelAssigner.reset();
297 const auto& trk = merger.OutputTracksTPCO2()[i];
298 for (int32_t j = 0; j < trk.getNClusters(); j++) {
299 uint8_t sectorIndex, rowIndex;
300 uint32_t clusterIndex;
301 trk.getClusterReference(clusRefs, j, sectorIndex, rowIndex, clusterIndex);
302 uint32_t clusterIdGlobal = clusters->clusterOffset[sectorIndex][rowIndex] + clusterIndex;
303 labelAssigner.addLabel(clusterIdGlobal);
304 }
305 merger.OutputTracksTPCO2MC()[i] = labelAssigner.computeLabel();
306 }
307#endif
308}
benchmark::State & state
A const (ready only) version of MCTruthContainer.
Helper class to access correction maps.
int32_t i
#define GPUsharedref()
#define get_global_size(dim)
#define GPUrestrict()
#define get_global_id(dim)
#define GPUCA_RTC_CONSTEXPR
#define GPUCA_DETERMINISTIC_CODE(det, indet)
#define GPUCA_TRACKLET_SELECTOR_MIN_HITS_B5(QPTB5)
#define GPUCA_GET_CONSTEXPR(obj, val)
GPUdii() void GPUTPCGMO2Output
void output(const std::map< std::string, ChannelStat > &channels)
Definition rawdump.cxx:197
uint32_t j
Definition RawData.h:0
uint16_t pid
Definition RawData.h:2
Definition of TPCFastTransform class.
PID response class.
Definition PIDResponse.h:39
GLboolean GLboolean GLboolean b
Definition glcorearb.h:1233
typedef void(APIENTRYP PFNGLCULLFACEPROC)(GLenum mode)
GLenum GLfloat param
Definition glcorearb.h:271
GLboolean GLboolean GLboolean GLboolean a
Definition glcorearb.h:1233
GLuint GLfloat GLfloat GLfloat GLfloat GLfloat GLfloat GLfloat GLfloat GLfloat t1
Definition glcorearb.h:5034
constexpr int MAXSECTOR
Definition Constants.h:28
Global TPC definitions and constants.
Definition SimTraits.h:167
GPUdi() T BetheBlochAleph(T bg
TrackParCovF TrackParCov
Definition Track.h:33
GPUReconstruction * rec
GPUd() bool operator()(const GPUTPCGMMerger
uint32_t x
uint32_t y
std::vector< Cluster > clusters