Newer
Older

Johannes Junggeburth
committed
/*
Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration

Johannes Junggeburth
committed
*/
#include "MuonReadoutGeomCnvAlg.h"
#include <GeoPrimitives/GeoPrimitivesHelpers.h>

Johannes Junggeburth
committed
#include <StoreGate/WriteCondHandle.h>
#include <StoreGate/ReadCondHandle.h>
#include <GeoModelKernel/GeoFullPhysVol.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometryR4/MdtReadoutElement.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometryR4/RpcReadoutElement.h>
#include <MuonReadoutGeometryR4/TgcReadoutElement.h>
#include <MuonReadoutGeometryR4/MmReadoutElement.h>
#include <MuonReadoutGeometryR4/sTgcReadoutElement.h>
#include <MuonReadoutGeometryR4/MuonChamber.h>

Johannes Junggeburth
committed
#include <MuonAlignmentDataR4/MdtAlignmentStore.h>

Johannes Junggeburth
committed
#include <MuonAlignmentDataR4/sTgcAlignmentStore.h>

Johannes Junggeburth
committed
#include <MuonAlignmentDataR4/MmAlignmentStore.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometry/MuonStation.h>
#include <MuonReadoutGeometry/MdtReadoutElement.h>
#include <MuonReadoutGeometry/RpcReadoutElement.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometry/sTgcReadoutElement.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometry/MMReadoutElement.h>

Johannes Junggeburth
committed
#include <MuonReadoutGeometry/MuonChannelDesign.h>

Johannes Junggeburth
committed
#include <AthenaKernel/IOVInfiniteRange.h>
#include <GaudiKernel/SystemOfUnits.h>

Johannes Junggeburth
committed
#include <GeoModelHelpers/defineWorld.h>
#include <GeoModelHelpers/cloneVolume.h>
#include <GeoModelHelpers/getChildNodesWithTrf.h>
#include <GeoModelHelpers/TransformToStringConverter.h>
#include <GeoModelHelpers/GeoShapeUtils.h>

Johannes Junggeburth
committed
#include <map>
#include <GaudiKernel/SystemOfUnits.h>

Johannes Junggeburth
committed
namespace {
Amg::Transform3D readOutToStation(const GeoVFullPhysVol* readOutVol) {
return readOutVol->getAbsoluteTransform().inverse() *
readOutVol->getParent()->getX();
}

Johannes Junggeburth
committed
using SubDetAlignment = ActsGeometryContext::AlignmentStorePtr;

Johannes Junggeburth
committed
}
MuonReadoutGeomCnvAlg::MuonReadoutGeomCnvAlg(const std::string& name, ISvcLocator* pSvcLocator):
AthReentrantAlgorithm{name, pSvcLocator} {}
StatusCode MuonReadoutGeomCnvAlg::initialize() {
ATH_CHECK(m_idHelperSvc.retrieve());
ATH_CHECK(m_writeKey.initialize());
ATH_CHECK(m_alignStoreKeys.initialize());

Johannes Junggeburth
committed
ATH_CHECK(detStore()->retrieve(m_detMgr));
return StatusCode::SUCCESS;
}
StatusCode MuonReadoutGeomCnvAlg::execute(const EventContext& ctx) const {
SG::WriteCondHandle<MuonGM::MuonDetectorManager> writeHandle{m_writeKey, ctx};
if (writeHandle.isValid()) {
ATH_MSG_DEBUG("The current readout geometry is still valid.");
return StatusCode::SUCCESS;
}
writeHandle.addDependency(IOVInfiniteRange::infiniteRunLB());
/// Prepare the Geometry context
ActsGeometryContext geoContext{};
using TrackingAlignment = ActsTrk::DetectorAlignStore::TrackingAlignStore;
for (const SG::ReadCondHandleKey<ActsTrk::DetectorAlignStore>& key : m_alignStoreKeys) {
SG::ReadCondHandle<ActsTrk::DetectorAlignStore> readHandle{key, ctx};
if (!readHandle.isValid()) {
ATH_MSG_FATAL("Failed to retrieve alignment store "<<key.fullKey());
return StatusCode::FAILURE;
}
writeHandle.addDependency(readHandle);
auto alignStore = std::make_unique<ActsTrk::DetectorAlignStore>(**readHandle);
/// Ensure that the position & tracking alignment caches are split from the conditions object
if (alignStore->geoModelAlignment) {
alignStore->geoModelAlignment->clearPosCache();
}
alignStore->trackingAlignment = std::make_unique<TrackingAlignment>(alignStore->detType);
geoContext.setStore(std::move(alignStore));
}
/// Check that for every detector technology there's an DetectorAlignStore in the geometry context
/// Otherwise create an empty one.
std::vector<ActsTrk::DetectorType> presentTechs = m_detMgr->getDetectorTypes();
for (const ActsTrk::DetectorType detType : presentTechs) {
if (geoContext.getStore(detType)) {
continue;
}
ATH_MSG_WARNING("No external detector alignment has been defined for technology "<<ActsTrk::to_string(detType));
geoContext.setStore(std::make_unique<ActsTrk::DetectorAlignStore>(detType));
}

Johannes Junggeburth
committed
std::unique_ptr<MuonGM::MuonDetectorManager> detMgr = std::make_unique<MuonGM::MuonDetectorManager>();
PVLink world{createGeoWorld()};

Johannes Junggeburth
committed
detMgr->addTreeTop(world);
ATH_CHECK(buildMdt(geoContext, detMgr.get(), world));

Johannes Junggeburth
committed
ATH_CHECK(buildSTGC(geoContext, detMgr.get(), world));

Johannes Junggeburth
committed
ATH_CHECK(buildMM(geoContext, detMgr.get(), world));
ATH_CHECK(buildRpc(geoContext, detMgr.get(), world));

Johannes Junggeburth
committed
ATH_CHECK(writeHandle.record(std::move(detMgr)));
return StatusCode::SUCCESS;
}
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
StatusCode MuonReadoutGeomCnvAlg::buildStation(const ActsGeometryContext& gctx,
MuonGM::MuonDetectorManager& mgr,
const Identifier& stationId,
PVLink world) const {
const std::string stName{m_idHelperSvc->stationNameString(stationId)};
const int stEta{m_idHelperSvc->stationEta(stationId)};
const int stPhi{m_idHelperSvc->stationPhi(stationId)};
MuonGM::MuonStation* station = mgr.getMuonStation(stName, stEta, stPhi);
if (station) {
ATH_MSG_DEBUG("Station "<<stName<<" "<<stEta<<" "<<stPhi<<" already exists.");
return StatusCode::SUCCESS;
}
/// Fetch the readout element to get its parent volume
const MuonGMR4::MuonReadoutElement* copyMe = m_detMgr->getReadoutElement(stationId);
/// Retrieve the full phyiscal volume
const GeoVFullPhysVol* readOutVol = copyMe->getMaterialGeom();
PVConstLink parentVolume = readOutVol->getParent();
/// Construct the aligned station transformation
const Amg::Transform3D stationTransform = copyMe->localToGlobalTrans(gctx) *
readOutToStation(readOutVol);
/// Copy the full physical volume of the muon station
PVLink parentPhysVol{make_intrusive<GeoFullPhysVol>(parentVolume->getLogVol())};
/// Make sure to copy all the children from the original tree that're not FullPhysVols -> represent
/// They represent the passive material inside the station and are needed for the TrackinGeometry building
const std::vector<GeoChildNodeWithTrf> children = getChildrenWithRef(parentVolume, false);
double minX{1.e9}, maxX{-1.e9}, minY1{1.e9}, maxY1{-1.e9}, minY2{1.e9}, maxY2{-1.e9}, minZ{1.e9}, maxZ{-1.e9};
for (const GeoChildNodeWithTrf& child : children) {
GeoVPhysVol* childVol = const_pointer_cast<GeoVPhysVol>(child.volume);
std::vector<Amg::Vector3D> edges = getPolyShapeEdges(childVol->getLogVol()->getShape(),
child.transform);
for (const Amg::Vector3D& edge : edges) {
minX = std::min(minX, edge.x());
maxX = std::max(maxX, edge.x());
minZ = std::min(minZ, edge.z());
maxZ = std::max(maxZ, edge.z());
if (edge.z() < 0) {
minY1 = std::min(minY1, edge.y());
maxY1 = std::max(maxY1, edge.y());
} else {
minY2 = std::min(minY2, edge.y());
maxY2 = std::max(maxY2, edge.y());
}
}
/// Skip the full physical volumes as they represent the readout elements
if (typeid(*childVol) == typeid(GeoFullPhysVol)) {
continue;
}
world->add(make_intrusive<GeoTransform>(child.transform));
world->add(cloneVolume(childVol));
}
/// Add the physical volume to the world
world->add(make_intrusive<GeoTransform>(stationTransform));
world->add(parentPhysVol);
/// To create the muon station, we need to extract the dimensions
/// --> Recieve the edge points from the shapes
const double shortS = (maxY1 - minY1);
const double longS = (maxY2 - minY2);
const double lengthR = (maxX - minX);
const double lengthZ = (maxZ - minZ);
auto newStation = std::make_unique<MuonGM::MuonStation>(stName,
shortS, lengthR, lengthZ, /// S / R / Z size
longS, lengthR, lengthZ, /// S / R / Z size (long)
stEta, stPhi, false);
newStation->setPhysVol(parentPhysVol);
mgr.addMuonStation(std::move(newStation));
return StatusCode::SUCCESS;
}
StatusCode MuonReadoutGeomCnvAlg::buildRpc(const ActsGeometryContext& gctx,
MuonGM::MuonDetectorManager* mgr,
PVLink world) const {
const std::vector<const MuonGMR4::RpcReadoutElement*> readoutEles = m_detMgr->getAllRpcReadoutElements();
ATH_MSG_DEBUG("Going to build "<<readoutEles.size()<<" Rpc readout elements.");
const RpcIdHelper& idHelper{m_idHelperSvc->rpcIdHelper()};
for (const MuonGMR4::RpcReadoutElement* copyMe : readoutEles) {
const Identifier reId = copyMe->identify();
const MuonGMR4::RpcReadoutElement::parameterBook& pars{copyMe->getParameters()};
/// Build the mother station if it's not already existing
ATH_CHECK(buildStation(gctx, *mgr, reId, world));
const std::string stName{m_idHelperSvc->stationNameString(reId)};
MuonGM::MuonStation* station = mgr->getMuonStation(stName,
m_idHelperSvc->stationEta(reId),
m_idHelperSvc->stationPhi(reId));
PVLink parentPhysVol{station->getPhysVol()};
GeoIntrusivePtr<const GeoVFullPhysVol> readOutVol{copyMe->getMaterialGeom()};
parentPhysVol->add(make_intrusive<GeoTransform>(readOutToStation(readOutVol).inverse()));
PVLink clonedVol{cloneVolume(const_pointer_cast<GeoVFullPhysVol>(readOutVol))};
GeoIntrusivePtr<GeoVFullPhysVol> physVol{dynamic_pointer_cast<GeoVFullPhysVol>(clonedVol)};
parentPhysVol->add(physVol);
std::unique_ptr<MuonGM::RpcReadoutElement> newElement = std::make_unique<MuonGM::RpcReadoutElement>(physVol, stName, 1, 1, false, mgr);
const bool aSide{copyMe->stationEta() > 0};
newElement->setDoubletPhi(copyMe->doubletPhi());
newElement->setDoubletR(copyMe->doubletR());
newElement->setDoubletZ(copyMe->doubletZ());
newElement->setIdentifier(reId);
newElement->setParentMuonStation(station);
/// Define the dimensions
newElement->setLongRsize(pars.halfLength);
newElement->setLongSsize(pars.halfWidth);
newElement->setLongZsize(pars.halfThickness);
newElement->setRsize(pars.halfLength);
newElement->setSsize(pars.halfWidth);
newElement->setZsize(pars.halfThickness);
newElement->m_nlayers = copyMe->nGasGaps();
newElement->m_phistripwidth = copyMe->stripPhiWidth();
newElement->m_etastripwidth = copyMe->stripEtaWidth();
newElement->m_phistrippitch = copyMe->stripPhiPitch();
newElement->m_etastrippitch = (aSide > 0 ? 1. : -1.) *copyMe->stripEtaPitch();
newElement->m_phistriplength = copyMe->stripPhiLength();
newElement->m_etastriplength = copyMe->stripEtaLength();
newElement->m_nphistripsperpanel = copyMe->nPhiStrips();
newElement->m_netastripsperpanel = copyMe->nEtaStrips();
newElement->m_nphistrippanels = copyMe->nPhiPanels();
newElement->m_hasDEDontop = true;
newElement->m_descratzneg = false;
for (unsigned int gasGap = 1; gasGap <= copyMe->nGasGaps(); ++gasGap) {
for (int doubPhi = copyMe->doubletPhiMax(); doubPhi >= copyMe->doubletPhi(); --doubPhi) {
for (bool measPhi : {false, true}) {
const int channel = 1;
const Identifier gapId = idHelper.channelID(copyMe->identify(),
copyMe->doubletZ(),
doubPhi, gasGap, measPhi,
channel);
const Amg::Vector3D locStripPos = copyMe->globalToLocalTrans(gctx) * copyMe->stripPosition(gctx, gapId);
ATH_MSG_VERBOSE("GasGap "<<m_idHelperSvc->toString(gapId)<<", local strip position: "<<Amg::toString(locStripPos));
newElement->m_gasGap_xPos[gasGap -1] = locStripPos.x();
/// Hack to assign the proper strip positions for REs having doubletPhi =2
/// in their Identifier
const int dbPIdx = copyMe->doubletPhi() == 2 ? 1 : doubPhi;
if (measPhi) {
newElement->m_first_phistrip_s[dbPIdx -1] = locStripPos.y();
newElement->m_phistrip_z = locStripPos.z();
} else{
newElement->m_first_etastrip_z = locStripPos.z();
newElement->m_etastrip_s[dbPIdx-1] = locStripPos.y();
}
}
}
newElement->m_mirrored = true;
newElement->fillCache();
newElement->m_mirrored = false;
ATH_CHECK(dumpAndCompare(gctx, *copyMe, *newElement));
mgr->addRpcReadoutElement(std::move(newElement));
}
return StatusCode::SUCCESS;
}

Johannes Junggeburth
committed
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
StatusCode MuonReadoutGeomCnvAlg::dumpAndCompare(const ActsGeometryContext& gctx,
const MuonGMR4::MmReadoutElement& refEle,
const MuonGM::MMReadoutElement& testEle) const {
if (!m_checkGeo) {
return StatusCode::SUCCESS;
}
ATH_MSG_VERBOSE("Compare basic readout transforms"<<std::endl
<<GeoTrf::toString(testEle.absTransform(),true)<<std::endl
<<GeoTrf::toString(refEle.localToGlobalTrans(gctx), true));
const MmIdHelper& idHelper{m_idHelperSvc->mmIdHelper()};
for (unsigned int gasGap = 1; gasGap <= refEle.nGasGaps(); ++ gasGap) {
const Identifier gapId = idHelper.channelID(refEle.identify(), refEle.multilayer(), gasGap, 1);
const Amg::Transform3D& refTrf{refEle.localToGlobalTrans(gctx, gapId)};
const Amg::Transform3D& testTrf{testEle.transform(gapId)};
if (!Amg::doesNotDeform(refTrf.inverse()*testTrf)) {
ATH_MSG_FATAL("The layer "<<m_idHelperSvc->toStringGasGap(gapId)<<" does not transform equally"
<<GeoTrf::toString(refTrf, true) <<" vs. "<<GeoTrf::toString(testTrf, true));
return StatusCode::FAILURE;
}
const MuonGMR4::StripDesign& stripDesign{refEle.stripLayer(gapId).design()};
for (int strip = stripDesign.firstStripNumber(); strip <= stripDesign.numStrips(); ++strip) {
const Identifier stripId = idHelper.channelID(refEle.identify(), refEle.multilayer(), gasGap, strip);
const Amg::Vector3D refStripPos{refEle.stripPosition(gctx, stripId)};
const Amg::Vector3D refStripDir{refEle.localToGlobalTrans(gctx, refEle.layerHash(stripId)).linear() * Amg::Vector3D::UnitX()};
Amg::Vector3D testStripPos{Amg::Vector3D::Zero()};
if (!testEle.stripGlobalPosition(stripId, testStripPos)) {
ATH_MSG_FATAL("Failed to retrieve strip position "<<m_idHelperSvc->toString(stripId));
return StatusCode::FAILURE;
}
const double dist = refStripDir.dot(refStripPos - testStripPos);
if (std::abs(dist) > 10. * Gaudi::Units::micrometer) {
ATH_MSG_FATAL("The strip "<<Amg::toString(testStripPos)<<" is not describing the same strip as "
<<Amg::toString(refStripPos)<<". Channel "<<m_idHelperSvc->toString(stripId)
<<" distance: "<<dist<<" "<<(dist / testEle.m_etaDesign[gasGap -1].inputWidth));
return StatusCode::FAILURE;
}
ATH_MSG_VERBOSE("Channel postion "<<m_idHelperSvc->toString(stripId)<<" match between legacy & new");
}
}
return StatusCode::SUCCESS;
}
StatusCode MuonReadoutGeomCnvAlg::dumpAndCompare(const ActsGeometryContext& gctx,
const MuonGMR4::RpcReadoutElement& refEle,
const MuonGM::RpcReadoutElement& testEle) const {

Johannes Junggeburth
committed
if (!m_checkGeo) {
return StatusCode::SUCCESS;
}
ATH_MSG_VERBOSE("Compare basic readout transforms"<<std::endl
<<GeoTrf::toString(testEle.absTransform(),true)<<std::endl
<<GeoTrf::toString(refEle.localToGlobalTrans(gctx), true));
const RpcIdHelper& idHelper{m_idHelperSvc->rpcIdHelper()};
for (unsigned int gasGap = 1; gasGap <= refEle.nGasGaps(); ++gasGap) {
for (int doubPhi = refEle.doubletPhi(); doubPhi <= refEle.doubletPhiMax(); ++doubPhi) {
for (bool measPhi : {false, true}) {
for (int strip = 1; strip <= testEle.Nstrips(measPhi); ++strip) {
const Identifier stripId = idHelper.channelID(refEle.identify(),
refEle.doubletZ(),
doubPhi, gasGap, measPhi, strip);
const Amg::Transform3D& refTrans{refEle.localToGlobalTrans(gctx, stripId)};
const Amg::Transform3D& testTrans{testEle.transform(stripId)};
if (strip == 1 && Amg::doesNotDeform(refTrans.inverse()*testTrans)) {
ATH_MSG_ERROR("Transformation for "<<m_idHelperSvc->toString(stripId)<<std::endl
<<" *** ref: "<<GeoTrf::toString(refTrans, true)<<std::endl
<<" *** test: "<<GeoTrf::toString(testTrans, true));
return StatusCode::FAILURE;
}
const Amg::Vector3D refStripPos = refEle.stripPosition(gctx, stripId);
const Amg::Vector3D testStripPos = testEle.stripPos(stripId);
if ((refStripPos - testStripPos).mag() > std::numeric_limits<float>::epsilon()){
ATH_MSG_ERROR("Mismatch in strip positions "<<m_idHelperSvc->toString(stripId)
<<" ref: "<<Amg::toString(refStripPos)<<" test: "<<Amg::toString(testStripPos)
<<" local coordinates -- ref: "<<Amg::toString(testEle.absTransform().inverse()*refStripPos)
<<" test: "<<Amg::toString(testEle.absTransform().inverse()*testStripPos));
return StatusCode::FAILURE;
}
ATH_MSG_VERBOSE("Agreement between new and old geometry for channel "<<m_idHelperSvc->toString(stripId)
<<" strip position "<<Amg::toString(refStripPos));
}
}
}
}
return StatusCode::SUCCESS;
}

Johannes Junggeburth
committed
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
StatusCode MuonReadoutGeomCnvAlg::buildMM(const ActsGeometryContext& gctx,
MuonGM::MuonDetectorManager* mgr,
PVLink world) const {
SubDetAlignment alignItr = gctx.getStore(ActsTrk::DetectorType::Mm);
const auto alignStore = alignItr ?
static_cast<const MmAlignmentStore*>(alignItr->internalAlignment.get()) : nullptr;
const std::vector<const MuonGMR4::MmReadoutElement*> mmReadouts{m_detMgr->getAllMmReadoutElements()};
for (const MuonGMR4::MmReadoutElement* copyMe : mmReadouts) {
const Identifier reId = copyMe->identify();
GeoIntrusivePtr<const GeoVFullPhysVol> readOutVol{copyMe->getMaterialGeom()};
PVLink clonedVol{cloneVolume(const_pointer_cast<GeoVFullPhysVol>(readOutVol))};
GeoIntrusivePtr<GeoFullPhysVol> physVol{dynamic_pointer_cast<GeoFullPhysVol>(clonedVol)};
world->add(make_intrusive<GeoTransform>(copyMe->localToGlobalTrans(gctx)));
world->add(physVol);
auto newRE = std::make_unique<MuonGM::MMReadoutElement>(physVol,
m_idHelperSvc->stationNameString(reId),
copyMe->stationEta(),
copyMe->stationPhi(),
copyMe->multilayer(), mgr,
alignStore ? alignStore->passivation : nullptr);
/// Loop over the gas gaps & efine the
for (unsigned int gasGap = 0; gasGap < copyMe->nGasGaps(); ++gasGap) {

Johannes Junggeburth
committed
const MuonGMR4::StripLayer& stripLayer{copyMe->stripLayer(MuonGMR4::MmReadoutElement::createHash(gasGap +1, 0))};

Johannes Junggeburth
committed
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
const MuonGMR4::StripDesign& designFrom{stripLayer.design()};
newRE->m_Xlg[gasGap] = stripLayer.toOrigin() *
Amg::getRotateZ3D(-designFrom.stereoAngle()) *
Amg::getRotateY3D(90. * Gaudi::Units::deg);
ATH_MSG_VERBOSE("Layer transform "<<gasGap<<" "<<GeoTrf::toString(newRE->m_Xlg[gasGap], true));
MuonGM::MuonChannelDesign& designTo{newRE->m_etaDesign[gasGap]};
designTo.defineTrapezoid(designFrom.shortHalfHeight(),
designFrom.longHalfHeight(),
designFrom.halfWidth(),
designFrom.stereoAngle());
designTo.type = MuonGM::MuonChannelDesign::ChannelType::etaStrip;
designTo.detType = MuonGM::MuonChannelDesign::DetType::MM;
designTo.inputPitch = designFrom.stripPitch();
designTo.inputWidth = designTo.inputPitch * std::cos(designTo.stereoAngle());
designTo.nMissedBottomEta = designTo.nMissedBottomStereo = designFrom.firstStripNumber() - 1;
designTo.totalStrips = designFrom.numStrips();
designTo.nch = designFrom.numStrips();
designTo.setFirstPos(designFrom.firstStripPos().x() + 0.5*designTo.inputPitch);
}
newRE->fillCache();
if (alignStore && alignStore->getBLine(reId)) {
newRE->setBLinePar(*alignStore->getBLine(reId));
}
ATH_CHECK(dumpAndCompare(gctx, *copyMe, *newRE));
mgr->addMMReadoutElement(std::move(newRE));
}
return StatusCode::SUCCESS;
}

Johannes Junggeburth
committed
StatusCode MuonReadoutGeomCnvAlg::buildSTGC(const ActsGeometryContext& gctx,
MuonGM::MuonDetectorManager* mgr,
PVLink world) const{
SubDetAlignment alignItr = gctx.getStore(ActsTrk::DetectorType::sTgc);
auto alignStore = alignItr ? static_cast<const sTgcAlignmentStore*>(alignItr->internalAlignment.get()) : nullptr;
const std::vector<const MuonGMR4::sTgcReadoutElement*> sTgcReadOuts{m_detMgr->getAllsTgcReadoutElements()};
for (const MuonGMR4::sTgcReadoutElement* copyMe : sTgcReadOuts) {
const Identifier reId = copyMe->identify();
GeoIntrusivePtr<const GeoVFullPhysVol> readOutVol{copyMe->getMaterialGeom()};
PVLink clonedVol{cloneVolume(const_pointer_cast<GeoVFullPhysVol>(readOutVol))};
GeoIntrusivePtr<GeoFullPhysVol> physVol{dynamic_pointer_cast<GeoFullPhysVol>(clonedVol)};
world->add(make_intrusive<GeoTransform>(copyMe->localToGlobalTrans(gctx)));
world->add(physVol);
auto newRE = std::make_unique<MuonGM::sTgcReadoutElement>(physVol,

Johannes Junggeburth
committed
m_idHelperSvc->stationNameString(reId).substr(1),

Johannes Junggeburth
committed
copyMe->stationEta(),
copyMe->stationPhi(),
copyMe->multilayer(), mgr);
if (alignStore && alignStore->getBLine(reId)) {
newRE->setBLinePar(*alignStore->getBLine(reId));
}
for (unsigned int layer = 1; layer <= copyMe->numLayers(); ++layer) {

Johannes Junggeburth
committed
using channelType = MuonGMR4::sTgcReadoutElement::ReadoutChannelType;
using ChannelDesign = MuonGM::MuonChannelDesign;
const IdentifierHash layerHash = MuonGMR4::sTgcReadoutElement::createHash(layer,channelType::Strip,0);
newRE->m_Xlg[layer -1] = Amg::getTranslate3D(copyMe->globalToLocalTrans(gctx) *
copyMe->center(gctx, layerHash));
const MuonGMR4::StripDesign& copyEtaDesign{copyMe->stripDesign(layerHash)};
/// Initialize the eta design
ChannelDesign& etaDesign{newRE->m_etaDesign[layer-1]};
etaDesign.type = ChannelDesign::ChannelType::etaStrip;
etaDesign.detType = ChannelDesign::DetType::STGC;

Johannes Junggeburth
committed
etaDesign.defineDiamond(copyEtaDesign.shortHalfHeight(),
copyEtaDesign.longHalfHeight(),
copyEtaDesign.halfWidth(),

Johannes Junggeburth
committed
} else {
etaDesign.defineTrapezoid(copyEtaDesign.shortHalfHeight(),
copyEtaDesign.longHalfHeight(),
copyEtaDesign.halfWidth());
}
etaDesign.inputPitch = copyEtaDesign.stripPitch();
etaDesign.inputWidth = copyEtaDesign.stripWidth();
etaDesign.nch = copyEtaDesign.numStrips();

Johannes Junggeburth
committed
etaDesign.setFirstPos(copyEtaDesign.firstStripPos().x());

Johannes Junggeburth
committed
/// Initialize the phi design
const MuonGMR4::WireGroupDesign& copyPhiDesign{copyMe->wireDesign(layerHash)};
ChannelDesign& phiDesign{newRE->m_phiDesign[layer-1]};
phiDesign.type = ChannelDesign::ChannelType::phiStrip;
phiDesign.detType = ChannelDesign::DetType::STGC;

Johannes Junggeburth
committed
phiDesign.defineTrapezoid(copyPhiDesign.shortHalfHeight(),
copyPhiDesign.longHalfHeight(),
copyPhiDesign.halfWidth());
} else {
phiDesign.defineDiamond(copyPhiDesign.shortHalfHeight(),
copyPhiDesign.longHalfHeight(),
copyPhiDesign.halfWidth(),

Johannes Junggeburth
committed
}
phiDesign.inputPitch = copyPhiDesign.stripPitch();
phiDesign.inputWidth = 0.015;

Johannes Junggeburth
committed
phiDesign.setFirstPos(copyPhiDesign.firstStripPos().x()); // Position of 1st wire, accounts for staggering

Johannes Junggeburth
committed
// phiDesign.firstPitch = firstWireGroup[il]; // Number of Wires in 1st group, group staggering
//phiDesign.groupWidth = wireGroupWidth; // Number of Wires normal group
phiDesign.nGroups = copyPhiDesign.numStrips(); // Number of Wire Groups
phiDesign.wireCutout = copyPhiDesign.wireCutout(); // Size of "active" wire region for digits

Johannes Junggeburth
committed
phiDesign.nch = copyPhiDesign.nAllWires();
const MuonGMR4::PadDesign& copyPadDesign{copyMe->padDesign(layerHash)};
MuonGM::MuonPadDesign& padDesign{newRE->m_padDesign[layer-1]};
padDesign.sPadWidth = 2.*copyPadDesign.shortHalfHeight();
padDesign.lPadWidth = 2.*copyPadDesign.longHalfHeight();
padDesign.Size = 2.*copyPadDesign.halfWidth();

Johannes Junggeburth
committed
}

Johannes Junggeburth
committed
mgr->addsTgcReadoutElement(std::move(newRE));
}
return StatusCode::SUCCESS;
}

Johannes Junggeburth
committed
StatusCode MuonReadoutGeomCnvAlg::buildMdt(const ActsGeometryContext& gctx,
MuonGM::MuonDetectorManager* mgr,
PVLink world) const {
/// Access the B-Line and As-built parameters
SubDetAlignment alignItr = gctx.getStore(ActsTrk::DetectorType::Mdt);
const MdtAlignmentStore* alignStore = alignItr ?
static_cast<const MdtAlignmentStore*>(alignItr->internalAlignment.get()) : nullptr;

Johannes Junggeburth
committed
const std::vector<const MuonGMR4::MdtReadoutElement*> mdtReadOuts{m_detMgr->getAllMdtReadoutElements()};
ATH_MSG_INFO("Copy "<<mdtReadOuts.size()<<" Mdt readout elements to the legacy system");
for (const MuonGMR4::MdtReadoutElement* copyMe : mdtReadOuts) {
const Identifier reId = copyMe->identify();
/// Build the mother station
ATH_CHECK(buildStation(gctx, *mgr, reId, world));

Johannes Junggeburth
committed
const std::string stName{m_idHelperSvc->stationNameString(reId)};
MuonGM::MuonStation* station = mgr->getMuonStation(stName, m_idHelperSvc->stationEta(reId), m_idHelperSvc->stationPhi(reId));
// cppcheck-suppress invalidLifetime; ok: mgr took ownership.

Johannes Junggeburth
committed
PVLink parentPhysVol{station->getPhysVol()};
GeoIntrusivePtr<const GeoVFullPhysVol> readOutVol{copyMe->getMaterialGeom()};
PVLink clonedVol{cloneVolume(const_pointer_cast<GeoVFullPhysVol>(readOutVol))};
GeoIntrusivePtr<GeoFullPhysVol> physVol{dynamic_pointer_cast<GeoFullPhysVol>(clonedVol)};
parentPhysVol->add(make_intrusive<GeoTransform>(readOutToStation(readOutVol).inverse()));

Johannes Junggeburth
committed
parentPhysVol->add(physVol);
const MuonGMR4::MdtReadoutElement::parameterBook& pars{copyMe->getParameters()};

Johannes Junggeburth
committed
std::unique_ptr<MuonGM::MdtReadoutElement> newElement = std::make_unique<MuonGM::MdtReadoutElement>(physVol, stName, mgr);
newElement->setIdentifier(reId);
newElement->setMultilayer(copyMe->multilayer());
// cppcheck-suppress invalidLifetime; ok: mgr took ownership.

Johannes Junggeburth
committed
newElement->setParentMuonStation(station);
/// Define the dimensions
newElement->setLongRsize(2*pars.halfY);
/// 1 cm is added as safety margin to the Mdt multilayer envelope
newElement->setLongSsize(2*pars.longHalfX - 1.*Gaudi::Units::cm);

Johannes Junggeburth
committed
newElement->setLongZsize(2*pars.halfHeight);
newElement->setRsize(2*pars.halfY);
newElement->setSsize(2*pars.shortHalfX - 1.*Gaudi::Units::cm);

Johannes Junggeburth
committed
newElement->setZsize(2*pars.halfHeight);
newElement->m_nlayers = copyMe->numLayers();
newElement->m_ntubesperlayer = copyMe->numTubesInLay();
newElement->m_deadlength = pars.deadLength;
newElement->m_endpluglength = pars.endPlugLength;

Johannes Junggeburth
committed
newElement->m_innerRadius = pars.tubeInnerRad;
newElement->m_tubeWallThickness = pars.tubeWall;
newElement->m_tubepitch = pars.tubePitch;
/// Need to check how to obtain this parameter from the new geometry
/// newElement->m_cutoutShift;
/// Determine the tube length's

Johannes Junggeburth
committed
const MuonGMR4::MdtTubeLayer& tubeLay{*pars.tubeLayers[0]};
unsigned int step{1};

Johannes Junggeburth
committed
double lastLength{2.*tubeLay.tubeHalfLength(1)};
for (unsigned tube = 0; tube < copyMe->numTubesInLay(); ++tube) {
const double currLength = 2.*tubeLay.tubeHalfLength(tube);
if (std::abs(lastLength - currLength) > std::numeric_limits<float>::epsilon() ||
tube == copyMe->numTubesInLay() -1) {
newElement->m_tubelength[step-1] = lastLength;
newElement->m_tubelength[step] = currLength;
if (step == 1) {
newElement->m_ntubesinastep = tube;
}

Johannes Junggeburth
committed
lastLength = currLength;
++step;
}
}
newElement->m_nsteps = step;
/// Define the tube staggering

Johannes Junggeburth
committed
const Amg::Transform3D globToLoc{copyMe->globalToLocalTrans(gctx)};

Johannes Junggeburth
committed
double xOffSet{pars.halfY}, yOffSet{pars.halfHeight};
if (newElement->barrel()) std::swap(xOffSet, yOffSet);
for (unsigned lay = 1; lay <= copyMe->numLayers(); ++lay) {
const IdentifierHash tubeHash{copyMe->measurementHash(lay, 1)};
const Amg::Vector3D locTube = globToLoc * copyMe->globalTubePos(gctx, tubeHash);
newElement->m_firstwire_x[lay-1] = locTube.z() + xOffSet;
newElement->m_firstwire_y[lay-1] = locTube.x() + yOffSet;
}
MdtAlignmentStore::chamberDistortions distort = alignStore ? alignStore->getDistortion(reId) :
MdtAlignmentStore::chamberDistortions{};
newElement->setBLinePar(distort.bLine);
station->setMdtAsBuiltParams(distort.asBuilt);
newElement->geoInitDone();
newElement->fillCache();
/// Add the readout element to the manager
ATH_CHECK(dumpAndCompare(gctx, *copyMe, *newElement));

Johannes Junggeburth
committed
mgr->addMdtReadoutElement(std::move(newElement));
}
return StatusCode::SUCCESS;
}
StatusCode MuonReadoutGeomCnvAlg::dumpAndCompare(const ActsGeometryContext& gctx,
const MuonGMR4::MdtReadoutElement& refEle,
const MuonGM::MdtReadoutElement& testEle) const {

Johannes Junggeburth
committed
if (!m_checkGeo) {
return StatusCode::SUCCESS;
}
ATH_MSG_VERBOSE("Detector element "<<m_idHelperSvc->toString(refEle.identify())
<<std::endl<<GeoTrf::toString(refEle.localToGlobalTrans(gctx))
<<std::endl<<GeoTrf::toString(testEle.getMaterialGeom()->getAbsoluteTransform())
<<std::endl<<"r-size: "<<testEle.getRsize()<<"/"<<testEle.getLongRsize()
<<" s-size: "<<testEle.getSsize()<<"/"<<testEle.getLongSsize()
<<" z-size: "<<testEle.getZsize()<<"/"<<testEle.getLongZsize());
for (unsigned int lay = 1; lay <= refEle.numLayers(); ++lay){
for (unsigned int tube = 1; tube <= refEle.numTubesInLay(); ++tube) {
const IdentifierHash tubeHash {refEle.measurementHash(lay,tube)};
if (!refEle.isValid(tubeHash)) {
ATH_MSG_VERBOSE("SKip layer / tube "<<lay <<","<<tube);
continue;
}
const Amg::Transform3D globToLocal = refEle.globalToLocalTrans(gctx, tubeHash);

Johannes Junggeburth
committed
const Amg::Vector3D refPos = refEle.globalTubePos(gctx, tubeHash);
const Amg::Vector3D tubePos = testEle.tubePos(lay, tube);

Johannes Junggeburth
committed
if ( (refPos - tubePos).mag() > Gaudi::Units::micrometer &&
(globToLocal*refPos - globToLocal * tubePos).perp() > Gaudi::Units::micrometer) {
ATH_MSG_ERROR("Tube positions differ for "<<m_idHelperSvc->toString(refEle.measurementId(tubeHash))
<<" reference: "<<GeoTrf::toString(refPos)<<" vs. test: "
<<GeoTrf::toString(tubePos) <<" delta: "<<(refPos - tubePos).mag());
return StatusCode::FAILURE;
}
ATH_MSG_VERBOSE("Tube positions layer: "<<lay<<", tube: "<<tube
<<std::endl<<"reference: "<<GeoTrf::toString(refPos)
<<std::endl<<"test: "<<GeoTrf::toString(tubePos)
<<std::endl<<testEle.tubeLength(lay, tube)<<"/"
<<testEle.getActiveTubeLength(lay, tube)<<"/"
<<testEle.getWireLength(lay,tube)
<<" vs. "<<refEle.tubeLength(tubeHash)<<"/"<<refEle.activeTubeLength(tubeHash)
<<"/"<<refEle.wireLength(tubeHash));
if (std::abs(testEle.tubeLength(lay,tube) - refEle.tubeLength(tubeHash)) >
std::numeric_limits<float>::epsilon() ) {
ATH_MSG_WARNING("Different tube length's detected for "<<m_idHelperSvc->toStringDetEl(refEle.identify())
<< " layer: "<<lay<<", tube: "<<tube<<" "<<testEle.tubeLength(lay,tube)<<" (new) vs. "
<<refEle.tubeLength(tubeHash)<<" (ref)");
}
}
}

Johannes Junggeburth
committed
return StatusCode::SUCCESS;
}