Plugins/EDM4hep/EDM4HepMeasurementTests.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include <boost/test/unit_test.hpp>
0010 #include <boost/test/unit_test_suite.hpp>
0011
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Definitions/TrackParametrization.hpp"
0014 #include "Acts/EventData/SubspaceHelpers.hpp"
0015 #include "Acts/EventData/detail/GenerateParameters.hpp"
0016 #include "Acts/Geometry/GeometryContext.hpp"
0017 #include "Acts/Surfaces/DiscSurface.hpp"
0018 #include "ActsPlugins/EDM4hep/EDM4hepUtil.hpp"
0019 #include "ActsPodioEdm/detail/SubspaceIndexHelpers.hpp"
0020 #include "ActsTests/CommonHelpers/FloatComparisons.hpp"
0021 #include <ActsPodioEdm/TrackerHitLocalCollection.h>
0022
0023 using namespace Acts;
0024 using namespace ActsPlugins;
0025 using namespace Acts::UnitLiterals;
0026 using namespace Acts::detail::Test;
0027
0028 namespace {
0029 std::default_random_engine rng(123);
0030 auto gctx = GeometryContext::dangerouslyDefaultConstruct();
0031 }  // namespace
0032
0033 BOOST_AUTO_TEST_SUITE(EDM4hepMeasurementTest)
0034
0035 BOOST_AUTO_TEST_CASE(WriteMeasurement) {
0036   auto [parameters, covariance] =
0037       generateParametersCovariance<double, eBoundSize>(rng);
0038
0039   std::uint64_t cellId = 1234;
0040
0041   ActsPodioEdm::TrackerHitLocalCollection hits;
0042   auto to = hits.create();
0043
0044   std::vector<std::uint8_t> indices = {eBoundLoc0, eBoundLoc1};
0045   FixedSubspaceHelper<eBoundSize, 2> helper(indices);
0046
0047   Vector2 measPos = helper.projectVector(parameters);
0048   SquareMatrix<2> measCov = helper.projectMatrix(covariance);
0049
0050   auto surface = Surface::makeShared<DiscSurface>(
0051       Transform3::Identity() * Translation3(Vector3{1, 2, 3}), 1_mm, 1_mm);
0052
0053   Vector3 global = surface->localToGlobal(gctx, measPos, Vector3::UnitZ());
0054
0055   EDM4hepUtil::writeMeasurement(gctx, {measPos.data(), 2},
0056                                 {measCov.data(), 2, 2}, indices, cellId,
0057                                 *surface, to);
0058
0059   BOOST_CHECK_EQUAL(to.getCellID(), cellId);
0060   BOOST_CHECK_EQUAL(to.getPosition()[0] * 1_mm, global.x());
0061   BOOST_CHECK_EQUAL(to.getPosition()[1] * 1_mm, global.y());
0062   BOOST_CHECK_EQUAL(to.getPosition()[2] * 1_mm, global.z());
0063
0064   // Time should be zero since we don't have time measurement
0065   BOOST_CHECK_EQUAL(to.getTime() * 1_ns, 0.0);
0066
0067   auto meas = to.getMeasurement();
0068   BOOST_CHECK_EQUAL(meas.size(), 2);
0069   CHECK_CLOSE_REL(meas[0], measPos.x(), 1e-6);
0070   CHECK_CLOSE_REL(meas[1], measPos.y(), 1e-6);
0071
0072   auto cov = to.getCovariance();
0073   BOOST_CHECK_EQUAL(cov.size(), 4);
0074   CHECK_CLOSE_REL(cov[0], measCov(0, 0), 1e-6);
0075   CHECK_CLOSE_REL(cov[1], measCov(1, 0), 1e-6);
0076   CHECK_CLOSE_REL(cov[2], measCov(0, 1), 1e-6);
0077   CHECK_CLOSE_REL(cov[3], measCov(1, 1), 1e-6);
0078
0079   auto unpackedIndices = ActsPodioEdm::detail::decodeIndices(
0080       static_cast<std::uint32_t>(to.getType()));
0081   BOOST_CHECK_EQUAL(unpackedIndices.size(), 2);
0082   BOOST_CHECK_EQUAL(unpackedIndices[0], eBoundLoc0);
0083   BOOST_CHECK_EQUAL(unpackedIndices[1], eBoundLoc1);
0084
0085   // Round-trip: read back and verify
0086   auto read = EDM4hepUtil::readMeasurement(to);
0087   BOOST_CHECK_EQUAL(read.cellId, cellId);
0088   BOOST_CHECK_EQUAL(read.indices.size(), 2);
0089   BOOST_CHECK_EQUAL(read.indices[0], eBoundLoc0);
0090   BOOST_CHECK_EQUAL(read.indices[1], eBoundLoc1);
0091   CHECK_CLOSE_REL(read.parameters(eBoundLoc0), measPos.x(), 1e-6);
0092   CHECK_CLOSE_REL(read.parameters(eBoundLoc1), measPos.y(), 1e-6);
0093   CHECK_CLOSE_REL(read.covariance(eBoundLoc0, eBoundLoc0), measCov(0, 0), 1e-6);
0094   CHECK_CLOSE_REL(read.covariance(eBoundLoc0, eBoundLoc1), measCov(0, 1), 1e-6);
0095   CHECK_CLOSE_REL(read.covariance(eBoundLoc1, eBoundLoc0), measCov(1, 0), 1e-6);
0096   CHECK_CLOSE_REL(read.covariance(eBoundLoc1, eBoundLoc1), measCov(1, 1), 1e-6);
0097 }
0098
0099 BOOST_AUTO_TEST_CASE(WriteMeasurementNoPosition) {
0100   auto [parameters, covariance] =
0101       generateParametersCovariance<double, eBoundSize>(rng);
0102
0103   std::uint64_t cellId = 1234;
0104
0105   ActsPodioEdm::TrackerHitLocalCollection hits;
0106   auto to = hits.create();
0107
0108   // Only measure phi and theta
0109   std::vector<std::uint8_t> indices = {eBoundPhi, eBoundTheta};
0110   FixedSubspaceHelper<eBoundSize, 2> helper(indices);
0111
0112   Vector2 measPos = helper.projectVector(parameters);
0113   SquareMatrix<2> measCov = helper.projectMatrix(covariance);
0114
0115   auto surface = Surface::makeShared<DiscSurface>(
0116       Transform3::Identity() * Translation3(Vector3{1, 2, 3}), 1_mm, 1_mm);
0117
0118   EDM4hepUtil::writeMeasurement(gctx, {measPos.data(), 2},
0119                                 {measCov.data(), 2, 2}, indices, cellId,
0120                                 *surface, to);
0121
0122   BOOST_CHECK_EQUAL(to.getCellID(), cellId);
0123   // Position should be zero since we don't have loc0/loc1
0124   BOOST_CHECK_EQUAL(to.getPosition()[0] * 1_mm, 0.0);
0125   BOOST_CHECK_EQUAL(to.getPosition()[1] * 1_mm, 0.0);
0126   BOOST_CHECK_EQUAL(to.getPosition()[2] * 1_mm, 0.0);
0127
0128   auto meas = to.getMeasurement();
0129   BOOST_CHECK_EQUAL(meas.size(), 2);
0130   CHECK_CLOSE_REL(meas[0], measPos.x(), 1e-6);
0131   CHECK_CLOSE_REL(meas[1], measPos.y(), 1e-6);
0132
0133   auto cov = to.getCovariance();
0134   BOOST_CHECK_EQUAL(cov.size(), 4);
0135   CHECK_CLOSE_REL(cov[0], measCov(0, 0), 1e-6);
0136   CHECK_CLOSE_REL(cov[1], measCov(1, 0), 1e-6);
0137   CHECK_CLOSE_REL(cov[2], measCov(0, 1), 1e-6);
0138   CHECK_CLOSE_REL(cov[3], measCov(1, 1), 1e-6);
0139
0140   auto unpackedIndices = ActsPodioEdm::detail::decodeIndices(
0141       static_cast<std::uint32_t>(to.getType()));
0142   BOOST_CHECK_EQUAL(unpackedIndices.size(), 2);
0143   BOOST_CHECK_EQUAL(unpackedIndices[0], eBoundPhi);
0144   BOOST_CHECK_EQUAL(unpackedIndices[1], eBoundTheta);
0145 }
0146
0147 BOOST_AUTO_TEST_CASE(WriteMeasurementWithTime) {
0148   auto [parameters, covariance] =
0149       generateParametersCovariance<double, eBoundSize>(rng);
0150
0151   std::uint64_t cellId = 1234;
0152
0153   ActsPodioEdm::TrackerHitLocalCollection hits;
0154   auto to = hits.create();
0155
0156   // Measure loc0, loc1, and time
0157   std::vector<std::uint8_t> indices = {eBoundLoc0, eBoundLoc1, eBoundTime};
0158   FixedSubspaceHelper<eBoundSize, 3> helper(indices);
0159
0160   Vector3 measPos = helper.projectVector(parameters);
0161   SquareMatrix<3> measCov = helper.projectMatrix(covariance);
0162
0163   auto surface = Surface::makeShared<DiscSurface>(
0164       Transform3::Identity() * Translation3(Vector3{1, 2, 3}), 1_mm, 1_mm);
0165
0166   Vector3 global =
0167       surface->localToGlobal(gctx, measPos.head<2>(), Vector3::UnitZ());
0168
0169   EDM4hepUtil::writeMeasurement(gctx, {measPos.data(), 3},
0170                                 {measCov.data(), 3, 3}, indices, cellId,
0171                                 *surface, to);
0172
0173   BOOST_CHECK_EQUAL(to.getCellID(), cellId);
0174   BOOST_CHECK_EQUAL(to.getPosition()[0] * 1_mm, global.x());
0175   BOOST_CHECK_EQUAL(to.getPosition()[1] * 1_mm, global.y());
0176   BOOST_CHECK_EQUAL(to.getPosition()[2] * 1_mm, global.z());
0177   // Time should be set since we have time measurement
0178   CHECK_CLOSE_REL(to.getTime() * 1_ns, parameters[eBoundTime], 1e-6);
0179
0180   auto meas = to.getMeasurement();
0181   BOOST_CHECK_EQUAL(meas.size(), 3);
0182   CHECK_CLOSE_REL(meas[0], measPos.x(), 1e-6);
0183   CHECK_CLOSE_REL(meas[1], measPos.y(), 1e-6);
0184   CHECK_CLOSE_REL(meas[2], measPos.z(), 1e-6);
0185
0186   auto cov = to.getCovariance();
0187   BOOST_CHECK_EQUAL(cov.size(), 9);
0188   for (int i = 0; i < 3; ++i) {
0189     for (int j = 0; j < 3; ++j) {
0190       CHECK_CLOSE_REL(cov[j * 3 + i], measCov(i, j), 1e-6);
0191     }
0192   }
0193
0194   auto unpackedIndices = ActsPodioEdm::detail::decodeIndices(
0195       static_cast<std::uint32_t>(to.getType()));
0196   BOOST_CHECK_EQUAL(unpackedIndices.size(), 3);
0197   BOOST_CHECK_EQUAL(unpackedIndices[0], eBoundLoc0);
0198   BOOST_CHECK_EQUAL(unpackedIndices[1], eBoundLoc1);
0199   BOOST_CHECK_EQUAL(unpackedIndices[2], eBoundTime);
0200 }
0201
0202 BOOST_AUTO_TEST_CASE(EncodeDecodeIndices) {
0203   // Test empty span
0204   {
0205     std::vector<std::uint8_t> indices = {};
0206     std::uint32_t encoded = ActsPodioEdm::detail::encodeIndices(indices);
0207     auto decoded = ActsPodioEdm::detail::decodeIndices(encoded);
0208     BOOST_CHECK_EQUAL(decoded.size(), 0);
0209   }
0210
0211   // Test single value
0212   {
0213     std::vector<std::uint8_t> indices = {3};
0214     std::uint32_t encoded = ActsPodioEdm::detail::encodeIndices(indices);
0215     auto decoded = ActsPodioEdm::detail::decodeIndices(encoded);
0216     BOOST_CHECK_EQUAL(decoded.size(), 1);
0217     BOOST_CHECK_EQUAL(decoded.at(0), 3);
0218   }
0219
0220   // Test maximum length (6)
0221   {
0222     std::vector<std::uint8_t> indices = {0, 1, 2, 3, 4, 5};
0223     auto encoded = ActsPodioEdm::detail::encodeIndices(indices);
0224     auto decoded = ActsPodioEdm::detail::decodeIndices(encoded);
0225     BOOST_CHECK_EQUAL(decoded.size(), 6);
0226     for (std::size_t i = 0; i < 6; ++i) {
0227       BOOST_CHECK_EQUAL(decoded.at(i), i);
0228     }
0229   }
0230
0231   // Test maximum value (6)
0232   {
0233     std::vector<std::uint8_t> indices = {6, 6, 6};
0234     auto encoded = ActsPodioEdm::detail::encodeIndices(indices);
0235     auto decoded = ActsPodioEdm::detail::decodeIndices(encoded);
0236     BOOST_CHECK_EQUAL(decoded.size(), 3);
0237     for (std::size_t i = 0; i < 3; ++i) {
0238       BOOST_CHECK_EQUAL(decoded.at(i), 6);
0239     }
0240   }
0241
0242   // Test mixed values
0243   {
0244     std::vector<std::uint8_t> indices = {2, 5, 1, 4};
0245     auto encoded = ActsPodioEdm::detail::encodeIndices(indices);
0246     auto decoded = ActsPodioEdm::detail::decodeIndices(encoded);
0247     BOOST_CHECK_EQUAL(decoded.size(), 4);
0248     BOOST_CHECK_EQUAL(decoded.at(0), 2);
0249     BOOST_CHECK_EQUAL(decoded.at(1), 5);
0250     BOOST_CHECK_EQUAL(decoded.at(2), 1);
0251     BOOST_CHECK_EQUAL(decoded.at(3), 4);
0252   }
0253 }
0254
0255 BOOST_AUTO_TEST_CASE(EncodeDecodeIndicesErrors) {
0256   // Test exceeding maximum length (7 values)
0257   {
0258     std::vector<std::uint8_t> indices = {0, 1, 2, 3, 4, 5, 6};
0259     BOOST_CHECK_THROW(ActsPodioEdm::detail::encodeIndices(indices),
0260                       std::runtime_error);
0261   }
0262
0263   // Test exceeding maximum value (7)
0264   {
0265     std::vector<std::uint8_t> indices = {7};
0266     BOOST_CHECK_THROW(ActsPodioEdm::detail::encodeIndices(indices),
0267                       std::runtime_error);
0268   }
0269
0270   // Test mixed valid/invalid values
0271   {
0272     std::vector<std::uint8_t> indices = {2, 7, 1};
0273     BOOST_CHECK_THROW(ActsPodioEdm::detail::encodeIndices(indices),
0274                       std::runtime_error);
0275   }
0276 }
0277
0278 BOOST_AUTO_TEST_CASE(TrackerHitLocalIndexApi) {
0279   ActsPodioEdm::TrackerHitLocalCollection hits;
0280   auto hit = hits.create();
0281
0282   // Initialize dimensionality via subspace indices (also sizes vectors)
0283   const std::array<std::uint8_t, 2> subspace{eBoundLoc0, eBoundLoc1};
0284   hit.setSubspaceIndices(subspace);
0285   BOOST_CHECK_EQUAL(hit.getSubspaceIndices().size(), 2u);
0286   BOOST_CHECK_EQUAL(hit.getSize(), 2u);
0287   hit.setValue(1.5f, 0);
0288   hit.setValue(-2.0f, 1);
0289   hit.setCov(4.0f, 0, 0);
0290   hit.setCov(0.25f, 0, 1);
0291   hit.setCov(0.25f, 1, 0);
0292   hit.setCov(9.0f, 1, 1);
0293
0294   CHECK_CLOSE_REL(hit.getValue(0), 1.5f, 1e-6f);
0295   CHECK_CLOSE_REL(hit.getValue(1), -2.0f, 1e-6f);
0296   CHECK_CLOSE_REL(hit.getCov(0, 0), 4.0f, 1e-6f);
0297   CHECK_CLOSE_REL(hit.getCov(0, 1), 0.25f, 1e-6f);
0298   CHECK_CLOSE_REL(hit.getCov(1, 0), 0.25f, 1e-6f);
0299   CHECK_CLOSE_REL(hit.getCov(1, 1), 9.0f, 1e-6f);
0300   BOOST_CHECK_THROW(hit.getCov(2, 0), std::runtime_error);
0301   BOOST_CHECK_THROW(hit.getCov(0, 2), std::runtime_error);
0302   BOOST_CHECK_THROW(hit.setCov(1.0f, 2, 0), std::runtime_error);
0303   BOOST_CHECK_THROW(hit.setCov(1.0f, 0, 2), std::runtime_error);
0304
0305   auto dims = hit.getSubspaceIndices();
0306   BOOST_REQUIRE_EQUAL(dims.size(), 2);
0307   BOOST_CHECK_EQUAL(dims[0], eBoundLoc0);
0308   BOOST_CHECK_EQUAL(dims[1], eBoundLoc1);
0309 }
0310
0311 BOOST_AUTO_TEST_CASE(TrackerHitLocalEnumApi) {
0312   ActsPodioEdm::TrackerHitLocalCollection hits;
0313   auto hit = hits.create();
0314
0315   const std::array<std::uint8_t, 2> subspace{eBoundLoc0, eBoundLoc1};
0316   hit.setSubspaceIndices(subspace);
0317
0318   // findSubspaceIndex: maps enum value to position in measurement vector
0319   BOOST_CHECK_EQUAL(hit.findSubspaceIndex(eBoundLoc0), 0u);
0320   BOOST_CHECK_EQUAL(hit.findSubspaceIndex(eBoundLoc1), 1u);
0321   BOOST_CHECK_THROW(hit.findSubspaceIndex(eBoundPhi), std::runtime_error);
0322
0323   // setValue/getCov via enum should round-trip identically to index-based API
0324   hit.setValue(1.5f, eBoundLoc0);
0325   hit.setValue(-2.0f, eBoundLoc1);
0326   hit.setCov(4.0f, eBoundLoc0, eBoundLoc0);
0327   hit.setCov(0.25f, eBoundLoc0, eBoundLoc1);
0328   hit.setCov(0.25f, eBoundLoc1, eBoundLoc0);
0329   hit.setCov(9.0f, eBoundLoc1, eBoundLoc1);
0330
0331   CHECK_CLOSE_REL(hit.getValue(eBoundLoc0), 1.5f, 1e-6f);
0332   CHECK_CLOSE_REL(hit.getValue(eBoundLoc1), -2.0f, 1e-6f);
0333   CHECK_CLOSE_REL(hit.getCov(eBoundLoc0, eBoundLoc0), 4.0f, 1e-6f);
0334   CHECK_CLOSE_REL(hit.getCov(eBoundLoc0, eBoundLoc1), 0.25f, 1e-6f);
0335   CHECK_CLOSE_REL(hit.getCov(eBoundLoc1, eBoundLoc0), 0.25f, 1e-6f);
0336   CHECK_CLOSE_REL(hit.getCov(eBoundLoc1, eBoundLoc1), 9.0f, 1e-6f);
0337
0338   // Enum access and index access return the same values
0339   CHECK_CLOSE_REL(hit.getValue(eBoundLoc0), hit.getValue(0u), 1e-6f);
0340   CHECK_CLOSE_REL(hit.getValue(eBoundLoc1), hit.getValue(1u), 1e-6f);
0341   CHECK_CLOSE_REL(hit.getCov(eBoundLoc0, eBoundLoc1), hit.getCov(0u, 1u),
0342                   1e-6f);
0343
0344   // Out-of-subspace enum values throw
0345   BOOST_CHECK_THROW(hit.getValue(eBoundPhi), std::runtime_error);
0346   BOOST_CHECK_THROW(hit.getCov(eBoundLoc0, eBoundPhi), std::runtime_error);
0347   BOOST_CHECK_THROW(hit.setValue(0.f, eBoundPhi), std::runtime_error);
0348   BOOST_CHECK_THROW(hit.setCov(0.f, eBoundLoc0, eBoundPhi), std::runtime_error);
0349
0350   // Read-only view has the same enum API
0351   ActsPodioEdm::TrackerHitLocal roHit = hit;
0352   CHECK_CLOSE_REL(roHit.getValue(eBoundLoc0), 1.5f, 1e-6f);
0353   CHECK_CLOSE_REL(roHit.getCov(eBoundLoc1, eBoundLoc1), 9.0f, 1e-6f);
0354   BOOST_CHECK_THROW(roHit.getValue(eBoundPhi), std::runtime_error);
0355 }
0356
0357 BOOST_AUTO_TEST_SUITE_END()