16 #include "Acts/Utilities/Result.hpp"
17 #include "Acts/Utilities/Delegate.hpp"
18 #include "Acts/EventData/SourceLink.hpp"
19 #include "Acts/TrackFinding/CombinatorialKalmanFilterError.hpp"
20 #include "Acts/Definitions/Algebra.hpp"
21 #include "Acts/Surfaces/Surface.hpp"
22 #include "Acts/Geometry/GeometryHierarchyMap.hpp"
23 #include "Acts/EventData/Types.hpp"
24 #include "Acts/EventData/TrackStatePropMask.hpp"
25 #include "boost/container/small_vector.hpp"
28 #include "Acts/EventData/TrackStateProxy.hpp"
29 #include "Acts/Utilities/CalibrationContext.hpp"
30 #include "Acts/EventData/TrackParameters.hpp"
33 #include "Acts/EventData/MultiTrajectory.hpp"
36 #include <type_traits>
42 template <
typename derived_t>
46 template <std::
size_t N>
50 template <std::
size_t N>
55 template <std::
size_t N>
60 template <std::
size_t N>
64 template <std::
size_t N>
65 using Predicted =
typename Acts::detail_lt::Types<N>::Coefficients;
68 template <std::
size_t N>
72 template <
typename T_Container>
87 using BoundState = std::tuple<BoundTrackParameters, BoundMatrix, double>;
101 template <
class T_Matrix>
102 static constexpr std::size_t
matrixColumns() {
return T_Matrix::ColsAtCompileTime;}
103 template <
class T_Matrix>
104 static constexpr std::size_t
matrixRows() {
return T_Matrix::RowsAtCompileTime;}
106 template <
typename T_Float,
class T_Matrix>
109 template <
class T_Matrix>
112 template <
class T_Matrix>
126 template <std::
size_t N>
127 using type = std::array<unsigned char, N>;
129 template <std::
size_t N>
132 for(
int i=0;
i<
N; ++
i) {
141 template <std::
size_t N,
class T_ResultType,
class T_Matrix>
144 using MatrixIndexMapType =
unsigned char;
146 using MatrixIndexType =
unsigned int;
153 if constexpr(MeasurementSelectorMatrixTraits::matrixColumns<T_Matrix>() == 1) {
155 for (MatrixIndexType meas_i=0; meas_i<
N; ++meas_i) {
156 assert( meas_i < parameter_map.size() );
157 ret(meas_i,0) =
matrix( parameter_map[meas_i], 0);
163 for (MatrixIndexType meas_i=0; meas_i<
N; ++meas_i) {
164 assert( meas_i < parameter_map.size());
165 MatrixIndexType param_i = parameter_map[meas_i];
166 for (MatrixIndexType meas_j=0; meas_j<
N; ++meas_j) {
167 assert( meas_j < parameter_map.size());
168 ret(meas_i,meas_j) =
matrix(param_i, parameter_map[meas_j]);
179 template <
typename measurement_vector_t,
typename measurement_cov_matrix_t,
180 typename predicted_vector_t,
typename predicted_cov_matrix_t>
182 const measurement_cov_matrix_t &a_cov,
183 const predicted_vector_t &
b,
184 const predicted_cov_matrix_t &b_cov) {
188 static_assert( MeasurementSelectorMatrixTraits::matrixColumns<measurement_vector_t>() == 1);
189 static_assert( MeasurementSelectorMatrixTraits::matrixColumns<predicted_vector_t>() == 1);
190 static_assert( MeasurementSelectorMatrixTraits::matrixRows<measurement_cov_matrix_t>()
191 == MeasurementSelectorMatrixTraits::matrixColumns<measurement_cov_matrix_t>() );
192 static_assert( MeasurementSelectorMatrixTraits::matrixRows<predicted_cov_matrix_t>()
193 == MeasurementSelectorMatrixTraits::matrixColumns<predicted_cov_matrix_t>() );
194 static_assert( MeasurementSelectorMatrixTraits::matrixRows<measurement_cov_matrix_t>()
195 == MeasurementSelectorMatrixTraits::matrixRows<measurement_vector_t>() );
196 static_assert( MeasurementSelectorMatrixTraits::matrixRows<predicted_cov_matrix_t>()
197 == MeasurementSelectorMatrixTraits::matrixRows<predicted_vector_t>() );
198 static_assert( MeasurementSelectorMatrixTraits::matrixRows<measurement_cov_matrix_t>()
199 == MeasurementSelectorMatrixTraits::matrixRows<predicted_cov_matrix_t>() );
211 template <std::
size_t N,
class PayloadType >
250 void acceptAndSort(std::function<
bool(
const PayloadType &
a,
const PayloadType &
b)> comparison) {
271 typename std::array<IndexType, N+1>::const_iterator
begin() {
return m_order.begin(); }
297 template <std::size_t NMeasMax,
302 using Config = Acts::GeometryHierarchyMap<AtlasMeasurementSelectorCuts>;
313 const derived_t &
derived()
const {
return *
static_cast<const derived_t *
>(
this); }
317 template <std::
size_t N>
320 constexpr std::size_t nrows = Acts::MultiTrajectoryTraits::MeasurementSizeMax;
321 constexpr std::size_t ncols = Acts::eBoundSize;
323 std::bitset<nrows * ncols> proj_bitset {};
325 for (
unsigned int col_i=0; col_i<
N; ++col_i) {
326 unsigned int row_i = parameter_map[col_i];
327 unsigned int idx = col_i *nrows + row_i;
328 proj_bitset[ (nrows * ncols - 1) -
idx ] = 1;
330 return proj_bitset.to_ullong();
343 const Acts::ProjectorBitset &projector_bitset,
344 boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface> &track_states,
345 const Acts::Logger&
logger,
346 bool outlier_states) {
348 const auto &boundParams = derived_t::boundParams(boundState);
349 const auto &pathLength = derived_t::pathLength(boundState);
351 using PM = Acts::TrackStatePropMask;
353 track_states.reserve( n_new_states );
354 std::optional<TrackStateProxy> firstTrackState{
356 for (
unsigned int state_i=0; state_i<n_new_states; ++state_i) {
357 PM
mask = PM::Predicted | PM::Filtered | PM::Jacobian | PM::Calibrated;
359 if (firstTrackState.has_value()) {
362 mask &= ~PM::Predicted & ~PM::Jacobian;
365 if (outlier_states) {
367 mask &= ~PM::Filtered;
371 ACTS_VERBOSE(
"Create SourceLink output track state #"
372 << trackState.index() <<
" with mask: " <<
mask);
374 if (firstTrackState.has_value()) {
375 trackState.shareFrom(*firstTrackState, PM::Predicted);
376 trackState.shareFrom(*firstTrackState, PM::Jacobian);
380 trackState.predicted() = boundParams.parameters();
381 if (boundParams.covariance()) {
382 trackState.predictedCovariance() = *boundParams.covariance();
384 trackState.jacobian() = derived_t::boundJacobiMatrix(boundState);
385 firstTrackState = trackState;
387 trackState.pathLength() = pathLength;
389 trackState.setReferenceSurface(boundParams.referenceSurface().getSharedPtr());
390 trackState.setProjectorBitset(projector_bitset);
393 if (trackState.referenceSurface().surfaceMaterial() !=
nullptr) {
394 typeFlags.set(Acts::TrackStateFlag::MaterialFlag);
396 typeFlags.set(Acts::TrackStateFlag::ParameterFlag);
400 track_states.push_back( trackState.index());
406 template <std::
size_t DIM,
typename T_SourceLink>
423 template <
typename T>
431 template <
typename Iterator>
446 template <std::
size_t DIM,
typename source_link_iterator_t,
typename T_Container>
447 Acts::Result<boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface> >
449 const Acts::CalibrationContext& calibrationContext,
450 const Acts::Surface& surface,
452 const source_link_iterator_t& sourceLinkBegin,
453 const source_link_iterator_t& sourceLinkEnd,
456 const Acts::Logger&
logger,
457 const std::size_t numMeasurementsCut,
458 const std::pair<float,float>& maxChi2Cut,
459 const T_Container &container)
const {
460 Acts::Result<boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface> >
461 result = boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface>{};
463 using BaseElementType = std::remove_cv_t<std::remove_pointer_t< container_value_t > >;
465 using TheMatchingMeasurement = MatchingMeasurement<DIM, container_value_t >;
466 using MeasCovPair = TheMatchingMeasurement::MeasCovPair;
478 = std::make_pair( project<DIM, Predicted>(parameter_map,
480 project<DIM, PredictedCovariance>(parameter_map,
481 derived().boundParams(boundState).covariance().
value()));
484 auto preCalibrator =
derived().template preCalibrator<DIM, BaseElementType>();
487 for (
const auto &measurement : MeasurementRange<T_Container>(container, sourceLinkBegin, sourceLinkEnd) ) {
488 TheMatchingMeasurement &matching_measurement=selected_measurements.
slot();
489 matching_measurement.m_measurement = preCalibrator(geometryContext,
491 derived().
template forwardToCalibrator(measurement),
492 derived().boundParams(boundState));
493 matching_measurement.m_chi2 =
computeChi2(matching_measurement.m_measurement.first,
494 matching_measurement.m_measurement.second,
498 if (matching_measurement.m_chi2<maxChi2Cut.second) {
499 matching_measurement.m_sourceLink=measurement;
500 selected_measurements.
acceptAndSort([](
const TheMatchingMeasurement &
a,
501 const TheMatchingMeasurement &
b) {
502 return a.m_chi2 <
b.m_chi2;
509 auto postCalibrator =
derived().template postCalibrator<DIM, BaseElementType>();
510 using post_calib_meas_cov_pair_t
516 static constexpr
bool pre_and_post_calib_types_agree
522 using Empty =
struct {};
523 typename std::conditional< !pre_and_post_calib_types_agree,
524 std::array< post_calib_meas_cov_pair_t, NMeasMax>,
526 if (postCalibrator) {
531 idx: selected_measurements) {
532 TheMatchingMeasurement &a_selected_measurement = selected_measurements.
getSlot(
idx);
536 post_calib_meas_cov_pair_t &calibrated_measurement
537 = [&calibrated, &a_selected_measurement, calibrated_meas_cov_i]() -> post_calib_meas_cov_pair_t & {
538 if constexpr(pre_and_post_calib_types_agree) {
540 (void) calibrated_meas_cov_i;
541 return a_selected_measurement.m_measurement;
544 (void) a_selected_measurement;
545 assert(calibrated_meas_cov_i < calibrated.size());
546 return calibrated[calibrated_meas_cov_i];
551 calibrated_measurement = postCalibrator(geometryContext,
553 derived().
template forwardToCalibrator(a_selected_measurement.m_sourceLink.value()),
554 derived().boundParams(boundState));
556 a_selected_measurement.m_chi2 =
computeChi2(calibrated_measurement.first,
557 calibrated_measurement.second,
561 a_selected_measurement.m_isOutLier = (a_selected_measurement.m_chi2 >= maxChi2Cut.first);
562 if constexpr(!pre_and_post_calib_types_agree) {
563 ++calibrated_meas_cov_i;
570 idx: selected_measurements) {
571 TheMatchingMeasurement &a_selected_measurement = selected_measurements.
getSlot(
idx);
572 a_selected_measurement.m_isOutLier = (a_selected_measurement.m_chi2 >= maxChi2Cut.first);
588 (!selected_measurements.
empty()
589 ? selected_measurements.
getSlot( *(selected_measurements.
begin())).m_isOutLier
591 assert(
result->size() == selected_measurements.
size() );
594 auto use_calibrated_storage = [&postCalibrator]() ->
bool {
595 if constexpr(pre_and_post_calib_types_agree) {
596 (void) postCalibrator;
602 return postCalibrator;
607 unsigned int state_i=0;
609 idx: selected_measurements) {
610 assert( state_i < result->
size());
612 TheMatchingMeasurement &a_selected_measurement = selected_measurements.
getSlot(
idx);
613 trackState.setUncalibratedSourceLink(
derived().makeSourceLink(std::move(a_selected_measurement.m_sourceLink.value())));
615 trackState.typeFlags().set( a_selected_measurement.m_isOutLier
616 ? Acts::TrackStateFlag::OutlierFlag
617 : Acts::TrackStateFlag::MeasurementFlag );
618 trackState.allocateCalibrated(DIM);
619 if (use_calibrated_storage()) {
622 assert( use_calibrated_storage() == !pre_and_post_calib_types_agree);
623 if constexpr(!pre_and_post_calib_types_agree) {
624 assert( state_i < calibrated.size());
625 trackState.template calibrated<DIM>()
626 = MeasurementSelectorMatrixTraits::matrixTypeCast<typename MeasurementSelectorTraits<derived_t>::MatrixFloatType>(calibrated[state_i].first);
627 trackState.template calibratedCovariance<DIM>()
628 = MeasurementSelectorMatrixTraits::matrixTypeCast<typename MeasurementSelectorTraits<derived_t>::MatrixFloatType>(calibrated[state_i].second);
629 trackState.chi2() = a_selected_measurement.m_chi2;
633 trackState.template calibrated<DIM>()
634 = MeasurementSelectorMatrixTraits::matrixTypeCast<typename MeasurementSelectorTraits<derived_t>::MatrixFloatType>(a_selected_measurement.m_measurement.first);
635 trackState.template calibratedCovariance<DIM>()
636 = MeasurementSelectorMatrixTraits::matrixTypeCast<typename MeasurementSelectorTraits<derived_t>::MatrixFloatType>(a_selected_measurement.m_measurement.second);
637 trackState.chi2() = a_selected_measurement.m_chi2;
644 template <
typename parameters_t>
645 static std::size_t
getEtaBin(
const parameters_t& boundParameters,
646 const std::vector<float> &
etaBins) {
650 const float eta = std::abs(std::atanh(
std::cos(boundParameters.parameters()[Acts::eBoundTheta])));
664 std::tuple< std::size_t, std::pair<float,float> >
getCuts(
const Acts::Surface& surface,
666 const Acts::Logger&
logger)
const {
667 std::tuple< std::size_t, std::pair<float,float> >
result;
668 std::size_t &numMeasurementsCut = std::get<0>(
result);
669 std::pair<float,float> &maxChi2Cut = std::get<1>(
result);
671 auto geoID = surface.geometryId();
676 numMeasurementsCut = 0;
682 numMeasurementsCut = (!
cuts->numMeasurementsCutOff.empty()
685 maxChi2Cut = !
cuts->chi2CutOff.empty()
689 ACTS_VERBOSE(
"Get cut for eta-bin="
691 <<
": chi2 (max,max-outlier) " << maxChi2Cut.first <<
", " << maxChi2Cut.second
692 <<
" max.meas." << numMeasurementsCut);
708 template <std::size_t NMeasMax,
710 typename measurement_container_variant_t >
724 template <
typename source_link_iterator_t>
725 Acts::Result<boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface> >
727 const Acts::CalibrationContext& calibrationContext,
728 const Acts::Surface& surface,
730 source_link_iterator_t sourceLinkBegin,
731 source_link_iterator_t sourceLinkEnd,
734 [[maybe_unused]] std::vector<TrackStateProxy> &trackStateCandidates,
736 const Acts::Logger&
logger)
const {
737 Acts::Result<boost::container::small_vector< typename TrackStateProxy::IndexType, s_maxBranchesPerSurface> >
738 result = Acts::CombinatorialKalmanFilterError::MeasurementSelectionFailed;
739 if (sourceLinkBegin != sourceLinkEnd) {
740 auto [numMeasurementsCut, maxChi2Cut] = this->
getCuts(surface,boundState,
logger);
743 if (numMeasurementsCut>0) {
744 const std::vector<measurement_container_variant_t> &
745 measurementContainer = sourceLinkBegin.m_iterator.measurementContainerList();
747 assert( sourceLinkBegin.m_iterator.containerIndex() == sourceLinkEnd.m_iterator.containerIndex() );
748 const measurement_container_variant_t &a_measurement_container_variant = measurementContainer.at(sourceLinkBegin.m_iterator.containerIndex());
749 result = std::visit( [
this,
760 &maxChi2Cut] (
const auto &measurement_container_with_dimension) {
761 using ArgType = std::remove_cv_t<std::remove_reference_t< decltype(measurement_container_with_dimension) > >;
762 constexpr std::size_t DIM = ArgType::dimension();
763 return this->
template selectMeasurementsCreateTrackStates<DIM>(geometryContext,
774 measurement_container_with_dimension.container());
776 a_measurement_container_variant);
783 template <std::size_t NMeasMax,
785 typename measurement_container_variant_t>
791 return std::get<0>(boundState);
796 return std::get<1>(boundState);
800 return std::get<2>(boundState);
804 template <
typename T_Value>
806 return Acts::SourceLink{
value};
812 template <
typename T>
814 if constexpr( std::is_same<
T, std::remove_pointer_t<T> >::
value ) {
826 template <std::
size_t DIM>
829 const Acts::CalibrationContext&,
830 const Acts::Surface&) {
831 return ParameterMapping::identity<DIM>();
835 template <std::
size_t DIM,
typename measurement_t>
839 (
const Acts::GeometryContext&,
840 const Acts::CalibrationContext&,
841 const measurement_t &,
846 template <std::
size_t DIM,
typename measurement_t>
850 (
const Acts::GeometryContext&,
851 const Acts::CalibrationContext&,
852 const measurement_t &,
858 template <std::
size_t DIM,
typename measurement_t>
866 template <std::
size_t DIM,
typename measurement_t>