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 244 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().

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

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