Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:15:42

0001 #include <cassert>
0002 
0003 #include "DataFormats/PortableTestObjects/interface/TestHostCollection.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/global/EDAnalyzer.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0012 #include "FWCore/Utilities/interface/EDGetToken.h"
0013 #include "FWCore/Utilities/interface/InputTag.h"
0014 
0015 namespace {
0016 
0017   template <typename T>
0018   class Column {
0019   public:
0020     Column(T const* data, size_t size) : data_(data), size_(size) {}
0021 
0022     void print(std::ostream& out) const {
0023       std::stringstream buffer;
0024       buffer << "{ ";
0025       if (size_ > 0) {
0026         buffer << data_[0];
0027       }
0028       if (size_ > 1) {
0029         buffer << ", " << data_[1];
0030       }
0031       if (size_ > 2) {
0032         buffer << ", " << data_[2];
0033       }
0034       if (size_ > 3) {
0035         buffer << ", ...";
0036       }
0037       buffer << '}';
0038       out << buffer.str();
0039     }
0040 
0041   private:
0042     T const* const data_;
0043     size_t const size_;
0044   };
0045 
0046   template <typename T>
0047   std::ostream& operator<<(std::ostream& out, Column<T> const& column) {
0048     column.print(out);
0049     return out;
0050   }
0051 
0052   template <typename T>
0053   void checkViewAddresses(T const& view) {
0054     // columns
0055     assert(view.metadata().addressOf_x() == view.x());
0056     assert(view.metadata().addressOf_x() == &view.x(0));
0057     assert(view.metadata().addressOf_x() == &view[0].x());
0058     assert(view.metadata().addressOf_y() == view.y());
0059     assert(view.metadata().addressOf_y() == &view.y(0));
0060     assert(view.metadata().addressOf_y() == &view[0].y());
0061     assert(view.metadata().addressOf_z() == view.z());
0062     assert(view.metadata().addressOf_z() == &view.z(0));
0063     assert(view.metadata().addressOf_z() == &view[0].z());
0064     assert(view.metadata().addressOf_id() == view.id());
0065     assert(view.metadata().addressOf_id() == &view.id(0));
0066     assert(view.metadata().addressOf_id() == &view[0].id());
0067     // scalars
0068     assert(view.metadata().addressOf_r() == &view.r());
0069     //assert(view.metadata().addressOf_r() == &view.r(0));                  // cannot access a scalar with an index
0070     //assert(view.metadata().addressOf_r() == &view[0].r());                // cannot access a scalar via a SoA row-like accessor
0071     // columns of arrays
0072     assert(view.metadata().addressOf_flags() == view.flags());
0073     assert(view.metadata().addressOf_flags() == &view.flags(0));
0074     assert(view.metadata().addressOf_flags() == &view[0].flags());
0075     // columns of Eigen matrices
0076     assert(view.metadata().addressOf_m() == view.m());
0077     assert(view.metadata().addressOf_m() == &view.m(0).coeffRef(0, 0));
0078     assert(view.metadata().addressOf_m() == &view[0].m().coeffRef(0, 0));
0079   }
0080 
0081   template <typename T>
0082   void checkViewAddresses2(T const& view) {
0083     assert(view.metadata().addressOf_x2() == view.x2());
0084     assert(view.metadata().addressOf_x2() == &view.x2(0));
0085     assert(view.metadata().addressOf_x2() == &view[0].x2());
0086     assert(view.metadata().addressOf_y2() == view.y2());
0087     assert(view.metadata().addressOf_y2() == &view.y2(0));
0088     assert(view.metadata().addressOf_y2() == &view[0].y2());
0089     assert(view.metadata().addressOf_z2() == view.z2());
0090     assert(view.metadata().addressOf_z2() == &view.z2(0));
0091     assert(view.metadata().addressOf_z2() == &view[0].z2());
0092     assert(view.metadata().addressOf_id2() == view.id2());
0093     assert(view.metadata().addressOf_id2() == &view.id2(0));
0094     assert(view.metadata().addressOf_id2() == &view[0].id2());
0095     assert(view.metadata().addressOf_m2() == view.m2());
0096     assert(view.metadata().addressOf_m2() == &view.m2(0).coeffRef(0, 0));
0097     assert(view.metadata().addressOf_m2() == &view[0].m2().coeffRef(0, 0));
0098     assert(view.metadata().addressOf_r2() == &view.r2());
0099     //assert(view.metadata().addressOf_r2() == &view.r2(0));                  // cannot access a scalar with an index
0100     //assert(view.metadata().addressOf_r2() == &view[0].r2());                // cannot access a scalar via a SoA row-like accessor
0101   }
0102 
0103   template <typename T>
0104   void checkViewAddresses3(T const& view) {
0105     assert(view.metadata().addressOf_x3() == view.x3());
0106     assert(view.metadata().addressOf_x3() == &view.x3(0));
0107     assert(view.metadata().addressOf_x3() == &view[0].x3());
0108     assert(view.metadata().addressOf_y3() == view.y3());
0109     assert(view.metadata().addressOf_y3() == &view.y3(0));
0110     assert(view.metadata().addressOf_y3() == &view[0].y3());
0111     assert(view.metadata().addressOf_z3() == view.z3());
0112     assert(view.metadata().addressOf_z3() == &view.z3(0));
0113     assert(view.metadata().addressOf_z3() == &view[0].z3());
0114     assert(view.metadata().addressOf_id3() == view.id3());
0115     assert(view.metadata().addressOf_id3() == &view.id3(0));
0116     assert(view.metadata().addressOf_id3() == &view[0].id3());
0117     assert(view.metadata().addressOf_m3() == view.m3());
0118     assert(view.metadata().addressOf_m3() == &view.m3(0).coeffRef(0, 0));
0119     assert(view.metadata().addressOf_m3() == &view[0].m3().coeffRef(0, 0));
0120     assert(view.metadata().addressOf_r3() == &view.r3());
0121     //assert(view.metadata().addressOf_r3() == &view.r3(0));                  // cannot access a scalar with an index
0122     //assert(view.metadata().addressOf_r3() == &view[0].r3());                // cannot access a scalar via a SoA row-like accessor
0123   }
0124 
0125 }  // namespace
0126 
0127 class TestAlpakaAnalyzer : public edm::global::EDAnalyzer<> {
0128 public:
0129   TestAlpakaAnalyzer(edm::ParameterSet const& config)
0130       : source_{config.getParameter<edm::InputTag>("source")},
0131         token_{consumes(source_)},
0132         //tokenMulti_{consumes(source_)},
0133         tokenMulti2_{consumes(source_)},
0134         tokenMulti3_{consumes(source_)},
0135         expectSize_{config.getParameter<int>("expectSize")},
0136         expectXvalues_{config.getParameter<std::vector<double>>("expectXvalues")} {
0137     if (std::string const& eb = config.getParameter<std::string>("expectBackend"); not eb.empty()) {
0138       expectBackend_ = cms::alpakatools::toBackend(eb);
0139       backendToken_ = consumes(edm::InputTag(source_.label(), "backend", source_.process()));
0140     }
0141   }
0142 
0143   void analyze(edm::StreamID sid, edm::Event const& event, edm::EventSetup const&) const override {
0144     portabletest::TestHostCollection const& product = event.get(token_);
0145     auto const& view = product.const_view();
0146     auto& mview = product.view();
0147     auto const& cmview = product.view();
0148 
0149     if (expectSize_ >= 0 and expectSize_ != view.metadata().size()) {
0150       throw cms::Exception("Assert") << "Expected input collection size " << expectSize_ << ", got "
0151                                      << view.metadata().size();
0152     }
0153 
0154     {
0155       edm::LogInfo msg("TestAlpakaAnalyzer");
0156       msg << source_.encode() << ".size() = " << view.metadata().size() << '\n';
0157       msg << "  data  @ " << product.buffer().data() << ",\n"
0158           << "  x     @ " << view.metadata().addressOf_x() << " = " << Column(view.x(), view.metadata().size()) << ",\n"
0159           << "  y     @ " << view.metadata().addressOf_y() << " = " << Column(view.y(), view.metadata().size()) << ",\n"
0160           << "  z     @ " << view.metadata().addressOf_z() << " = " << Column(view.z(), view.metadata().size()) << ",\n"
0161           << "  id    @ " << view.metadata().addressOf_id() << " = " << Column(view.id(), view.metadata().size())
0162           << ",\n"
0163           << "  r     @ " << view.metadata().addressOf_r() << " = " << view.r() << '\n'
0164           << "  flags @ " << view.metadata().addressOf_flags() << " = " << Column(view.flags(), view.metadata().size())
0165           << ",\n"
0166           << "  m     @ " << view.metadata().addressOf_m() << " = { ... {" << view[1].m()(1, Eigen::indexing::all)
0167           << " } ... } \n";
0168       msg << std::hex << "  [y - x] = 0x"
0169           << reinterpret_cast<intptr_t>(view.metadata().addressOf_y()) -
0170                  reinterpret_cast<intptr_t>(view.metadata().addressOf_x())
0171           << "  [z - y] = 0x"
0172           << reinterpret_cast<intptr_t>(view.metadata().addressOf_z()) -
0173                  reinterpret_cast<intptr_t>(view.metadata().addressOf_y())
0174           << "  [id - z] = 0x"
0175           << reinterpret_cast<intptr_t>(view.metadata().addressOf_id()) -
0176                  reinterpret_cast<intptr_t>(view.metadata().addressOf_z())
0177           << "  [r - id] = 0x"
0178           << reinterpret_cast<intptr_t>(view.metadata().addressOf_r()) -
0179                  reinterpret_cast<intptr_t>(view.metadata().addressOf_id())
0180           << "  [flags - r] = 0x"
0181           << reinterpret_cast<intptr_t>(view.metadata().addressOf_flags()) -
0182                  reinterpret_cast<intptr_t>(view.metadata().addressOf_r())
0183           << "  [m - flags] = 0x"
0184           << reinterpret_cast<intptr_t>(view.metadata().addressOf_m()) -
0185                  reinterpret_cast<intptr_t>(view.metadata().addressOf_flags());
0186     }
0187 
0188     checkViewAddresses(view);
0189     checkViewAddresses(mview);
0190     checkViewAddresses(cmview);
0191 
0192     const portabletest::Matrix matrix{{1, 2, 3, 4, 5, 6}, {2, 4, 6, 8, 10, 12}, {3, 6, 9, 12, 15, 18}};
0193     const portabletest::Array flags{{6, 4, 2, 0}};
0194     assert(view.r() == 1.);
0195     for (int32_t i = 0; i < view.metadata().size(); ++i) {
0196       auto vi = view[i];
0197       if (not expectXvalues_.empty() and vi.x() != expectXvalues_[i % expectXvalues_.size()]) {
0198         throw cms::Exception("Assert") << "Index " << i << " expected value "
0199                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x();
0200       }
0201       assert(vi.y() == 0.);
0202       assert(vi.z() == 0.);
0203       assert(vi.id() == i);
0204       assert(vi.flags() == flags);
0205       assert(vi.m() == matrix * i);
0206     }
0207 
0208     if (expectBackend_) {
0209       auto backend = static_cast<cms::alpakatools::Backend>(event.get(backendToken_));
0210       if (expectBackend_ != backend) {
0211         throw cms::Exception("Assert") << "Expected input backend " << cms::alpakatools::toString(*expectBackend_)
0212                                        << ", got " << cms::alpakatools::toString(backend);
0213       }
0214     }
0215 
0216     //    portabletest::TestHostMultiCollection const& productMulti = event.get(tokenMulti_);
0217     //    auto const& viewMulti0 = productMulti.const_view<0>();
0218     //    auto& mviewMulti0 = productMulti.view<0>();
0219     //    auto const& cmviewMulti0 = productMulti.view<0>();
0220     //    auto const& viewMulti1 = productMulti.const_view<1>();
0221     //    auto& mviewMulti1 = productMulti.view<1>();
0222     //    auto const& cmviewMulti1 = productMulti.view<1>();
0223 
0224     portabletest::TestHostMultiCollection2 const& productMulti2 = event.get(tokenMulti2_);
0225     auto const& viewMulti2_0 = productMulti2.const_view<0>();
0226     auto& mviewMulti2_0 = productMulti2.view<0>();
0227     auto const& cmviewMulti2_0 = productMulti2.view<0>();
0228     auto const& viewMulti2_1 = productMulti2.const_view<1>();
0229     auto& mviewMulti2_1 = productMulti2.view<1>();
0230     auto const& cmviewMulti2_1 = productMulti2.view<1>();
0231 
0232     checkViewAddresses(viewMulti2_0);
0233     checkViewAddresses(mviewMulti2_0);
0234     checkViewAddresses(cmviewMulti2_0);
0235     checkViewAddresses2(viewMulti2_1);
0236     checkViewAddresses2(mviewMulti2_1);
0237     checkViewAddresses2(cmviewMulti2_1);
0238 
0239     assert(viewMulti2_0.r() == 1.);
0240     for (int32_t i = 0; i < viewMulti2_0.metadata().size(); ++i) {
0241       auto vi = viewMulti2_0[i];
0242       //      std::stringstream s;
0243       //      s << "i=" << i << " x=" << vi.x() << " y=" << vi.y() << " z=" << vi.z() << " id=" << vi.id() << "'\nm=" << vi.m();
0244       //      std::cout << s.str() << std::endl;
0245       if (not expectXvalues_.empty() and vi.x() != expectXvalues_[i % expectXvalues_.size()]) {
0246         throw cms::Exception("Assert") << "Index " << i << " expected value "
0247                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x();
0248       }
0249       //assert(vi.x() == 0.);
0250       assert(vi.y() == 0.);
0251       assert(vi.z() == 0.);
0252       assert(vi.id() == i);
0253       assert(vi.m() == matrix * i);
0254     }
0255     assert(viewMulti2_1.r2() == 2.);
0256     for (int32_t i = 0; i < viewMulti2_1.metadata().size(); ++i) {
0257       auto vi = viewMulti2_1[i];
0258       if (not expectXvalues_.empty() and vi.x2() != expectXvalues_[i % expectXvalues_.size()]) {
0259         throw cms::Exception("Assert") << "Index " << i << " expected value "
0260                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x2();
0261       }
0262       assert(vi.y2() == 0.);
0263       assert(vi.z2() == 0.);
0264       assert(vi.id2() == i);
0265       assert(vi.m2() == matrix * i);
0266     }
0267 
0268     portabletest::TestHostMultiCollection3 const& productMulti3 = event.get(tokenMulti3_);
0269     auto const& viewMulti3_0 = productMulti3.const_view<0>();
0270     auto& mviewMulti3_0 = productMulti3.view<0>();
0271     auto const& cmviewMulti3_0 = productMulti3.view<0>();
0272     auto const& viewMulti3_1 = productMulti3.const_view<1>();
0273     auto& mviewMulti3_1 = productMulti3.view<1>();
0274     auto const& cmviewMulti3_1 = productMulti3.view<1>();
0275     auto const& viewMulti3_2 = productMulti3.const_view<2>();
0276     auto& mviewMulti3_2 = productMulti3.view<2>();
0277     auto const& cmviewMulti3_2 = productMulti3.view<2>();
0278 
0279     checkViewAddresses(viewMulti3_0);
0280     checkViewAddresses(mviewMulti3_0);
0281     checkViewAddresses(cmviewMulti3_0);
0282     checkViewAddresses2(viewMulti3_1);
0283     checkViewAddresses2(mviewMulti3_1);
0284     checkViewAddresses2(cmviewMulti3_1);
0285     checkViewAddresses3(viewMulti3_2);
0286     checkViewAddresses3(mviewMulti3_2);
0287     checkViewAddresses3(cmviewMulti3_2);
0288 
0289     assert(viewMulti3_0.r() == 1.);
0290     for (int32_t i = 0; i < viewMulti3_0.metadata().size(); ++i) {
0291       auto vi = viewMulti3_0[i];
0292       if (not expectXvalues_.empty() and vi.x() != expectXvalues_[i % expectXvalues_.size()]) {
0293         throw cms::Exception("Assert") << "Index " << i << " expected value "
0294                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x();
0295       }
0296       assert(vi.y() == 0.);
0297       assert(vi.z() == 0.);
0298       assert(vi.id() == i);
0299       assert(vi.m() == matrix * i);
0300     }
0301     assert(viewMulti3_1.r2() == 2.);
0302     for (int32_t i = 0; i < viewMulti3_1.metadata().size(); ++i) {
0303       auto vi = viewMulti3_1[i];
0304       if (not expectXvalues_.empty() and vi.x2() != expectXvalues_[i % expectXvalues_.size()]) {
0305         throw cms::Exception("Assert") << "Index " << i << " expected value "
0306                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x2();
0307       }
0308       assert(vi.y2() == 0.);
0309       assert(vi.z2() == 0.);
0310       assert(vi.id2() == i);
0311       assert(vi.m2() == matrix * i);
0312     }
0313 
0314     assert(viewMulti3_2.r3() == 3.);
0315     for (int32_t i = 0; i < viewMulti3_2.metadata().size(); ++i) {
0316       auto vi = viewMulti3_2[i];
0317       if (not expectXvalues_.empty() and vi.x3() != expectXvalues_[i % expectXvalues_.size()]) {
0318         throw cms::Exception("Assert") << "Index " << i << " expected value "
0319                                        << expectXvalues_[i % expectXvalues_.size()] << ", got " << vi.x3();
0320       }
0321       assert(vi.y3() == 0.);
0322       assert(vi.z3() == 0.);
0323       assert(vi.id3() == i);
0324       assert(vi.m3() == matrix * i);
0325     }
0326   }
0327 
0328   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0329     edm::ParameterSetDescription desc;
0330     desc.add<edm::InputTag>("source");
0331     desc.add<int>("expectSize", -1)
0332         ->setComment("Expected size of the input collection. Values < 0 mean the check is not performed. Default: -1");
0333     desc.add<std::vector<double>>("expectXvalues", std::vector<double>(0.))
0334         ->setComment(
0335             "Expected values of the 'x' field in the input collection. Empty value means to not perform the check. If "
0336             "input collection has more elements than this parameter, the parameter values are looped over. Default: "
0337             "{0.}");
0338     desc.add<std::string>("expectBackend", "")
0339         ->setComment(
0340             "Expected backend of the input collection. Empty value means to not perform the check. Default: empty "
0341             "string");
0342     descriptions.addWithDefaultLabel(desc);
0343   }
0344 
0345 private:
0346   const edm::InputTag source_;
0347   const edm::EDGetTokenT<portabletest::TestHostCollection> token_;
0348   edm::EDGetTokenT<unsigned short> backendToken_;
0349   std::optional<cms::alpakatools::Backend> expectBackend_;
0350   //const edm::EDGetTokenT<portabletest::TestHostMultiCollection> tokenMulti_;
0351   const edm::EDGetTokenT<portabletest::TestHostMultiCollection2> tokenMulti2_;
0352   const edm::EDGetTokenT<portabletest::TestHostMultiCollection3> tokenMulti3_;
0353   const int expectSize_;
0354   const std::vector<double> expectXvalues_;
0355 };
0356 
0357 #include "FWCore/Framework/interface/MakerMacros.h"
0358 DEFINE_FWK_MODULE(TestAlpakaAnalyzer);