File indexing completed on 2023-02-02 16:37:49
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/stream/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 assert(view.metadata().addressOf_x() == view.x());
0055 assert(view.metadata().addressOf_x() == &view.x(0));
0056 assert(view.metadata().addressOf_x() == &view[0].x());
0057 assert(view.metadata().addressOf_y() == view.y());
0058 assert(view.metadata().addressOf_y() == &view.y(0));
0059 assert(view.metadata().addressOf_y() == &view[0].y());
0060 assert(view.metadata().addressOf_z() == view.z());
0061 assert(view.metadata().addressOf_z() == &view.z(0));
0062 assert(view.metadata().addressOf_z() == &view[0].z());
0063 assert(view.metadata().addressOf_id() == view.id());
0064 assert(view.metadata().addressOf_id() == &view.id(0));
0065 assert(view.metadata().addressOf_id() == &view[0].id());
0066 assert(view.metadata().addressOf_m() == view.m());
0067 assert(view.metadata().addressOf_m() == &view.m(0).coeffRef(0, 0));
0068 assert(view.metadata().addressOf_m() == &view[0].m().coeffRef(0, 0));
0069 assert(view.metadata().addressOf_r() == &view.r());
0070
0071
0072 }
0073
0074 }
0075
0076 class TestAlpakaAnalyzer : public edm::stream::EDAnalyzer<> {
0077 public:
0078 TestAlpakaAnalyzer(edm::ParameterSet const& config)
0079 : source_{config.getParameter<edm::InputTag>("source")},
0080 token_{consumes(source_)},
0081 expectSize_(config.getParameter<int>("expectSize")) {}
0082
0083 void analyze(edm::Event const& event, edm::EventSetup const&) override {
0084 portabletest::TestHostCollection const& product = event.get(token_);
0085 auto const& view = product.const_view();
0086 auto& mview = product.view();
0087 auto const& cmview = product.view();
0088
0089 if (expectSize_ >= 0 and expectSize_ != view.metadata().size()) {
0090 throw cms::Exception("Assert") << "Expected input collection size " << expectSize_ << ", got "
0091 << view.metadata().size();
0092 }
0093
0094 {
0095 edm::LogInfo msg("TestAlpakaAnalyzer");
0096 msg << source_.encode() << ".size() = " << view.metadata().size() << '\n';
0097 msg << " data @ " << product.buffer().data() << ",\n"
0098 << " x @ " << view.metadata().addressOf_x() << " = " << Column(view.x(), view.metadata().size()) << ",\n"
0099 << " y @ " << view.metadata().addressOf_y() << " = " << Column(view.y(), view.metadata().size()) << ",\n"
0100 << " z @ " << view.metadata().addressOf_z() << " = " << Column(view.z(), view.metadata().size()) << ",\n"
0101 << " id @ " << view.metadata().addressOf_id() << " = " << Column(view.id(), view.metadata().size())
0102 << ",\n"
0103 << " r @ " << view.metadata().addressOf_r() << " = " << view.r() << '\n'
0104 << " m @ " << view.metadata().addressOf_m() << " = { ... {" << view[1].m()(1, Eigen::all)
0105 << " } ... } \n";
0106 msg << std::hex << " [y - x] = 0x"
0107 << reinterpret_cast<intptr_t>(view.metadata().addressOf_y()) -
0108 reinterpret_cast<intptr_t>(view.metadata().addressOf_x())
0109 << " [z - y] = 0x"
0110 << reinterpret_cast<intptr_t>(view.metadata().addressOf_z()) -
0111 reinterpret_cast<intptr_t>(view.metadata().addressOf_y())
0112 << " [id - z] = 0x"
0113 << reinterpret_cast<intptr_t>(view.metadata().addressOf_id()) -
0114 reinterpret_cast<intptr_t>(view.metadata().addressOf_z())
0115 << " [r - id] = 0x"
0116 << reinterpret_cast<intptr_t>(view.metadata().addressOf_r()) -
0117 reinterpret_cast<intptr_t>(view.metadata().addressOf_id())
0118 << " [m - r] = 0x"
0119 << reinterpret_cast<intptr_t>(view.metadata().addressOf_m()) -
0120 reinterpret_cast<intptr_t>(view.metadata().addressOf_r());
0121 }
0122
0123 checkViewAddresses(view);
0124 checkViewAddresses(mview);
0125 checkViewAddresses(cmview);
0126
0127 const portabletest::Matrix matrix{{1, 2, 3, 4, 5, 6}, {2, 4, 6, 8, 10, 12}, {3, 6, 9, 12, 15, 18}};
0128 assert(view.r() == 1.);
0129 for (int32_t i = 0; i < view.metadata().size(); ++i) {
0130 auto vi = view[i];
0131 assert(vi.x() == 0.);
0132 assert(vi.y() == 0.);
0133 assert(vi.z() == 0.);
0134 assert(vi.id() == i);
0135 assert(vi.m() == matrix * i);
0136 }
0137 }
0138
0139 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0140 edm::ParameterSetDescription desc;
0141 desc.add<edm::InputTag>("source");
0142 desc.add<int>("expectSize", -1)
0143 ->setComment("Expected size of the input collection. Values < 0 mean the check is not performed. Default: -1");
0144 descriptions.addWithDefaultLabel(desc);
0145 }
0146
0147 private:
0148 const edm::InputTag source_;
0149 const edm::EDGetTokenT<portabletest::TestHostCollection> token_;
0150 const int expectSize_;
0151 };
0152
0153 #include "FWCore/Framework/interface/MakerMacros.h"
0154 DEFINE_FWK_MODULE(TestAlpakaAnalyzer);