CMS 3D CMS Logo

List of all members | Public Types | Public Member Functions | Private Attributes | Static Private Attributes
ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits > Class Template Reference

#include <HelixFit.h>

Public Types

using HitConstView = TrackingRecHitSoAConstView< TrackerTraits >
 
using HitView = TrackingRecHitSoAView< TrackerTraits >
 
using OutputSoAView = reco::TrackSoAView< TrackerTraits >
 
using ParamsOnDevice = pixelCPEforDevice::ParamsOnDeviceT< TrackerTraits >
 
using TrackingRecHitSoAs = TrackingRecHitSoA< TrackerTraits >
 
using TupleMultiplicity = caStructures::TupleMultiplicityT< TrackerTraits >
 
using Tuples = typename reco::TrackSoA< TrackerTraits >::HitContainer
 

Public Member Functions

void allocate (TupleMultiplicity const *tupleMultiplicity, OutputSoAView &helix_fit_results)
 
void deallocate ()
 
 HelixFit (float bf, bool fitNas4)
 
void launchBrokenLineKernels (const HitConstView &hv, ParamsOnDevice const *cpeParams, uint32_t nhits, uint32_t maxNumberOfTuples, Queue &queue)
 
void launchRiemannKernels (const HitConstView &hv, ParamsOnDevice const *cpeParams, uint32_t nhits, uint32_t maxNumberOfTuples, Queue &queue)
 
void setBField (double bField)
 
 ~HelixFit ()
 

Private Attributes

float bField_
 
const bool fitNas4_
 
OutputSoAView outputSoa_
 
TupleMultiplicity const * tupleMultiplicity_ = nullptr
 
Tuples const * tuples_ = nullptr
 

Static Private Attributes

static constexpr uint32_t maxNumberOfConcurrentFits_ = riemannFit::maxNumberOfConcurrentFits
 

Detailed Description

template<typename TrackerTraits>
class ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >

Definition at line 54 of file HelixFit.h.

Member Typedef Documentation

◆ HitConstView

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::HitConstView = TrackingRecHitSoAConstView<TrackerTraits>

Definition at line 59 of file HelixFit.h.

◆ HitView

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::HitView = TrackingRecHitSoAView<TrackerTraits>

Definition at line 58 of file HelixFit.h.

◆ OutputSoAView

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::OutputSoAView = reco::TrackSoAView<TrackerTraits>

Definition at line 62 of file HelixFit.h.

◆ ParamsOnDevice

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::ParamsOnDevice = pixelCPEforDevice::ParamsOnDeviceT<TrackerTraits>

Definition at line 66 of file HelixFit.h.

◆ TrackingRecHitSoAs

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::TrackingRecHitSoAs = TrackingRecHitSoA<TrackerTraits>

Definition at line 56 of file HelixFit.h.

◆ TupleMultiplicity

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::TupleMultiplicity = caStructures::TupleMultiplicityT<TrackerTraits>

Definition at line 64 of file HelixFit.h.

◆ Tuples

template<typename TrackerTraits >
using ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::Tuples = typename reco::TrackSoA<TrackerTraits>::HitContainer

Definition at line 61 of file HelixFit.h.

Constructor & Destructor Documentation

◆ HelixFit()

template<typename TrackerTraits >
ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::HelixFit ( float  bf,
bool  fitNas4 
)
inlineexplicit

◆ ~HelixFit()

template<typename TrackerTraits >
ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::~HelixFit ( )
inline

Member Function Documentation

◆ allocate()

template<typename TrackerTraits >
void ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::allocate ( TupleMultiplicity const *  tupleMultiplicity,
OutputSoAView helix_fit_results 
)

Definition at line 6 of file HelixFit.cc.

References ALPAKA_ACCELERATOR_NAMESPACE::caPixelDoublets::ALPAKA_ASSERT_ACC(), and tupleMultiplicity.

Referenced by ALPAKA_ACCELERATOR_NAMESPACE::CAHitNtupletGenerator< TrackerTraits >::makeTuplesAsync().

6  {
7  tuples_ = &helix_fit_results.hitIndices();
9  outputSoa_ = helix_fit_results;
10 
13  }
TupleMultiplicity const * tupleMultiplicity_
Definition: HelixFit.h:91
TupleMultiplicity< TrackerTraits > const *__restrict__ tupleMultiplicity

◆ deallocate()

template<typename TrackerTraits >
void ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::deallocate ( )

Definition at line 16 of file HelixFit.cc.

Referenced by ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::~HelixFit().

16 {}

◆ launchBrokenLineKernels()

template<typename TrackerTraits >
void ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::launchBrokenLineKernels ( const HitConstView hv,
ParamsOnDevice const *  cpeParams,
uint32_t  nhits,
uint32_t  maxNumberOfTuples,
Queue &  queue 
)

Definition at line 246 of file BrokenLineFit.dev.cc.

References ALPAKA_ACCELERATOR_NAMESPACE::caPixelDoublets::ALPAKA_ASSERT_ACC(), cms::alpakatools::divide_up_by(), mps_fire::i, hltrates_dqm_sourceclient-live_cfg::offset, and createBeamHaloJobs::queue.

Referenced by ALPAKA_ACCELERATOR_NAMESPACE::CAHitNtupletGenerator< TrackerTraits >::makeTuplesAsync().

251  {
253 
254  uint32_t blockSize = 64;
255  uint32_t numberOfBlocks = cms::alpakatools::divide_up_by(maxNumberOfConcurrentFits_, blockSize);
256  const WorkDiv1D workDivTriplets = cms::alpakatools::make_workdiv<Acc1D>(numberOfBlocks, blockSize);
257  const WorkDiv1D workDivQuadsPenta = cms::alpakatools::make_workdiv<Acc1D>(numberOfBlocks / 4, blockSize);
258 
259  // Fit internals
260  auto tkidDevice =
261  cms::alpakatools::make_device_buffer<typename TrackerTraits::tindex_type[]>(queue, maxNumberOfConcurrentFits_);
262  auto hitsDevice = cms::alpakatools::make_device_buffer<double[]>(
263  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Matrix3xNd<6>) / sizeof(double));
264  auto hits_geDevice = cms::alpakatools::make_device_buffer<float[]>(
265  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Matrix6xNf<6>) / sizeof(float));
266  auto fast_fit_resultsDevice = cms::alpakatools::make_device_buffer<double[]>(
267  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Vector4d) / sizeof(double));
268 
269  for (uint32_t offset = 0; offset < maxNumberOfTuples; offset += maxNumberOfConcurrentFits_) {
270  // fit triplets
271 
272  alpaka::exec<Acc1D>(queue,
273  workDivTriplets,
274  Kernel_BLFastFit<3, TrackerTraits>{},
275  tuples_,
277  hv,
278  cpeParams,
279  tkidDevice.data(),
280  hitsDevice.data(),
281  hits_geDevice.data(),
282  fast_fit_resultsDevice.data(),
283  3,
284  3,
285  offset);
286 
287  alpaka::exec<Acc1D>(queue,
288  workDivTriplets,
289  Kernel_BLFit<3, TrackerTraits>{},
291  bField_,
292  outputSoa_,
293  tkidDevice.data(),
294  hitsDevice.data(),
295  hits_geDevice.data(),
296  fast_fit_resultsDevice.data());
297 
298  if (fitNas4_) {
299  // fit all as 4
300  riemannFit::rolling_fits<4, TrackerTraits::maxHitsOnTrack, 1>([this,
301  &hv,
302  &cpeParams,
303  &tkidDevice,
304  &hitsDevice,
305  &hits_geDevice,
306  &fast_fit_resultsDevice,
307  &offset,
308  &queue,
309  &workDivQuadsPenta](auto i) {
310  alpaka::exec<Acc1D>(queue,
311  workDivQuadsPenta,
312  Kernel_BLFastFit<4, TrackerTraits>{},
313  tuples_,
315  hv,
316  cpeParams,
317  tkidDevice.data(),
318  hitsDevice.data(),
319  hits_geDevice.data(),
320  fast_fit_resultsDevice.data(),
321  4,
322  4,
323  offset);
324 
325  alpaka::exec<Acc1D>(queue,
326  workDivQuadsPenta,
327  Kernel_BLFit<4, TrackerTraits>{},
329  bField_,
330  outputSoa_,
331  tkidDevice.data(),
332  hitsDevice.data(),
333  hits_geDevice.data(),
334  fast_fit_resultsDevice.data());
335  });
336 
337  } else {
338  riemannFit::rolling_fits<4, TrackerTraits::maxHitsOnTrackForFullFit, 1>([this,
339  &hv,
340  &cpeParams,
341  &tkidDevice,
342  &hitsDevice,
343  &hits_geDevice,
344  &fast_fit_resultsDevice,
345  &offset,
346  &queue,
347  &workDivQuadsPenta](auto i) {
348  alpaka::exec<Acc1D>(queue,
349  workDivQuadsPenta,
350  Kernel_BLFastFit<i, TrackerTraits>{},
351  tuples_,
353  hv,
354  cpeParams,
355  tkidDevice.data(),
356  hitsDevice.data(),
357  hits_geDevice.data(),
358  fast_fit_resultsDevice.data(),
359  i,
360  i,
361  offset);
362 
363  alpaka::exec<Acc1D>(queue,
364  workDivQuadsPenta,
365  Kernel_BLFit<i, TrackerTraits>{},
367  bField_,
368  outputSoa_,
369  tkidDevice.data(),
370  hitsDevice.data(),
371  hits_geDevice.data(),
372  fast_fit_resultsDevice.data());
373  });
374 
375  static_assert(TrackerTraits::maxHitsOnTrackForFullFit < TrackerTraits::maxHitsOnTrack);
376 
377  //Fit all the rest using the maximum from previous call
378  alpaka::exec<Acc1D>(queue,
379  workDivQuadsPenta,
380  Kernel_BLFastFit<TrackerTraits::maxHitsOnTrackForFullFit, TrackerTraits>{},
381  tuples_,
383  hv,
384  cpeParams,
385  tkidDevice.data(),
386  hitsDevice.data(),
387  hits_geDevice.data(),
388  fast_fit_resultsDevice.data(),
389  TrackerTraits::maxHitsOnTrackForFullFit,
390  TrackerTraits::maxHitsOnTrack - 1,
391  offset);
392 
393  alpaka::exec<Acc1D>(queue,
394  workDivQuadsPenta,
395  Kernel_BLFit<TrackerTraits::maxHitsOnTrackForFullFit, TrackerTraits>{},
397  bField_,
398  outputSoa_,
399  tkidDevice.data(),
400  hitsDevice.data(),
401  hits_geDevice.data(),
402  fast_fit_resultsDevice.data());
403  }
404 
405  } // loop on concurrent fits
406  }
Eigen::Vector4d Vector4d
Definition: FitResult.h:12
constexpr Idx divide_up_by(Idx value, Idx divisor)
Definition: workdivision.h:20
TupleMultiplicity const * tupleMultiplicity_
Definition: HelixFit.h:91
Eigen::Matrix< float, 6, N > Matrix6xNf
Definition: HelixFit.h:35
WorkDiv< Dim1D > WorkDiv1D
Definition: config.h:32
Eigen::Matrix< double, 3, N > Matrix3xNd
Definition: HelixFit.h:30
static constexpr uint32_t maxNumberOfConcurrentFits_
Definition: HelixFit.h:87

◆ launchRiemannKernels()

template<typename TrackerTraits >
void ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::launchRiemannKernels ( const HitConstView hv,
ParamsOnDevice const *  cpeParams,
uint32_t  nhits,
uint32_t  maxNumberOfTuples,
Queue &  queue 
)

Definition at line 213 of file RiemannFit.dev.cc.

References cms::cuda::assert(), hltrates_dqm_sourceclient-live_cfg::offset, and createBeamHaloJobs::queue.

Referenced by ALPAKA_ACCELERATOR_NAMESPACE::CAHitNtupletGenerator< TrackerTraits >::makeTuplesAsync().

217  {
218  assert(tuples_);
219 
220  auto blockSize = 64;
221  auto numberOfBlocks = (maxNumberOfConcurrentFits_ + blockSize - 1) / blockSize;
222  const auto workDivTriplets = cms::alpakatools::make_workdiv<Acc1D>(numberOfBlocks, blockSize);
223  const auto workDivQuadsPenta = cms::alpakatools::make_workdiv<Acc1D>(numberOfBlocks / 4, blockSize);
224 
225  // Fit internals
226  auto hitsDevice = cms::alpakatools::make_device_buffer<double[]>(
227  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Matrix3xNd<4>) / sizeof(double));
228  auto hits_geDevice = cms::alpakatools::make_device_buffer<float[]>(
229  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Matrix6x4f) / sizeof(float));
230  auto fast_fit_resultsDevice = cms::alpakatools::make_device_buffer<double[]>(
231  queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::Vector4d) / sizeof(double));
232  auto circle_fit_resultsDevice_holder =
233  cms::alpakatools::make_device_buffer<char[]>(queue, maxNumberOfConcurrentFits_ * sizeof(riemannFit::CircleFit));
234  riemannFit::CircleFit *circle_fit_resultsDevice_ =
235  (riemannFit::CircleFit *)(circle_fit_resultsDevice_holder.data());
236 
237  for (uint32_t offset = 0; offset < maxNumberOfTuples; offset += maxNumberOfConcurrentFits_) {
238  // triplets
239  alpaka::exec<Acc1D>(queue,
240  workDivTriplets,
241  Kernel_FastFit<3, TrackerTraits>{},
242  tuples_,
244  3,
245  hv,
246  cpeParams,
247  hitsDevice.data(),
248  hits_geDevice.data(),
249  fast_fit_resultsDevice.data(),
250  offset);
251 
252  alpaka::exec<Acc1D>(queue,
253  workDivTriplets,
254  Kernel_CircleFit<3, TrackerTraits>{},
256  3,
257  bField_,
258  hitsDevice.data(),
259  hits_geDevice.data(),
260  fast_fit_resultsDevice.data(),
261  circle_fit_resultsDevice_,
262  offset);
263 
264  alpaka::exec<Acc1D>(queue,
265  workDivTriplets,
266  Kernel_LineFit<3, TrackerTraits>{},
268  3,
269  bField_,
270  outputSoa_,
271  hitsDevice.data(),
272  hits_geDevice.data(),
273  fast_fit_resultsDevice.data(),
274  circle_fit_resultsDevice_,
275  offset);
276 
277  // quads
278  alpaka::exec<Acc1D>(queue,
279  workDivQuadsPenta,
280  Kernel_FastFit<4, TrackerTraits>{},
281  tuples_,
283  4,
284  hv,
285  cpeParams,
286  hitsDevice.data(),
287  hits_geDevice.data(),
288  fast_fit_resultsDevice.data(),
289  offset);
290 
291  alpaka::exec<Acc1D>(queue,
292  workDivQuadsPenta,
293  Kernel_CircleFit<4, TrackerTraits>{},
295  4,
296  bField_,
297  hitsDevice.data(),
298  hits_geDevice.data(),
299  fast_fit_resultsDevice.data(),
300  circle_fit_resultsDevice_,
301  offset);
302 
303  alpaka::exec<Acc1D>(queue,
304  workDivQuadsPenta,
305  Kernel_LineFit<4, TrackerTraits>{},
307  4,
308  bField_,
309  outputSoa_,
310  hitsDevice.data(),
311  hits_geDevice.data(),
312  fast_fit_resultsDevice.data(),
313  circle_fit_resultsDevice_,
314  offset);
315 
316  if (fitNas4_) {
317  // penta
318  alpaka::exec<Acc1D>(queue,
319  workDivQuadsPenta,
320  Kernel_FastFit<4, TrackerTraits>{},
321  tuples_,
323  5,
324  hv,
325  cpeParams,
326  hitsDevice.data(),
327  hits_geDevice.data(),
328  fast_fit_resultsDevice.data(),
329  offset);
330 
331  alpaka::exec<Acc1D>(queue,
332  workDivQuadsPenta,
333  Kernel_CircleFit<4, TrackerTraits>{},
335  5,
336  bField_,
337  hitsDevice.data(),
338  hits_geDevice.data(),
339  fast_fit_resultsDevice.data(),
340  circle_fit_resultsDevice_,
341  offset);
342 
343  alpaka::exec<Acc1D>(queue,
344  workDivQuadsPenta,
345  Kernel_LineFit<4, TrackerTraits>{},
347  5,
348  bField_,
349  outputSoa_,
350  hitsDevice.data(),
351  hits_geDevice.data(),
352  fast_fit_resultsDevice.data(),
353  circle_fit_resultsDevice_,
354  offset);
355  } else {
356  // penta all 5
357  alpaka::exec<Acc1D>(queue,
358  workDivQuadsPenta,
359  Kernel_FastFit<5, TrackerTraits>{},
360  tuples_,
362  5,
363  hv,
364  cpeParams,
365  hitsDevice.data(),
366  hits_geDevice.data(),
367  fast_fit_resultsDevice.data(),
368  offset);
369 
370  alpaka::exec<Acc1D>(queue,
371  workDivQuadsPenta,
372  Kernel_CircleFit<5, TrackerTraits>{},
374  5,
375  bField_,
376  hitsDevice.data(),
377  hits_geDevice.data(),
378  fast_fit_resultsDevice.data(),
379  circle_fit_resultsDevice_,
380  offset);
381 
382  alpaka::exec<Acc1D>(queue,
383  workDivQuadsPenta,
384  Kernel_LineFit<5, TrackerTraits>{},
386  5,
387  bField_,
388  outputSoa_,
389  hitsDevice.data(),
390  hits_geDevice.data(),
391  fast_fit_resultsDevice.data(),
392  circle_fit_resultsDevice_,
393  offset);
394  }
395  }
396  }
Eigen::Vector4d Vector4d
Definition: FitResult.h:12
TupleMultiplicity const * tupleMultiplicity_
Definition: HelixFit.h:91
Eigen::Matrix< float, 6, 4 > Matrix6x4f
Definition: HelixFit.h:25
assert(be >=bs)
Eigen::Matrix< double, 3, N > Matrix3xNd
Definition: HelixFit.h:30
static constexpr uint32_t maxNumberOfConcurrentFits_
Definition: HelixFit.h:87

◆ setBField()

template<typename TrackerTraits >
void ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::setBField ( double  bField)
inline

Member Data Documentation

◆ bField_

template<typename TrackerTraits >
float ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::bField_
private

◆ fitNas4_

template<typename TrackerTraits >
const bool ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::fitNas4_
private

Definition at line 95 of file HelixFit.h.

◆ maxNumberOfConcurrentFits_

template<typename TrackerTraits >
constexpr uint32_t ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::maxNumberOfConcurrentFits_ = riemannFit::maxNumberOfConcurrentFits
staticprivate

Definition at line 87 of file HelixFit.h.

◆ outputSoa_

template<typename TrackerTraits >
OutputSoAView ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::outputSoa_
private

Definition at line 92 of file HelixFit.h.

◆ tupleMultiplicity_

template<typename TrackerTraits >
TupleMultiplicity const* ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::tupleMultiplicity_ = nullptr
private

Definition at line 91 of file HelixFit.h.

◆ tuples_

template<typename TrackerTraits >
Tuples const* ALPAKA_ACCELERATOR_NAMESPACE::HelixFit< TrackerTraits >::tuples_ = nullptr
private

Definition at line 90 of file HelixFit.h.