Project
Loading...
Searching...
No Matches
GPUTPCNNClusterizerKernels.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 "clusterFinderDefs.h"
16#include "PackedCharge.h"
18#include "GPUConstantMem.h"
19#include "GPUTPCClusterFinder.h"
20#include "GPUTPCCFClusterizer.h"
21#include "GPUTPCGeometry.h"
22
23using namespace o2::gpu;
24using namespace o2::gpu::tpccf;
25
26#include "CfConsts.h"
27#include "CfUtils.h"
28#include "ClusterAccumulator.h"
30
31#if !defined(GPUCA_GPUCODE)
32#include "GPUHostDataTypes.h"
33#include "MCLabelAccumulator.h"
34#endif
35
36#ifdef GPUCA_GPUCODE
37#include "GPUTPCCFClusterizer.inc"
38#endif
39
40// Defining individual thread functions for data filling, determining the class label and running the CF clusterizer
41template <>
42GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::runCfClusterizer>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
43{
44 uint32_t glo_idx = get_global_id(0);
45 auto& clusterer = processors.tpcClusterer[sector];
46 auto& clustererNN = processors.tpcNNClusterer[sector];
47 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
48 CPU_ONLY(MCLabelAccumulator labelAcc(clusterer));
49 tpc::ClusterNative* clusterOut = (withMC) ? nullptr : clusterer.mPclusterByRow;
50 int8_t isAccepted = (clustererNN.mNnClusterizerUseClassification ? (clustererNN.mOutputDataClass[CAMath::Min(glo_idx, (uint32_t)clusterer.mPmemory->counters.nClusters - 1)] > 0) : 1);
51 GPUTPCCFClusterizer::computeClustersImpl(get_num_groups(0), get_local_size(0), get_group_id(0), get_local_id(0), clusterer, clusterer.mPmemory->fragment, reinterpret_cast<GPUTPCCFClusterizer::GPUSharedMemory&>(smem), chargeMap, clusterer.mPfilteredPeakPositions, clusterer.Param().rec, CPU_PTR(&labelAcc), clusterer.mPmemory->counters.nClusters, clusterer.mNMaxClusterPerRow, clusterer.mPclusterInRow, clusterOut, clusterer.mPclusterPosInRow, isAccepted);
52}
53
54template <>
55GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::fillInputNNCPU>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
56{
57 auto& clusterer = processors.tpcClusterer[sector];
58 auto& clustererNN = processors.tpcNNClusterer[sector];
59
60 uint32_t glo_idx = get_global_id(0);
61 if (glo_idx + batchStart >= clusterer.mPmemory->counters.nClusters || glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
62 return;
63 }
64
65 uint32_t write_idx = glo_idx * clustererNN.mNnClusterizerElementSize;
66
67 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
68 CfArray2D<uint8_t> isPeakMap(clusterer.mPpeakMap);
69 CfChargePos peak = clusterer.mPfilteredPeakPositions[CAMath::Min(glo_idx + batchStart, (uint32_t)(clusterer.mPmemory->counters.nClusters - 1))];
70 int32_t row = static_cast<int>(peak.row());
71 int32_t pad = static_cast<int>(peak.pad());
72 int32_t time = static_cast<int>(peak.time());
73 float central_charge = static_cast<float>(chargeMap[peak].unpack());
74 int32_t row_offset = GPUTPCNNClusterizerKernels::rowOffset(row, clustererNN.mNnClusterizerSizeInputRow);
75
76 for (int32_t r = -clustererNN.mNnClusterizerSizeInputRow; r <= clustererNN.mNnClusterizerSizeInputRow; ++r) {
77 int32_t target_row = row + r;
78 bool is_row_boundary = (target_row < 0) || (target_row >= o2::tpc::constants::MAXGLOBALPADROW);
79 int32_t pad_offset = is_row_boundary ? 0 : GPUTPCNNClusterizerKernels::padOffset(row, target_row);
80
81 for (int32_t p = -clustererNN.mNnClusterizerSizeInputPad + pad_offset; p <= clustererNN.mNnClusterizerSizeInputPad + pad_offset; ++p) {
82 int32_t target_pad = pad + p;
83 bool is_boundary = is_row_boundary || GPUTPCNNClusterizerKernels::isBoundary(target_row + row_offset, target_pad, clustererNN.mNnClusterizerSizeInputRow);
84
85 for (int32_t t = -clustererNN.mNnClusterizerSizeInputTime; t <= clustererNN.mNnClusterizerSizeInputTime; ++t) {
86 int32_t target_time = time + t;
87
88 if (is_boundary || target_time < 0 || target_time >= TPC_MAX_FRAGMENT_LEN_GPU) {
89 // Fill boundary value
90 float boundary_value = static_cast<float>(clustererNN.mNnClusterizerBoundaryFillValue);
91 if (dtype == 0) {
92 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)boundary_value;
93 } else {
94 clustererNN.mInputData_32[write_idx] = boundary_value;
95 }
96 } else {
97 CfChargePos tmp_pos(target_row, target_pad, target_time);
98 float normalized_charge = static_cast<float>(chargeMap[tmp_pos].unpack()) / central_charge;
99 if (dtype == 0) {
100 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)normalized_charge;
101 } else {
102 clustererNN.mInputData_32[write_idx] = normalized_charge;
103 }
104 }
105 // if((CAMath::Abs(static_cast<float>(clustererNN.mInputData_16_Test[write_idx]) - static_cast<float>(clustererNN.mInputData_16[write_idx])) > 1e-4) && ((glo_idx + batchStart) < clusterer.mPmemory->counters.nClusters)) {
106 // printf("Warning: Input data mismatch at index %d, %d - row, pad, time: %d, %d, %d : %f -> %f\n", glo_idx, glo_idx + batchStart, r, p, t,
107 // static_cast<float>(clustererNN.mInputData_16_Test[write_idx]), static_cast<float>(clustererNN.mInputData_16[write_idx]));
108 // }
109 write_idx++;
110 }
111 }
112 }
113
114 if (clustererNN.mNnClusterizerAddIndexData) {
115 if (dtype == 0) {
116 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)(static_cast<float>(sector) / o2::tpc::constants::MAXSECTOR);
117 clustererNN.mInputData_16[write_idx + 1] = (OrtDataType::Float16_t)(static_cast<float>(row) / o2::tpc::constants::MAXGLOBALPADROW);
118 clustererNN.mInputData_16[write_idx + 2] = (OrtDataType::Float16_t)(static_cast<float>(pad) / GPUTPCGeometry::NPads(row));
119 } else {
120 clustererNN.mInputData_32[write_idx] = static_cast<float>(sector) / o2::tpc::constants::MAXSECTOR;
121 clustererNN.mInputData_32[write_idx + 1] = static_cast<float>(row) / o2::tpc::constants::MAXGLOBALPADROW;
122 clustererNN.mInputData_32[write_idx + 2] = static_cast<float>(pad) / GPUTPCGeometry::NPads(row);
123 }
124 }
125
126 if (!clustererNN.mNnClusterizerSetDeconvolutionFlags) {
127 clustererNN.mClusterFlags[2 * glo_idx] = 0;
128 clustererNN.mClusterFlags[2 * glo_idx + 1] = 0;
129
130 for (uint16_t i = 0; i < 8; ++i) {
131 Delta2 d = cfconsts::InnerNeighbors[i];
132 CfChargePos tmp_pos = peak.delta(d);
133 clustererNN.mClusterFlags[2 * glo_idx] += CfUtils::isPeak(isPeakMap[tmp_pos]);
134 }
135 clustererNN.mClusterFlags[2 * glo_idx + 1] = clustererNN.mClusterFlags[2 * glo_idx];
136 }
137}
138
139template <>
140GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::fillInputNNGPU>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
141{
142 uint32_t glo_idx = get_global_id(0);
143 auto& clusterer = processors.tpcClusterer[sector];
144 auto& clustererNN = processors.tpcNNClusterer[sector];
145
146 if (glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode * clustererNN.mNnClusterizerRowTimeSizeFull) {
147 return;
148 }
149
150 uint32_t base_idx = glo_idx / clustererNN.mNnClusterizerRowTimeSizeFull;
151 uint32_t transient_index = glo_idx - (base_idx * clustererNN.mNnClusterizerRowTimeSizeFull);
152
153 // Early exit for out-of-bounds threads
154 if (base_idx + batchStart >= clusterer.mPmemory->counters.nClusters) {
155 return;
156 }
157 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
158 CfArray2D<uint8_t> isPeakMap(clusterer.mPpeakMap);
159
160 // Use dedicated neural network shared memory arrays for warp-level caching
161 // First thread in each warp loads shared data
162 CfChargePos peak = clusterer.mPfilteredPeakPositions[CAMath::Min(base_idx + batchStart, (uint32_t)(clusterer.mPmemory->counters.nClusters - 1))];
163 float central_charge = static_cast<float>(chargeMap[peak].unpack());
164 int32_t row = static_cast<int>(peak.row());
165 int32_t pad = static_cast<int>(peak.pad());
166 int32_t time = static_cast<int>(peak.time());
167
168 // Handle index data with fewer branches
169 if (clustererNN.mNnClusterizerAddIndexData && transient_index >= clustererNN.mNnClusterizerRowTimeSize) {
170 int32_t data_idx = transient_index - clustererNN.mNnClusterizerRowTimeSize;
171 uint32_t write_idx = base_idx * clustererNN.mNnClusterizerElementSize + clustererNN.mNnClusterizerChargeArraySize + data_idx;
172
173 float index_values[3] = {
174 static_cast<float>(sector) / o2::tpc::constants::MAXSECTOR,
175 static_cast<float>(row) / o2::tpc::constants::MAXGLOBALPADROW,
176 static_cast<float>(pad) / GPUTPCGeometry::NPads(row)};
177
178 if (dtype == 0) {
179 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)index_values[data_idx];
180 } else {
181 clustererNN.mInputData_32[write_idx] = index_values[data_idx];
182 }
183
184 // Handle deconvolution flags only once per cluster (last thread in element)
185 if (!clustererNN.mNnClusterizerSetDeconvolutionFlags && data_idx == 2) {
186 uint8_t cluster_flags = 0;
187 for (uint16_t i = 0; i < 8; i++) {
188 Delta2 d = cfconsts::InnerNeighbors[i];
189 CfChargePos tmp_pos = peak.delta(d);
190 cluster_flags += CfUtils::isPeak(isPeakMap[tmp_pos]);
191 }
192 clustererNN.mClusterFlags[2 * base_idx] = cluster_flags;
193 clustererNN.mClusterFlags[2 * base_idx + 1] = cluster_flags;
194 }
195 return;
196 }
197
198 // Main data processing - optimize index calculations
199 if (transient_index < clustererNN.mNnClusterizerRowTimeSize) {
200 // Optimize 3D index calculation
201 int32_t row_idx = transient_index / clustererNN.mNnClusterizerFullTimeSize;
202 int32_t r_local = row_idx - clustererNN.mNnClusterizerSizeInputRow;
203 int32_t time_idx = transient_index - row_idx * clustererNN.mNnClusterizerFullTimeSize;
204 int32_t t_local = time_idx - clustererNN.mNnClusterizerSizeInputTime;
205 int32_t write_idx = base_idx * clustererNN.mNnClusterizerElementSize + row_idx * clustererNN.mNnClusterizerPadTimeSize + time_idx;
206
207 // Early boundary check for row
208 int32_t target_row = row + r_local;
209 int8_t is_row_boundary = (target_row < 0) || (target_row > (o2::tpc::constants::MAXGLOBALPADROW - 1));
210
211 // Calculate offsets
212 int32_t row_offset = GPUTPCNNClusterizerKernels::rowOffset(row, clustererNN.mNnClusterizerSizeInputRow);
213 int32_t pad_offset = GPUTPCNNClusterizerKernels::padOffset(row, target_row);
214 for (int32_t p_local = -clustererNN.mNnClusterizerSizeInputPad + pad_offset; p_local <= clustererNN.mNnClusterizerSizeInputPad + pad_offset; p_local++) {
215 if (is_row_boundary) {
216 // Use boundary fill value
217 float boundary_val = static_cast<float>(clustererNN.mNnClusterizerBoundaryFillValue);
218 if (dtype == 0) {
219 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)boundary_val;
220 } else {
221 clustererNN.mInputData_32[write_idx] = boundary_val;
222 }
223 write_idx += clustererNN.mNnClusterizerFullTimeSize; // Move to next pad position
224 continue;
225 }
226
227 // Calculate target pad and time
228 int32_t target_pad = pad + p_local;
229 int32_t target_time = time + t_local;
230
231 // Optimized boundary check
232 int8_t is_boundary = GPUTPCNNClusterizerKernels::isBoundary(target_row + row_offset, target_pad, clustererNN.mNnClusterizerSizeInputRow) || (target_time < 0) || (target_time >= TPC_MAX_FRAGMENT_LEN_GPU);
233
234 float output_value;
235 if (is_boundary) {
236 output_value = static_cast<float>(clustererNN.mNnClusterizerBoundaryFillValue);
237 } else {
238 // Coalesced memory access - create position and read charge
239 CfChargePos tmp_pos(target_row, target_pad, target_time);
240 output_value = static_cast<float>(chargeMap[tmp_pos].unpack()) / central_charge; // Normalize by central charge
241 }
242
243 // Write output with reduced branching
244 if (dtype == 0) {
245 clustererNN.mInputData_16[write_idx] = (OrtDataType::Float16_t)output_value;
246 } else {
247 clustererNN.mInputData_32[write_idx] = output_value;
248 }
249 // if (write_idx >= clustererNN.mNnClusterizerElementSize * clustererNN.mNnClusterizerBatchedMode) {
250 // printf("Error: Write index out of bounds (central array)! %d >= %d (write_idx: %d, base_idx: %d, transient_index: %d, row_idx: %d, time_idx: %d, r_local: %d, t_local: %d)\n",
251 // write_idx, (int)(clustererNN.mNnClusterizerElementSize * clustererNN.mNnClusterizerBatchedMode), write_idx, base_idx, transient_index, row_idx, time_idx, r_local, t_local);
252 // }
253 // if ((clusterer.mPmemory->counters.nClusters - batchStart) < clustererNN.mNnClusterizerBatchedMode) {
254 // if (write_idx >= ((clusterer.mPmemory->counters.nClusters - batchStart) * clustererNN.mNnClusterizerElementSize)) {
255 // printf("Error: Write index out of bounds (end of array)! %d >= %d (write_idx: %d, base_idx: %d, transient_index: %d, row_idx: %d, time_idx: %d, r_local: %d, t_local: %d)\n",
256 // write_idx, (int)((clusterer.mPmemory->counters.nClusters - batchStart) * clustererNN.mNnClusterizerElementSize), write_idx, base_idx, transient_index, row_idx, time_idx, r_local, t_local);
257 // }
258 // if (write_idx > ((clusterer.mPmemory->counters.nClusters - batchStart) * clustererNN.mNnClusterizerElementSize - 5)) {
259 // printf("Sanity check (should appear only once) %d == %d (write_idx: %d, base_idx: %d, transient_index: %d, row_idx: %d, time_idx: %d, r_local: %d, t_local: %d)\n",
260 // write_idx, (int)((clusterer.mPmemory->counters.nClusters - batchStart) * clustererNN.mNnClusterizerElementSize - 4), write_idx, base_idx, transient_index, row_idx, time_idx, r_local, t_local);
261 // }
262 // }
263
264 write_idx += clustererNN.mNnClusterizerFullTimeSize; // Move to next pad position
265 }
266 }
267}
268
269template <>
270GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::determineClass1Labels>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
271{
272 uint32_t glo_idx = get_global_id(0);
273 auto& clusterer = processors.tpcClusterer[sector];
274 auto& clustererNN = processors.tpcNNClusterer[sector];
275 if (glo_idx + batchStart >= clusterer.mPmemory->counters.nClusters || glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
276 return;
277 }
278 if (clustererNN.mNnClusterizerUseClassification) {
279 if (dtype == 0) {
280 clustererNN.mOutputDataClass[glo_idx + batchStart] = (int32_t)((clustererNN.mModelProbabilities_16[glo_idx]).ToFloat() > clustererNN.mNnClassThreshold);
281 } else if (dtype == 1) {
282 clustererNN.mOutputDataClass[glo_idx + batchStart] = (int32_t)(clustererNN.mModelProbabilities_32[glo_idx] > clustererNN.mNnClassThreshold);
283 }
284 } else {
285 clustererNN.mOutputDataClass[glo_idx + batchStart] = 1;
286 }
287}
288
289template <>
290GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::determineClass2Labels>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
291{
292 uint32_t glo_idx = get_global_id(0);
293 auto& clusterer = processors.tpcClusterer[sector];
294 auto& clustererNN = processors.tpcNNClusterer[sector];
295 if (glo_idx + batchStart >= clusterer.mPmemory->counters.nClusters || glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
296 return;
297 }
298 if (clustererNN.mNnClusterizerUseClassification) {
299 uint32_t elem_iterator = glo_idx * clustererNN.mNnClusterizerModelClassNumOutputNodes;
300 float current_max_prob = 0.f; // If the neural network doesn't contain the softmax as a last layer, the outputs can range in [-infty, infty]
301 uint32_t class_label = 0;
302 for (uint32_t pIdx = elem_iterator; pIdx < elem_iterator + clustererNN.mNnClusterizerModelClassNumOutputNodes; pIdx++) {
303 if (pIdx == elem_iterator) {
304 if (dtype == 0) {
305 current_max_prob = static_cast<float>(clustererNN.mModelProbabilities_16[pIdx]);
306 } else if (dtype == 1) {
307 current_max_prob = clustererNN.mModelProbabilities_32[pIdx];
308 }
309 } else {
310 if (dtype == 0) {
311 current_max_prob = CAMath::Max(current_max_prob, clustererNN.mModelProbabilities_16[pIdx].ToFloat());
312 } else if (dtype == 1) {
313 current_max_prob = CAMath::Max(current_max_prob, clustererNN.mModelProbabilities_32[pIdx]);
314 }
315 }
316 }
317 // uint32_t class_label = std::distance(elem_iterator, std::max_element(elem_iterator, elem_iterator + clustererNN.mNnClusterizerModelClassNumOutputNodes)); // Multiple outputs of the class network are the probabilities for each class. The highest one "wins"
318 clustererNN.mOutputDataClass[glo_idx + batchStart] = class_label;
319 if (class_label > 1) {
320 clustererNN.mClusterFlags[2 * glo_idx] = 1;
321 clustererNN.mClusterFlags[2 * glo_idx + 1] = 1;
322 }
323 } else {
324 clustererNN.mOutputDataClass[glo_idx + batchStart] = 1;
325 }
326}
327
328template <>
329GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::publishClass1Regression>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
330{
331 uint32_t glo_idx = get_global_id(0);
332 auto& clusterer = processors.tpcClusterer[sector];
333 auto& clustererNN = processors.tpcNNClusterer[sector];
334 if (glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
335 return;
336 }
337
338 uint32_t maxClusterNum = clusterer.mPmemory->counters.nClusters;
339 uint32_t full_glo_idx = glo_idx + batchStart;
340 int32_t model_output_index = glo_idx * clustererNN.mNnClusterizerModelReg1NumOutputNodes;
341
342 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
343 CfChargePos peak = clusterer.mPfilteredPeakPositions[CAMath::Min(full_glo_idx, maxClusterNum - 1)];
344 float central_charge = static_cast<float>(chargeMap[peak].unpack());
345
346 CPU_ONLY(MCLabelAccumulator labelAccElem(clusterer));
347 MCLabelAccumulator* labelAcc = CPU_PTR(&labelAccElem);
348
349 if (full_glo_idx >= maxClusterNum) {
350 if (withMC) {
351 ClusterAccumulator dummy_pc;
352 CPU_ONLY(labelAcc->collect(peak, central_charge));
353 GPUTPCCFClusterizer::buildCluster(
354 clusterer.Param().rec,
355 chargeMap,
356 peak,
357 smem.posBcast,
358 smem.buf,
359 smem.innerAboveThreshold,
360 &dummy_pc,
361 labelAcc);
362 }
363 return;
364 }
365
366 tpc::ClusterNative* clusterOut = clusterer.mPclusterByRow;
367
368 // LOG(info) << glo_idx << " -- " << model_output_index << " / " << clustererNN.outputDataReg1.size() << " / " << clustererNN.mNnClusterizerModelReg1NumOutputNodes << " -- " << clusterer.peakPositions.size() << " -- " << clusterer.centralCharges.size();
369
371
372 // Publishing logic is taken from default clusterizer
373 if (withMC) {
374 ClusterAccumulator dummy_pc;
375 CPU_ONLY(labelAcc->collect(peak, central_charge));
376 GPUTPCCFClusterizer::buildCluster(
377 clusterer.Param().rec,
378 chargeMap,
379 peak,
380 smem.posBcast,
381 smem.buf,
382 smem.innerAboveThreshold,
383 &dummy_pc,
384 labelAcc);
385 }
386 if ((clusterer.mPmemory->fragment).isOverlap(peak.time())) {
387 if (clusterer.mPclusterPosInRow) {
388 clusterer.mPclusterPosInRow[full_glo_idx] = clusterer.mNMaxClusterPerRow;
389 }
390 return;
391 }
392
393 // For flag influence on cluster error setting: O2/GPU/GPUTracking/Base/GPUParam.inc -> UpdateClusterError2ByState
394 bool notSinglePad = false, notSingleTime = false;
395 for (uint16_t i = 0; i < 8; i++) {
396 Delta2 d = cfconsts::InnerNeighbors[i];
397 CfChargePos tmp_pos = peak.delta(d);
398 notSinglePad |= (d.x != 0) && (static_cast<float>(chargeMap[tmp_pos].unpack()) > 0);
399 notSingleTime |= (d.y != 0) && (static_cast<float>(chargeMap[tmp_pos].unpack()) > 0);
400 }
401
402 if (dtype == 0) {
403 pc.setFull(central_charge * clustererNN.mOutputDataReg1_16[model_output_index + 4].ToFloat(),
404 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg1_16[model_output_index].ToFloat(),
405 notSinglePad ? clustererNN.mOutputDataReg1_16[model_output_index + 2].ToFloat() : 0.f,
406 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg1_16[model_output_index + 1].ToFloat(),
407 notSingleTime ? clustererNN.mOutputDataReg1_16[model_output_index + 3].ToFloat() : 0.f,
408 clustererNN.mClusterFlags[2 * glo_idx],
409 clustererNN.mClusterFlags[2 * glo_idx + 1]);
410 } else if (dtype == 1) {
411 pc.setFull(central_charge * clustererNN.mOutputDataReg1_32[model_output_index + 4],
412 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg1_32[model_output_index],
413 notSinglePad ? clustererNN.mOutputDataReg1_32[model_output_index + 2] : 0.f,
414 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg1_32[model_output_index + 1],
415 notSingleTime ? clustererNN.mOutputDataReg1_32[model_output_index + 3] : 0.f,
416 clustererNN.mClusterFlags[2 * glo_idx],
417 clustererNN.mClusterFlags[2 * glo_idx + 1]);
418 }
419
420 tpc::ClusterNative myCluster;
421 bool rejectCluster = !pc.toNative(peak, central_charge, myCluster, clusterer.Param(), chargeMap);
422 if (clustererNN.mNnClusterizerUseClassification) {
423 rejectCluster |= (clustererNN.mOutputDataClass[CAMath::Min(full_glo_idx, (uint32_t)clusterer.mPmemory->counters.nClusters - 1)] <= 0);
424 }
425 if (rejectCluster) {
426 if (clusterer.mPclusterPosInRow) {
427 clusterer.mPclusterPosInRow[full_glo_idx] = clusterer.mNMaxClusterPerRow;
428 }
429 return;
430 }
431
432 uint32_t rowIndex = 0;
433 if (clusterOut != nullptr) {
434 rowIndex = GPUTPCCFClusterizer::sortIntoBuckets(
435 clusterer,
436 myCluster,
437 peak.row(),
438 clusterer.mNMaxClusterPerRow,
439 clusterer.mPclusterInRow,
440 clusterOut);
441 if (clusterer.mPclusterPosInRow != nullptr) {
442 clusterer.mPclusterPosInRow[full_glo_idx] = rowIndex;
443 }
444 } else if (clusterer.mPclusterPosInRow) {
445 rowIndex = clusterer.mPclusterPosInRow[full_glo_idx];
446 }
447 CPU_ONLY(labelAcc->commit(peak.row(), rowIndex, clusterer.mNMaxClusterPerRow));
448}
449
450template <>
451GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::publishClass2Regression>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint32_t batchStart)
452{
453 uint32_t glo_idx = get_global_id(0);
454 auto& clusterer = processors.tpcClusterer[sector];
455 auto& clustererNN = processors.tpcNNClusterer[sector];
456 if (glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
457 return;
458 }
459
460 uint32_t maxClusterNum = clusterer.mPmemory->counters.nClusters;
461 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
462 CfChargePos peak = clusterer.mPfilteredPeakPositions[CAMath::Min(glo_idx + batchStart, (uint32_t)(clusterer.mPmemory->counters.nClusters - 1))];
463 float central_charge = static_cast<float>(chargeMap[peak].unpack());
464
465 CPU_ONLY(MCLabelAccumulator labelAccElem(clusterer));
466 MCLabelAccumulator* labelAcc = CPU_PTR(&labelAccElem);
467 tpc::ClusterNative* clusterOut = (withMC) ? nullptr : clusterer.mPclusterByRow;
468 uint32_t full_glo_idx = glo_idx + batchStart;
469
470 if (full_glo_idx >= maxClusterNum) {
471 if (withMC) {
472 ClusterAccumulator dummy_pc;
473 CPU_ONLY(labelAcc->collect(peak, central_charge));
474 GPUTPCCFClusterizer::buildCluster(
475 clusterer.Param().rec,
476 chargeMap,
477 peak,
478 smem.posBcast,
479 smem.buf,
480 smem.innerAboveThreshold,
481 &dummy_pc,
482 labelAcc);
483 }
484 return;
485 }
486
487 uint32_t model_output_index = glo_idx * clustererNN.mNnClusterizerModelReg2NumOutputNodes;
488
490
491 if (withMC) {
492 ClusterAccumulator dummy_pc;
493 CPU_ONLY(labelAcc->collect(peak, central_charge));
494 GPUTPCCFClusterizer::buildCluster(
495 clusterer.Param().rec,
496 chargeMap,
497 peak,
498 smem.posBcast,
499 smem.buf,
500 smem.innerAboveThreshold,
501 &dummy_pc,
502 labelAcc);
503 }
504 if ((clusterer.mPmemory->fragment).isOverlap(peak.time())) {
505 if (clusterer.mPclusterPosInRow) {
506 clusterer.mPclusterPosInRow[full_glo_idx] = clusterer.mNMaxClusterPerRow;
507 }
508 return;
509 }
510
511 // Cluster 1
512 if (dtype == 0) {
513 pc.setFull(central_charge * clustererNN.mOutputDataReg2_16[model_output_index + 8].ToFloat(),
514 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg2_16[model_output_index].ToFloat(),
515 clustererNN.mOutputDataReg2_16[model_output_index + 4].ToFloat(),
516 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg2_16[model_output_index + 2].ToFloat(),
517 clustererNN.mOutputDataReg2_16[model_output_index + 6].ToFloat(),
518 clustererNN.mClusterFlags[2 * glo_idx],
519 clustererNN.mClusterFlags[2 * glo_idx + 1]);
520 } else if (dtype == 1) {
521 pc.setFull(central_charge * clustererNN.mOutputDataReg2_32[model_output_index + 8],
522 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg2_32[model_output_index],
523 clustererNN.mOutputDataReg2_32[model_output_index + 4],
524 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg2_32[model_output_index + 2],
525 clustererNN.mOutputDataReg2_32[model_output_index + 6],
526 clustererNN.mClusterFlags[2 * glo_idx],
527 clustererNN.mClusterFlags[2 * glo_idx + 1]);
528 }
529
530 tpc::ClusterNative myCluster;
531 bool rejectCluster = !pc.toNative(peak, central_charge, myCluster, clusterer.Param(), chargeMap);
532 if (clustererNN.mNnClusterizerUseClassification) {
533 rejectCluster |= (clustererNN.mOutputDataClass[CAMath::Min(full_glo_idx, (uint32_t)clusterer.mPmemory->counters.nClusters - 1)] <= 0);
534 }
535 if (rejectCluster) {
536 if (clusterer.mPclusterPosInRow) {
537 clusterer.mPclusterPosInRow[full_glo_idx] = clusterer.mNMaxClusterPerRow;
538 }
539 return;
540 }
541
542 uint32_t rowIndex = 0;
543 if (clusterOut != nullptr) {
544 rowIndex = GPUTPCCFClusterizer::sortIntoBuckets(
545 clusterer,
546 myCluster,
547 peak.row(),
548 clusterer.mNMaxClusterPerRow,
549 clusterer.mPclusterInRow,
550 clusterOut);
551 if (clusterer.mPclusterPosInRow != nullptr) {
552 clusterer.mPclusterPosInRow[full_glo_idx] = rowIndex;
553 }
554 } else if (clusterer.mPclusterPosInRow) {
555 rowIndex = clusterer.mPclusterPosInRow[full_glo_idx];
556 }
557 CPU_ONLY(labelAcc->commit(peak.row(), rowIndex, clusterer.mNMaxClusterPerRow));
558
559 // Cluster 2
560 if (dtype == 0) {
561 pc.setFull(central_charge * clustererNN.mOutputDataReg2_16[model_output_index + 9].ToFloat(),
562 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg2_16[model_output_index + 1].ToFloat(),
563 clustererNN.mOutputDataReg2_16[model_output_index + 5].ToFloat(),
564 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg2_16[model_output_index + 3].ToFloat(),
565 clustererNN.mOutputDataReg2_16[model_output_index + 7].ToFloat(),
566 clustererNN.mClusterFlags[2 * glo_idx],
567 clustererNN.mClusterFlags[2 * glo_idx + 1]);
568 } else if (dtype == 1) {
569 pc.setFull(central_charge * clustererNN.mOutputDataReg2_32[model_output_index + 9],
570 static_cast<float>(peak.pad()) + clustererNN.mOutputDataReg2_32[model_output_index + 1],
571 clustererNN.mOutputDataReg2_32[model_output_index + 5],
572 (clusterer.mPmemory->fragment).start + static_cast<float>(peak.time()) + clustererNN.mOutputDataReg2_32[model_output_index + 3],
573 clustererNN.mOutputDataReg2_32[model_output_index + 7],
574 clustererNN.mClusterFlags[2 * glo_idx],
575 clustererNN.mClusterFlags[2 * glo_idx + 1]);
576 }
577
578 rejectCluster = !pc.toNative(peak, central_charge, myCluster, clusterer.Param(), chargeMap);
579 if (clustererNN.mNnClusterizerUseClassification) {
580 rejectCluster |= (clustererNN.mOutputDataClass[CAMath::Min(full_glo_idx, (uint32_t)clusterer.mPmemory->counters.nClusters - 1)] <= 0);
581 }
582 if (rejectCluster) {
583 if (clusterer.mPclusterPosInRow) {
584 clusterer.mPclusterPosInRow[full_glo_idx] = clusterer.mNMaxClusterPerRow;
585 }
586 return;
587 }
588
589 if (clusterOut != nullptr) {
590 rowIndex = GPUTPCCFClusterizer::sortIntoBuckets(
591 clusterer,
592 myCluster,
593 peak.row(),
594 clusterer.mNMaxClusterPerRow,
595 clusterer.mPclusterInRow,
596 clusterOut);
597 if (clusterer.mPclusterPosInRow != nullptr) {
598 clusterer.mPclusterPosInRow[full_glo_idx] = rowIndex;
599 }
600 } else if (clusterer.mPclusterPosInRow) {
601 rowIndex = clusterer.mPclusterPosInRow[full_glo_idx];
602 }
603 // CPU_ONLY(labelAcc->commit(peak.row(), rowIndex, clusterer.mNMaxClusterPerRow)); // -> Is this needed? How to handle MC labels for split clusters?
604}
605
606// ---------------------------------
607template <>
608GPUdii() void GPUTPCNNClusterizerKernels::Thread<GPUTPCNNClusterizerKernels::publishDeconvolutionFlags>(int32_t nBlocks, int32_t nThreads, int32_t iBlock, int32_t iThread, GPUSharedMemory& smem, processorType& processors, uint8_t sector, int8_t dtype, int8_t withMC, uint batchStart)
609{
610 // Implements identical publishing logic as the heuristic clusterizer and deconvolution kernel
611 uint32_t glo_idx = get_global_id(0);
612 auto& clusterer = processors.tpcClusterer[sector];
613 auto& clustererNN = processors.tpcNNClusterer[sector];
614 if (glo_idx + batchStart >= clusterer.mPmemory->counters.nClusters || glo_idx >= (uint32_t)clustererNN.mNnClusterizerBatchedMode) {
615 return;
616 }
617 CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
618 CfChargePos peak = clusterer.mPfilteredPeakPositions[glo_idx + batchStart];
619
620 clustererNN.mClusterFlags[2 * glo_idx] = 0;
621 clustererNN.mClusterFlags[2 * glo_idx + 1] = 0;
622 for (int i = 0; i < 8; i++) {
623 Delta2 d = cfconsts::InnerNeighbors[i];
624 CfChargePos tmp_pos = peak.delta(d);
625 PackedCharge charge = chargeMap[tmp_pos];
626 clustererNN.mClusterFlags[2 * glo_idx] += (d.y != 0 && charge.isSplit());
627 clustererNN.mClusterFlags[2 * glo_idx + 1] += (d.x != 0 && charge.isSplit());
628 }
629 for (int i = 0; i < 16; i++) {
630 Delta2 d = cfconsts::OuterNeighbors[i];
631 CfChargePos tmp_pos = peak.delta(d);
632 PackedCharge charge = chargeMap[tmp_pos];
633 clustererNN.mClusterFlags[2 * glo_idx] += (d.y != 0 && charge.isSplit() && !charge.has3x3Peak());
634 clustererNN.mClusterFlags[2 * glo_idx + 1] += (d.x != 0 && charge.isSplit() && !charge.has3x3Peak());
635 }
636}
637
638// THe following arithmetic is done because the network is trained with a split between IROC and OROC boundary
639GPUd() int32_t GPUTPCNNClusterizerKernels::padOffset(int32_t row_ref, int32_t row_current)
640{
641 if (row_current < 0 || row_current >= o2::tpc::constants::MAXGLOBALPADROW) {
642 return 0; // Short-circuit for out-of-bound rows
643 } else {
644 return (int)((GPUTPCGeometry::NPads(row_current) - GPUTPCGeometry::NPads(row_ref)) / 2);
645 }
646}
647
648GPUd() int32_t GPUTPCNNClusterizerKernels::rowOffset(int32_t row, int32_t offset)
649{
650 return (row > 62 ? offset : 0);
651}
652
653GPUd() bool GPUTPCNNClusterizerKernels::isBoundary(int32_t row, int32_t pad, int32_t offset)
654{
655 if (pad < 0 || row < 0) { // Faster short-circuit
656 return true;
657 } else if (row < 63) {
658 return (pad >= static_cast<int>(GPUTPCGeometry::NPads(row)));
659 } else if (row < (63 + offset)) { // to account for the gap between IROC and OROC. Charge will be set to the boundary fill value in order to signal boundaries to the neural network
660 return true;
662 return (pad >= static_cast<int>(GPUTPCGeometry::NPads(row - offset)));
663 } else {
664 return true;
665 }
666}
int16_t charge
Definition RawEventData.h:5
int16_t time
Definition RawEventData.h:4
int32_t i
#define get_local_size(dim)
#define get_local_id(dim)
#define get_num_groups(dim)
#define get_global_id(dim)
#define get_group_id(dim)
GPUdii() void GPUTPCNNClusterizerKernels
GPUd() int32_t GPUTPCNNClusterizerKernels
void collect(const CfChargePos &, tpccf::Charge)
void commit(tpccf::Row, uint32_t, uint32_t)
#define TPC_MAX_FRAGMENT_LEN_GPU
#define CPU_ONLY(x)
#define CPU_PTR(x)
GLintptr offset
Definition glcorearb.h:660
typedef void(APIENTRYP PFNGLCULLFACEPROC)(GLenum mode)
GLboolean r
Definition glcorearb.h:1233
GLuint start
Definition glcorearb.h:469
constexpr int MAXSECTOR
Definition Constants.h:28
constexpr int MAXGLOBALPADROW
Definition Constants.h:34
int16_t y
int16_t x
std::vector< int > row