diff --git a/CI/physmon/reference/performance_gx2f.root b/CI/physmon/reference/performance_gx2f.root index 3f9c4482d04..0c1c3facef3 100644 Binary files a/CI/physmon/reference/performance_gx2f.root and b/CI/physmon/reference/performance_gx2f.root differ diff --git a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp index 51be82c6885..f6acef88ddd 100644 --- a/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp +++ b/Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp @@ -399,6 +399,11 @@ class Gx2Fitter { // End the propagation and return to the fitter if (result.finished) { + // Remove the missing surfaces that occur after the last measurement + if (result.measurementStates > 0) { + result.missedActiveSurfaces.resize(result.measurementHoles); + } + return; } @@ -415,15 +420,13 @@ class Gx2Fitter { // Update: // - Waiting for a current surface auto surface = navigator.currentSurface(state.navigation); - if (surface != nullptr && - surface->associatedDetectorElement() != nullptr) { + if (surface != nullptr) { ++result.surfaceCount; ACTS_VERBOSE("Surface " << surface->geometryId() << " detected."); - // check if measurement surface - auto sourcelink_it = inputMeasurements->find(surface->geometryId()); - - if (sourcelink_it != inputMeasurements->end()) { + // Check if we have a measurement surface + if (auto sourcelink_it = inputMeasurements->find(surface->geometryId()); + sourcelink_it != inputMeasurements->end()) { ACTS_VERBOSE("Measurement surface " << surface->geometryId() << " detected."); @@ -528,6 +531,107 @@ class Gx2Fitter { // We count the processed state ++result.processedStates; + + // Update the number of holes count only when encountering a + // measurement + result.measurementHoles = result.missedActiveSurfaces.size(); + } else if (surface->associatedDetectorElement() != nullptr || + surface->surfaceMaterial() != nullptr) { + // Here we handle material and holes + // TODO add material handling + ACTS_VERBOSE("Non-Measurement surface " << surface->geometryId() + << " detected."); + + // We only create track states here if there is already a measurement + // detected or if the surface has material (no holes before the first + // measurement) + if (result.measurementStates > 0 + // || surface->surfaceMaterial() != nullptr + ) { + ACTS_VERBOSE("Handle hole."); + + auto& fittedStates = *result.fittedStates; + + // Mask for the track states. We don't need Smoothed and Filtered + TrackStatePropMask mask = TrackStatePropMask::Predicted | + TrackStatePropMask::Jacobian | + TrackStatePropMask::Calibrated; + + ACTS_VERBOSE(" processSurface: addTrackState"); + + // Add a TrackState entry multi trajectory. This allocates + // storage for all components, which we will set later. + typename traj_t::TrackStateProxy trackStateProxy = + fittedStates.makeTrackState(mask, result.lastTrackIndex); + std::size_t currentTrackIndex = trackStateProxy.index(); + { + // Set the trackStateProxy components with the state from the + // ongoing propagation + { + trackStateProxy.setReferenceSurface(surface->getSharedPtr()); + // Bind the transported state to the current surface + auto res = stepper.boundState(state.stepping, *surface, false, + freeToBoundCorrection); + if (!res.ok()) { + result.result = res.error(); + return; + } + const auto& [boundParams, jacobian, pathLength] = *res; + + // Fill the track state + trackStateProxy.predicted() = boundParams.parameters(); + trackStateProxy.predictedCovariance() = state.stepping.cov; + + trackStateProxy.jacobian() = jacobian; + trackStateProxy.pathLength() = pathLength; + } + + // Get and set the type flags + auto typeFlags = trackStateProxy.typeFlags(); + typeFlags.set(TrackStateFlag::ParameterFlag); + if (surface->surfaceMaterial() != nullptr) { + typeFlags.set(TrackStateFlag::MaterialFlag); + } + + // Set hole only, if we are on a sensitive surface + if (surface->associatedDetectorElement() != nullptr) { + ACTS_VERBOSE("Detected hole on " << surface->geometryId()); + // If the surface is sensitive, set the hole type flag + typeFlags.set(TrackStateFlag::HoleFlag); + } else { + ACTS_VERBOSE("Detected in-sensitive surface " + << surface->geometryId()); + } + } + + ACTS_VERBOSE( + "Actor - indices after processing, before over writing:" + << "\n " + << "result.lastMeasurementIndex: " + << result.lastMeasurementIndex << "\n " + << "trackStateProxy.index(): " << trackStateProxy.index() + << "\n " + << "result.lastTrackIndex: " << result.lastTrackIndex + << "\n " + << "currentTrackIndex: " << currentTrackIndex) + result.lastTrackIndex = currentTrackIndex; + + if (trackStateProxy.typeFlags().test(TrackStateFlag::HoleFlag)) { + // Count the missed surface + result.missedActiveSurfaces.push_back(surface); + } + + ++result.processedStates; + } else { + ACTS_VERBOSE("Ignoring hole, because no preceding measurements."); + } + + if (surface->surfaceMaterial() != nullptr) { + // TODO write similar to KF? + // Update state and stepper with material effects + // materialInteractor(surface, state, stepper, navigator, + // MaterialUpdateStage::FullUpdate); + } } else { ACTS_INFO("Actor: This case is not implemented yet") } diff --git a/Tests/UnitTests/Core/TrackFitting/Gx2fTests.cpp b/Tests/UnitTests/Core/TrackFitting/Gx2fTests.cpp index 52b836d0fc7..a896d7b96a2 100644 --- a/Tests/UnitTests/Core/TrackFitting/Gx2fTests.cpp +++ b/Tests/UnitTests/Core/TrackFitting/Gx2fTests.cpp @@ -802,5 +802,116 @@ BOOST_AUTO_TEST_CASE(NotEnoughMeasurements) { ACTS_INFO("*** Test: NotEnoughMeasurements -- Finish"); } + +BOOST_AUTO_TEST_CASE(FindHoles) { + ACTS_INFO("*** Test: FindHoles -- Start"); + + std::default_random_engine rng(42); + + ACTS_DEBUG("Create the detector"); + // const std::size_t nSurfaces = 7; + const std::size_t nSurfaces = 8; + Detector detector; + detector.geometry = makeToyDetector(geoCtx, nSurfaces); + + ACTS_DEBUG("Set the start parameters for measurement creation and fit"); + const auto parametersMeasurements = makeParameters(); + const auto startParametersFit = makeParameters( + 7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e); + + ACTS_DEBUG("Create the measurements"); + using SimPropagator = + Acts::Propagator; + const SimPropagator simPropagator = makeStraightPropagator(detector.geometry); + const auto measurements = + createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements, + resMapAllPixel, rng); + auto sourceLinks = prepareSourceLinks(measurements.sourceLinks); + ACTS_VERBOSE("sourceLinks.size() [before] = " << sourceLinks.size()); + + // We remove the first measurement in the list. This does not create a hole. + sourceLinks.erase(std::next(sourceLinks.begin(), 0)); + ACTS_VERBOSE( + "sourceLinks.size() [after first erase] = " << sourceLinks.size()); + + // We remove the last measurement in the list. This does not create a hole. + sourceLinks.pop_back(); + ACTS_VERBOSE("sourceLinks.size() [after pop] = " << sourceLinks.size()); + + // We remove the second to last measurement in the list. This effectivly + // creates a hole on that surface. + const std::size_t indexHole = sourceLinks.size() - 2; + ACTS_VERBOSE("Remove measurement " << indexHole); + sourceLinks.erase(std::next(sourceLinks.begin(), indexHole)); + ACTS_VERBOSE("sourceLinks.size() [after second-to-last erase]= " + << sourceLinks.size()); + + // We removed 3 measurements + // const std::size_t nMeasurements = nSurfaces - 2; + const std::size_t nMeasurements = nSurfaces - 3; + BOOST_REQUIRE_EQUAL(sourceLinks.size(), nMeasurements); + + ACTS_DEBUG("Set up the fitter"); + const Surface* rSurface = ¶metersMeasurements.referenceSurface(); + + using RecoStepper = EigenStepper<>; + const auto recoPropagator = + makeConstantFieldPropagator(detector.geometry, 0_T); + + using RecoPropagator = decltype(recoPropagator); + using Gx2Fitter = + Experimental::Gx2Fitter; + const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone()); + + Experimental::Gx2FitterExtensions extensions; + extensions.calibrator + .connect<&testSourceLinkCalibrator>(); + TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry}; + extensions.surfaceAccessor + .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor); + + const Experimental::Gx2FitterOptions gx2fOptions( + geoCtx, magCtx, calCtx, extensions, PropagatorPlainOptions(), rSurface, + false, false, FreeToBoundCorrection(false), 20, true, 1e-5); + + Acts::TrackContainer tracks{Acts::VectorTrackContainer{}, + Acts::VectorMultiTrajectory{}}; + + ACTS_DEBUG("Fit the track"); + ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements); + ACTS_VERBOSE("startParameter fit:\n" << startParametersFit); + const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(), + startParametersFit, gx2fOptions, tracks); + + BOOST_REQUIRE(res.ok()); + + const auto& track = *res; + + // -1, because the index starts at 0 + // -2, because the first and the last surface are not part of the track + BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1 - 2); + BOOST_CHECK(track.hasReferenceSurface()); + + // Track quantities + CHECK_CLOSE_ABS(track.chi2(), 6.5, 2.); + BOOST_CHECK_EQUAL(track.nDoF(), 10u); + BOOST_CHECK_EQUAL(track.nHoles(), 1u); + BOOST_CHECK_EQUAL(track.nMeasurements(), nMeasurements); + BOOST_CHECK_EQUAL(track.nSharedHits(), 0u); + BOOST_CHECK_EQUAL(track.nOutliers(), 0u); + + // Parameters + // We need quite coarse checks here, since on different builds + // the created measurements differ in the randomness + BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 7e0); + BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0); + BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1e3); + BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], M_PI / 2, 1e-3); + BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1); + BOOST_CHECK_CLOSE(track.parameters()[eBoundTime], 12591.2832360000, 1e-6); + BOOST_CHECK_CLOSE(track.covariance().determinant(), 4.7e-28, 2e0); + + ACTS_INFO("*** Test: FindHoles -- Finish"); +} BOOST_AUTO_TEST_SUITE_END() } // namespace Acts::Test