Skip to content

Commit

Permalink
add Trackstate filling
Browse files Browse the repository at this point in the history
  • Loading branch information
@AJPfleger
AJPfleger committed Oct 4, 2023
1 parent 691f1d5 commit c6c17e7
Show file tree
Hide file tree
Showing 2 changed files with 134 additions and 36 deletions.
156 changes: 127 additions & 29 deletions Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,6 @@ struct Gx2FitterExtensions {
/// Retrieves the associated surface from a source link
SourceLinkSurfaceAccessor surfaceAccessor;

// @TODO get an own Calibrator and Updater
/// Default constructor which connects the default void components
Gx2FitterExtensions() {
calibrator.template connect<&detail::voidFitterCalibrator<traj_t>>();
Expand Down Expand Up @@ -241,6 +240,8 @@ void collector(typename traj_t::TrackStateProxy& trackStateProxy,
auto projPredicted =
(trackStateProxy.effectiveProjector() * predicted).eval();

/// TODO change effectiveProjector because Eigen puts memory on the heap (ask
/// Paul)
ACTS_VERBOSE("Processing and collecting measurements in Actor:\n"
<< "\tMeasurement:\t" << measurement.transpose()
<< "\n\tPredicted:\t" << predicted.transpose()
Expand Down Expand Up @@ -345,6 +346,9 @@ class Gx2Fitter {
/// Calibration context for the fit
const CalibrationContext* calibrationContext{nullptr};

/// The number of the current iteration of the update
size_t nUpdate = Acts::MultiTrajectoryTraits::kInvalid;

/// @brief Gx2f actor operation
///
/// @tparam propagator_state_t is the type of Propagator state
Expand Down Expand Up @@ -391,38 +395,122 @@ class Gx2Fitter {
auto sourcelink_it = inputMeasurements->find(surface->geometryId());

if (sourcelink_it != inputMeasurements->end()) {
ACTS_VERBOSE("Measurement surface " << surface->geometryId()
<< " detected.");

if (surface and surface->surfaceMaterial()) {
ACTS_VERBOSE(
"dbgActor: materialInteractor PREupdate could be needed.")
}

// Transport the covariance to the surface
stepper.transportCovarianceToBound(state.stepping, *surface,
freeToBoundCorrection);
auto res = stepper.boundState(state.stepping, *surface, false,
freeToBoundCorrection);
if (!res.ok()) {
std::cout << "dbgActor: res = stepper.boundState res not ok"
<< std::endl;
return;

std::cout << "actorPRE - result.lastMeasurementIndex: "
<< result.lastMeasurementIndex << std::endl;
std::cout << "SIZE - result.fittedStates: "
<< result.fittedStates->size() << std::endl;
std::cout << "actorPRE - result.lastTrackIndex: "
<< result.lastTrackIndex << std::endl;

size_t currentTrackIndex = Acts::MultiTrajectoryTraits::kInvalid;

std::cout << " processSurface: result.lastTrackIndex = "
<< result.lastTrackIndex << std::endl;
// Checks if during the update an existing surface is found.
// If not, there will be a new index generated afterwards
if (nUpdate == 0) {
ACTS_VERBOSE(" processSurface: nUpdate == 0 decision");

// add a full TrackState entry multi trajectory
// (this allocates storage for all components, we will set them
// later)
currentTrackIndex = result.fittedStates->addTrackState(
~(TrackStatePropMask::Smoothed | TrackStatePropMask::Filtered),
result.lastTrackIndex);
} else {
ACTS_VERBOSE(" processSurface: nUpdate > 0 decision");

if (result.lastTrackIndex ==
Acts::MultiTrajectoryTraits::kInvalid) {
currentTrackIndex = 0;
ACTS_VERBOSE(" processSurface: currentTrackIndex (kInv->0) = "
<< currentTrackIndex);
} else if (result.lastTrackIndex <
result.fittedStates->size() - 1) {
currentTrackIndex = result.lastTrackIndex + 1;
ACTS_VERBOSE(" processSurface: currentTrackIndex (n+1) = "
<< currentTrackIndex);
} else {
// add a full TrackState entry multi trajectory
// (this allocates storage for all components, we will set them
// later)
currentTrackIndex = result.fittedStates->addTrackState(
~(TrackStatePropMask::Smoothed |
TrackStatePropMask::Filtered),
result.lastTrackIndex);
ACTS_VERBOSE(" processSurface: currentTrackIndex (ADD NEW)= "
<< currentTrackIndex);
}
}
auto& [boundParams, jacobian, pathLength] = *res;
result.jacobianFromStart = jacobian * result.jacobianFromStart;

// add a full TrackState entry multi trajectory
// (this allocates storage for all components, we will set them later)
auto fittedStates = *result.fittedStates;
const auto newTrackIndex = fittedStates.addTrackState(
TrackStatePropMask::All, result.lastTrackIndex);
auto& fittedStates = *result.fittedStates;
/// use trackstate API

/// START of getTrackStateProxy()
// add a full TrackState entry multi trajectory
// (this allocates storage for all components, we will set them later)
/// const auto newTrackIndex = fittedStates.addTrackState(mask,
/// lastTrackIndex);

// now get track state proxy back
auto trackStateProxy = fittedStates.getTrackState(newTrackIndex);
auto trackStateProxy = fittedStates.getTrackState(currentTrackIndex);

trackStateProxy.setReferenceSurface(surface->getSharedPtr());
// assign the source link to the track state
trackStateProxy.setUncalibratedSourceLink(sourcelink_it->second);

// Bind the transported state to the current surface
auto res = stepper.boundState(state.stepping, *surface, false,
freeToBoundCorrection);
if (!res.ok()) {
ACTS_ERROR("Propagate to surface " << surface->geometryId()
<< " failed: " << res.error());
return;
}
auto& [boundParams, jacobian, pathLength] = *res;

// Fill the track state
trackStateProxy.predicted() = std::move(boundParams.parameters());
/// Do we need this block? it was a conflict
if (boundParams.covariance().has_value()) {
trackStateProxy.predictedCovariance() =
std::move(*boundParams.covariance());
}
trackStateProxy.jacobian() = std::move(jacobian);
trackStateProxy.pathLength() = std::move(pathLength);
/// END of getTrackStateProxy()

// We have predicted parameters, so calibrate the uncalibrated input
// measurement
extensions.calibrator(state.geoContext, *calibrationContext,
sourcelink_it->second, trackStateProxy);

// Get and set the type flags
auto typeFlags = trackStateProxy.typeFlags();
typeFlags.set(TrackStateFlag::ParameterFlag);
if (surface->surfaceMaterial() != nullptr) {
typeFlags.set(TrackStateFlag::MaterialFlag);
}

result.jacobianFromStart =
trackStateProxy.jacobian() * result.jacobianFromStart;

// Collect:
// - Residuals
// - Covariances
// - ProjectedJacobians
if (trackStateProxy.calibratedSize() == 1) {
collector<1>(trackStateProxy, result, *actorLogger);
} else if (trackStateProxy.calibratedSize() == 2) {
Expand All @@ -432,17 +520,26 @@ class Gx2Fitter {
"Only measurements of 1 and 2 dimensions are implemented yet.");
}

if (boundParams.covariance().has_value()) {
trackStateProxy.predictedCovariance() =
std::move(*boundParams.covariance());
}

trackStateProxy.jacobian() = std::move(jacobian);
trackStateProxy.pathLength() = std::move(pathLength);
// Set the measurement type flag
typeFlags.set(TrackStateFlag::MeasurementFlag);
// We count the processed state
++result.processedStates;
std::cout << "actor - result.lastMeasurementIndex: "
<< result.lastMeasurementIndex << std::endl;
std::cout << "actor - trackStateProxy.index(): "
<< trackStateProxy.index() << std::endl;
std::cout << "actor - result.lastTrackIndex: "
<< result.lastTrackIndex << std::endl;
std::cout << "actor - currentTrackIndex: " << currentTrackIndex
<< std::endl;
result.lastMeasurementIndex = currentTrackIndex;
result.lastTrackIndex = currentTrackIndex;
} else {
ACTS_INFO("Actor: This case is not implemented yet")
}
}

if (result.surfaceCount > 17) {
if (result.surfaceCount > 900) {
ACTS_INFO("Actor: finish due to limit. Result might be garbage.");
result.finished = true;
}
Expand Down Expand Up @@ -524,12 +621,14 @@ class Gx2Fitter {
using PropagatorOptions = Acts::PropagatorOptions<Actors, Aborters>;

const size_t reducedMatrixSize = 4;
Acts::CurvilinearTrackParameters params = sParameters;
start_parameters_t params = sParameters;
BoundVector deltaParams = BoundVector::Zero();
double chi2sum = 0;
BoundMatrix aMatrix = BoundMatrix::Zero();
BoundVector bVector = BoundVector::Zero();

GX2FResult gx2fResult;

ACTS_VERBOSE("params:\n" << params);

/// Actual Fitting /////////////////////////////////////////////////////////
Expand All @@ -555,6 +654,7 @@ class Gx2Fitter {
gx2fActor.extensions = gx2fOptions.extensions;
gx2fActor.calibrationContext = &gx2fOptions.calibrationContext.get();
gx2fActor.actorLogger = m_actorLogger.get();
gx2fActor.nUpdate = nUpdate;

typename propagator_t::template action_list_t_result_t<
CurvilinearTrackParameters, Actors>
Expand All @@ -570,7 +670,7 @@ class Gx2Fitter {
// TODO Improve Propagator + Actor [allocate before loop], rewrite
// makeMeasurements
auto& propRes = *result;
auto gx2fResult = std::move(propRes.template get<GX2FResult>());
gx2fResult = std::move(propRes.template get<GX2FResult>());

ACTS_VERBOSE("gx2fResult.collectorResiduals.size() = "
<< gx2fResult.collectorResiduals.size());
Expand Down Expand Up @@ -642,12 +742,10 @@ class Gx2Fitter {

// Prepare track for return
auto track = trackContainer.getTrack(trackContainer.addTrack());
track.tipIndex() = gx2fResult.lastMeasurementIndex;
track.parameters() = params.parameters();
track.covariance() = fullCovariancePredicted;
// TODO track.tipIndex() = gx2fResult.lastMeasurementIndex;
// TODO track.setReferenceSurface(params.referenceSurface().getSharedPtr());
// TODO track.nMeasurements() = gx2fResult.measurementStates;
// TODO track.nHoles() = gx2fResult.measurementHoles;
track.setReferenceSurface(params.referenceSurface().getSharedPtr());
calculateTrackQuantities(
track); // TODO write test for calculateTrackQuantities

Expand Down
14 changes: 7 additions & 7 deletions Tests/UnitTests/Core/TrackFitting/Gx2fTests.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -246,7 +246,7 @@ BOOST_AUTO_TEST_CASE(NoFit) {

auto& track = *res;
BOOST_CHECK_EQUAL(track.tipIndex(), Acts::MultiTrajectoryTraits::kInvalid);
BOOST_CHECK(!track.hasReferenceSurface());
BOOST_CHECK(track.hasReferenceSurface());
BOOST_CHECK_EQUAL(track.nMeasurements(), 0u);
BOOST_CHECK_EQUAL(track.nHoles(), 0u);
BOOST_CHECK_EQUAL(track.parameters(), startParametersFit.parameters());
Expand Down Expand Up @@ -341,9 +341,9 @@ BOOST_AUTO_TEST_CASE(Fit5Iterations) {
BOOST_REQUIRE(res.ok());

auto& track = *res;
BOOST_CHECK_EQUAL(track.tipIndex(), Acts::MultiTrajectoryTraits::kInvalid);
BOOST_CHECK(!track.hasReferenceSurface());
BOOST_CHECK_EQUAL(track.nMeasurements(), 0u);
BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
BOOST_CHECK(track.hasReferenceSurface());
BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
BOOST_CHECK_EQUAL(track.nHoles(), 0u);
// We need quite coarse checks here, since on different builds
// the created measurements differ in the randomness
Expand Down Expand Up @@ -451,9 +451,9 @@ BOOST_AUTO_TEST_CASE(MixedDetector) {
BOOST_REQUIRE(res.ok());

auto& track = *res;
BOOST_CHECK_EQUAL(track.tipIndex(), Acts::MultiTrajectoryTraits::kInvalid);
BOOST_CHECK(!track.hasReferenceSurface());
BOOST_CHECK_EQUAL(track.nMeasurements(), 0u);
BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
BOOST_CHECK(track.hasReferenceSurface());
BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
BOOST_CHECK_EQUAL(track.nHoles(), 0u);
// We need quite coarse checks here, since on different builds
// the created measurements differ in the randomness
Expand Down

0 comments on commit c6c17e7

Please sign in to comment.