diff --git a/ACTSTracking/ACTSSeededCKFTrackingProc.hxx b/ACTSTracking/ACTSSeededCKFTrackingProc.hxx index d7b7e35..e123f3e 100644 --- a/ACTSTracking/ACTSSeededCKFTrackingProc.hxx +++ b/ACTSTracking/ACTSSeededCKFTrackingProc.hxx @@ -64,6 +64,10 @@ class ACTSSeededCKFTrackingProc : public ACTSProcBase { float _seedFinding_rMax = 150; float _seedFinding_deltaRMin = 5; float _seedFinding_deltaRMax = 80; + float _seedFinding_deltaRMinTop = 0; + float _seedFinding_deltaRMaxTop = 0; + float _seedFinding_deltaRMinBottom = 0; + float _seedFinding_deltaRMaxBottom = 0; float _seedFinding_collisionRegion = 75; float _seedFinding_zMax = 600; float _seedFinding_sigmaScattering = 50; @@ -71,6 +75,12 @@ class ACTSSeededCKFTrackingProc : public ACTSProcBase { float _seedFinding_minPt = 500; float _seedFinding_impactMax = 3 * Acts::UnitConstants::mm; + StringVec _seedFinding_zBinEdges; + int _zTopBinLen = 1; + int _zBottomBinLen = 1; + int _phiTopBinLen = 1; + int _phiBottomBinLen = 1; + // Track fit parameters double _initialTrackError_pos; double _initialTrackError_phi; diff --git a/ACTSTracking/ACTSTruthTrackingProc.hxx b/ACTSTracking/ACTSTruthTrackingProc.hxx index 9dd97b6..d495467 100644 --- a/ACTSTracking/ACTSTruthTrackingProc.hxx +++ b/ACTSTracking/ACTSTruthTrackingProc.hxx @@ -1,11 +1,13 @@ #ifndef ACTSTruthTrackingProc_h #define ACTSTruthTrackingProc_h 1 +#include #include #include #include +#include "Acts/EventData/ParticleHypothesis.hpp" #include "ACTSProcBase.hxx" diff --git a/ACTSTracking/GeometryContainers.hxx b/ACTSTracking/GeometryContainers.hxx index fb7f149..73ddf29 100644 --- a/ACTSTracking/GeometryContainers.hxx +++ b/ACTSTracking/GeometryContainers.hxx @@ -101,6 +101,14 @@ struct GeometryIdMultisetAccessor { return container->equal_range(geoId); } + // get the range of elements with requested geoId + std::pair range( + const Acts::Surface& surface) const { + assert(container != nullptr); + auto [begin, end] = container->equal_range(surface.geometryId()); + return {Iterator{begin}, Iterator{end}}; + } + // get the element using the iterator const Value& at(const Iterator& it) const { return *it; } }; diff --git a/ACTSTracking/Helpers.hxx b/ACTSTracking/Helpers.hxx index 406e2b0..3cdd4cd 100644 --- a/ACTSTracking/Helpers.hxx +++ b/ACTSTracking/Helpers.hxx @@ -1,18 +1,26 @@ #pragma once #include +#include #include #include #include +#include "Acts/EventData/ParticleHypothesis.hpp" #include #include #include +#include +#include #include "SourceLink.hxx" namespace ACTSTracking { +using TrackResult = Acts::TrackContainer::TrackProxy; + //! Get path to a resource file /** * Get absolute file of a file `inpath` by looking in the following places: @@ -43,13 +51,14 @@ std::string findFile(const std::string& inpath); * * \return Track with equivalent parameters of the ACTS track */ +/* EVENT::Track* ACTS2Marlin_track( const Acts::CombinatorialKalmanFilterResult& fitOutput, std::size_t trackTip, std::shared_ptr magneticField, Acts::MagneticFieldProvider::Cache& magCache); - +*/ //! Convert ACTS KF result to LCIO track class /** * Converted propertie are: @@ -63,10 +72,16 @@ EVENT::Track* ACTS2Marlin_track( * * \return Track with equivalent parameters of the ACTS track */ +/* EVENT::Track* ACTS2Marlin_track( const Acts::KalmanFitterResult& fitOutput, std::shared_ptr magneticField, Acts::MagneticFieldProvider::Cache& magCache); +*/ +EVENT::Track* ACTS2Marlin_track( + const TrackResult& fitter_res, + std::shared_ptr magneticField, + Acts::MagneticFieldProvider::Cache& magCache); //! Convert ACTS track state class to Marlin class /** @@ -93,4 +108,8 @@ EVENT::TrackState* ACTS2Marlin_trackState(int location, */ EVENT::LCCollection* getCollection(EVENT::LCEvent* evt, const std::string& name); + +Acts::ParticleHypothesis convertParticle(const EVENT::MCParticle* mcParticle); + } // namespace ACTSTracking + diff --git a/ACTSTracking/MeasurementCalibrator.hxx b/ACTSTracking/MeasurementCalibrator.hxx index 345c9b1..51be357 100644 --- a/ACTSTracking/MeasurementCalibrator.hxx +++ b/ACTSTracking/MeasurementCalibrator.hxx @@ -1,12 +1,16 @@ #pragma once #include +#include "Acts/EventData/SourceLink.hpp" +#include "Acts/EventData/VectorMultiTrajectory.hpp" +#include "Acts/Geometry/GeometryContext.hpp" +#include "Acts/Utilities/CalibrationContext.hpp" #include "SourceLink.hxx" namespace ACTSTracking { //! Hit stored as an measurement -using Measurement = Acts::BoundVariantMeasurement; +using Measurement = Acts::BoundVariantMeasurement; //! Collection of measurements using MeasurementContainer = std::vector; @@ -17,7 +21,7 @@ class MeasurementCalibrator { MeasurementCalibrator() = default; /// Construct using a user-provided container to chose measurements from. MeasurementCalibrator(const MeasurementContainer& measurements) - : m_measurements(&measurements) {} + : m_measurements(measurements) {} //! Find the measurement corresponding to the source link. /** @@ -28,16 +32,28 @@ class MeasurementCalibrator { template const Measurement& operator()(const SourceLink& sourceLink, const parameters_t& /* parameters */) const { - assert(m_measurements and - "Undefined measurement container in DigitizedCalibrator"); - assert((sourceLink.index() < m_measurements->size()) and + assert((sourceLink.index() < m_measurements.size()) and "Source link index is outside the container bounds"); - return (*m_measurements)[sourceLink.index()]; + return m_measurements[sourceLink.index()]; + } + + void calibrate(const Acts::GeometryContext& gctx, + const Acts::CalibrationContext& cctx, + const Acts::SourceLink& sourceLink, + Acts::VectorMultiTrajectory::TrackStateProxy trackState) const { + trackState.setUncalibratedSourceLink(sourceLink); + const auto& idxSourceLink = sourceLink.get(); + + assert((idxSourceLink.index() < m_measurements.size()) and + "Source link index is outside the container bounds"); + + const auto& meas = std::get<1>(m_measurements[idxSourceLink.index()]); //TODO workaround + trackState.allocateCalibrated(meas.size()); + trackState.setCalibrated(meas); } private: - // use pointer so the calibrator is copyable and default constructible. - const MeasurementContainer* m_measurements = nullptr; + const MeasurementContainer& m_measurements; }; } // namespace ACTSTracking diff --git a/ACTSTracking/SeedSpacePoint.hxx b/ACTSTracking/SeedSpacePoint.hxx index 330c077..53af863 100644 --- a/ACTSTracking/SeedSpacePoint.hxx +++ b/ACTSTracking/SeedSpacePoint.hxx @@ -42,6 +42,12 @@ class SeedSpacePoint { constexpr SourceLink sourceLink() const { return m_sourceLink; } + const std::optional t() const + { + return m_sourceLink.lciohit()->getTime(); + } + // TODO missing: const std::optional varianceT() const + private: // Global position float m_x; diff --git a/ACTSTracking/SourceLink.hxx b/ACTSTracking/SourceLink.hxx index 5301918..ed69002 100644 --- a/ACTSTracking/SourceLink.hxx +++ b/ACTSTracking/SourceLink.hxx @@ -2,6 +2,8 @@ #include +#include "Acts/Surfaces/Surface.hpp" + #include "GeometryContainers.hxx" namespace ACTSTracking { @@ -46,6 +48,22 @@ class SourceLink final { /// Container of index source links using SourceLinkContainer = GeometryIdMultiset; -// Wrapper for SourceLinkContainer for use with CKF -using SourceLinkAccessor = GeometryIdMultisetAccessor; +/// Accessor for the above source link container +/// +/// It wraps up a few lookup methods to be used in the Combinatorial Kalman +/// Filter +struct SourceLinkAccessor : GeometryIdMultisetAccessor { + using BaseIterator = GeometryIdMultisetAccessor::Iterator; + + using Iterator = Acts::SourceLinkAdapterIterator; + + // get the range of elements with requested geoId + std::pair range(const Acts::Surface& surface) const { + assert(container != nullptr); + auto [begin, end] = container->equal_range(surface.geometryId()); + return {Iterator{begin}, Iterator{end}}; + } +}; + } // namespace ACTSTracking + diff --git a/src/ACTSProcBase.cxx b/src/ACTSProcBase.cxx index 0cba549..2985ab3 100644 --- a/src/ACTSProcBase.cxx +++ b/src/ACTSProcBase.cxx @@ -596,7 +596,6 @@ void ACTSProcBase::buildDetector() { ringLayoutConfiguration(lbc.layerConfigurations[0]); ringLayoutConfiguration(lbc.layerConfigurations[2]); volumeConfig.layerBuilder = layerBuilder; - volumeConfig.volumeSignature = 0; auto volumeBuilder = std::make_shared( volumeConfig, Acts::getDefaultLogger(lbc.configurationName + "VolumeBuilder", diff --git a/src/ACTSSeededCKFTrackingProc.cxx b/src/ACTSSeededCKFTrackingProc.cxx index 6b15f67..788ca75 100644 --- a/src/ACTSSeededCKFTrackingProc.cxx +++ b/src/ACTSSeededCKFTrackingProc.cxx @@ -16,10 +16,10 @@ #include #include #include -#include -#include +// #include +// #include #include -#include +#include #include #include #include @@ -36,14 +36,9 @@ using namespace Acts::UnitLiterals; // Track fitting definitions using TrackFinderOptions = - Acts::CombinatorialKalmanFilterOptions; - -using TrackFinderResult = Acts::Result< - Acts::CombinatorialKalmanFilterResult>; - -using TrackFinderResultContainer = std::vector; + Acts::CombinatorialKalmanFilterOptions; +using SSPoint = ACTSTracking::SeedSpacePoint; ACTSSeededCKFTrackingProc aACTSSeededCKFTrackingProc; @@ -96,6 +91,42 @@ ACTSSeededCKFTrackingProc::ACTSSeededCKFTrackingProc() "Maximum dR between hits in a seed.", _seedFinding_deltaRMax, _seedFinding_deltaRMax); + registerProcessorParameter("SeedFinding_DeltaRMinTop", + "Minimum dR between the reference hit and outer ones in a seed.", + _seedFinding_deltaRMinTop, 0.f); + + registerProcessorParameter("SeedFinding_DeltaRMaxTop", + "Maximum dR between the reference hit and outer ones in a seed.", + _seedFinding_deltaRMaxTop, 0.f); + + registerProcessorParameter("SeedFinding_DeltaRMinBottom", + "Minimum dR between the reference hit and inner ones in a seed.", + _seedFinding_deltaRMinBottom, 0.f); + + registerProcessorParameter("SeedFinding_DeltaRMaxBottom", + "Maximum dR between the reference hit and inner ones in a seed.", + _seedFinding_deltaRMaxBottom, 0.f); + + registerProcessorParameter("SeedFinding_zTopBinLen", + "Number of top bins along Z for seeding", + _zTopBinLen, 1); + + registerProcessorParameter("SeedFinding_zBottomBinLen", + "Number of bottom bins along Z for seeding", + _zBottomBinLen, 1); + + registerProcessorParameter("SeedFinding_phiTopBinLen", + "Number of top bins along phi for seeding", + _phiTopBinLen, 1); + + registerProcessorParameter("SeedFinding_phiBottomBinLen", + "Number of bottom bins along phi for seeding", + _phiBottomBinLen, 1); + + registerProcessorParameter("SeedFinding_zBinEdges", + "Bins placement along Z for seeding.", + _seedFinding_zBinEdges, StringVec(0)); + registerProcessorParameter( "SeedFinding_CollisionRegion", "Size of the collision region in one direction (assumed symmetric).", @@ -180,6 +211,11 @@ void ACTSSeededCKFTrackingProc::init() { } _seedGeometrySelection = ACTSTracking::GeometryIdSelector(geoSelection); + + if (_seedFinding_deltaRMinTop == 0.f) _seedFinding_deltaRMinTop = _seedFinding_deltaRMin; + if (_seedFinding_deltaRMaxTop == 0.f) _seedFinding_deltaRMaxTop = _seedFinding_deltaRMax; + if (_seedFinding_deltaRMinBottom == 0.f) _seedFinding_deltaRMinBottom = _seedFinding_deltaRMin; + if (_seedFinding_deltaRMaxBottom == 0.f) _seedFinding_deltaRMaxBottom = _seedFinding_deltaRMax; } void ACTSSeededCKFTrackingProc::processRunHeader(LCRunHeader *) {} @@ -238,6 +274,7 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { sortedHits) { // Convert to Acts hit const Acts::Surface *surface = trackingGeometry()->findSurface(hit.first); + if (surface == nullptr) throw std::runtime_error("Surface not found"); const double *lcioglobalpos = hit.second->getPosition(); Acts::Vector3 globalPos = {lcioglobalpos[0], lcioglobalpos[1], @@ -250,7 +287,7 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { Acts::Vector2 loc = lpResult.value(); - Acts::SymMatrix2 localCov = Acts::SymMatrix2::Zero(); + Acts::SquareMatrix2 localCov = Acts::SquareMatrix2::Zero(); const EVENT::TrackerHitPlane *hitplane = dynamic_cast(hit.second); if (hitplane) { @@ -262,8 +299,9 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { ACTSTracking::SourceLink sourceLink(surface->geometryId(), measurements.size(), hit.second); - ACTSTracking::Measurement meas = Acts::makeMeasurement( - sourceLink, loc, localCov, Acts::eBoundLoc0, Acts::eBoundLoc1); + Acts::SourceLink src_wrap { sourceLink }; + Acts::Measurement meas = Acts::makeMeasurement( + src_wrap, loc, localCov, Acts::eBoundLoc0, Acts::eBoundLoc1); measurements.push_back(meas); sourceLinks.emplace_hint(sourceLinks.end(), sourceLink); @@ -325,7 +363,7 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { using Stepper = Acts::EigenStepper<>; using Navigator = Acts::Navigator; using Propagator = Acts::Propagator; - using CKF = Acts::CombinatorialKalmanFilter; + using CKF = Acts::CombinatorialKalmanFilter; // Configurations Navigator::Config navigatorCfg{trackingGeometry()}; @@ -342,12 +380,12 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { // Set the options Acts::MeasurementSelector::Config measurementSelectorCfg = { {Acts::GeometryIdentifier(), - {_CKF_chi2CutOff, (std::size_t)(_CKF_numMeasurementsCutOff)}}}; + { {}, { _CKF_chi2CutOff }, { (std::size_t)(_CKF_numMeasurementsCutOff) }}}}; Acts::PropagatorPlainOptions pOptions; pOptions.maxSteps = 10000; if (_propagateBackward) { - pOptions.direction = Acts::backward; + pOptions.direction = Acts::Direction::Backward; } // Construct a perigee surface as the target surface @@ -355,26 +393,53 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { Acts::Surface::makeShared( Acts::Vector3{0., 0., 0.}); + Acts::GainMatrixUpdater kfUpdater; + Acts::GainMatrixSmoother kfSmoother; + + Acts::MeasurementSelector measSel { measurementSelectorCfg }; + ACTSTracking::MeasurementCalibrator measCal { measurements }; + Acts::CombinatorialKalmanFilterExtensions + extensions; + extensions.calibrator.connect< + &ACTSTracking::MeasurementCalibrator::calibrate>( + &measCal); + extensions.updater.connect< + &Acts::GainMatrixUpdater::operator()>( + &kfUpdater); + extensions.smoother.connect< + &Acts::GainMatrixSmoother::operator()>( + &kfSmoother); + extensions.measurementSelector + .connect<&Acts::MeasurementSelector::select>( + &measSel); + + using ACTSTracking::SourceLinkAccessor; + SourceLinkAccessor slAccessor; + slAccessor.container = &sourceLinks; + Acts::SourceLinkAccessorDelegate slAccessorDelegate; + slAccessorDelegate.connect<&SourceLinkAccessor::range>(&slAccessor); + // std::unique_ptr // logger=Acts::getDefaultLogger("TrackFitting", // Acts::Logging::Level::VERBOSE); TrackFinderOptions ckfOptions = TrackFinderOptions( geometryContext(), magneticFieldContext(), calibrationContext(), - ACTSTracking::SourceLinkAccessor(), - ACTSTracking::MeasurementCalibrator(std::move(measurements)), - Acts::MeasurementSelector(measurementSelectorCfg), - // Acts::LoggerWrapper{*logger}, pOptions, - Acts::getDummyLogger(), pOptions, &(*perigeeSurface)); + slAccessorDelegate, + extensions, pOptions, perigeeSurface.get()); // // Finder configuration static const Acts::Vector3 zeropos(0, 0, 0); - Acts::SeedfinderConfig finderCfg; + Acts::SeedFinderConfig finderCfg; finderCfg.rMax = _seedFinding_rMax; finderCfg.deltaRMin = _seedFinding_deltaRMin; finderCfg.deltaRMax = _seedFinding_deltaRMax; + finderCfg.deltaRMinTopSP = _seedFinding_deltaRMinTop; + finderCfg.deltaRMaxTopSP = _seedFinding_deltaRMaxTop; + finderCfg.deltaRMinBottomSP = _seedFinding_deltaRMinBottom; + finderCfg.deltaRMaxBottomSP = _seedFinding_deltaRMaxBottom; finderCfg.collisionRegionMin = -_seedFinding_collisionRegion; finderCfg.collisionRegionMax = _seedFinding_collisionRegion; finderCfg.zMin = -_seedFinding_zMax; @@ -384,35 +449,60 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { finderCfg.sigmaScattering = _seedFinding_sigmaScattering; finderCfg.radLengthPerSeed = _seedFinding_radLengthPerSeed; finderCfg.minPt = _seedFinding_minPt * Acts::UnitConstants::MeV; - finderCfg.bFieldInZ = (*magneticField()->getField(zeropos, magCache))[2]; - finderCfg.beamPos = {0, 0}; finderCfg.impactMax = _seedFinding_impactMax * Acts::UnitConstants::mm; + finderCfg.useVariableMiddleSPRange = true; + + Acts::SeedFilterConfig filterCfg; + filterCfg.maxSeedsPerSpM = finderCfg.maxSeedsPerSpM; - Acts::SpacePointGridConfig gridCfg; - gridCfg.bFieldInZ = finderCfg.bFieldInZ; + finderCfg.seedFilter = + std::make_unique>( + Acts::SeedFilter(filterCfg.toInternalUnits())); + finderCfg = finderCfg.toInternalUnits().calculateDerivedQuantities(); + + Acts::SeedFinderOptions finderOpts; + finderOpts.bFieldInZ = (*magneticField()->getField(zeropos, magCache))[2]; // TODO investigate + finderOpts.beamPos = {0, 0}; + finderOpts = finderOpts.toInternalUnits(); + finderOpts = finderOpts.calculateDerivedQuantities(finderCfg); + + Acts::CylindricalSpacePointGridConfig gridCfg; gridCfg.cotThetaMax = finderCfg.cotThetaMax; gridCfg.deltaRMax = finderCfg.deltaRMax; gridCfg.minPt = finderCfg.minPt; gridCfg.rMax = finderCfg.rMax; gridCfg.zMax = finderCfg.zMax; gridCfg.zMin = finderCfg.zMin; + gridCfg.impactMax = finderCfg.impactMax; + if (_seedFinding_zBinEdges.size() > 0) + { + gridCfg.zBinEdges.resize(_seedFinding_zBinEdges.size()); + for (int k = 0; k < _seedFinding_zBinEdges.size(); k++) + { + float pos = std::atof(_seedFinding_zBinEdges[k].c_str()); + if (pos >= finderCfg.zMin && pos < finderCfg.zMax) + { + gridCfg.zBinEdges[k] = pos; + } + else + { + streamlog_out(WARNING) << "Wrong parameter SeedFinding_zBinEdges; " + << "default used" << std::endl; + gridCfg.zBinEdges.clear(); + break; + } + } + } - Acts::SeedFilterConfig filterCfg; - filterCfg.maxSeedsPerSpM = finderCfg.maxSeedsPerSpM; - - finderCfg.seedFilter = - std::make_unique>( - Acts::SeedFilter(filterCfg)); + Acts::CylindricalSpacePointGridOptions gridOpts; + gridOpts.bFieldInZ = (*magneticField()->getField(zeropos, magCache))[2]; // Create tools - std::function( - const ACTSTracking::SeedSpacePoint &, float, float, float)> - extractGlobalQuantities = - [](const ACTSTracking::SeedSpacePoint &sp, float, float, - float) -> std::pair { - Acts::Vector3 position{sp.x(), sp.y(), sp.z()}; - Acts::Vector2 covariance{sp.varianceR(), sp.varianceZ()}; - return std::make_pair(position, covariance); + auto extractGlobalQuantities = [](const SSPoint& sp, float, float, + float) { + Acts::Vector3 position { sp.x(), sp.y(), sp.z() }; + Acts::Vector2 covariance { sp.varianceR(), sp.varianceZ() }; + return std::make_tuple(position, covariance, sp.t()); }; std::vector spacePointPtrs( @@ -420,48 +510,55 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { std::transform(spacePoints.begin(), spacePoints.end(), spacePointPtrs.begin(), [](const ACTSTracking::SeedSpacePoint &sp) { return &sp; }); - std::shared_ptr> - bottomBinFinder = - std::make_shared>(); - std::shared_ptr> topBinFinder = - std::make_shared>(); - - std::unique_ptr> grid = - Acts::SpacePointGridCreator::createGrid( - gridCfg); + Acts::Extent rRangeSPExtent; - Acts::BinnedSPGroup spacePointsGrouping( + Acts::CylindricalSpacePointGrid grid = + Acts::CylindricalSpacePointGridCreator::createGrid( + gridCfg.toInternalUnits(), gridOpts.toInternalUnits()); + Acts::CylindricalSpacePointGridCreator::fillGrid(finderCfg, finderOpts, grid, spacePointPtrs.begin(), spacePointPtrs.end(), extractGlobalQuantities, - bottomBinFinder, topBinFinder, std::move(grid), finderCfg); + rRangeSPExtent); + + const Acts::GridBinFinder<2ul> bottomBinFinder(_phiBottomBinLen, _zBottomBinLen); + const Acts::GridBinFinder<2ul> topBinFinder(_phiTopBinLen, _zTopBinLen); + + auto spacePointsGrouping = Acts::CylindricalBinnedGroup( + std::move(grid), bottomBinFinder, topBinFinder); + + Acts::SeedFinder finder(finderCfg); + decltype(finder)::SeedingState state; + std::vector> seeds; + + state.spacePointData.resize(spacePointPtrs.size(), + finderCfg.useDetailedDoubleMeasurementInfo); + + float up = Acts::clampValue( + std::floor(rRangeSPExtent.max(Acts::binR) / 2) * 2); + const Acts::Range1D rMiddleSPRange( + std::floor(rRangeSPExtent.min(Acts::binR) / 2) * 2 + + finderCfg.deltaRMiddleMinSPRange, + up - finderCfg.deltaRMiddleMaxSPRange); // TODO investigate - Acts::Seedfinder finder(finderCfg); - decltype(finder)::State state; - std::vector> seeds; std::vector paramseeds; - Acts::BinnedSPGroupIterator group = - spacePointsGrouping.begin(); - Acts::BinnedSPGroupIterator groupEnd = - spacePointsGrouping.end(); - for (; !(group == groupEnd); ++group) { - // - // Run seeding + for (const auto [bottom, middle, top] : spacePointsGrouping) + { seeds.clear(); - finder.createSeedsForGroup(state, std::back_inserter(seeds), group.bottom(), - group.middle(), group.top()); + finder.createSeedsForGroup( + finderOpts, state, spacePointsGrouping.grid(), + std::back_inserter(seeds), bottom, middle, top, rMiddleSPRange); // // Loop over seeds and get track parameters paramseeds.clear(); - for (const Acts::Seed &seed : seeds) { - // Get the bottom space point and its reference surface - // @todo do we need to sort the sps first - const ACTSTracking::SeedSpacePoint *bottomSP = seed.sp().front(); - const ACTSTracking::SourceLink &sourceLink = bottomSP->sourceLink(); - const Acts::GeometryIdentifier &geoId = sourceLink.geometryId(); - - const Acts::Surface *surface = trackingGeometry()->findSurface(geoId); + for (const Acts::Seed &seed : seeds) + { + const SSPoint* bottomSP = seed.sp().front(); + + const auto& sourceLink = bottomSP->sourceLink(); + const Acts::GeometryIdentifier& geoId = sourceLink.geometryId(); + const Acts::Surface* surface = trackingGeometry()->findSurface(geoId); if (surface == nullptr) { std::cout << "surface with geoID " << geoId << " is not found in the tracking gemetry"; @@ -493,7 +590,7 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { float p = std::abs(1 / params[Acts::eBoundQOverP]); // build the track covariance matrix using the smearing sigmas - Acts::BoundSymMatrix cov = Acts::BoundSymMatrix::Zero(); + Acts::BoundSquareMatrix cov = Acts::BoundSquareMatrix::Zero(); cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = std::pow(_initialTrackError_pos, 2); cov(Acts::eBoundLoc1, Acts::eBoundLoc1) = @@ -508,7 +605,7 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { std::pow(_initialTrackError_relP * p / (p * p), 2); Acts::BoundTrackParameters paramseed(surface->getSharedPtr(), params, - charge, cov); + cov, Acts::ParticleHypothesis::pion()); paramseeds.push_back(paramseed); // @@ -549,49 +646,43 @@ void ACTSSeededCKFTrackingProc::processEvent(LCEvent *evt) { // // Find the tracks - TrackFinderResultContainer results; - if (_runCKF) - results = trackFinder.findTracks(sourceLinks, paramseeds, ckfOptions); + if (!_runCKF) continue; - for (TrackFinderResult &result : results) { - if (result.ok()) { - const Acts::CombinatorialKalmanFilterResult - &fitOutput = result.value(); - for (const size_t &trackTip : fitOutput.lastMeasurementIndices) { - if (fitOutput.fittedParameters.count(trackTip) == 0) { - streamlog_out(WARNING) << "No fitted track parameters for " - "trajectory with entry index = " - << trackTip << std::endl; - continue; - } + using TrackContainer = Acts::TrackContainer; + auto trackContainer = std::make_shared(); + auto trackStateContainer = std::make_shared(); + TrackContainer tracks(trackContainer, trackStateContainer); + + for (std::size_t iseed = 0; iseed < paramseeds.size(); ++iseed) { + + tracks.clear(); + auto result = trackFinder.findTracks(paramseeds.at(iseed), ckfOptions, tracks); + if (result.ok()) { + const auto& fitOutput = result.value(); + for (const TrackContainer::TrackProxy& trackTip : fitOutput) + { // // Helpful debug output - Acts::MultiTrajectoryHelpers::TrajectoryState trajState = - Acts::MultiTrajectoryHelpers::trajectoryState( - fitOutput.fittedStates, trackTip); streamlog_out(DEBUG) << "Trajectory Summary" << std::endl; streamlog_out(DEBUG) - << "\tchi2Sum " << trajState.chi2Sum << std::endl; + << "\tchi2Sum " << trackTip.chi2() << std::endl; streamlog_out(DEBUG) - << "\tNDF " << trajState.NDF << std::endl; + << "\tNDF " << trackTip.nDoF() << std::endl; streamlog_out(DEBUG) - << "\tnHoles " << trajState.nHoles << std::endl; + << "\tnHoles " << trackTip.nHoles() << std::endl; streamlog_out(DEBUG) - << "\tnMeasurements " << trajState.nMeasurements << std::endl; + << "\tnMeasurements " << trackTip.nMeasurements() << std::endl; streamlog_out(DEBUG) - << "\tnOutliers " << trajState.nOutliers << std::endl; + << "\tnOutliers " << trackTip.nOutliers() << std::endl; streamlog_out(DEBUG) - << "\tnStates " << trajState.nStates << std::endl; - - const Acts::BoundTrackParameters ¶ms = - fitOutput.fittedParameters.at(trackTip); - streamlog_out(DEBUG) << "Fitted Paramemeters" << std::endl - << params << std::endl; + << "\tnStates " << trackTip.nTrackStates() << std::endl; // Make track object EVENT::Track *track = ACTSTracking::ACTS2Marlin_track( - fitOutput, trackTip, magneticField(), magCache); + trackTip, magneticField(), magCache); // Save results trackCollection->addElement(track); diff --git a/src/ACTSTruthCKFTrackingProc.cxx b/src/ACTSTruthCKFTrackingProc.cxx index a1f6e3b..7baaf88 100644 --- a/src/ACTSTruthCKFTrackingProc.cxx +++ b/src/ACTSTruthCKFTrackingProc.cxx @@ -21,6 +21,9 @@ #include #include +#include "Acts/EventData/VectorTrackContainer.hpp" +#include "Acts/EventData/VectorMultiTrajectory.hpp" + using namespace Acts::UnitLiterals; #include "Helpers.hxx" @@ -28,15 +31,14 @@ using namespace Acts::UnitLiterals; #include "SourceLink.hxx" using TrackFinderOptions = - Acts::CombinatorialKalmanFilterOptions; - + Acts::CombinatorialKalmanFilterOptions; +/* using TrackFinderResult = Acts::Result< Acts::CombinatorialKalmanFilterResult>; using TrackFinderResultContainer = std::vector; - +*/ ACTSTruthCKFTrackingProc aACTSTruthCKFTrackingProc; ACTSTruthCKFTrackingProc::ACTSTruthCKFTrackingProc() @@ -128,7 +130,7 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { params[Acts::eBoundQOverP] = mcParticle->getCharge() / p; // build the track covariance matrix using the smearing sigmas - Acts::BoundSymMatrix cov = Acts::BoundSymMatrix::Zero(); + Acts::BoundSquareMatrix cov = Acts::BoundSquareMatrix::Zero(); cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = std::pow(_initialTrackError_d0, 2); cov(Acts::eBoundLoc1, Acts::eBoundLoc1) = @@ -140,8 +142,8 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { cov(Acts::eBoundQOverP, Acts::eBoundQOverP) = std::pow(_initialTrackError_relP * p / (p * p), 2); - Acts::BoundTrackParameters seed(perigeeSurface, params, - mcParticle->getCharge(), cov); + Acts::BoundTrackParameters seed(perigeeSurface, params, cov, + ACTSTracking::convertParticle(mcParticle)); seeds.push_back(seed); } @@ -198,7 +200,7 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { Acts::Vector2 loc = lpResult.value(); - Acts::SymMatrix2 localCov = Acts::SymMatrix2::Zero(); + Acts::SquareMatrix2 localCov = Acts::SquareMatrix2::Zero(); const EVENT::TrackerHitPlane* hitplane = dynamic_cast(hit.second); if (hitplane) { @@ -210,8 +212,9 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { ACTSTracking::SourceLink sourceLink(surface->geometryId(), measurements.size(), hit.second); - ACTSTracking::Measurement meas = Acts::makeMeasurement( - sourceLink, loc, localCov, Acts::eBoundLoc0, Acts::eBoundLoc1); + Acts::SourceLink src_wrap { sourceLink }; + Acts::Measurement meas = Acts::makeMeasurement( + src_wrap, loc, localCov, Acts::eBoundLoc0, Acts::eBoundLoc1); measurements.push_back(meas); sourceLinks.emplace_hint(sourceLinks.end(), sourceLink); @@ -230,7 +233,7 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { using Stepper = Acts::EigenStepper<>; using Navigator = Acts::Navigator; using Propagator = Acts::Propagator; - using CKF = Acts::CombinatorialKalmanFilter; + using CKF = Acts::CombinatorialKalmanFilter; // Configurations Navigator::Config navigatorCfg{trackingGeometry()}; @@ -247,21 +250,44 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { // Set the options Acts::MeasurementSelector::Config measurementSelectorCfg = { {Acts::GeometryIdentifier(), - {_CKF_chi2CutOff, (std::size_t)(_CKF_numMeasurementsCutOff)}}}; + { {}, { _CKF_chi2CutOff }, { (std::size_t)(_CKF_numMeasurementsCutOff) }}}}; Acts::PropagatorPlainOptions pOptions; pOptions.maxSteps = 10000; + Acts::GainMatrixUpdater kfUpdater; + Acts::GainMatrixSmoother kfSmoother; + + Acts::MeasurementSelector measSel { measurementSelectorCfg }; + ACTSTracking::MeasurementCalibrator measCal { measurements }; + Acts::CombinatorialKalmanFilterExtensions + extensions; + extensions.calibrator.connect< + &ACTSTracking::MeasurementCalibrator::calibrate>( + &measCal); + extensions.updater.connect< + &Acts::GainMatrixUpdater::operator()>( + &kfUpdater); + extensions.smoother.connect< + &Acts::GainMatrixSmoother::operator()>( + &kfSmoother); + extensions.measurementSelector + .connect<&Acts::MeasurementSelector::select>( + &measSel); + + using ACTSTracking::SourceLinkAccessor; + SourceLinkAccessor slAccessor; + slAccessor.container = &sourceLinks; + Acts::SourceLinkAccessorDelegate slAccessorDelegate; + slAccessorDelegate.connect<&SourceLinkAccessor::range>(&slAccessor); + // std::unique_ptr // logger=Acts::getDefaultLogger("TrackFitting", // Acts::Logging::Level::VERBOSE); TrackFinderOptions ckfOptions = TrackFinderOptions( geometryContext(), magneticFieldContext(), calibrationContext(), - ACTSTracking::SourceLinkAccessor(), - ACTSTracking::MeasurementCalibrator(std::move(measurements)), - Acts::MeasurementSelector(measurementSelectorCfg), - // Acts::LoggerWrapper{*logger}, pOptions, - Acts::getDummyLogger(), pOptions, &(*perigeeSurface)); + slAccessorDelegate, + extensions, pOptions, perigeeSurface.get()); // // Output @@ -276,6 +302,31 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { // // Find the tracks + using TrackContainer = + Acts::TrackContainer; + auto trackContainer = std::make_shared(); + auto trackStateContainer = std::make_shared(); + TrackContainer tracks(trackContainer, trackStateContainer); + + for (std::size_t iseed = 0; iseed < seeds.size(); ++iseed) { + auto result = trackFinder.findTracks(seeds.at(iseed), ckfOptions, tracks); + if (result.ok()) { + const auto& fitOutput = result.value(); + for (const auto& trackTip : fitOutput) + { + std::cout << trackTip.chi2() << std::endl; + EVENT::Track* track = ACTSTracking::ACTS2Marlin_track( + trackTip, magneticField(), magCache); + trackCollection->addElement(track); + } + } else { + streamlog_out(WARNING) << "Track fit error: " << result.error() << std::endl; + _fitFails++; + } + + } +/* auto t1 = std::chrono::high_resolution_clock::now(); TrackFinderResultContainer results = trackFinder.findTracks(sourceLinks, seeds, ckfOptions); @@ -336,6 +387,7 @@ void ACTSTruthCKFTrackingProc::processEvent(LCEvent* evt) { _fitFails++; } } +*/ // Save the output track collection evt->addCollection(trackCollection, _outputTrackCollection); diff --git a/src/ACTSTruthTrackingProc.cxx b/src/ACTSTruthTrackingProc.cxx index 1b7273f..beb5837 100644 --- a/src/ACTSTruthTrackingProc.cxx +++ b/src/ACTSTruthTrackingProc.cxx @@ -7,6 +7,7 @@ #include #include #include +#include #include #include @@ -14,6 +15,7 @@ #include #include +#include #include #include #include @@ -25,8 +27,10 @@ using namespace Acts::UnitLiterals; #include "MeasurementCalibrator.hxx" #include "SourceLink.hxx" -using TrackFitterResult = - Acts::Result>; +using TrackContainer = + Acts::TrackContainer; +using TrackFitterResult = Acts::Result; // sorting by value of R(=x^2+y^2) in global coordinated so the hits are always // sorted from close to the IP outward @@ -96,7 +100,7 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { using Stepper = Acts::EigenStepper<>; using Navigator = Acts::Navigator; using Propagator = Acts::Propagator; - using Fitter = Acts::KalmanFitter; + using Fitter = Acts::KalmanFitter; // Configurations Navigator::Config navigatorCfg{trackingGeometry()}; @@ -218,7 +222,7 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { // Make container // MeasurementContainer track; - std::vector trackSourceLinks; + std::vector trackSourceLinks; ACTSTracking::MeasurementContainer track; for (EVENT::TrackerHit* hit : trackFilteredByRHits) { // Convert to Acts hit @@ -234,7 +238,7 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { Acts::Vector2 loc = lpResult.value(); - Acts::SymMatrix2 cov = Acts::SymMatrix2::Zero(); + Acts::SquareMatrix2 cov = Acts::SquareMatrix2::Zero(); const EVENT::TrackerHitPlane* hitplane = dynamic_cast(hit); if (hitplane) { @@ -244,8 +248,8 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { throw std::runtime_error("Currently only support TrackerHitPlane."); } - ACTSTracking::SourceLink sourceLink(surface->geometryId(), track.size(), - hit); + ACTSTracking::SourceLink s_link(surface->geometryId(), track.size(), hit); + Acts::SourceLink sourceLink { std::move(s_link) }; ACTSTracking::Measurement meas = Acts::makeMeasurement( sourceLink, loc, cov, Acts::eBoundLoc0, Acts::eBoundLoc1); @@ -260,21 +264,29 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { Acts::Surface::makeShared( Acts::Vector3{0., 0., 0.}); + Updater kfUpdater; + Smoother kfSmoother; + + Acts::KalmanFitterExtensions extensions; + extensions.updater.connect< + &Acts::GainMatrixUpdater::operator()>( + &kfUpdater); + extensions.smoother.connect< + &Acts::GainMatrixSmoother::operator()>( + &kfSmoother); + // Set the KalmanFitter options // std::unique_ptr // logger=Acts::getDefaultLogger("TrackFitting", // Acts::Logging::Level::VERBOSE); - Acts::KalmanFitterOptions - kfOptions = - Acts::KalmanFitterOptions( - geometryContext(), magneticFieldContext(), calibrationContext(), - ACTSTracking::MeasurementCalibrator(track), - Acts::VoidOutlierFinder(), - // Acts::LoggerWrapper{*logger}, Acts::PropagatorPlainOptions(), - Acts::getDummyLogger(), Acts::PropagatorPlainOptions(), - &(*perigeeSurface)); + Acts::KalmanFitterOptions kfOptions { + geometryContext(), + magneticFieldContext(), + calibrationContext(), + extensions, + Acts::PropagatorPlainOptions(), + &(*perigeeSurface) + }; double px = mcParticle->getMomentum()[0]; double py = mcParticle->getMomentum()[1]; @@ -294,7 +306,7 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { params[Acts::eBoundQOverP] = mcParticle->getCharge() / p; // build the track covariance matrix using the smearing sigmas - Acts::BoundSymMatrix cov = Acts::BoundSymMatrix::Zero(); + Acts::BoundSquareMatrix cov = Acts::BoundSquareMatrix::Zero(); cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = std::pow(_initialTrackError_d0, 2); cov(Acts::eBoundLoc1, Acts::eBoundLoc1) = @@ -311,16 +323,23 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { Acts::Vector3(mcParticle->getVertex())); Acts::BoundTrackParameters initialparams(perigeeSurface, params, - mcParticle->getCharge(), cov); + cov, ACTSTracking::convertParticle(mcParticle)); + // reference Examples TruthTracking/ParticleSmearing.cpp streamlog_out(DEBUG) << "Initial Paramemeters" << std::endl << initialparams << std::endl; - TrackFitterResult result = - trackFitter.fit(trackSourceLinks, initialparams, kfOptions); + auto trackContainer = std::make_shared(); + auto trackStateContainer = std::make_shared(); + TrackContainer tracks(trackContainer, trackStateContainer); + + TrackFitterResult result = trackFitter.fit(trackSourceLinks.begin(), + trackSourceLinks.end(), + initialparams, kfOptions, + tracks); if (result.ok()) { - const Acts::KalmanFitterResult& fitOutput = - result.value(); + const auto& fitOutput = result.value(); +/* if (fitOutput.fittedParameters) { // Make the track object and relations object IMPL::LCRelationImpl* relationTrack = new IMPL::LCRelationImpl; @@ -366,6 +385,25 @@ void ACTSTruthTrackingProc::processEvent(LCEvent* evt) { << "No fitted paramemeters for track" << std::endl; _fitFails++; } +*/ + if (fitOutput.hasReferenceSurface()) + { + EVENT::Track* track = ACTSTracking::ACTS2Marlin_track( + fitOutput, magneticField(), magCache); + trackCollection->addElement(track); + + IMPL::LCRelationImpl* relationTrack = new IMPL::LCRelationImpl; + relationTrack->setFrom(track); + relationTrack->setTo(mcParticle); + relationTrack->setWeight(1.0); + trackRelationCollection->addElement(relationTrack); + } + else + { + streamlog_out(WARNING) + << "No fitted paramemeters for track" << std::endl; + _fitFails++; + } } else { streamlog_out(WARNING) << "Track fit error: " << result.error() << std::endl; diff --git a/src/GeometryIdMappingTool.cxx b/src/GeometryIdMappingTool.cxx index 7f9287c..9c76aa3 100644 --- a/src/GeometryIdMappingTool.cxx +++ b/src/GeometryIdMappingTool.cxx @@ -215,7 +215,7 @@ uint64_t GeometryIdMappingTool::getGeometryID(uint32_t systemID, sensitive_id = ladderID; break; } - geometry_id |= sensitive_id << (0 * 4); + geometry_id |= sensitive_id << (2 * 4); return geometry_id; } diff --git a/src/Helpers.cxx b/src/Helpers.cxx index 41e40a0..2eae839 100644 --- a/src/Helpers.cxx +++ b/src/Helpers.cxx @@ -35,7 +35,7 @@ std::string findFile(const std::string& inpath) { // nothing was found :( return inpath; } - +/* EVENT::Track* ACTS2Marlin_track( const Acts::CombinatorialKalmanFilterResult& fitOutput, @@ -225,6 +225,84 @@ EVENT::Track* ACTS2Marlin_track( return track; } +*/ + +EVENT::Track* ACTS2Marlin_track( + const TrackResult& fitter_res, + std::shared_ptr magneticField, + Acts::MagneticFieldProvider::Cache& magCache) +{ + IMPL::TrackImpl* track = new IMPL::TrackImpl; + + track->setChi2(fitter_res.chi2()); + track->setNdf(fitter_res.nDoF()); + + const Acts::Vector3 zeroPos(0, 0, 0); + Acts::Result fieldRes = magneticField->getField(zeroPos, magCache); + if (!fieldRes.ok()) { + throw std::runtime_error("Field lookup error: " + fieldRes.error().value()); + } + Acts::Vector3 field = *fieldRes; + + const Acts::BoundVector& params = fitter_res.parameters(); + const Acts::BoundMatrix& covariance = fitter_res.covariance(); + EVENT::TrackState* trackStateAtIP = ACTSTracking::ACTS2Marlin_trackState( + EVENT::TrackState::AtIP, params, covariance, field[2] / Acts::UnitConstants::T); + track->trackStates().push_back(trackStateAtIP); + + EVENT::TrackerHitVec hitsOnTrack; + EVENT::TrackStateVec statesOnTrack; + + for (const auto& trk_state : fitter_res.trackStatesReversed()) + { + if (!trk_state.hasUncalibratedSourceLink()) continue; + + auto sl = trk_state.getUncalibratedSourceLink() + .get(); + EVENT::TrackerHit* curr_hit = sl.lciohit(); + hitsOnTrack.push_back(curr_hit); + + const Acts::Vector3 hitPos(curr_hit->getPosition()[0], + curr_hit->getPosition()[1], + curr_hit->getPosition()[2]); + + Acts::Result fieldRes = + magneticField->getField(hitPos, magCache); + if (!fieldRes.ok()) { + throw std::runtime_error("Field lookup error: " + + fieldRes.error().value()); + } + Acts::Vector3 field = *fieldRes; + + EVENT::TrackState* trackState = ACTSTracking::ACTS2Marlin_trackState( + EVENT::TrackState::AtOther, trk_state.smoothed(), + trk_state.smoothedCovariance(), field[2] / Acts::UnitConstants::T); + statesOnTrack.push_back(trackState); + } + + std::reverse(hitsOnTrack.begin(), hitsOnTrack.end()); + std::reverse(statesOnTrack.begin(), statesOnTrack.end()); + + for (EVENT::TrackerHit* hit : hitsOnTrack) { + track->addHit(hit); + } + + if (statesOnTrack.size() > 0) { + dynamic_cast(statesOnTrack.back()) + ->setLocation(EVENT::TrackState::AtLastHit); + dynamic_cast(statesOnTrack.front()) + ->setLocation(EVENT::TrackState::AtFirstHit); + } + + EVENT::TrackStateVec& myTrackStates = track->trackStates(); + myTrackStates.insert(myTrackStates.end(), statesOnTrack.begin(), + statesOnTrack.end()); + + return track; +} + + + EVENT::TrackState* ACTS2Marlin_trackState( int location, const Acts::BoundTrackParameters& params, double Bz) { return ACTS2Marlin_trackState(location, params.parameters(), @@ -313,4 +391,45 @@ EVENT::LCCollection* getCollection(EVENT::LCEvent* evt, return nullptr; } } + +Acts::ParticleHypothesis convertParticle(const EVENT::MCParticle* mcParticle) +{ + switch (mcParticle->getPDG()) { + case 11: + return Acts::ParticleHypothesis {Acts::PdgParticle::eElectron}; + case -11: + return Acts::ParticleHypothesis {Acts::PdgParticle::ePositron}; + case 13: + return Acts::ParticleHypothesis {Acts::PdgParticle::eMuon}; + case -13: + return Acts::ParticleHypothesis {Acts::PdgParticle::eAntiMuon}; + case 15: + return Acts::ParticleHypothesis {Acts::PdgParticle::eTau}; + case -15: + return Acts::ParticleHypothesis {Acts::PdgParticle::eAntiTau}; + case 22: + return Acts::ParticleHypothesis {Acts::PdgParticle::eGamma}; + case 111: + return Acts::ParticleHypothesis {Acts::PdgParticle::ePionZero}; + case 211: + return Acts::ParticleHypothesis {Acts::PdgParticle::ePionPlus}; + case -211: + return Acts::ParticleHypothesis {Acts::PdgParticle::ePionMinus}; + case 2112: + return Acts::ParticleHypothesis {Acts::PdgParticle::eNeutron}; + case -2112: + return Acts::ParticleHypothesis {Acts::PdgParticle::eAntiNeutron}; + case 2212: + return Acts::ParticleHypothesis {Acts::PdgParticle::eProton}; + case -2212: + return Acts::ParticleHypothesis {Acts::PdgParticle::eAntiProton}; + } + + Acts::PdgParticle pdg = Acts::PdgParticle::eInvalid; + float mass = 0.0f; + Acts::AnyCharge charge_type { 0.0f }; + return Acts::ParticleHypothesis { pdg, mass, charge_type }; +} + + } // namespace ACTSTracking