File indexing completed on 2025-01-22 07:33:56
0001 #include <alpaka/alpaka.hpp>
0002
0003 #include "DataFormats/VertexSoA/interface/ZVertexDevice.h"
0004 #include "DataFormats/VertexSoA/interface/ZVertexHost.h"
0005 #include "DataFormats/VertexSoA/interface/alpaka/ZVertexSoACollection.h"
0006 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0007 #include "HeterogeneousCore/AlpakaInterface/interface/workdivision.h"
0008
0009 namespace ALPAKA_ACCELERATOR_NAMESPACE::testZVertexSoAT {
0010
0011 class TestFillKernel {
0012 public:
0013 ALPAKA_FN_ACC void operator()(Acc1D const& acc,
0014 reco::ZVertexSoAView zvertex_view,
0015 reco::ZVertexTracksSoAView ztracks_view) const {
0016 if (cms::alpakatools::once_per_grid(acc)) {
0017 zvertex_view.nvFinal() = 420;
0018 }
0019
0020 for (int32_t j : cms::alpakatools::uniform_elements(acc, zvertex_view.metadata().size())) {
0021 zvertex_view[j].zv() = (float)j;
0022 zvertex_view[j].wv() = (float)j;
0023 zvertex_view[j].chi2() = (float)j;
0024 zvertex_view[j].ptv2() = (float)j;
0025 zvertex_view[j].sortInd() = (uint16_t)j;
0026 }
0027 for (int32_t j : cms::alpakatools::uniform_elements(acc, ztracks_view.metadata().size())) {
0028 ztracks_view[j].idv() = (int16_t)j;
0029 ztracks_view[j].ndof() = (int32_t)j;
0030 }
0031 }
0032 };
0033
0034 class TestVerifyKernel {
0035 public:
0036 ALPAKA_FN_ACC void operator()(Acc1D const& acc,
0037 reco::ZVertexSoAView zvertex_view,
0038 reco::ZVertexTracksSoAView ztracks_view) const {
0039 if (cms::alpakatools::once_per_grid(acc)) {
0040 ALPAKA_ASSERT_ACC(zvertex_view.nvFinal() == 420);
0041 }
0042
0043 for (int32_t j : cms::alpakatools::uniform_elements(acc, zvertex_view.nvFinal())) {
0044 ALPAKA_ASSERT(zvertex_view[j].zv() - (float)j < 0.0001);
0045 ALPAKA_ASSERT(zvertex_view[j].wv() - (float)j < 0.0001);
0046 ALPAKA_ASSERT(zvertex_view[j].chi2() - (float)j < 0.0001);
0047 ALPAKA_ASSERT(zvertex_view[j].ptv2() - (float)j < 0.0001);
0048 ALPAKA_ASSERT(zvertex_view[j].sortInd() == uint32_t(j));
0049 }
0050 for (int32_t j : cms::alpakatools::uniform_elements(acc, ztracks_view.metadata().size())) {
0051 ALPAKA_ASSERT(ztracks_view[j].idv() == j);
0052 ALPAKA_ASSERT(ztracks_view[j].ndof() == j);
0053 }
0054 }
0055 };
0056
0057 void runKernels(reco::ZVertexSoAView zvertex_view, reco::ZVertexTracksSoAView ztracks_view, Queue& queue) {
0058 uint32_t items = 64;
0059 uint32_t groups = cms::alpakatools::divide_up_by(zvertex_view.metadata().size(), items);
0060 auto workDiv = cms::alpakatools::make_workdiv<Acc1D>(groups, items);
0061 alpaka::exec<Acc1D>(queue, workDiv, TestFillKernel{}, zvertex_view, ztracks_view);
0062 alpaka::exec<Acc1D>(queue, workDiv, TestVerifyKernel{}, zvertex_view, ztracks_view);
0063 }
0064
0065 }