diff --git a/Tracking/Acts/FaserActsKalmanFilter/CMakeLists.txt b/Tracking/Acts/FaserActsKalmanFilter/CMakeLists.txt index 7bb3abe1e6a5f703af8a5aebd053354462d6db6f..34d9584f9f6b30b1683c962d6baa59037bc8d6f3 100755 --- a/Tracking/Acts/FaserActsKalmanFilter/CMakeLists.txt +++ b/Tracking/Acts/FaserActsKalmanFilter/CMakeLists.txt @@ -36,6 +36,7 @@ atlas_add_component(FaserActsKalmanFilter # FaserActsKalmanFilter/TruthTrackFinderTool.h FaserActsKalmanFilter/Measurement.h # FaserActsKalmanFilter/MultiTrackFinderTool.h + FaserActsKalmanFilter/MyAmbiguitySolver.h FaserActsKalmanFilter/PerformanceWriterTool.h FaserActsKalmanFilter/PlotHelpers.h FaserActsKalmanFilter/ResPlotTool.h diff --git a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CombinatorialKalmanFilterAlg.h b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CombinatorialKalmanFilterAlg.h index f98dc333b971dd92d2dd91d7fab3ee90e8e16199..9a8a08e2ba7018a3fcb066afbeb58d3644958599 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CombinatorialKalmanFilterAlg.h +++ b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/CombinatorialKalmanFilterAlg.h @@ -7,8 +7,8 @@ #include "TrackerPrepRawData/FaserSCT_ClusterContainer.h" #include "FaserActsGeometryInterfaces/IFaserActsTrackingGeometryTool.h" #include "FaserActsKalmanFilter/TruthBasedInitialParameterTool.h" -#include "FaserActsKalmanFilter/SPSimpleInitialParameterTool.h" -#include "FaserActsKalmanFilter/SPSeedBasedInitialParameterTool.h" +//#include "FaserActsKalmanFilter/SPSimpleInitialParameterTool.h" +//#include "FaserActsKalmanFilter/SPSeedBasedInitialParameterTool.h" #include "Acts/TrackFinding/CombinatorialKalmanFilter.hpp" #include "Acts/TrackFinding/MeasurementSelector.hpp" #include "FaserActsKalmanFilter/Measurement.h" @@ -70,6 +70,7 @@ public: const TrackerDD::SCT_DetectorManager* m_detManager {nullptr}; Gaudi::Property<std::string> m_actsLogging {this, "ActsLogging", "VERBOSE"}; + Gaudi::Property<int> m_minNumberMeasurements {this, "MinNumberMeasurements", 12}; Gaudi::Property<bool> m_backwardPropagation {this, "BackwardPropagation", false}; Gaudi::Property<bool> m_performanceWriter {this, "PerformanceWriter", true}; Gaudi::Property<bool> m_summaryWriter {this, "SummaryWriter", true}; diff --git a/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/MyAmbiguitySolver.h b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/MyAmbiguitySolver.h new file mode 100644 index 0000000000000000000000000000000000000000..47257befa2fcbd70da95d1ee88bd41d0fdabeb2a --- /dev/null +++ b/Tracking/Acts/FaserActsKalmanFilter/FaserActsKalmanFilter/MyAmbiguitySolver.h @@ -0,0 +1,80 @@ +#ifndef FASERACTSKALMANFILTER_AMBIGUITYSOLVER_H +#define FASERACTSKALMANFILTER_AMBIGUITYSOLVER_H + +#include "Acts/TrackFinding/CombinatorialKalmanFilter.hpp" +#include "FaserActsKalmanFilter/FaserActsRecMultiTrajectory.h" + +using CombinatorialKalmanFilterResult = Acts::CombinatorialKalmanFilterResult<IndexSourceLink>; +using TrackFitterResult = Acts::Result<CombinatorialKalmanFilterResult>; +using TrackFinderResult = std::vector<TrackFitterResult>; + + +size_t numberMeasurements(const CombinatorialKalmanFilterResult& ckfResult) { + auto traj = FaserActsRecMultiTrajectory(ckfResult.fittedStates, ckfResult.lastMeasurementIndices, ckfResult.fittedParameters); + const auto& mj = traj.multiTrajectory(); + const auto& trackTips = traj.tips(); + size_t maxMeasurements = 0; + for (const auto& trackTip : trackTips) { + auto trajState = Acts::MultiTrajectoryHelpers::trajectoryState(mj, trackTip); + size_t nMeasurements = trajState.nMeasurements; + if (nMeasurements > maxMeasurements) { + maxMeasurements = nMeasurements; + } + std::cout << "# measurements: " << trajState.nMeasurements << std::endl; + } + return maxMeasurements; +} + +int countSharedHits(const CombinatorialKalmanFilterResult& result1, const CombinatorialKalmanFilterResult& result2) { + int count = 0; + std::vector<size_t> hitIndices {}; + + for (auto measIndex : result1.lastMeasurementIndices) { + result1.fittedStates.visitBackwards(measIndex, [&](const auto& state) { + if (not state.typeFlags().test(Acts::TrackStateFlag::MeasurementFlag)) + return; + size_t hitIndex = state.uncalibrated().index(); + hitIndices.emplace_back(hitIndex); + }); + } + + for (auto measIndex : result2.lastMeasurementIndices) { + result2.fittedStates.visitBackwards(measIndex, [&](const auto& state) { + if (not state.typeFlags().test(Acts::TrackStateFlag::MeasurementFlag)) + return; + size_t hitIndex = state.uncalibrated().index(); + if (std::find(hitIndices.begin(), hitIndices.end(), hitIndex) != hitIndices.end()) { + count += 1; + } + }); + } + return count; +} + + +std::pair<int, int> solveAmbiguity(TrackFinderResult& results, size_t minMeasurements = 13) { + std::map<std::pair<size_t, size_t>, size_t> trackPairs {}; + for (size_t i = 0; i < results.size(); ++i) { + // if (not results.at(i).ok()) continue; + // if (numberMeasurements(results.at(i).value()) < minMeasurements) continue; + for (size_t j = i+1; j < results.size(); ++j) { + // if (not results.at(j).ok()) continue; + // if (numberMeasurements(results.at(j).value()) < minMeasurements) continue; + int n = countSharedHits(results.at(i).value(), results.at(j).value()); + trackPairs[std::make_pair(i, j)] = n; + } + } + + std::pair<size_t, size_t> bestTrackPair; + size_t minSharedHits = std::numeric_limits<size_t>::max(); + for (const auto& trackPair : trackPairs) { + if (trackPair.second < minSharedHits) { + minSharedHits = trackPair.second; + bestTrackPair = trackPair.first; + } + } + + return bestTrackPair; +} + +#endif //FASERACTSKALMANFILTER_AMBIGUITYSOLVER_H diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/CombinatorialKalmanFilterAlg.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/CombinatorialKalmanFilterAlg.cxx index 1c7047960d5f98a78a9bd5dbe91eb15ed49e74dd..cc70de2becc3b4832c9abdbecf5813aa2f7c56e9 100644 --- a/Tracking/Acts/FaserActsKalmanFilter/src/CombinatorialKalmanFilterAlg.cxx +++ b/Tracking/Acts/FaserActsKalmanFilter/src/CombinatorialKalmanFilterAlg.cxx @@ -19,6 +19,7 @@ #include "Acts/Surfaces/PerigeeSurface.hpp" #include "Acts/MagneticField/MagneticFieldContext.hpp" #include "FaserActsKalmanFilter/TrackSelection.h" +#include "FaserActsKalmanFilter/MyAmbiguitySolver.h" #include <algorithm> using TrajectoriesContainer = std::vector<FaserActsRecMultiTrajectory>; @@ -123,6 +124,43 @@ StatusCode CombinatorialKalmanFilterAlg::execute() { for (size_t i = 0; i < std::min(trackQuality.size(), (size_t)m_nTrajectories); ++i) { selectedResults.push_back(Acts::Result<Acts::CombinatorialKalmanFilterResult<IndexSourceLink>>(trackQuality[i].track)); } + /* + TrackFinderResult minTrackRequirements {}; + for (auto& result : results) { + if (not result.ok()) { + continue; + } + auto& ckfResult = result.value(); + auto traj = FaserActsRecMultiTrajectory(ckfResult.fittedStates, ckfResult.lastMeasurementIndices, ckfResult.fittedParameters); + const auto& mj = traj.multiTrajectory(); + const auto& trackTips = traj.tips(); + size_t maxMeasurements = 0; + for (const auto& trackTip : trackTips) { + auto trajState = Acts::MultiTrajectoryHelpers::trajectoryState(mj, trackTip); + size_t nMeasurements = trajState.nMeasurements; + if (nMeasurements > maxMeasurements) { + maxMeasurements = nMeasurements; + } + } + if (maxMeasurements >= m_minNumberMeasurements) { + minTrackRequirements.push_back(Acts::Result<Acts::CombinatorialKalmanFilterResult<IndexSourceLink>>(ckfResult)); + } + } + + TrackFinderResult selectedResults {}; + if (minTrackRequirements.size() > 2) { + std::pair<std::size_t, std::size_t> trackIndices = solveAmbiguity(results); + selectedResults.push_back( + Acts::Result<Acts::CombinatorialKalmanFilterResult<IndexSourceLink>>(results.at(trackIndices.first).value())); + selectedResults.push_back( + Acts::Result<Acts::CombinatorialKalmanFilterResult<IndexSourceLink>>(results.at(trackIndices.second).value())); + } else { + for (auto& result : minTrackRequirements) { + selectedResults.push_back( + Acts::Result<Acts::CombinatorialKalmanFilterResult<IndexSourceLink>>(result.value())); + } + } + */ computeSharedHits(sourceLinks.get(), selectedResults); // Loop over the track finding results for all initial parameters diff --git a/Tracking/Acts/FaserActsKalmanFilter/src/MyAmbiguitySolver.cxx b/Tracking/Acts/FaserActsKalmanFilter/src/MyAmbiguitySolver.cxx new file mode 100644 index 0000000000000000000000000000000000000000..998269a83cf6228092dd3b31987696f929820e4e --- /dev/null +++ b/Tracking/Acts/FaserActsKalmanFilter/src/MyAmbiguitySolver.cxx @@ -0,0 +1,3 @@ +#include "FaserActsKalmanFilter/MyAmbiguitySolver.h" + +