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

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

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