17 #include "GaudiKernel/PhysicalConstants.h" 
   19 #include "CLHEP/Random/RandGaussZiggurat.h" 
   24 #include "Acts/Surfaces/PlaneSurface.hpp" 
   27 using namespace Acts::UnitLiterals;
 
   29     constexpr 
double straightQoverP = 1. / (100._TeV);
 
   36     using namespace SegmentFit;
 
   45         ATH_CHECK(m_auxMeasProv.initialize(m_writeKey.key()));
 
   49         ATH_CHECK(m_trackingGeometryTool.retrieve());
 
   50         ATH_CHECK(m_extrapolationTool.retrieve());
 
   54         auto trkGeo = m_trackingGeometryTool->trackingGeometry();
 
   55         trkGeo->visitSurfaces([
this](
const Acts::Surface* surf){
 
   56             const auto detEl = surf->associatedDetectorElement();
 
   62             if (!m_idHelperSvc->isMuon(
det->identify())){
 
   65             ATH_MSG_ALWAYS(
"Sensitive muon surface "<<m_idHelperSvc->toString(
det->identify())<<
" -> geoId: " 
   66                         <<surf->geometryId());
 
   69         return StatusCode::SUCCESS;
 
   71     std::tuple<Amg::Vector3D, Amg::Vector3D> 
 
   74                                           CLHEP::HepRandomEngine* engine)
 const{
 
   77         auto smearedPars = segPars;
 
   84             const unsigned idx = Acts::toUnderlying(precPar);
 
   86             smearedPars[
idx] = CLHEP::RandGaussZiggurat::shoot(engine, segPars[
idx], uncert);
 
   88         <<
" parameter -- cov: "<<uncert
 
   89     <<
", original: "<<segPars[
idx]<<
", smeared: "<<smearedPars[
idx]<<
", deviation: " 
   90             <<(smearedPars[
idx] - segPars[
idx]) / uncert );
 
   93         auto [smearLocPos, smearLocDir] = 
makeLine(smearedPars);
 
   97             SeedingAux::strawSigns(smearLocPos, smearLocDir, segment.
measurements())) {
 
   99                         <<
" changes the L/R ambiguity -> avoid for this test");
 
  100             return smearSegment(gctx, segment, engine);
 
  104         if (smearLocDir.z() < 0) {
 
  105             smearLocDir = -smearLocDir;
 
  107         return std::make_tuple(locToGlob * smearLocPos, locToGlob.linear() * smearLocDir);
 
  115         const Acts::GeometryContext tgContext = gctx.
context();
 
  116         const Acts::MagneticFieldContext mfContext = m_extrapolationTool->getMagneticFieldContext(ctx);
 
  122         CLHEP::HepRandomEngine* randEngine = rngWrapper->
getEngine(ctx);
 
  125         ATH_CHECK(outHandle.record(std::make_unique<xAOD::MuonSegmentContainer>(),
 
  126                                    std::make_unique<xAOD::MuonSegmentAuxContainer>()));
 
  129         ATH_CHECK(surfaceHandle.record(std::make_unique<xAOD::TrackSurfaceContainer>(),
 
  130                                        std::make_unique<xAOD::TrackStateAuxContainer>()));
 
  133         auto auxMeasHandle = m_auxMeasProv.makeHandle(ctx, *surfaceHandle);
 
  140         ParDecor_t dec_locPars{m_localParsKey, ctx};
 
  141         ParDecor_t dec_seedPars{m_seedParsKey, ctx};
 
  144             const auto msSector = m_detMgr->getSectorEnvelope(reFitMe->chamberIndex(), 
 
  146                                                               reFitMe->etaIndex());
 
  153                 auto invTrf = sectorTrf.inverse();
 
  156                 auto& seedPars = dec_seedPars(*reFitMe);
 
  157                 seedPars[Acts::toUnderlying(ParamDefs::x0)] = locSeedPos.x();
 
  158                 seedPars[Acts::toUnderlying(ParamDefs::y0)] = locSeedPos.y();
 
  163             const GeoTrf::CoordEulerAngles sectorAngles = GeoTrf::getCoordRotationAngles(sectorTrf);
 
  167             const auto* refMeas = startMeas.front();
 
  169             const Amg::Vector3D firstSurfPos{surfAcc.get(refMeas)->transform(tgContext).translation()};
 
  171             if (reFitMe->nPhiLayers() < 1) {
 
  176                 const Amg::Vector3D lastSurfPos = surfAcc.get(startMeas.back())->transform(tgContext).translation();
 
  179                 const Amg::Transform3D trfBeneath = GeoTrf::GeoTransformRT{sectorAngles, firstSurfPos - pseudoSurfDist * planeNormal}; 
 
  180                 const Amg::Transform3D trfAbove   = GeoTrf::GeoTransformRT{sectorAngles, lastSurfPos + pseudoSurfDist * planeNormal}; 
 
  182                 auto surfBeneath = Acts::Surface::makeShared<Acts::PlaneSurface>(trfBeneath);
 
  183                 auto surfAbove   = Acts::Surface::makeShared<Acts::PlaneSurface>(trfAbove);
 
  185                 startMeas.insert(startMeas.begin(), auxMeasHandle.newMeasurement<1>(surfBeneath, ProjectorType::e1DimNoTime, 
AmgSymMatrix(1){covVal}));
 
  186                 startMeas.insert(startMeas.end(), auxMeasHandle.newMeasurement<1>(surfAbove, ProjectorType::e1DimNoTime, 
AmgSymMatrix(1){covVal}));
 
  197             const double extDist = Amg::intersect<3>(seedPos, seedDir, trfZ, 
 
  206                 std::stringstream sstr{};
 
  208                     <<reFitMe->chiSquared() / reFitMe->numberDoF()<<
", nDoF: "<<reFitMe->numberDoF()<<
", " 
  209                     <<reFitMe->nPrecisionHits()<<
", "<<reFitMe->nPhiLayers()<<std::endl;
 
  211                     sstr<<
"  **** "<<(*meas)<<
", chi2: "<<SeedingAux::chi2Term(locPos, 
locDir, *meas)
 
  212                         <<
", sign: "<<(meas->isStraw() ? 
 
  213                                 (SeedingAux::strawSign(locPos,
locDir, *meas) == 1 ? 
"R" : 
"L") : 
"-")
 
  214                         <<
", geoId: "<<(meas->spacePoint() ? surfAcc.get(meas->spacePoint()->primaryMeasurement())->geometryId()
 
  215                                                  : Acts::GeometryIdentifier{})<<std::endl;
 
  220                 ATH_MSG_VERBOSE(
"Run G2F fit on "<<msSector->identString()<<std::endl<<sstr.str());
 
  223             auto target = Acts::Surface::makeShared<Acts::PlaneSurface>(
trf);
 
  226             Acts::BoundMatrix initialCov{Acts::BoundMatrix::Identity()};
 
  228             auto initialPars = Acts::BoundTrackParameters::create(tgContext, 
target, fourPos, seedDir, straightQoverP,
 
  230             if (!initialPars.ok()) {
 
  235             auto fitTraject = m_trackFitTool->fit(startMeas, *initialPars, 
 
  236                                                   tgContext, mfContext, calContext, 
target.get());
 
  243             auto track = fitTraject->getTrack(0);
 
  248             std::vector<const xAOD::UncalibratedMeasurement*> goodMeas{};
 
  250             fitTraject->trackStateContainer().visitBackwards(
track.tipIndex(),[&](
const auto& state){
 
  251                 if (!state.hasUncalibratedSourceLink()){
 
  254                 goodMeas.insert(goodMeas.begin(),
 
  257                                 <<
", id: "<<surfAcc.get(goodMeas.front())->geometryId());
 
  266                     summary.nEtaTrigHits += (m->type() != xAOD::UncalibMeasType::MdtDriftCircleType);
 
  270             Acts::BoundTrackParameters 
parameters = 
track.createParametersAtReference();
 
  279             const Amg::Vector3D refitSeg = refitPos + Amg::intersect<3>(refitPos, refitDir, Amg::Vector3D::UnitZ(), 0).value_or(0.) * refitDir;
 
  280             const Amg::Vector3D globPos{msSector->localToGlobalTrans(gctx) * refitSeg};
 
  282             auto newSegment = outHandle->push_back(std::make_unique<xAOD::MuonSegment>());
 
  283             dec_segLink(*newSegment) = Link_t{*segments, reFitMe->index(), ctx};
 
  285             newSegment->setDirection(globDir.x(), globDir.y(), globDir.z());
 
  286             newSegment->setPosition(globPos.x(), globPos.y(), globPos.z());
 
  288             newSegment->setFitQuality(
track.chi2(), 
track.nDoF());
 
  290             auto& locFitPars = dec_locPars(*newSegment);
 
  291             locFitPars[Acts::toUnderlying(ParamDefs::x0)] = refitSeg.x();
 
  292             locFitPars[Acts::toUnderlying(ParamDefs::y0)] = refitSeg.y();
 
  296         return StatusCode::SUCCESS;