Skip to content

Commit 2d959e0

Browse files
committed
Try general radius criterion
1 parent 01d32b0 commit 2d959e0

File tree

8 files changed

+304
-118
lines changed

8 files changed

+304
-118
lines changed

core/CMakeLists.txt

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -92,6 +92,7 @@ traccc_add_library( traccc_core core TYPE SHARED
9292
"include/traccc/fitting/kalman_filter/kalman_actor.hpp"
9393
"include/traccc/fitting/kalman_filter/kalman_fitter.hpp"
9494
"include/traccc/fitting/kalman_filter/kalman_step_aborter.hpp"
95+
"include/traccc/fitting/kalman_filter/measurement_selector.hpp"
9596
"include/traccc/fitting/kalman_filter/statistics_updater.hpp"
9697
"include/traccc/fitting/kalman_filter/two_filters_smoother.hpp"
9798
"include/traccc/fitting/details/kalman_fitting_types.hpp"

core/include/traccc/finding/details/combinatorial_kalman_filter.hpp

Lines changed: 29 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@
1717
#include "traccc/finding/finding_config.hpp"
1818
#include "traccc/fitting/kalman_filter/gain_matrix_updater.hpp"
1919
#include "traccc/fitting/kalman_filter/is_line_visitor.hpp"
20+
#include "traccc/fitting/kalman_filter/measurement_selector.hpp"
2021
#include "traccc/fitting/status_codes.hpp"
2122
#include "traccc/sanity/contiguous_on.hpp"
2223
#include "traccc/utils/logging.hpp"
@@ -248,6 +249,9 @@ combinatorial_kalman_filter(
248249
bound_track_parameters<algebra_type>>>
249250
best_links;
250251

252+
const bool is_line = detail::is_line(sf);
253+
measurement_selector::config calib_cfg{};
254+
251255
// Iterate over the measurements
252256
TRACCC_VERBOSE_HOST("No. measurements: " << (up - lo));
253257
for (unsigned int meas_id = lo; meas_id < up; meas_id++) {
@@ -256,63 +260,53 @@ combinatorial_kalman_filter(
256260
// The measurement on surface to handle.
257261
const auto meas = measurements.at(meas_id);
258262

263+
const scalar_type chi2 = measurement_selector::predicted_chi2(
264+
meas, in_param, calib_cfg, is_line);
265+
266+
// If the measurement is outside the chi2 cut, skip it
267+
if (chi2 > config.chi2_max) {
268+
continue;
269+
}
270+
259271
// Create a standalone track state object.
260272
auto trk_state =
261273
edm::make_track_state<algebra_type>(measurements, meas_id);
262274

275+
// Kalman filter status code
263276
kalman_fitter_status res{kalman_fitter_status::ERROR_OTHER};
264-
// Check if the measurement is even remotely close to the track
277+
278+
// Don't run the filter on the first measurement
265279
if (step == 0 && !sf.has_material()) {
266-
detray::dmatrix<algebra_type, 2, 1> meas_local;
267-
edm::get_measurement_local<algebra_type>(meas, meas_local);
268-
269-
const subspace<algebra_type, e_bound_size> subs(
270-
meas.subspace());
271-
detray::dmatrix<algebra_type, 2, e_bound_size> H =
272-
subs.template projector<2>();
273-
const auto residual = meas_local - H * in_param.vector();
274-
275-
const scalar_type res_0{getter::element(residual, 0, 0)};
276-
const scalar_type res_1{
277-
meas.dimensions() == 1
278-
? 0.f
279-
: getter::element(residual, 1, 0)};
280-
const scalar_type search_rad =
281-
math::sqrt(res_0 * res_0 + res_1 * res_1);
282-
if (search_rad > 0.f) {
283-
TRACCC_VERBOSE_HOST(
284-
"Measurement outside search radius: "
285-
<< search_rad << " mm for meas. " << meas_id);
286-
continue;
287-
}
288280
res = kalman_fitter_status::SUCCESS;
281+
282+
// Copy the full track parameters
283+
// TODO: Apply calibration ?
289284
trk_state.filtered_params() = in_param;
290-
trk_state.filtered_chi2() = 0.f;
291285

292-
// Set the covariance to the measurement variance
293-
detray::dmatrix<algebra_type, 2, 2> V;
294-
edm::get_measurement_covariance<algebra_type>(meas, V);
286+
// Update measurement covariance
287+
const auto V =
288+
measurement_selector::calibrated_measurement_covariance<
289+
algebra_type, 2>(meas, calib_cfg);
290+
295291
auto& filtered_cov =
296292
trk_state.filtered_params().covariance();
297293
getter::element(filtered_cov, e_bound_loc0, e_bound_loc0) =
298294
getter::element(V, 0, 0);
299295
getter::element(filtered_cov, e_bound_loc1, e_bound_loc1) =
300296
getter::element(V, 1, 1);
301297
} else {
302-
const bool is_line = detail::is_line(sf);
303-
304-
// Run the Kalman update on a copy of the track parameters
305-
res = gain_matrix_updater<algebra_type>{}(
306-
trk_state, measurements, in_param, is_line);
298+
// Run the Kalman update on the track state
299+
constexpr gain_matrix_updater<algebra_type>
300+
kalman_updater{};
301+
res = kalman_updater(trk_state, meas, in_param, is_line);
307302
}
308303

309-
const traccc::scalar chi2 = trk_state.filtered_chi2();
304+
trk_state.filtered_chi2() = chi2;
310305

311306
TRACCC_DEBUG_HOST("KF status: " << fitter_debug_msg{res}());
312307

313308
// The chi2 from Kalman update should be less than chi2_max
314-
if (res == kalman_fitter_status::SUCCESS &&
315-
(chi2 < config.chi2_max)) {
309+
if (res == kalman_fitter_status::SUCCESS) {
316310

317311
TRACCC_VERBOSE_HOST("Found measurement: " << meas_id);
318312

core/include/traccc/fitting/fitting_config.hpp

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,10 @@ struct fitting_config {
2929
float min_p = 100.f * traccc::unit<float>::MeV;
3030
float min_pT = 600.f * traccc::unit<float>::MeV;
3131

32+
// Maximal local radius a measurement is allowed to be away from the
33+
// predicted track position in order to be considered as a condidate
34+
float max_measurement_radius{5.f * traccc::unit<float>::mm};
35+
3236
/// Particle hypothesis
3337
traccc::pdg_particle<traccc::scalar> ptc_hypothesis =
3438
traccc::muon<traccc::scalar>();

core/include/traccc/fitting/kalman_filter/gain_matrix_updater.hpp

Lines changed: 7 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -11,12 +11,10 @@
1111
#include "traccc/definitions/qualifiers.hpp"
1212
#include "traccc/definitions/track_parametrization.hpp"
1313
#include "traccc/edm/measurement_collection.hpp"
14-
#include "traccc/edm/measurement_helpers.hpp"
15-
#include "traccc/edm/track_state_collection.hpp"
1614
#include "traccc/fitting/details/regularize_covariance.hpp"
15+
#include "traccc/fitting/kalman_filter/measurement_selector.hpp"
1716
#include "traccc/fitting/status_codes.hpp"
1817
#include "traccc/utils/logging.hpp"
19-
#include "traccc/utils/subspace.hpp"
2018

2119
namespace traccc {
2220

@@ -42,18 +40,16 @@ struct gain_matrix_updater {
4240
/// @param bound_params bound parameter
4341
///
4442
/// @return true if the update succeeds
45-
template <typename track_state_backend_t>
43+
template <typename track_state_backend_t, typename measurement_backend_t>
4644
[[nodiscard]] TRACCC_HOST_DEVICE inline kalman_fitter_status operator()(
4745
typename edm::track_state<track_state_backend_t>& trk_state,
48-
const edm::measurement_collection<default_algebra>::const_device&
49-
measurements,
46+
const edm::measurement<measurement_backend_t>& measurement,
5047
const bound_track_parameters<algebra_t>& bound_params,
5148
const bool is_line) const {
5249

5350
static constexpr unsigned int D = 2;
5451

55-
[[maybe_unused]] const unsigned int dim{
56-
measurements.at(trk_state.measurement_index()).dimensions()};
52+
[[maybe_unused]] const unsigned int dim{measurement.dimensions()};
5753

5854
TRACCC_VERBOSE_HOST_DEVICE("In gain-matrix-updater...");
5955
TRACCC_VERBOSE_HOST_DEVICE("Measurement dim: %d", dim);
@@ -70,8 +66,7 @@ struct gain_matrix_updater {
7066

7167
// Measurement data on surface
7268
matrix_type<D, 1> meas_local;
73-
edm::get_measurement_local<algebra_t>(
74-
measurements.at(trk_state.measurement_index()), meas_local);
69+
edm::get_measurement_local<algebra_t>(measurement, meas_local);
7570

7671
assert((dim > 1) || (getter::element(meas_local, 1u, 0u) == 0.f));
7772

@@ -83,8 +78,7 @@ struct gain_matrix_updater {
8378
// Predicted covaraince of bound track parameters
8479
const bound_matrix_type& predicted_cov = bound_params.covariance();
8580

86-
const subspace<algebra_t, e_bound_size> subs(
87-
measurements.at(trk_state.measurement_index()).subspace());
81+
const subspace<algebra_t, e_bound_size> subs(measurement.subspace());
8882
matrix_type<D, e_bound_size> H = subs.template projector<D>();
8983

9084
// Flip the sign of projector matrix element in case the first element
@@ -101,8 +95,7 @@ struct gain_matrix_updater {
10195

10296
// Spatial resolution (Measurement covariance)
10397
matrix_type<D, D> V;
104-
edm::get_measurement_covariance<algebra_t>(
105-
measurements.at(trk_state.measurement_index()), V);
98+
edm::get_measurement_covariance<algebra_t>(measurement, V);
10699
// @TODO: Fix properly
107100
if (/*dim == 1*/ getter::element(meas_local, 1u, 0u) == 0.f) {
108101
getter::element(V, 1u, 1u) = 1000.f;

core/include/traccc/fitting/kalman_filter/kalman_actor.hpp

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -406,6 +406,8 @@ struct kalman_actor : detray::actor {
406406
const auto sf = navigation.current_surface();
407407
const bool is_line = detail::is_line(sf);
408408

409+
const auto measurement =
410+
actor_state.m_measurements.at(trk_state.measurement_index());
409411
if (!actor_state.backward_mode) {
410412
if constexpr (direction_e ==
411413
kalman_actor_direction::FORWARD_ONLY ||
@@ -417,8 +419,7 @@ struct kalman_actor : detray::actor {
417419
// Forward filter
418420
TRACCC_DEBUG_HOST_DEVICE("Run filtering...");
419421
actor_state.fit_result = gain_matrix_updater<algebra_t>{}(
420-
trk_state, actor_state.m_measurements, bound_param,
421-
is_line);
422+
trk_state, measurement, bound_param, is_line);
422423

423424
// Update the propagation flow
424425
bound_param = trk_state.filtered_params();
@@ -447,8 +448,7 @@ struct kalman_actor : detray::actor {
447448
} else {
448449
actor_state.fit_result =
449450
two_filters_smoother<algebra_t>{}(
450-
trk_state, actor_state.m_measurements,
451-
bound_param, is_line);
451+
trk_state, measurement, bound_param, is_line);
452452
}
453453
} else {
454454
assert(false);

0 commit comments

Comments
 (0)