diff --git a/algo/ca/core/CMakeLists.txt b/algo/ca/core/CMakeLists.txt index eb0d7c8d51895945b0aac9f9793a6b07a7e35744..566fbed38071a4db66aa23f6268d27ff8193e949 100644 --- a/algo/ca/core/CMakeLists.txt +++ b/algo/ca/core/CMakeLists.txt @@ -18,6 +18,7 @@ set(SRCS ${CMAKE_CURRENT_SOURCE_DIR}/data/CaMeasurementU.cxx ${CMAKE_CURRENT_SOURCE_DIR}/data/CaMeasurementXy.cxx ${CMAKE_CURRENT_SOURCE_DIR}/data/CaMeasurementTime.cxx + ${CMAKE_CURRENT_SOURCE_DIR}/data/CaWindowData.cxx ${CMAKE_CURRENT_SOURCE_DIR}/pars/CaConfigReader.cxx ${CMAKE_CURRENT_SOURCE_DIR}/pars/CaField.cxx ${CMAKE_CURRENT_SOURCE_DIR}/pars/CaInitManager.cxx @@ -84,12 +85,14 @@ install( data/CaMeasurementU.h data/CaMeasurementXy.h data/CaMeasurementTime.h + data/CaWindowData.h pars/CaConfigReader.h data/CaGridEntry.h data/CaGrid.h data/CaGridArea.h data/CaTriplet.h data/CaBranch.h + data/CaWindowData.h pars/CaConstants.h pars/CaField.h pars/CaInitManager.h diff --git a/algo/ca/core/data/CaWindowData.cxx b/algo/ca/core/data/CaWindowData.cxx new file mode 100644 index 0000000000000000000000000000000000000000..cff42c32f476b896a27f536d83daeaab9d35b37b --- /dev/null +++ b/algo/ca/core/data/CaWindowData.cxx @@ -0,0 +1,20 @@ +/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Sergei Zharko [committer] */ + +/// \file CaWindowData.h +/// \author Sergei Zharko <s.zharko@gsi.de> +/// \brief Container for all data, which are processed within a single sub-timeslice (implementation) +/// \since 29.01.2024 + +#include "CaWindowData.h" + +using cbm::algo::ca::WindowData; + +// --------------------------------------------------------------------------------------------------------------------- +// +void WindowData::ResetHitData(int nHits) +{ + fvHits.reset(nHits); + fvbHitSuppressed.reset(nHits, 0); +} diff --git a/algo/ca/core/data/CaWindowData.h b/algo/ca/core/data/CaWindowData.h new file mode 100644 index 0000000000000000000000000000000000000000..a0d76d345d5320d79bfa38f73837dea4f1d046a0 --- /dev/null +++ b/algo/ca/core/data/CaWindowData.h @@ -0,0 +1,150 @@ +/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Sergei Zharko [committer] */ + +/// \file CaWindowData.h +/// \author Sergei Zharko <s.zharko@gsi.de> +/// \brief Container for all data, which are processed within a single sub-timeslice +/// \since 29.01.2024 + +#pragma once + +#include "CaConstants.h" +#include "CaGrid.h" +#include "CaHit.h" +#include "CaTrack.h" + +#include <array> + +namespace cbm::algo::ca +{ + /// \class WindowData + /// \brief Container for internal data, processed on a single sub-TS window + class WindowData { + public: + /// \brief Constructor + WindowData() = default; + + /// \brief Destructor + ~WindowData() = default; + + /// \brief Gets grid for station index + /// \param iStation Index of the tracking station + [[gnu::always_inline]] ca::Grid& Grid(int iStation) { return fvGrid[iStation]; } + + /// \brief Gets grid for station index + /// \param iStation Index of the tracking station + [[gnu::always_inline]] const ca::Grid& Grid(int iStation) const { return fvGrid[iStation]; } + + /// \brief Gets hit by index + /// \param iHit Hit index + [[gnu::always_inline]] ca::Hit& Hit(int iHit) { return fvHits[iHit]; } + + /// \brief Gets hit by index + /// \param iHit Hit index + [[gnu::always_inline]] const ca::Hit& Hit(int iHit) const { return fvHits[iHit]; } + + /// \brief Gets hit vector + [[gnu::always_inline]] Vector<ca::Hit>& Hits() { return fvHits; } + + /// \brief Gets hit vector + [[gnu::always_inline]] const Vector<ca::Hit>& Hits() const { return fvHits; } + + /// \brief Gets hit suppression flag + /// \param iHit Hit index + [[gnu::always_inline]] uint8_t IsHitSuppressed(int iHit) const { return fvbHitSuppressed[iHit]; } + + /// \brief Resets hit data + /// \param nHits Number of hits in the sample + void ResetHitData(int nHits); + + /// \brief Reset suppressed hit flags + [[gnu::always_inline]] void ResetHitSuppressionFlags() { fvbHitSuppressed.reset(fvHits.size()); } + + /// \brief Set hit suppression flag + [[gnu::always_inline]] void SuppressHit(int iHit) { fvbHitSuppressed[iHit] = 1; } + + /// \brief Index of the first hit on the station + /// \param iStation Index of the tracking station + [[gnu::always_inline]] ca::HitIndex_t& HitStartIndexOnStation(int iStation) + { + return fvHitStartIndexOnStation[iStation]; + } + + /// \brief Index of the first hit on the station + /// \param iStation Index of the tracking station + [[gnu::always_inline]] ca::HitIndex_t HitStartIndexOnStation(int iStation) const + { + return fvHitStartIndexOnStation[iStation]; + } + + /// \brief Number of hits on station + /// \param iStation Index of the tracking station + [[gnu::always_inline]] ca::HitIndex_t& NofHitsOnStation(int iStation) { return fvNofHitsOnStation[iStation]; } + + /// \brief Number of hits on station + /// \param iStation Index of the tracking station + [[gnu::always_inline]] ca::HitIndex_t NofHitsOnStation(int iStation) const { return fvNofHitsOnStation[iStation]; } + + /// \brief Accesses indices of hits, used by reconstructed tracks + [[gnu::always_inline]] Vector<ca::HitIndex_t>& RecoHitIndices() { return fvRecoHitIndices; } + + /// \brief Accesses indices of hits + [[gnu::always_inline]] const Vector<ca::HitIndex_t>& RecoHitIndices() const { return fvRecoHitIndices; } + + /// \brief Accesses index of hit in the input data + /// \param iHit Index of reconstructed hit, used in reconstructed tracks + /// \return Index of reconstructed hit in the algorithm input data object + [[gnu::always_inline]] ca::HitIndex_t& RecoHitIndex(int iHit) { return fvRecoHitIndices[iHit]; } + + /// \brief Accesses index of hit in the input data + /// \param iHit Index of reconstructed hit, used in reconstructed tracks + /// \return Index of reconstructed hit in the algorithm input data object + [[gnu::always_inline]] ca::HitIndex_t RecoHitIndex(int iHit) const { return fvRecoHitIndices[iHit]; } + + /// \brief Accesses reconstructed track by index + /// \param iTrack Track index + [[gnu::always_inline]] ca::Track& RecoTrack(int iTrack) { return fvRecoTracks[iTrack]; } + + /// \brief Accesses reconstructed track by index + /// \param iTrack Track index + [[gnu::always_inline]] const ca::Track& RecoTrack(int iTrack) const { return fvRecoTracks[iTrack]; } + + /// \brief Accesses reconstructed track container + [[gnu::always_inline]] Vector<ca::Track>& RecoTracks() { return fvRecoTracks; } + + /// \brief Accesses reconstructed track container + [[gnu::always_inline]] const Vector<ca::Track>& RecoTracks() const { return fvRecoTracks; } + + private: + static constexpr int kMaxNofStations = constants::size::MaxNstations; ///< Alias to max number of stations + + /// \brief Grid vs. station index + std::array<ca::Grid, kMaxNofStations> fvGrid; + + /// \brief Hits of the current time window + /// + /// It is a portion of fInputData to process in the current time window + /// hit.Id is replaced by the hit index in fInputData + Vector<ca::Hit> fvHits{"WindowData::fHits"}; + + /// \brief Flag, if the hit is suppressed for tracking + Vector<unsigned char> fvbHitSuppressed{"WindowData::fbHitSuppressed"}; + + /// \brief First hit index of the station + std::array<ca::HitIndex_t, kMaxNofStations + 1> fvHitStartIndexOnStation = {0}; + + /// \brief Number of hits on the station + std::array<ca::HitIndex_t, kMaxNofStations + 1> fvNofHitsOnStation = {0}; + + /// \brief Sample of reconstructed tracks + Vector<ca::Track> fvRecoTracks{"WindowData::fvRecoTracks"}; + + /// \brief Sample of reconstructed hit indices + Vector<ca::HitIndex_t> fvRecoHitIndices{"WindowData::fvRecoHitIndices"}; + + /// \brief Sample of track candidates + Vector<ca::Track> fvTrackCandidates{"WindowData::fvTrackCandidates"}; + + } _fvecalignment; +} // namespace cbm::algo::ca diff --git a/algo/ca/core/tracking/CaFramework.h b/algo/ca/core/tracking/CaFramework.h index 66ea7f788e4497849a793fa60a3b0a1684027aea..c7f1a3322dcb56c6fad8bc90da69ea8ea01df1e8 100644 --- a/algo/ca/core/tracking/CaFramework.h +++ b/algo/ca/core/tracking/CaFramework.h @@ -25,6 +25,7 @@ #include "CaTrackingMonitor.h" #include "CaTriplet.h" #include "CaVector.h" +#include "CaWindowData.h" #include <array> #include <iomanip> @@ -221,23 +222,11 @@ namespace cbm::algo::ca public: Vector<CaHitTimeInfo> fHitTimeInfo; - ca::Grid vGrid[constants::size::MaxNstations]; ///< + ca::WindowData fWData{}; /// ----- Data used during track finding ----- /// - ///\brief hits of the current time window - /// it is a portion of fInputData to process in the current time window - /// hit.Id is replaced by the hit index in fInputData - Vector<ca::Hit> fWindowHits{"Framework::fWindowHits"}; - - Vector<unsigned char> fIsWindowHitSuppressed{"Framework::fIsWindowHitSuppressed"}; - - ///\brief first index of the station hits in the time window - ca::HitIndex_t fStationHitsStartIndex[constants::size::MaxNstations + 1]{0}; - - ///\brief number of station hits in the time window - ca::HitIndex_t fStationNhits[constants::size::MaxNstations + 1]{0}; /// ---------- @@ -246,9 +235,6 @@ namespace cbm::algo::ca Vector<Track> fRecoTracks{"Framework::fRecoTracks"}; ///< reconstructed tracks Vector<ca::HitIndex_t> fRecoHits{"Framework::fRecoHits"}; ///< packed hits of reconstructed tracks - Vector<Track> fSliceRecoTracks{"Framework::fSliceRecoTracks"}; ///< reconstructed tracks in sub-timeslice - Vector<ca::HitIndex_t> fSliceRecoHits{"Framework::fSliceRecoHits"}; ///< packed hits of reconstructed tracks - /// Created triplets vs station index Vector<ca::Triplet> fTriplets[constants::size::MaxNstations]{{"Framework::fTriplets"}}; diff --git a/algo/ca/core/tracking/CaTrackExtender.cxx b/algo/ca/core/tracking/CaTrackExtender.cxx index 901c43a6d8216bec5429336917f5e6f81d86c8dc..b44b8b28e6e824d9cb1d19a879cf2d25f8d33a0f 100644 --- a/algo/ca/core/tracking/CaTrackExtender.cxx +++ b/algo/ca/core/tracking/CaTrackExtender.cxx @@ -207,10 +207,10 @@ void TrackExtender::FindMoreHits(ca::Branch& t, TrackParamV& Tout, const bool up TrackParamV& tr = fit.Tr(); - ca::GridArea area( - frAlgo.vGrid[ista], tr.X()[0], tr.Y()[0], - (sqrt(frAlgo.fPickGather * tr.C00()) + frAlgo.vGrid[ista].GetMaxRangeX() + frAlgo.fMaxDZ * abs(tr.Tx()))[0], - (sqrt(frAlgo.fPickGather * tr.C11()) + frAlgo.vGrid[ista].GetMaxRangeY() + frAlgo.fMaxDZ * abs(tr.Ty()))[0]); + const auto& grid = frAlgo.fWData.Grid(ista); + ca::GridArea area(grid, tr.X()[0], tr.Y()[0], + (sqrt(frAlgo.fPickGather * tr.C00()) + grid.GetMaxRangeX() + frAlgo.fMaxDZ * abs(tr.Tx()))[0], + (sqrt(frAlgo.fPickGather * tr.C11()) + grid.GetMaxRangeY() + frAlgo.fMaxDZ * abs(tr.Ty()))[0]); if (frAlgo.GetParameters().DevIsIgnoreHitSearchAreas()) { area.DoLoopOverEntireGrid(); @@ -220,10 +220,10 @@ void TrackExtender::FindMoreHits(ca::Branch& t, TrackParamV& Tout, const bool up while (area.GetNextObjectId(ih)) { // loop over the hits in the area - if (frAlgo.fIsWindowHitSuppressed[ih]) { + if (frAlgo.fWData.IsHitSuppressed(ih)) { continue; } - const ca::Hit& hit = frAlgo.fWindowHits[ih]; + const ca::Hit& hit = frAlgo.fWData.Hit(ih); if (sta.timeInfo && tr.NdfTime()[0] > -2.) { fscal dt = hit.T() - tr.Time()[0]; @@ -248,11 +248,10 @@ void TrackExtender::FindMoreHits(ca::Branch& t, TrackParamV& Tout, const bool up fscal d_y = hit.Y() - y[0]; fscal d2 = d_x * d_x + d_y * d_y; if (d2 > r2_best) continue; - - fscal dxm_est = sqrt(pickGather2 * C00)[0] + frAlgo.vGrid[ista].GetMaxRangeX(); + fscal dxm_est = sqrt(pickGather2 * C00)[0] + grid.GetMaxRangeX(); if (fabs(d_x) > dxm_est) continue; - fscal dym_est = sqrt(pickGather2 * C11)[0] + frAlgo.vGrid[ista].GetMaxRangeY(); + fscal dym_est = sqrt(pickGather2 * C11)[0] + grid.GetMaxRangeY(); if (fabs(d_y) > dym_est) continue; r2_best = d2; @@ -261,9 +260,9 @@ void TrackExtender::FindMoreHits(ca::Branch& t, TrackParamV& Tout, const bool up if (iHit_best < 0) break; - newHits.push_back(frAlgo.fWindowHits[iHit_best].Id()); - const ca::Hit& hit = frAlgo.fWindowHits[iHit_best]; + const ca::Hit& hit = frAlgo.fWData.Hit(iHit_best); + newHits.push_back(hit.Id()); fit.Extrapolate(hit.Z(), fld); fit.FilterHit(sta, hit); diff --git a/algo/ca/core/tracking/CaTrackFinder.cxx b/algo/ca/core/tracking/CaTrackFinder.cxx index 7ee4c05c8ca1061dfbb4b6f925a5d293c5b7a863..98c47fe74ebfa361a3de8de77f31b636fe3b1350 100644 --- a/algo/ca/core/tracking/CaTrackFinder.cxx +++ b/algo/ca/core/tracking/CaTrackFinder.cxx @@ -116,6 +116,7 @@ void TrackFinder::FindTracks() 1.5 * l * sqrt(1. + constants::phys::ProtonMassD * constants::phys::ProtonMassD / minProtonMomentum / minProtonMomentum) * constants::phys::SpeedOfLightInvD; + // TODO: Is it possible, that the proton mass selection affects the search of heavier particles? fscal dt = h.RangeT(); CaHitTimeInfo& info = frAlgo.fHitTimeInfo[caHitId]; @@ -130,7 +131,7 @@ void TrackFinder::FindTracks() if (info.fEventTimeMin > 500.e6 || info.fEventTimeMax < -500.) { // cut hits with bogus start time > 500 ms frAlgo.fvHitKeyFlags[h.FrontKey()] = 1; frAlgo.fvHitKeyFlags[h.BackKey()] = 1; - LOG(error) << "CA: skip bogus hit " << h.ToString(); + LOG(warning) << "CA: skip bogus hit " << h.ToString(); continue; } @@ -305,11 +306,10 @@ void TrackFinder::FindTracks() // TODO: only add those hits from the region before tsStartNew that belong to the not stored tracks int trackFirstHit = 0; - for (auto it = frAlgo.fSliceRecoTracks.begin(); it != frAlgo.fSliceRecoTracks.end(); - trackFirstHit += it->fNofHits, it++) { + for (const auto& track : frAlgo.fWData.RecoTracks()) { bool isTrackCompletelyInOverlap = true; - for (int ih = 0; ih < it->fNofHits; ih++) { - int caHitId = frAlgo.fSliceRecoHits[trackFirstHit + ih]; + for (int ih = 0; ih < track.fNofHits; ih++) { + int caHitId = frAlgo.fWData.RecoHitIndex(trackFirstHit + ih); CaHitTimeInfo& info = frAlgo.fHitTimeInfo[caHitId]; if (info.fEventTimeMax < tsStart) { // this hit is before the overlap isTrackCompletelyInOverlap = false; @@ -327,24 +327,25 @@ void TrackFinder::FindTracks() // // release the track hits - for (int i = 0; i < it->fNofHits; i++) { - int caHitId = frAlgo.fSliceRecoHits[trackFirstHit + i]; + for (int i = 0; i < track.fNofHits; i++) { + int caHitId = frAlgo.fWData.RecoHitIndex(trackFirstHit + i); const auto& h = frAlgo.fInputData.GetHit(caHitId); frAlgo.fvHitKeyFlags[h.FrontKey()] = 0; frAlgo.fvHitKeyFlags[h.BackKey()] = 0; } } else { // save the track - frAlgo.fRecoTracks.push_back(*it); + frAlgo.fRecoTracks.push_back(track); // mark the track hits as used - for (int i = 0; i < it->fNofHits; i++) { - int caHitId = frAlgo.fSliceRecoHits[trackFirstHit + i]; + for (int i = 0; i < track.fNofHits; i++) { + int caHitId = frAlgo.fWData.RecoHitIndex(trackFirstHit + i); const auto& h = frAlgo.fInputData.GetHit(caHitId); frAlgo.fvHitKeyFlags[h.FrontKey()] = 1; frAlgo.fvHitKeyFlags[h.BackKey()] = 1; frAlgo.fRecoHits.push_back(caHitId); } } + trackFirstHit += track.fNofHits; } // sub-timeslice tracks tsStart -= 5; // do 5 ns margin diff --git a/algo/ca/core/tracking/CaTrackFinderWindow.cxx b/algo/ca/core/tracking/CaTrackFinderWindow.cxx index 97aff8d80804b70dabef090ff6fe3b1b73a8d6f0..a6eba0e1db35d697e98aa64c393145a76098a050 100644 --- a/algo/ca/core/tracking/CaTrackFinderWindow.cxx +++ b/algo/ca/core/tracking/CaTrackFinderWindow.cxx @@ -120,27 +120,26 @@ void TrackFinderWindow::ReadWindowData() { int nHits = 0; for (int iS = 0; iS < frAlgo.GetParameters().GetNstationsActive(); iS++) { - frAlgo.fStationHitsStartIndex[iS] = nHits; - frAlgo.fStationNhits[iS] = frAlgo.fSliceHitIds[iS].size(); - nHits += frAlgo.fStationNhits[iS]; + frAlgo.fWData.HitStartIndexOnStation(iS) = nHits; + frAlgo.fWData.NofHitsOnStation(iS) = frAlgo.fSliceHitIds[iS].size(); + nHits += frAlgo.fWData.NofHitsOnStation(iS); } - frAlgo.fStationHitsStartIndex[frAlgo.GetParameters().GetNstationsActive()] = nHits; - - frAlgo.fWindowHits.reset(nHits); - frAlgo.fIsWindowHitSuppressed.reset(nHits, 0); + frAlgo.fWData.HitStartIndexOnStation(frAlgo.GetParameters().GetNstationsActive()) = nHits; + frAlgo.fWData.ResetHitData(nHits); for (int iS = 0; iS < frAlgo.GetParameters().GetNstationsActive(); iS++) { + int iFstHit = frAlgo.fWData.HitStartIndexOnStation(iS); for (ca::HitIndex_t ih = 0; ih < frAlgo.fSliceHitIds[iS].size(); ++ih) { ca::Hit h = frAlgo.fInputData.GetHit(frAlgo.fSliceHitIds[iS][ih]); h.SetId(frAlgo.fSliceHitIds[iS][ih]); - frAlgo.fWindowHits[frAlgo.fStationHitsStartIndex[iS] + ih] = h; + frAlgo.fWData.Hit(iFstHit + ih) = h; } } if constexpr (fDebug) { LOG(info) << "===== Sliding Window hits: "; for (int i = 0; i < nHits; ++i) { - LOG(info) << " " << frAlgo.fWindowHits[i].ToString(); + LOG(info) << " " << frAlgo.fWData.Hit(i).ToString(); } LOG(info) << "===== "; } @@ -149,11 +148,11 @@ void TrackFinderWindow::ReadWindowData() frAlgo.fHitNtriplets.reset(nHits, 0); - frAlgo.fSliceRecoTracks.clear(); - frAlgo.fSliceRecoTracks.reserve(2 * nHits / frAlgo.GetParameters().GetNstationsActive()); + frAlgo.fWData.RecoTracks().clear(); + frAlgo.fWData.RecoTracks().reserve(2 * nHits / frAlgo.GetParameters().GetNstationsActive()); - frAlgo.fSliceRecoHits.clear(); - frAlgo.fSliceRecoHits.reserve(2 * nHits); + frAlgo.fWData.RecoHitIndices().clear(); + frAlgo.fWData.RecoHitIndices().reserve(2 * nHits); frAlgo.fTrackCandidates.clear(); frAlgo.fTrackCandidates.reserve(nHits / 10); @@ -232,15 +231,20 @@ void TrackFinderWindow::CaTrackFinderSlice() if (xStep > 0.3 * scale) xStep = 0.3 * scale; if (xStep < 0.01 * scale) xStep = 0.01 * scale; - frAlgo.vGrid[iS].BuildBins(gridMinX, gridMaxX, gridMinY, gridMaxY, xStep, yStep); - frAlgo.vGrid[iS].StoreHits(frAlgo.fWindowHits, frAlgo.fStationHitsStartIndex[iS], frAlgo.fStationNhits[iS], - frAlgo.fvHitKeyFlags); + auto& grid = frAlgo.fWData.Grid(iS); + grid.BuildBins(gridMinX, gridMaxX, gridMinY, gridMaxY, xStep, yStep); + /* clang-format off */ + grid.StoreHits(frAlgo.fWData.Hits(), + frAlgo.fWData.HitStartIndexOnStation(iS), + frAlgo.fWData.NofHitsOnStation(iS), + frAlgo.fvHitKeyFlags); + /* clang-format on */ } // ---- Loop over Track Finder iterations ----------------------------------------------------------------// - + // TODO: replace with exception assert(frAlgo.fNFindIterations == (int) frAlgo.GetParameters().GetCAIterations().size()); for (frAlgo.fCurrentIterationIndex = 0; @@ -253,18 +257,19 @@ void TrackFinderWindow::CaTrackFinderSlice() if (frAlgo.fCurrentIterationIndex > 0) { for (int ista = 0; ista < frAlgo.GetParameters().GetNstationsActive(); ++ista) { - frAlgo.vGrid[ista].RemoveUsedHits(frAlgo.fWindowHits, frAlgo.fvHitKeyFlags); + frAlgo.fWData.Grid(ista).RemoveUsedHits(frAlgo.fWData.Hits(), frAlgo.fvHitKeyFlags); } } - frAlgo.fIsWindowHitSuppressed.reset(frAlgo.fWindowHits.size(), 0); - + int nHits = frAlgo.fWData.Hits().size(); + // --> frAlgo.fIsWindowHitSuppressed.reset(frAlgo.fWindowHits.size(), 0); + frAlgo.fWData.ResetHitSuppressionFlags(); // TODO: ??? No effect? for (int j = 0; j < frAlgo.GetParameters().GetNstationsActive(); j++) { frAlgo.fTriplets[j].clear(); } - frAlgo.fHitFirstTriplet.reset(frAlgo.fWindowHits.size(), 0); - frAlgo.fHitNtriplets.reset(frAlgo.fWindowHits.size(), 0); + frAlgo.fHitFirstTriplet.reset(nHits, 0); + frAlgo.fHitNtriplets.reset(nHits, 0); { // #pragma omp task @@ -332,9 +337,10 @@ void TrackFinderWindow::CaTrackFinderSlice() for (int istal = frAlgo.GetParameters().GetNstationsActive() - 2; istal >= frAlgo.fFirstCAstation; istal--) { // start with downstream chambers - Tindex nGridEntriesL = frAlgo.vGrid[istal].GetEntries().size(); + const auto& grid = frAlgo.fWData.Grid(istal); + Tindex nGridEntriesL = grid.GetEntries().size(); for (Tindex iel = 0; iel < nGridEntriesL; ++iel) { - ca::HitIndex_t ihitl = frAlgo.vGrid[istal].GetEntries()[iel].GetObjectId(); + ca::HitIndex_t ihitl = grid.GetEntries()[iel].GetObjectId(); auto oldSize = frAlgo.fTriplets[istal].size(); for (auto& pattern : staPattern) { constructor.CreateTripletsForHit(istal, istal + pattern.first, istal + pattern.second, ihitl); @@ -436,7 +442,7 @@ void TrackFinderWindow::CaTrackFinderSlice() frAlgo.fTrackCandidates.clear(); - for (const auto& h : frAlgo.fWindowHits) { + for (const auto& h : frAlgo.fWData.Hits()) { frAlgo.fHitKeyToTrack[h.FrontKey()] = -1; frAlgo.fHitKeyToTrack[h.BackKey()] = -1; } @@ -451,8 +457,8 @@ void TrackFinderWindow::CaTrackFinderSlice() for (Tindex itrip = 0; itrip < (Tindex) frAlgo.fTriplets[istaF].size(); ++itrip) { ca::Triplet& first_trip = (frAlgo.fTriplets[istaF][itrip]); - if (frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[first_trip.GetLHit()].FrontKey()] - || frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[first_trip.GetLHit()].BackKey()]) { + const auto& fstTripLHit = frAlgo.fWData.Hit(first_trip.GetLHit()); + if (frAlgo.fvHitKeyFlags[fstTripLHit.FrontKey()] || frAlgo.fvHitKeyFlags[fstTripLHit.BackKey()]) { continue; } @@ -460,7 +466,7 @@ void TrackFinderWindow::CaTrackFinderSlice() int minNhits = frAlgo.fpCurrentIteration->GetMinNhits(); - if (frAlgo.fWindowHits[first_trip.GetLHit()].Station() == 0) { + if (fstTripLHit.Station() == 0) { minNhits = frAlgo.fpCurrentIteration->GetMinNhitsStation0(); } if (frAlgo.fpCurrentIteration->GetTrackFromTripletsFlag()) { @@ -480,7 +486,7 @@ void TrackFinderWindow::CaTrackFinderSlice() curr_tr.SetChi2(0.); curr_tr.SetStation(0); curr_tr.ResetHits(); - curr_tr.AddHit(frAlgo.fWindowHits[first_trip.GetLHit()].Id()); + curr_tr.AddHit(fstTripLHit.Id()); curr_tr.SetChi2(first_trip.GetChi2()); best_tr = curr_tr; @@ -650,11 +656,11 @@ void TrackFinderWindow::CaTrackFinderSlice() frAlgo.fvHitKeyFlags[hit.FrontKey()] = 1; frAlgo.fvHitKeyFlags[hit.BackKey()] = 1; - frAlgo.fSliceRecoHits.push_back(iHit); + frAlgo.fWData.RecoHitIndices().push_back(iHit); } Track t; t.fNofHits = tr.NofHits(); - frAlgo.fSliceRecoTracks.push_back(t); + frAlgo.fWData.RecoTracks().push_back(t); if (0) { // SG debug std::stringstream s; s << "store track " << iCandidate << " chi2= " << tr.Chi2() << "\n"; @@ -671,9 +677,9 @@ void TrackFinderWindow::CaTrackFinderSlice() // suppress strips of suppressed hits - for (unsigned int ih = 0; ih < frAlgo.fWindowHits.size(); ih++) { - if (frAlgo.fIsWindowHitSuppressed[ih]) { - const ca::Hit& hit = frAlgo.fWindowHits[ih]; + for (unsigned int ih = 0; ih < frAlgo.fWData.Hits().size(); ih++) { + if (frAlgo.fWData.IsHitSuppressed(ih)) { + const ca::Hit& hit = frAlgo.fWData.Hit(ih); frAlgo.fvHitKeyFlags[hit.FrontKey()] = 1; frAlgo.fvHitKeyFlags[hit.BackKey()] = 1; } @@ -684,7 +690,7 @@ void TrackFinderWindow::CaTrackFinderSlice() frAlgo.fTrackFitter.FitCaTracks(); // Merge clones - frAlgo.fCloneMerger.Exec(frAlgo.fSliceRecoTracks, frAlgo.fSliceRecoHits); + frAlgo.fCloneMerger.Exec(frAlgo.fWData.RecoTracks(), frAlgo.fWData.RecoHitIndices()); } @@ -708,22 +714,15 @@ void TrackFinderWindow::CAFindTrack(int ista, ca::Branch& best_tr, const ca::Tri { // -- finish with current track // add rest of hits - const ca::HitIndex_t& ihitm = curr_trip->GetMHit(); - const ca::HitIndex_t& ihitr = curr_trip->GetRHit(); - - if (!(frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[ihitm].FrontKey()] - || frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[ihitm].BackKey()])) { + const auto& hitM = frAlgo.fWData.Hit(curr_trip->GetMHit()); + const auto& hitR = frAlgo.fWData.Hit(curr_trip->GetRHit()); - // curr_tr.Hits.push_back(fGridHitIds[ihitm]); - - curr_tr.AddHit(frAlgo.fWindowHits[ihitm].Id()); + if (!(frAlgo.fvHitKeyFlags[hitM.FrontKey()] || frAlgo.fvHitKeyFlags[hitM.BackKey()])) { + curr_tr.AddHit(hitM.Id()); } - if (!(frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[ihitr].FrontKey()] - || frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[ihitr].BackKey()])) { - - //curr_tr.Hits.push_back(fGridHitIds[ihitr]); - curr_tr.AddHit(frAlgo.fWindowHits[ihitr].Id()); + if (!(frAlgo.fvHitKeyFlags[hitR.FrontKey()] || frAlgo.fvHitKeyFlags[hitR.BackKey()])) { + curr_tr.AddHit(hitR.Id()); } //if( curr_tr.NofHits() < min_best_l - 1 ) return; // suppouse that only one hit can be added by extender @@ -745,7 +744,6 @@ void TrackFinderWindow::CAFindTrack(int ista, ca::Branch& best_tr, const ca::Tri || ((curr_tr.NofHits() == best_tr.NofHits()) && (curr_tr.Chi2() < best_tr.Chi2()))) { best_tr = curr_tr; } - return; } else //MEANS level ! = 0 @@ -764,8 +762,8 @@ void TrackFinderWindow::CAFindTrack(int ista, ca::Branch& best_tr, const ca::Tri fscal dchi2 = 0.; if (!checkTripletMatch(*curr_trip, new_trip, dchi2)) continue; - if (frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[new_trip.GetLHit()].FrontKey()] - || frAlgo.fvHitKeyFlags[frAlgo.fWindowHits[new_trip.GetLHit()].BackKey()]) { //hits are used + const auto& hitL = frAlgo.fWData.Hit(new_trip.GetLHit()); // left hit of new triplet + if (frAlgo.fvHitKeyFlags[hitL.FrontKey()] || frAlgo.fvHitKeyFlags[hitL.BackKey()]) { //hits are used // no used hits allowed -> compare and store track if ((curr_tr.NofHits() > best_tr.NofHits()) || ((curr_tr.NofHits() == best_tr.NofHits()) && (curr_tr.Chi2() < best_tr.Chi2()))) { @@ -810,7 +808,7 @@ void TrackFinderWindow::CAFindTrack(int ista, ca::Branch& best_tr, const ca::Tri // add new hit new_tr[ista] = curr_tr; - new_tr[ista].AddHit(frAlgo.fWindowHits[new_trip.GetLHit()].Id()); + new_tr[ista].AddHit(hitL.Id()); const int new_ista = ista + new_trip.GetMSta() - new_trip.GetLSta(); new_tr[ista].SetChi2(new_chi2); diff --git a/algo/ca/core/tracking/CaTrackFinderWindow.h b/algo/ca/core/tracking/CaTrackFinderWindow.h index a3fdc588900dfd5f49dd8f7d5b9f2c3e18f122ba..ee5d3d717db29de89350639ebf100b669c976875 100644 --- a/algo/ca/core/tracking/CaTrackFinderWindow.h +++ b/algo/ca/core/tracking/CaTrackFinderWindow.h @@ -10,12 +10,12 @@ #pragma once // include this header only once per compilation unit #include "CaBranch.h" +#include "CaGrid.h" #include "CaHit.h" #include "CaSimd.h" #include "CaTrackParam.h" #include "CaVector.h" - namespace cbm::algo::ca { @@ -68,8 +68,6 @@ namespace cbm::algo::ca /// Data members ca::Framework& frAlgo; ///< Reference to the main track finder algorithm class - - private: static constexpr bool fDebug = false; // print debug info }; diff --git a/algo/ca/core/tracking/CaTrackFitter.cxx b/algo/ca/core/tracking/CaTrackFitter.cxx index c5a8e6386a49af27a060e22e1e214f743b78c1de..7167e5a0950907e4b90456a8f8623473c6358453 100644 --- a/algo/ca/core/tracking/CaTrackFitter.cxx +++ b/algo/ca/core/tracking/CaTrackFitter.cxx @@ -31,7 +31,7 @@ TrackFitter::~TrackFitter() {} void TrackFitter::FitCaTracks() { // LOG(info) << " Start CA Track Fitter "; - int start_hit = 0; // for interation in frAlgo.fSliceRecoHits[] + int start_hit = 0; // for interation in frAlgo.fWData.RecoHitIndices() // static ca::FieldValue fldB0, fldB1, fldB2 _fvecalignment; // static ca::FieldRegion fld _fvecalignment; @@ -102,19 +102,19 @@ void TrackFitter::FitCaTracks() mxy[ista].SetCov(1., 0., 1.); } - unsigned short N_vTracks = frAlgo.fSliceRecoTracks.size(); // number of tracks processed per one SSE register + unsigned short N_vTracks = frAlgo.fWData.RecoTracks().size(); // number of tracks processed per one SSE register for (unsigned short itrack = 0; itrack < N_vTracks; itrack += fvec::size()) { if (N_vTracks - itrack < static_cast<unsigned short>(fvec::size())) nTracks_SIMD = N_vTracks - itrack; for (int i = 0; i < nTracks_SIMD; i++) { - t[i] = &frAlgo.fSliceRecoTracks[itrack + i]; + t[i] = &frAlgo.fWData.RecoTrack(itrack + i); } // fill the rest of the SIMD vectors with something reasonable for (uint i = nTracks_SIMD; i < fvec::size(); i++) { - t[i] = &frAlgo.fSliceRecoTracks[itrack]; + t[i] = &frAlgo.fWData.RecoTrack(itrack); } // get hits of current track @@ -133,7 +133,7 @@ void TrackFitter::FitCaTracks() for (int ih = 0; ih < nHitsTrack; ih++) { - const ca::Hit& hit = frAlgo.GetInputData().GetHit(frAlgo.fSliceRecoHits[start_hit++]); + const ca::Hit& hit = frAlgo.GetInputData().GetHit(frAlgo.fWData.RecoHitIndex(start_hit++)); const int ista = hit.Station(); //if (sta[ista].fieldStatus) { isFieldPresent[iVec] = true; } diff --git a/algo/ca/core/tracking/CaTripletConstructor.cxx b/algo/ca/core/tracking/CaTripletConstructor.cxx index 529a9c448038bbc142075ee42c5d39a431c67333..a3535c533f38ca846b6f1c38d157dd1de15b029e 100644 --- a/algo/ca/core/tracking/CaTripletConstructor.cxx +++ b/algo/ca/core/tracking/CaTripletConstructor.cxx @@ -96,7 +96,7 @@ void TripletConstructor::CreateTripletsForHit(int istal, int istam, int istar, c if (!fIsInitialized) return; fIhitL = iHitL; - ca::Hit& hitl = fAlgo->fWindowHits[fIhitL]; + const auto& hitl = fAlgo->fWData.Hit(fIhitL); // fit a straight line through the target and the left hit. { @@ -217,9 +217,11 @@ void TripletConstructor::FitDoublets() for (unsigned int i2 = 0; i2 < hitsMtmp.size(); i2++) { const ca::HitIndex_t indM = hitsMtmp[i2]; - const ca::Hit& hitm = fAlgo->fWindowHits[indM]; + const ca::Hit& hitm = fAlgo->fWData.Hit(indM); - if (fAlgo->fIsWindowHitSuppressed[indM]) continue; + if (fAlgo->fWData.IsHitSuppressed(indM)) { + continue; + } TrackParamV& T2 = fFit.Tr(); T2 = fTrL; @@ -255,7 +257,7 @@ void TripletConstructor::FitDoublets() for (unsigned int iClone = i2 + 1; iClone < hitsMtmp.size(); iClone++) { int indClone = hitsMtmp[iClone]; - const ca::Hit& hitClone = fAlgo->fWindowHits[indClone]; + const ca::Hit& hitClone = fAlgo->fWData.Hit(indClone); fscal dz = hitClone.Z() - T2.Z()[0]; @@ -282,7 +284,7 @@ void TripletConstructor::FitDoublets() continue; } } - fAlgo->fIsWindowHitSuppressed[indClone] = 1; + fAlgo->fWData.SuppressHit(indClone); } } @@ -333,7 +335,7 @@ void TripletConstructor::FindRightHit() if constexpr (fDebugDublets) { ca::HitIndex_t iwhit[2] = {fIhitL, fHitsM_2[i2]}; - ca::HitIndex_t ihit[2] = {fAlgo->fWindowHits[iwhit[0]].Id(), fAlgo->fWindowHits[iwhit[1]].Id()}; + ca::HitIndex_t ihit[2] = {fAlgo->fWData.Hit(iwhit[0]).Id(), fAlgo->fWData.Hit(iwhit[1]).Id()}; int ih0 = ihit[0]; int ih1 = ihit[1]; @@ -403,11 +405,11 @@ void TripletConstructor::FindRightHit() for (unsigned int ih = 0; ih < collectedHits.size(); ih++) { ca::HitIndex_t irh = collectedHits[ih]; if constexpr (fDebugDublets) { - ca::HitIndex_t ihit = fAlgo->fWindowHits[irh].Id(); + ca::HitIndex_t ihit = fAlgo->fWData.Hit(irh).Id(); const ca::Hit& h = fAlgo->fInputData.GetHit(ihit); LOG(info) << " hit " << ihit << " " << h.ToString(); } - if (fAlgo->fIsWindowHitSuppressed[irh]) { + if (fAlgo->fWData.IsHitSuppressed(irh)) { if constexpr (fDebugDublets) { LOG(info) << " the hit is suppressed"; } @@ -469,8 +471,8 @@ void TripletConstructor::FitTriplets() // prepare data ca::HitIndex_t iwhit[NHits] = {fIhitL, fHitsM_3[i3], fHitsR_3[i3]}; - ca::HitIndex_t ihit[NHits] = {fAlgo->fWindowHits[iwhit[0]].Id(), fAlgo->fWindowHits[iwhit[1]].Id(), - fAlgo->fWindowHits[iwhit[2]].Id()}; + ca::HitIndex_t ihit[NHits] = {fAlgo->fWData.Hit(iwhit[0]).Id(), fAlgo->fWData.Hit(iwhit[1]).Id(), + fAlgo->fWData.Hit(iwhit[2]).Id()}; if (fAlgo->fParameters.DevIsMatchTripletsViaMc()) { int mc1 = fAlgo->GetMcTrackIdForCaHit(ihit[0]); @@ -652,9 +654,9 @@ void TripletConstructor::StoreTriplets() const ca::HitIndex_t ihitm = fHitsM_3[i3]; const ca::HitIndex_t ihitr = fHitsR_3[i3]; - CBMCA_DEBUG_ASSERT(ihitl < fAlgo->fStationHitsStartIndex[fIstaL + 1]); - CBMCA_DEBUG_ASSERT(ihitm < fAlgo->fStationHitsStartIndex[fIstaM + 1]); - CBMCA_DEBUG_ASSERT(ihitr < fAlgo->fStationHitsStartIndex[fIstaR + 1]); + CBMCA_DEBUG_ASSERT(ihitl < fAlgo->fWData.HitStartIndexOnStation(fIstaL + 1)); + CBMCA_DEBUG_ASSERT(ihitm < fAlgo->fWData.HitStartIndexOnStation(fIstaM + 1)); + CBMCA_DEBUG_ASSERT(ihitr < fAlgo->fWData.HitStartIndexOnStation(fIstaR + 1)); if (!fAlgo->fpCurrentIteration->GetTrackFromTripletsFlag()) { if (chi2 > fAlgo->fTripletFinalChi2Cut) { @@ -703,7 +705,7 @@ void TripletConstructor::CollectHits(const TrackParamV& Tr, const int iSta, cons const fscal timeError2 = T.C55()[0]; const fscal time = T.Time()[0]; - const auto& grid = fAlgo->vGrid[iSta]; + const auto& grid = fAlgo->fWData.Grid(iSta); ca::GridArea area(grid, T.X()[0], T.Y()[0], (sqrt(Pick_m22 * T.C00()) + grid.GetMaxRangeX() + fAlgo->fMaxDZ * abs(T.Tx()))[0], (sqrt(Pick_m22 * T.C11()) + grid.GetMaxRangeY() + fAlgo->fMaxDZ * abs(T.Ty()))[0]); @@ -721,9 +723,11 @@ void TripletConstructor::CollectHits(const TrackParamV& Tr, const int iSta, cons while (area.GetNextObjectId(ih) && ((int) collectedHits.size() < maxNhits)) { // loop over station hits - if (fAlgo->fIsWindowHitSuppressed[ih]) continue; + if (fAlgo->fWData.IsHitSuppressed(ih)) { + continue; + } - const ca::Hit& hit = fAlgo->fWindowHits[ih]; + const ca::Hit& hit = fAlgo->fWData.Hit(ih); if constexpr (fDebugCollectHits) { LOG(info) << "hit in the search area: " << hit.ToString(); LOG(info) << " check the hit.. "; diff --git a/reco/L1/L1Algo/utils/L1AlgoDraw.cxx b/reco/L1/L1Algo/utils/L1AlgoDraw.cxx index fc384e6de3b1737fed82126ece530ea470c191d5..5cf56fb5f334952e173cd2df5a1706a07bc224d2 100644 --- a/reco/L1/L1Algo/utils/L1AlgoDraw.cxx +++ b/reco/L1/L1Algo/utils/L1AlgoDraw.cxx @@ -79,17 +79,18 @@ void L1AlgoDraw::InitL1Draw(ca::Framework* algo_) vHits.clear(); vHitsQa.clear(); - vHits.reserve(algo->fWindowHits.size()); - vHitsQa.reserve(algo->fWindowHits.size()); - for (unsigned int i = 0; i < algo->fWindowHits.size(); i++) { - vHits.push_back(algo->fWindowHits[i]); - int iQaHit = algo->fWindowHits[i].Id(); + auto nHits = algo->fWData.Hits().size(); + vHits.reserve(nHits); + vHitsQa.reserve(nHits); + for (const auto& hit : algo->fWData.Hits()) { + vHits.push_back(hit); + int iQaHit = hit.Id(); vHitsQa.push_back(CbmL1::Instance()->GetQaHits()[iQaHit]); } NStations = algo->GetParameters().GetNstationsActive(); for (int i = 0; i < NStations; i++) { - HitsStartIndex[i] = algo->fStationHitsStartIndex[i]; - HitsStopIndex[i] = algo->fStationHitsStartIndex[i] + algo->fStationNhits[i]; + HitsStartIndex[i] = algo->fWData.HitStartIndexOnStation(i); + HitsStopIndex[i] = HitsStartIndex[i] + algo->fWData.NofHitsOnStation(i); vStations[i] = algo->GetParameters().GetStation(i); } } @@ -244,8 +245,8 @@ void L1AlgoDraw::DrawRecoTracks() // CbmL1 &L1 = *CbmL1::Instance(); int curRecoHit = 0; - cbm::algo::ca::Vector<ca::HitIndex_t>& recoHits = algo->fSliceRecoHits; - for (vector<Track>::iterator it = algo->fSliceRecoTracks.begin(); it != algo->fSliceRecoTracks.end(); ++it) { + cbm::algo::ca::Vector<ca::HitIndex_t>& recoHits = algo->fWData.RecoHitIndices(); + for (vector<Track>::iterator it = algo->fWData.RecoTracks().begin(); it != algo->fWData.RecoTracks().end(); ++it) { Track& T = *it; int nHits = T.fNofHits; // if (nHits > 5) continue; // draw clones @@ -483,8 +484,8 @@ void L1AlgoDraw::DrawDoublet(int il, int ir) void L1AlgoDraw::DrawInfo() { cout << " vHits.size = " << algo->GetInputData().GetNhits() << endl; - cout << " vRecoHits.size = " << algo->fSliceRecoHits.size() << endl; - cout << " vTracks.size = " << algo->fSliceRecoTracks.size() << endl; + cout << " vRecoHits.size = " << algo->fWData.RecoHitIndices().size() << endl; + cout << " vTracks.size = " << algo->fWData.RecoTracks().size() << endl; } void L1AlgoDraw::DrawTarget()