Skip to content

Commit

Permalink
update
Browse files Browse the repository at this point in the history
  • Loading branch information
Benjamin Huth committed Aug 29, 2024
1 parent 5e6e1cb commit 592399a
Show file tree
Hide file tree
Showing 3 changed files with 40 additions and 29 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -13,17 +13,12 @@

namespace Acts {

/// Specialization of the GeoModelDetectorElement for the ITk. This allows
/// mapping of Acts::GeometryIdentifiers to ITk modules in a straight-forward
/// way.
class GeoModelDetectorElementITk : public GeoModelDetectorElement {
class ITkIdentifier {
Acts::MultiIndex<std::size_t, 1, 2, 20, 20, 20, 1> m_identifier;

public:
GeoModelDetectorElementITk(const PVConstLink& geoPhysVol,
std::shared_ptr<Surface> surface,
const Transform3& sfTransform,
ActsScalar thickness, int hardware,
int barrelEndcap, int layerWheel, int etaModule,
int phiModule, int side);
ITkIdentifier(int hardware, int barrelEndcap, int layerWheel, int etaModule,
int phiModule, int side);

/// Access the hardware specifier (pixel=0, strip=1)
int hardware() const;
Expand All @@ -45,6 +40,24 @@ class GeoModelDetectorElementITk : public GeoModelDetectorElement {

/// A unique identifier that represents the combination of specifiers
std::size_t value() const;
};

/// Specialization of the GeoModelDetectorElement for the ITk. This allows
/// mapping of Acts::GeometryIdentifiers to ITk modules in a straight-forward
/// way.
class GeoModelDetectorElementITk : public GeoModelDetectorElement {
public:
GeoModelDetectorElementITk(const PVConstLink& geoPhysVol,
std::shared_ptr<Surface> surface,
const Transform3& sfTransform,
ActsScalar thickness, int hardware,
int barrelEndcap, int layerWheel, int etaModule,
int phiModule, int side)
: GeoModelDetectorElement(geoPhysVol, surface, sfTransform, thickness),
m_identifier(hardware, barrelEndcap, layerWheel, etaModule, phiModule,
side) {}

ITkIdentifier identifier() const { return m_identifier; }

/// Convert a GeoModelDetectorElement to a GeoModelDetectorElementITk
/// A new surface is constructed.
Expand All @@ -57,7 +70,7 @@ class GeoModelDetectorElementITk : public GeoModelDetectorElement {
int etaModule, int phiModule, int side);

private:
Acts::MultiIndex<std::size_t, 1, 2, 20, 20, 20, 1> m_identifier;
ITkIdentifier m_identifier;
};

} // namespace Acts
22 changes: 10 additions & 12 deletions Plugins/GeoModel/src/GeoModelDetectorElementITk.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -21,11 +21,9 @@ namespace Acts {
constexpr static std::array<std::pair<unsigned, int>, 3> s_barrelEndcapMap{
{{0, 0}, {1, 2}, {2, -2}}};

Acts::GeoModelDetectorElementITk::GeoModelDetectorElementITk(
const PVConstLink& geoPhysVol, std::shared_ptr<Surface> surface,
const Transform3& sfTransform, ActsScalar thickness, int hardware,
int barrelEndcap, int layerWheel, int etaModule, int phiModule, int side)
: GeoModelDetectorElement(geoPhysVol, surface, sfTransform, thickness) {
Acts::ITkIdentifier::ITkIdentifier(int hardware, int barrelEndcap,
int layerWheel, int etaModule, int phiModule,
int side) {
m_identifier.set(0, hardware);

auto found = std::ranges::find(s_barrelEndcapMap, barrelEndcap,
Expand All @@ -40,11 +38,11 @@ Acts::GeoModelDetectorElementITk::GeoModelDetectorElementITk(
m_identifier.set(5, side);
}

int Acts::GeoModelDetectorElementITk::hardware() const {
int Acts::ITkIdentifier::hardware() const {
return m_identifier.level(0);
}

int Acts::GeoModelDetectorElementITk::barrelEndcap() const {
int Acts::ITkIdentifier::barrelEndcap() const {
auto found = std::ranges::find(s_barrelEndcapMap, m_identifier.level(1),
&std::pair<unsigned, int>::first);
if (found == s_barrelEndcapMap.end()) {
Expand All @@ -53,23 +51,23 @@ int Acts::GeoModelDetectorElementITk::barrelEndcap() const {
return found->second;
}

int Acts::GeoModelDetectorElementITk::layerWheel() const {
int Acts::ITkIdentifier::layerWheel() const {
return m_identifier.level(2);
}

int Acts::GeoModelDetectorElementITk::phiModule() const {
int Acts::ITkIdentifier::phiModule() const {
return m_identifier.level(3);
}

int Acts::GeoModelDetectorElementITk::etaModule() const {
int Acts::ITkIdentifier::etaModule() const {
return m_identifier.level(4);
}

int Acts::GeoModelDetectorElementITk::side() const {
int Acts::ITkIdentifier::side() const {
return m_identifier.level(5);
}

std::size_t Acts::GeoModelDetectorElementITk::value() const {
std::size_t Acts::ITkIdentifier::value() const {
return m_identifier.value();
}

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -45,12 +45,12 @@ BOOST_AUTO_TEST_CASE(GeoModelDetectorElementConstruction) {
itkElement->surface().bounds().type());
BOOST_CHECK_NE(element->surface().associatedDetectorElement(),
itkElement->surface().associatedDetectorElement());
BOOST_CHECK_EQUAL(itkElement->barrelEndcap(), barrelEndcap);
BOOST_CHECK_EQUAL(itkElement->hardware(), hardware);
BOOST_CHECK_EQUAL(itkElement->layerWheel(), layerWheel);
BOOST_CHECK_EQUAL(itkElement->phiModule(), phiModule);
BOOST_CHECK_EQUAL(itkElement->etaModule(), etaModule);
BOOST_CHECK_EQUAL(itkElement->side(), side);
BOOST_CHECK_EQUAL(itkElement->identifier().barrelEndcap(), barrelEndcap);
BOOST_CHECK_EQUAL(itkElement->identifier().hardware(), hardware);
BOOST_CHECK_EQUAL(itkElement->identifier().layerWheel(), layerWheel);
BOOST_CHECK_EQUAL(itkElement->identifier().phiModule(), phiModule);
BOOST_CHECK_EQUAL(itkElement->identifier().etaModule(), etaModule);
BOOST_CHECK_EQUAL(itkElement->identifier().side(), side);
}

BOOST_AUTO_TEST_SUITE_END()

0 comments on commit 592399a

Please sign in to comment.