Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-12-15 10:35:54 

0001 #include <cmath>
0002 #include <algorithm>
0003 #include <unordered_map>
0004 
0005 #include "Acts/ActsVersion.hpp"
0006 #include "Acts/Definitions/Units.hpp"
0007 #include "Acts/Definitions/Common.hpp"
0008 #include "Acts/Geometry/TrackingGeometry.hpp"
0009 
0010 #include "Acts/Seeding/BinFinder.hpp"
0011 #include "Acts/Seeding/BinnedSPGroup.hpp"
0012 #include "Acts/Seeding/SpacePointGrid.hpp"
0013 #include "Acts/Seeding/SeedFilter.hpp"
0014 #include "Acts/Seeding/EstimateTrackParamsFromSeed.hpp"
0015 #include "Acts/Surfaces/PerigeeSurface.hpp"
0016 #include "Acts/Seeding/SeedFinderConfig.hpp"
0017 #include "Acts/Seeding/SeedFinder.hpp"
0018 
0019 // Gaudi
0020 #include "Gaudi/Algorithm.h"
0021 #include "GaudiKernel/ToolHandle.h"
0022 #include "GaudiKernel/RndmGenerators.h"
0023 #include "Gaudi/Property.h"
0024 
0025 #include <k4FWCore/DataHandle.h>
0026 #include <k4Interface/IGeoSvc.h>
0027 #include "JugTrack/IActsGeoSvc.h"
0028 #include "JugTrack/DD4hepBField.h"
0029 #include "ActsExamples/EventData/Index.hpp"
0030 #include "ActsExamples/EventData/Measurement.hpp"
0031 #include "ActsExamples/EventData/Track.hpp"
0032 
0033 #include "DDRec/CellIDPositionConverter.h"
0034 
0035 #include "edm4eic/TrackerHitCollection.h"
0036 #include "Math/Vector3D.h"
0037 
0038 
0039   ///// (Reconstructed) track parameters e.g. close to the vertex.
0040   //using TrackParameters = Acts::CurvilinearTrackParameters;
0041 
0042   ///// Container of reconstructed track states for multiple tracks.
0043   //using TrackParametersContainer = std::vector<TrackParameters>;
0044 
0045   ///// MultiTrajectory definition
0046   //using Trajectory = Acts::MultiTrajectory<SourceLink>;
0047 
0048   ///// Container for the truth fitting/finding track(s)
0049   //using TrajectoryContainer = std::vector<SimMultiTrajectory>;
0050 
0051 namespace Jug::Reco {
0052 
0053     /** Initial Track parameters from MC truth.
0054      *
0055      *  TrackParmetersContainer
0056      *  \ingroup tracking
0057      */
0058     class TrackParamACTSSeeding : public Gaudi::Algorithm {
0059     public:
0060         DataHandle<edm4eic::TrackerHitCollection>
0061         m_inputHitCollection { "inputHitCollection",
0062             Gaudi::DataHandle::Reader, this };
0063         DataHandle<ActsExamples::TrackParametersContainer>
0064         m_outputInitialTrackParameters {
0065             "outputInitialTrackParameters",
0066             Gaudi::DataHandle::Writer, this };
0067 
0068         SmartIF<IGeoSvc> m_geoSvc;
0069         SmartIF<IActsGeoSvc> m_actsGeoSvc;
0070         std::shared_ptr<const dd4hep::rec::CellIDPositionConverter> m_converter;
0071 
0072         // Acts::GeometryContext m_geoContext;
0073         std::shared_ptr<const Jug::BField::DD4hepBField> m_BField =
0074             nullptr;
0075         Acts::MagneticFieldContext m_fieldContext;
0076 
0077         /// Index type to reference elements in a container.
0078         ///
0079         /// We do not expect to have more than 2^32 elements in any
0080         /// given container so a fixed sized integer type is
0081         /// sufficient.
0082         //using Index = eic::Index;
0083 
0084         /// Space point representation of eic::TrackerHitData suitable
0085         /// for ACTS track seeding.
0086         class SpacePoint
0087         // : edm4eic::TrackerHitData
0088         {
0089         public:
0090             // Acts::Vector3 _dummy[16];
0091             Acts::Vector3 _position;
0092             Acts::Vector3 _positionError;
0093             int32_t _measurementIndex;
0094             const Acts::Surface *_surface;
0095             // Constructor to circumvent the fact that eic::TrackerHit
0096             // and associated classes are all non-polymorphic
0097             SpacePoint(const edm4eic::TrackerHit h,
0098                        const int32_t measurementIndex,
0099                        SmartIF<IActsGeoSvc> m_actsGeoSvc,
0100                        std::shared_ptr<const dd4hep::rec::CellIDPositionConverter> m_converter)
0101                 : _measurementIndex(measurementIndex)
0102             {
0103                 _position[0] = h.getPosition().x;
0104                 _position[1] = h.getPosition().y;
0105                 _position[2] = h.getPosition().z;
0106                 _positionError[0] = h.getPositionError().xx;
0107                 _positionError[1] = h.getPositionError().yy;
0108                 _positionError[2] = h.getPositionError().zz;
0109                 const auto volumeId =
0110                     m_converter->findContext(h.getCellID())->identifier;
0111                 const auto its = m_actsGeoSvc->surfaceMap().find(volumeId);
0112                 if (its == m_actsGeoSvc->surfaceMap().end()) {
0113                     _surface = nullptr;
0114                 }
0115                 else {
0116                     _surface = its->second;
0117                 }
0118             }
0119             SpacePoint(const SpacePoint &sp)
0120                 : _position(sp._position),
0121                   _positionError(sp._positionError),
0122                   _measurementIndex(sp._measurementIndex),
0123                   _surface(sp._surface)
0124             {
0125             }
0126             float x() const { return _position[0]; }
0127             float y() const { return _position[1]; }
0128             float z() const { return _position[2]; }
0129             float r() const { return std::hypot(x(), y()); }
0130             float varianceR() const
0131             {
0132                 return (std::pow(x(), 2) * _positionError[0] +
0133                         std::pow(y(), 2) * _positionError[1]) /
0134                     (std::pow(x(), 2) + std::pow(y(), 2));
0135             }
0136             float varianceZ() const { return _positionError[2]; }
0137             constexpr uint32_t measurementIndex() const {
0138                 return _measurementIndex; }
0139             bool isOnSurface() const {
0140                 if (_surface == nullptr) {
0141                     return false;
0142                 }
0143                 return _surface->isOnSurface(
0144                      Acts::GeometryContext(), {x(), y(), z()},
0145                      {0, 0, 0});
0146             }
0147         };
0148 
0149         static bool spCompare(SpacePoint r, SpacePoint s)
0150         {
0151             return
0152                 std::hypot(r.x(), r.y(), r.z()) <
0153                 std::hypot(s.x(), s.y(), s.z());
0154         }
0155 
0156         /// Container of sim seed
0157         using SeedContainer = std::vector<Acts::Seed<SpacePoint>>;
0158 
0159         /// A proto track is a collection of hits identified by their
0160         /// indices.
0161         using ProtoTrack = std::vector<ActsExamples::Index>;
0162         /// Container of proto tracks. Each proto track is identified by
0163         /// its index.
0164         using ProtoTrackContainer = std::vector<ProtoTrack>;
0165 
0166         using SpacePointContainer = std::vector<SpacePoint>;
0167 
0168         struct Config {
0169             /// Input space point collections.
0170             ///
0171             /// We allow multiple space point collections to allow
0172             /// different parts of the detector to use different
0173             /// algorithms for space point construction, e.g.
0174             /// single-hit space points for pixel-like detectors or
0175             /// double-hit space points for strip-like detectors.
0176             std::vector<std::string> inputSpacePoints;
0177             /// Output track seed collection.
0178             std::string outputSeeds;
0179             /// Output proto track collection.
0180             std::string outputProtoTracks;
0181 
0182             float cotThetaMax = std::sinh(4.01);
0183             float minPt = 100 * Acts::UnitConstants::MeV / cotThetaMax;
0184             float rMax = 440 * Acts::UnitConstants::mm;
0185             float zMin = -1500 * Acts::UnitConstants::mm;
0186             float zMax = 1700 * Acts::UnitConstants::mm;
0187             float deltaRMin = 0 * Acts::UnitConstants::mm;
0188             float deltaRMax = 600 * Acts::UnitConstants::mm;
0189             //
0190             float collisionRegionMin = -250 * Acts::UnitConstants::mm;
0191             float collisionRegionMax = 250 * Acts::UnitConstants::mm;
0192             float maxSeedsPerSpM = 2;
0193             float sigmaScattering = 5;
0194             float radLengthPerSeed = 0.1;
0195             float impactMax = 3 * Acts::UnitConstants::mm;
0196             //
0197             float bFieldInZ = 1.7 * Acts::UnitConstants::T;
0198             float beamPosX = 0 * Acts::UnitConstants::mm;
0199             float beamPosY = 0 * Acts::UnitConstants::mm;
0200 
0201             /// The minimum magnetic field to trigger the track
0202             /// parameters estimation
0203             double bFieldMin = 0.1 * Acts::UnitConstants::T;
0204 
0205             /// Constant term of the loc0 resolution.
0206             double sigmaLoc0 = 25 * Acts::UnitConstants::um;
0207             /// Constant term of the loc1 resolution.
0208             double sigmaLoc1 = 100 * Acts::UnitConstants::um;
0209             /// Phi angular resolution.
0210             double sigmaPhi = 0.02 * Acts::UnitConstants::degree;
0211             /// Theta angular resolution.
0212             double sigmaTheta = 0.02 * Acts::UnitConstants::degree;
0213             /// q/p resolution.
0214             double sigmaQOverP = 0.1 / Acts::UnitConstants::GeV;
0215             /// Time resolution.
0216             double sigmaT0 = 1400 * Acts::UnitConstants::s;
0217 
0218         int numPhiNeighbors = 1;
0219 
0220             float deltaRMiddleMinSPRange = 10. * Acts::UnitConstants::mm;
0221             float deltaRMiddleMaxSPRange = 10. * Acts::UnitConstants::mm;
0222 
0223             // vector containing the map of z bins in the top and bottom layers
0224             std::vector<std::pair<int, int> > zBinNeighborsTop;
0225             std::vector<std::pair<int, int> > zBinNeighborsBottom;
0226         } m_cfg;
0227         Acts::SpacePointGridConfig m_gridCfg;
0228         Acts::SpacePointGridOptions m_gridOpt;
0229         Acts::SeedFinderConfig<SpacePoint> m_finderCfg;
0230         Acts::SeedFinderOptions m_finderOpt;
0231         /// The track parameters covariance (assumed to be the same
0232         /// for all estimated track parameters for the moment)
0233         Acts::BoundSquareMatrix m_covariance =
0234             Acts::BoundSquareMatrix::Zero();
0235 
0236     public:
0237         TrackParamACTSSeeding(const std::string &name,
0238                               ISvcLocator* svcLoc)
0239             : Gaudi::Algorithm(name, svcLoc)
0240         {
0241             declareProperty("inputHitCollection",
0242                             m_inputHitCollection, "");
0243             declareProperty("outputInitialTrackParameters",
0244                             m_outputInitialTrackParameters, "");
0245         }
0246 
0247         StatusCode initialize() override;
0248 
0249         StatusCode execute(const EventContext&) const override;
0250     };
0251 
0252 
0253     StatusCode TrackParamACTSSeeding::initialize()
0254     {
0255         if (Gaudi::Algorithm::initialize().isFailure()) {
0256             return StatusCode::FAILURE;
0257         }
0258 
0259         m_geoSvc = service("GeoSvc");
0260         if (m_geoSvc == nullptr) {
0261             error() << "Unable to locate Geometry Service. " << endmsg;
0262             return StatusCode::FAILURE;
0263         }
0264         m_actsGeoSvc = service("ActsGeoSvc");
0265         if (m_actsGeoSvc == nullptr) {
0266             error() << "Unable to locate ACTS Geometry Service. " << endmsg;
0267             return StatusCode::FAILURE;
0268         }              
0269         m_converter = std::make_shared<const dd4hep::rec::CellIDPositionConverter>(*(m_geoSvc->getDetector()));
0270 
0271         m_BField = std::dynamic_pointer_cast<
0272             const Jug::BField::DD4hepBField>(
0273                 m_actsGeoSvc->getFieldProvider());
0274         m_fieldContext = Jug::BField::BFieldVariant(m_BField);
0275 
0276         m_gridCfg.minPt = m_cfg.minPt;
0277         m_gridCfg.rMax = m_cfg.rMax;
0278         m_gridCfg.zMax = m_cfg.zMax;
0279         m_gridCfg.zMin = m_cfg.zMin;
0280         m_gridCfg.cotThetaMax = m_cfg.cotThetaMax;
0281 
0282         m_gridCfg.deltaRMax = m_cfg.deltaRMax;
0283 
0284         // Construct seed filter
0285         Acts::SeedFilterConfig filterCfg;
0286         filterCfg.maxSeedsPerSpM = m_cfg.maxSeedsPerSpM;
0287         m_finderCfg.seedFilter =
0288             std::make_unique<Acts::SeedFilter<SpacePoint>>(
0289                 Acts::SeedFilter<SpacePoint>(filterCfg));
0290 
0291         m_finderCfg.rMax = m_cfg.rMax;
0292         m_finderCfg.deltaRMin = m_cfg.deltaRMin;
0293         m_finderCfg.deltaRMax = m_cfg.deltaRMax;
0294         m_finderCfg.collisionRegionMin = m_cfg.collisionRegionMin;
0295         m_finderCfg.collisionRegionMax = m_cfg.collisionRegionMax;
0296         m_finderCfg.zMin = m_cfg.zMin;
0297         m_finderCfg.zMax = m_cfg.zMax;
0298         m_finderCfg.maxSeedsPerSpM = m_cfg.maxSeedsPerSpM;
0299         m_finderCfg.cotThetaMax = m_cfg.cotThetaMax;
0300         m_finderCfg.sigmaScattering = m_cfg.sigmaScattering;
0301         m_finderCfg.radLengthPerSeed = m_cfg.radLengthPerSeed;
0302         m_finderCfg.minPt = m_cfg.minPt;
0303         m_finderCfg.impactMax = m_cfg.impactMax;
0304 
0305         m_finderOpt.bFieldInZ = m_cfg.bFieldInZ;
0306         m_finderOpt.beamPos =
0307             Acts::Vector2(m_cfg.beamPosX, m_cfg.beamPosY);
0308 
0309         m_finderCfg =
0310           m_finderCfg.toInternalUnits().calculateDerivedQuantities();
0311         m_finderOpt =
0312           m_finderOpt.toInternalUnits().calculateDerivedQuantities(
0313             m_finderCfg);
0314 
0315         // Set up the track parameters covariance (the same for all
0316         // tracks)
0317         m_covariance(Acts::eBoundLoc0, Acts::eBoundLoc0) =
0318             std::pow(m_cfg.sigmaLoc0, 2);
0319         m_covariance(Acts::eBoundLoc1, Acts::eBoundLoc1) =
0320             std::pow(m_cfg.sigmaLoc1, 2);
0321         m_covariance(Acts::eBoundPhi, Acts::eBoundPhi) =
0322             std::pow(m_cfg.sigmaPhi, 2);
0323         m_covariance(Acts::eBoundTheta, Acts::eBoundTheta) =
0324             std::pow(m_cfg.sigmaTheta, 2);
0325         m_covariance(Acts::eBoundQOverP, Acts::eBoundQOverP) =
0326             std::pow(m_cfg.sigmaQOverP, 2);
0327         m_covariance(Acts::eBoundTime, Acts::eBoundTime) =
0328             std::pow(m_cfg.sigmaT0, 2);
0329 
0330         return StatusCode::SUCCESS;
0331     }
0332 
0333     StatusCode TrackParamACTSSeeding::execute(const EventContext&) const
0334     {
0335         const edm4eic::TrackerHitCollection *hits =
0336             m_inputHitCollection.get();
0337         // Create output collections
0338         auto initTrackParameters =
0339             m_outputInitialTrackParameters.createAndPut();
0340 
0341         static SeedContainer seeds;
0342         static Acts::SeedFinder<SpacePoint>::SeedingState state;
0343 
0344         // Sadly, eic::TrackerHit and eic::TrackerHitData are
0345     // non-polymorphic
0346         std::vector<SpacePoint> spacePoint;
0347         std::vector<const SpacePoint *> spacePointPtrs;
0348 
0349         std::shared_ptr<const Acts::TrackingGeometry>
0350             trackingGeometry = m_actsGeoSvc->trackingGeometry();
0351 
0352 #ifdef USE_LOCAL_COORD
0353         // Currently broken, possibly geometry issues
0354 #error Do not use, broken
0355         if (msgLevel(MSG::DEBUG)) {
0356             debug() << __FILE__ << ':' << __LINE__ << ": "
0357                     << sourceLinks->size() << ' '
0358                     << hits->size() << endmsg;
0359         }
0360         auto its = sourceLinks->begin();
0361         auto itm = measurements->begin();
0362         for (; its != sourceLinks->end() &&
0363                  itm != measurements->end();
0364              its++, itm++) {
0365             const Acts::Surface *surface =
0366                 trackingGeometry->findSurface(its->get().geometryId());
0367             if (surface != nullptr) {
0368                 Acts::Vector3 v =
0369                     surface->localToGlobal(
0370                         Acts::GeometryContext(),
0371                         {std::get<Acts::Measurement<
0372                          Acts::BoundIndices, 2>>(*itm).
0373                          parameters()[0],
0374                          std::get<Acts::Measurement<
0375                          Acts::BoundIndices, 2>>(*itm).
0376                          parameters()[1]},
0377                         {0, 0, 0});
0378                 if (msgLevel(MSG::DEBUG)) {
0379                     debug() << __FILE__ << ':' << __LINE__ << ": "
0380                             << its - sourceLinks->begin() << ' '
0381                         // << itm - measurements->begin() << ' '
0382                             << v[0] << ' ' << v[1] << ' ' << v[2]
0383                             << endmsg;
0384                 }
0385                 spacePoint.push_back(
0386                     SpacePoint(
0387                         edm4eic::TrackerHit(
0388                             static_cast<uint64_t>(spacePoint.size() + 1),
0389                             edm4eic::Vector3f(v[0], v[1], v[2]),
0390                             // Dummy uncertainties
0391                             edm4eic::CovDiag3f(25.0e-6 / 3.0,
0392                                                25.0e-6 / 3.0, 0.0),
0393                             0.0, 10.0, 0.05, 0.0),
0394                         static_cast<int32_t>(spacePoint.size())));
0395                 spacePointPtrs.push_back(&spacePoint.back());
0396         }
0397 #else // USE_LOCAL_COORD
0398         for(const auto &h : *hits) {
0399             spacePoint.push_back(SpacePoint(
0400                 h, static_cast<int32_t>(spacePoint.size() + 1),
0401                 m_actsGeoSvc, m_converter));
0402             spacePointPtrs.push_back(&spacePoint.back());
0403             if (msgLevel(MSG::DEBUG)) {
0404                 debug() << __FILE__ << ':' << __LINE__ << ": "
0405                         << ' ' << h.getPosition().x
0406                         << ' ' << h.getPosition().y
0407                         << ' ' << h.getPosition().z
0408                         << ' ' << h.getPositionError().xx
0409                         << ' ' << h.getPositionError().yy
0410                         << ' ' << h.getPositionError().zz
0411                         << ' ' << h.getTime()
0412                         << ' ' << h.getTimeError()
0413                         << ' ' << h.getEdep()
0414                         << ' ' << h.getEdepError()
0415                         << ' ' << spacePointPtrs.back()
0416                         << ' ' << spacePointPtrs.back()->measurementIndex()
0417                         << ' ' << spacePointPtrs.back()->isOnSurface()
0418                         << endmsg;
0419             }
0420         }
0421 #endif // USE_LOCAL_COORD
0422         if (msgLevel(MSG::DEBUG)) {
0423             debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg;
0424         }
0425 
0426         auto extractGlobalQuantities =
0427             [=](const SpacePoint& sp, float, float, float) ->
0428             std::pair<Acts::Vector3, Acts::Vector2> {
0429                 Acts::Vector3 position { sp.x(), sp.y(), sp.z() };
0430                 Acts::Vector2 covariance {
0431                     sp.varianceR(), sp.varianceZ() };
0432                 return std::make_pair(position, covariance);
0433         };
0434         if (msgLevel(MSG::DEBUG)) {
0435             debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg;
0436         }
0437 
0438         // extent used to store r range for middle spacepoint
0439         Acts::Extent rRangeSPExtent;
0440 
0441         auto bottomBinFinder =
0442             std::make_shared<Acts::BinFinder<SpacePoint>>(
0443                 Acts::BinFinder<SpacePoint>(m_cfg.zBinNeighborsBottom,
0444                                             m_cfg.numPhiNeighbors));
0445         if (msgLevel(MSG::DEBUG)) {
0446             debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg;
0447         }
0448         auto topBinFinder =
0449             std::make_shared<Acts::BinFinder<SpacePoint>>(
0450                 Acts::BinFinder<SpacePoint>(m_cfg.zBinNeighborsTop,
0451                                             m_cfg.numPhiNeighbors));
0452         if (msgLevel(MSG::DEBUG)) {
0453             debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg;
0454         }
0455 
0456         const float bFieldInZSave = m_gridOpt.bFieldInZ;
0457         const float minPtSave = m_gridCfg.minPt;
0458         m_gridCfg.minPt = 400 * Acts::UnitConstants::MeV;
0459         m_gridOpt.bFieldInZ =
0460             (m_gridCfg.minPt / Acts::UnitConstants::MeV) /
0461             (150.0 * (1.0 + FLT_EPSILON) *
0462              (m_gridCfg.rMax / Acts::UnitConstants::mm)) *
0463             1000.0 * Acts::UnitConstants::T;
0464         if (msgLevel(MSG::DEBUG)) {
0465             debug() << "createGrid() with temporary minPt = "
0466                     << m_gridCfg.minPt / Acts::UnitConstants::MeV
0467                     << " MeV, bFieldInZ = "
0468                     << m_gridOpt.bFieldInZ / (1000 * Acts::UnitConstants::T)
0469                     << " kT" << endmsg;
0470         }
0471         auto grid =
0472             Acts::SpacePointGridCreator::createGrid<SpacePoint>(
0473                 m_gridCfg, m_gridOpt);
0474         m_gridOpt.bFieldInZ = bFieldInZSave;
0475         m_gridCfg.minPt = minPtSave;
0476         if (msgLevel(MSG::DEBUG)) {
0477             debug() << "phiBins = "
0478                     << grid->axes().front()->getNBins()
0479                     << ", zBins = "
0480                     << grid->axes().back()->getNBins() << endmsg;
0481         }
0482 
0483         auto spacePointsGrouping =
0484             Acts::BinnedSPGroup<SpacePoint>(
0485                 spacePointPtrs.begin(), spacePointPtrs.end(),
0486                 extractGlobalQuantities, bottomBinFinder,
0487                 topBinFinder, std::move(grid),
0488                 rRangeSPExtent,
0489                 m_finderCfg, m_finderOpt);
0490         auto finder = Acts::SeedFinder<SpacePoint>(m_finderCfg);
0491 
0492         if (msgLevel(MSG::DEBUG)) {
0493             debug() << __FILE__ << ':' << __LINE__
0494                     << ": spacePointsGrouping.size() = "
0495                     << spacePointsGrouping.size() << endmsg;
0496         }
0497 
0498         const Acts::Range1D<float> rMiddleSPRange(
0499             std::floor(rRangeSPExtent.min(Acts::binR) / 2) * 2 +
0500             m_cfg.deltaRMiddleMinSPRange,
0501             std::floor(rRangeSPExtent.max(Acts::binR) / 2) * 2 -
0502             m_cfg.deltaRMiddleMaxSPRange);
0503 
0504         // Run the seeding
0505         seeds.clear();
0506 
0507         for (const auto [bottom, middle, top] : spacePointsGrouping) {
0508             finder.createSeedsForGroup(
0509                 m_finderOpt, state, spacePointsGrouping.grid(),
0510                 std::back_inserter(seeds), bottom, middle, top, rMiddleSPRange
0511             );
0512         }
0513 
0514         if (msgLevel(MSG::DEBUG)) {
0515             debug() << "seeds.size() = " << seeds.size() << endmsg;
0516         }
0517 
0518         ActsExamples::TrackParametersContainer trackParameters;
0519         ProtoTrackContainer tracks;
0520         trackParameters.reserve(seeds.size());
0521         tracks.reserve(seeds.size());
0522 
0523         std::shared_ptr<const Acts::MagneticFieldProvider>
0524             magneticField = m_actsGeoSvc->getFieldProvider();
0525 
0526         // if (msgLevel(MSG::DEBUG)) { debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg; }
0527         auto bCache = magneticField->makeCache(m_fieldContext);
0528 
0529         std::unordered_map<size_t, bool> spTaken;
0530 
0531         for (size_t iseed = 0; iseed < seeds.size(); iseed++) {
0532             const auto &seed = seeds[iseed];
0533             // Get the bottom space point and its reference surface
0534             const auto bottomSP = seed.sp().front();
0535             auto hitIdx = bottomSP->measurementIndex();
0536             // const Acts::Surface *surface = nullptr;
0537             const Acts::Surface *surface = bottomSP->_surface;
0538             if (surface == nullptr) {
0539                 if (msgLevel(MSG::DEBUG)) {
0540                     debug() << "hit " << hitIdx << " ("
0541                             << bottomSP->x() << ", " << bottomSP->y()
0542                             << ", " << bottomSP->z()
0543                             << ") lost its surface" << endmsg;
0544                 }
0545             }
0546             if (std::find_if(spacePoint.begin(), spacePoint.end(),
0547                              [&surface](const SpacePoint sp) {
0548                                  return surface == sp._surface;
0549                              }) == spacePoint.end()) {
0550                 if (msgLevel(MSG::DEBUG)) {
0551                     debug() << "hit " << hitIdx
0552                             << " has a surface that was never "
0553                             << "associated with a hit" << endmsg;
0554                 }
0555                 surface = nullptr;
0556             }
0557             if (surface == nullptr && msgLevel(MSG::DEBUG)) {
0558                 debug() << "hit " << hitIdx
0559                         << " is not found in the tracking gemetry"
0560                         << endmsg;
0561                 continue;
0562             }
0563 
0564             // Get the magnetic field at the bottom space point
0565             auto fieldRes = magneticField->getField(
0566                 {bottomSP->x(), bottomSP->y(), bottomSP->z()},
0567                 bCache);
0568             // if (msgLevel(MSG::DEBUG)) { debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg; }
0569             // Estimate the track parameters from seed
0570             auto optParams = Acts::estimateTrackParamsFromSeed(
0571                 Acts::GeometryContext(),
0572                 seed.sp().begin(), seed.sp().end(),
0573                 *surface, *fieldRes, m_cfg.bFieldMin);
0574             // if (msgLevel(MSG::DEBUG)) { debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg; }
0575             if (not optParams.has_value()) {
0576                 debug() << "Estimation of track parameters for seed "
0577                         << iseed << " failed." << endmsg;
0578                 continue;
0579             }
0580             else if (!(spTaken[seed.sp()[0]->measurementIndex()] ||
0581                        spTaken[seed.sp()[1]->measurementIndex()] ||
0582                        spTaken[seed.sp()[2]->measurementIndex()])) {
0583                 const auto& params = optParams.value();
0584                 initTrackParameters->emplace_back(
0585                     surface->getSharedPtr(), params,
0586                     m_covariance, Acts::ParticleHypothesis::pion());
0587                 // Activate/deactivate for unique seed filtering
0588 #if 0
0589                 spTaken[seed.sp()[0]->measurementIndex()] = true;
0590                 spTaken[seed.sp()[1]->measurementIndex()] = true;
0591                 spTaken[seed.sp()[2]->measurementIndex()] = true;
0592 #endif
0593             }
0594             // if (msgLevel(MSG::DEBUG)) { debug() << __FILE__ << ':' << __LINE__ << ": " << endmsg; }
0595         }
0596 
0597         return StatusCode::SUCCESS;
0598     }
0599     
0600     DECLARE_COMPONENT(TrackParamACTSSeeding)
0601 } // namespace Jug::reco
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!