diff --git a/algo/ca/core/CMakeLists.txt b/algo/ca/core/CMakeLists.txt index c65a8e0584ebbc61063bea2187babb1755e876da..26fe08f9084fa4acf26806b0d6a620d113277e7f 100644 --- a/algo/ca/core/CMakeLists.txt +++ b/algo/ca/core/CMakeLists.txt @@ -22,6 +22,7 @@ target_compile_definitions(CaCore PUBLIC NO_ROOT) target_link_libraries(CaCore Vc::Vc Boost::serialization + #external::fles_logging # conflicting with FairLogger in libL1.so ) install(TARGETS CaCore DESTINATION lib) @@ -36,6 +37,7 @@ install( utils/CaSimd.h utils/CaSimdVc.h utils/CaSimdPseudo.h + utils/CaVector.h DESTINATION include/ ) diff --git a/algo/ca/core/utils/CaSimdPseudo.h b/algo/ca/core/utils/CaSimdPseudo.h index 74f64da5afcd7d01aba7fcd850d8fa2c1b22a0e1..f1cdd98c3132ccd59b9b42a8b1aa031a037222eb 100644 --- a/algo/ca/core/utils/CaSimdPseudo.h +++ b/algo/ca/core/utils/CaSimdPseudo.h @@ -5,9 +5,10 @@ #ifndef Ca_CaSimdPseudo_H #define Ca_CaSimdPseudo_H +#include <boost/serialization/access.hpp> + #include <iomanip> #include <iostream> -#include <boost/serialization/access.hpp> #include <cmath> diff --git a/algo/ca/core/utils/CaSimdVc.h b/algo/ca/core/utils/CaSimdVc.h index 668622c655b05d69cb67846380d0b8c70c618d4d..f8b22b5d1cdafe63a7a9ec68b31766f1232fd5bc 100644 --- a/algo/ca/core/utils/CaSimdVc.h +++ b/algo/ca/core/utils/CaSimdVc.h @@ -5,11 +5,12 @@ #ifndef Ca_CaSimdVc_H #define Ca_CaSimdVc_H -#include "Vc/Vc" #include <boost/serialization/access.hpp> #include <boost/serialization/array.hpp> #include <boost/serialization/split_free.hpp> +#include "Vc/Vc" + namespace cbm::algo::ca { typedef Vc::float_v fvec; @@ -30,7 +31,9 @@ namespace boost::serialization void save(Archive& ar, const cbm::algo::ca::fvec& vect, unsigned int) { std::array<cbm::algo::ca::fscal, cbm::algo::ca::fvec::size()> buffer; - for (size_t i = 0; i < cbm::algo::ca::fvec::size(); ++i) { buffer[i] = vect[i]; } + for (size_t i = 0; i < cbm::algo::ca::fvec::size(); ++i) { + buffer[i] = vect[i]; + } ar << buffer; } @@ -49,7 +52,7 @@ namespace boost::serialization { split_free(ar, vect, version); } -} +} // namespace boost::serialization #endif diff --git a/algo/ca/core/utils/CaVector.h b/algo/ca/core/utils/CaVector.h new file mode 100644 index 0000000000000000000000000000000000000000..0746c0d5b7aae7337112ebb6e082784984fd9cd7 --- /dev/null +++ b/algo/ca/core/utils/CaVector.h @@ -0,0 +1,296 @@ +/* Copyright (C) 2021-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Sergey Gorbunov [committer], Sergei Zharko */ + +/// @file CaVector.h +/// @author Sergey Gorbunov +/// @date 2021-06-16 + +#ifndef Ca_CaVector_H +#define Ca_CaVector_H + + +#ifndef FAST_CODE +#include "Logger.h" +#endif +#include <boost/serialization/access.hpp> +#include <boost/serialization/base_object.hpp> +#include <boost/serialization/string.hpp> +#include <boost/serialization/vector.hpp> + +#include <sstream> + +namespace cbm::algo::ca +{ + /// \class Vector + /// + /// ca::Vector class is a wrapper around std::vector. + /// It does the following: + /// 1. gives names to vectors for better debugging + /// 2. controls the out-of-range access in debug mode + /// 3. supresses methods that are currently not controlled + /// 4. warns when slow memory operations are called, + /// i.e. when the preallocated capacity is reached and the entire vector should be copied to a new place + /// 5. blocks usage of boolean vectors, as they have a special + /// space-optimized but slow implementation in std. (Use ca::Vector<char> instead). + /// + template<class T> + class Vector : private std::vector<T> { + friend class boost::serialization::access; + + public: + typedef std::vector<T> Tbase; + + /// \brief Generic constructor from vairadic parameter list + template<typename... Tinput> + Vector(Tinput... value) : Tbase(value...) + { + } + + /// \brief Generic constructor from vairadic parameter list including the name of the vector + template<typename... Tinput> + Vector(const char* name, Tinput... value) : Tbase(value...) + , fName(name) + { + } + + /// \brief Constructor to make initializations like ca::Vector<int> myVector {"MyVector", {1, 2, 3}} + Vector(const std::string& name, std::initializer_list<T> init) : Tbase(init), fName(name) {} + + /// \brief Copy constructor + Vector(const Vector& v) : Tbase() { *this = v; } + + /// \brief Move constructor + Vector(Vector&& v) noexcept : Tbase(std::move(v)), fName(std::move(v.fName)) {} + + /// \brief Copy assignment operator + Vector& operator=(const Vector& v) + { + if (this != &v) { + fName = v.fName; + Tbase::reserve(v.capacity()); // make sure that the capacity is transmitted + Tbase::assign(v.begin(), v.end()); + } + return *this; + } + + /// \brief Move assignment operator + Vector& operator=(Vector&& v) noexcept + { + if (this != &v) { + std::swap(fName, v.fName); + Tbase::swap(v); + } + return *this; + } + + /// \brief Swap operator + void swap(Vector& v) noexcept + { + if (this != &v) { + std::swap(fName, v.fName); + Tbase::swap(v); + } + } + + /// \brief Sets the name of the vector + /// \param s Name of the vector + void SetName(const std::string& s) { fName = s; } + + /// \brief Sets the name of the vector + /// \param s Name of the vector (string stream) + void SetName(const std::basic_ostream<char>& s) + { + // helps to set a composed name in a single line via: + // SetName(std::stringstream()<<"my name "<<..whatever..); + fName = dynamic_cast<const std::stringstream&>(s).str(); + } + + /// \brief Gets name of the vector + std::string GetName() const + { + std::string s = "ca::Vector<"; + s += fName + "> "; + return s; + } + + /// \brief Clears vector and resizes it to the selected size with selected values + /// \param count New size of the vector + /// \param value Variadic list of the parameters to pass to the base std::vector::resize function + template<typename... Tinput> + void reset(std::size_t count, Tinput... value) + { + // does the same as Tbase::assign(), but works with the default T constructor too + // (no second parameter) + Tbase::clear(); + Tbase::resize(count, value...); + } + + /// \brief Enlarges the vector to the new size + /// \param count New size of the vector + /// \param value Value of the new elements + template<typename... Tinput> + void enlarge(std::size_t count, Tinput... value) + { + if (count < Tbase::size()) { + LOG(fatal) << "ca::Vector \"" << fName << "\"::enlarge(" << count + << "): the new size is smaller than the current one " << Tbase::size() << ", something goes wrong."; + assert(count >= Tbase::size()); + } + if ((!Tbase::empty()) && (count > Tbase::capacity())) { + LOG(warning) << "ca::Vector \"" << fName << "\"::enlarge(" << count << "): allocated capacity of " + << Tbase::capacity() << " is reached, the vector of size " << Tbase::size() + << " will be copied to the new place."; + } + Tbase::resize(count, value...); + } + + /// \brief Reduces the vector to a given size + /// \param count Size of the new vector + void reduce(std::size_t count) + { + if (count > Tbase::size()) { + LOG(fatal) << "ca::Vector \"" << fName << "\"::reduce(" << count + << "): the new size is bigger than the current one " << Tbase::size() << ", something goes wrong."; + assert(count < Tbase::size()); + } + Tbase::resize(count); + } + + /// \brief Reserves a new size for the vector + /// \param count New size of the vector + void reserve(std::size_t count) + { + if (!Tbase::empty()) { + LOG(fatal) << "ca::Vector \"" << fName << "\"::reserve(" << count << "): the vector is not empty; " + << " it will be copied to the new place."; + assert(Tbase::empty()); + } + Tbase::reserve(count); + } + + /// \brief Pushes back a value to the vector + /// \param value New value + /// \note Raises a warning, if the vector re-alocates memory + template<typename Tinput> + void push_back(Tinput value) + { +#ifndef FAST_CODE + if (Tbase::size() >= Tbase::capacity()) { + LOG(warning) << "ca::Vector \"" << fName << "\"::push_back(): allocated capacity of " << Tbase::capacity() + << " is reached, re-allocate and copy."; + } +#endif + Tbase::push_back(value); + } + + /// \brief Pushes back a value to the vector without testing for the memory re-alocation + /// \param value New value + template<typename Tinput> + void push_back_no_warning(Tinput value) + { + Tbase::push_back(value); + } + + /// \brief Creates a parameter in the end of the vector + /// \param value Variadic list of the parameters, which are passed to the element constructor + template<typename... Tinput> + void emplace_back(Tinput&&... value) + { +#ifndef FAST_CODE + if (Tbase::size() >= Tbase::capacity()) { + LOG(warning) << "ca::Vector \"" << fName << "\"::emplace_back(): allocated capacity of " << Tbase::capacity() + << " is reached, re-allocate and copy."; + } +#endif + Tbase::emplace_back(value...); + } + + /// \brief Mutable access to the element by its index + /// \param pos Index of the element + T& operator[](std::size_t pos) + { +#ifndef FAST_CODE + if (pos >= Tbase::size()) { + LOG(fatal) << "ca::Vector \"" << fName << "\": trying to access element " << pos + << " outside of the vector of the size of " << Tbase::size(); + assert(pos < Tbase::size()); + } +#endif + return Tbase::operator[](pos); + } + + /// \brief Constant access to the element by its index + /// \param pos Index of the element + const T& operator[](std::size_t pos) const + { +#ifndef FAST_CODE + if (pos >= Tbase::size()) { + LOG(fatal) << "ca::Vector \"" << fName << "\": trying to access element " << pos + << " outside of the vector of the size of " << Tbase::size(); + assert(pos < Tbase::size()); + } +#endif + return Tbase::operator[](pos); + } + + /// \brief Mutable access to the last element of the vector + T& back() + { +#ifndef FAST_CODE + if (Tbase::size() == 0) { + LOG(fatal) << "ca::Vector \"" << fName << "\": trying to access element of an empty vector"; + assert(Tbase::size() > 0); + } +#endif + return Tbase::back(); + } + + /// \brief Constant access to the last element of the vector + const T& back() const + { +#ifndef FAST_CODE + if (Tbase::size() == 0) { + LOG(fatal) << "ca::Vector \"" << fName << "\": trying to access element of an empty vector"; + assert(Tbase::size() > 0); + } +#endif + return Tbase::back(); + } + + using Tbase::begin; + using Tbase::capacity; + using Tbase::clear; + using Tbase::end; + using Tbase::insert; //TODO:: make it private + using Tbase::pop_back; + using Tbase::rbegin; + using Tbase::reserve; + using Tbase::size; + using typename Tbase::iterator; + + private: + std::string fName {"no name"}; ///< Name of the vector + using Tbase::assign; // use reset() instead + using Tbase::at; + using Tbase::resize; + + /// \brief Serialization function for the vector + template<class Archive> + void serialize(Archive& ar, const unsigned int /*version*/) + { + ar& boost::serialization::base_object<Tbase>(*this); + ar& fName; + } + }; + + /// + /// std::vector<bool> has a special implementation that is space-optimized + /// and therefore slow and not thread-safe. + /// That is why one should use ca::Vector<char> instead. + /// + template<> + class Vector<bool> { + }; +} // namespace cbm::algo::ca +#endif diff --git a/reco/KF/KF.cmake b/reco/KF/KF.cmake index d4aa313d967a136a3fa4e8e617205bd2edd7b440..273f26b0666816d20f457fa011c1f3dcb3c9dc65 100644 --- a/reco/KF/KF.cmake +++ b/reco/KF/KF.cmake @@ -91,6 +91,7 @@ ENDIF (SSE_FOUND) set(LIBRARY_NAME KF) set(LINKDEF ${LIBRARY_NAME}LinkDef.h) set(PUBLIC_DEPENDENCIES + CaCore CbmBase CbmData CbmRecoBase diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt index d30f029ba0d8dc843e2213ac31ee9741731f00bf..1ad447ddd20eb68ea5437e2cbc166ed5ecd75ebe 100644 --- a/reco/L1/CMakeLists.txt +++ b/reco/L1/CMakeLists.txt @@ -72,7 +72,6 @@ set(SRCS catools/CaToolsMCData.cxx catools/CaToolsMCPoint.cxx catools/CaToolsMCTrack.cxx - catools/CaToolsQa.cxx catools/CaToolsDebugger.cxx catools/CaToolsWindowFinder.cxx catools/CaToolsWFExpression.cxx @@ -110,13 +109,13 @@ set(HEADERS CbmL1TrackPar.h CbmL1Vtx.h L1Algo/L1Def.h - L1Algo/L1Vector.h L1Algo/L1EArray.h L1Algo/utils/CaUvConverter.h L1Algo/utils/CaMonitor.h catools/CaToolsWindowFinder.h catools/CaToolsLinkKey.h catools/CaToolsHitRecord.h + catools/CaToolsDef.h #utils/CbmCaIdealHitProducer.h catools/CaToolsMaterialHelper.h qa/CbmCaInputQaBase.h @@ -226,7 +225,6 @@ install(FILES L1Algo/L1Algo.h L1Algo/L1Branch.h L1Algo/L1Field.h L1Algo/L1Hit.h - L1Algo/L1Vector.h L1Algo/L1EArray.h DESTINATION include/L1Algo ) diff --git a/reco/L1/CbmCaMCModule.cxx b/reco/L1/CbmCaMCModule.cxx index 869ec615c1cfb49cc8e9c57fe18e85a5b8a01f9b..5cd92cd0e9078671fa2ff1217fbe919703875318 100644 --- a/reco/L1/CbmCaMCModule.cxx +++ b/reco/L1/CbmCaMCModule.cxx @@ -48,6 +48,7 @@ using cbm::algo::ca::constants::clrs::CL; // clear log using cbm::algo::ca::constants::clrs::RDb; // red bold log using cbm::ca::MCModule; +using cbm::ca::tools::MCTrack; // --------------------------------------------------------------------------------------------------------------------- // @@ -270,7 +271,7 @@ void MCModule::MatchRecoAndMCTracks() if (double(nHitsTrkMc) > double(nHitsTrkRe) * maxPurity) { maxPurity = double(nHitsTrkMc) / double(nHitsTrkRe); } - ::ca::tools::MCTrack& trkMc = fpMCData->GetTrack(iTmc); + auto& trkMc = fpMCData->GetTrack(iTmc); // Match reconstructed and MC tracks, if purity with this MC track passes the threshold if (double(nHitsTrkMc) >= CbmL1Constants::MinPurity * double(nHitsTrkRe)) { @@ -470,7 +471,7 @@ void MCModule::ReadMCTracks() << " not found"; } // Create a CbmL1MCTrack - ::ca::tools::MCTrack aTrk; + auto aTrk = MCTrack {}; aTrk.SetId(fpMCData->GetNofTracks()); // assign current number of tracks read so far as an ID of a new track aTrk.SetExternalId(iTrkExt); // external index of track is its index from CbmMCTrack objects container diff --git a/reco/L1/CbmCaMCModule.h b/reco/L1/CbmCaMCModule.h index 6cb47fec936f5de7a3d6eab1b291c047c1b09366..9e1090c2e0e7ac6e4c7b8ea7a182ed2153313845 100644 --- a/reco/L1/CbmCaMCModule.h +++ b/reco/L1/CbmCaMCModule.h @@ -42,6 +42,7 @@ #include "CaMonitor.h" #include "CaToolsMCData.h" #include "CaToolsMCPoint.h" +#include "CaVector.h" #include "L1InputData.h" #include "L1Parameters.h" @@ -50,10 +51,9 @@ class CbmMCDataObject; class CbmL1HitDebugInfo; class CbmL1Track; -enum class L1DetectorID; - namespace cbm::ca { + namespace ca = cbm::algo::ca; /// @brief Class CbmCaPerformance is an interface to communicate between /// @@ -89,7 +89,7 @@ namespace cbm::ca void Finish(); /// @brief Gets a pointer to MC data object - const ::ca::tools::MCData* GetMCData() const { return fpMCData; } + const tools::MCData* GetMCData() const { return fpMCData; } /// @brief Defines performance action in the beginning of each event or time slice /// @note This function should be called before hits initialization @@ -108,7 +108,7 @@ namespace cbm::ca /// This method finds a match for a given hit or matches for hits clusters (in case of STS), finds the best /// link in the match and returns the corresponding global ID of the MC points. template<L1DetectorID DetId> - std::tuple<int, vector<int>> MatchHitWithMc(int iHitExt); + std::tuple<int, std::vector<int>> MatchHitWithMc(int iHitExt); /// @brief Match reconstructed and MC data /// @@ -142,15 +142,15 @@ namespace cbm::ca /// @brief Registers MC data object /// @param mcData Instance of MC data - void RegisterMCData(::ca::tools::MCData& mcData) { fpMCData = &mcData; } + void RegisterMCData(tools::MCData& mcData) { fpMCData = &mcData; } /// @brief Registers reconstructed track container /// @param vRecoTrack Reference to reconstructed track container - void RegisterRecoTrackContainer(L1Vector<CbmL1Track>& vRecoTracks) { fpvRecoTracks = &vRecoTracks; } + void RegisterRecoTrackContainer(ca::Vector<CbmL1Track>& vRecoTracks) { fpvRecoTracks = &vRecoTracks; } /// @brief Registers hit index container /// @param vHitIds Reference to hit index container - void RegisterHitIndexContainer(L1Vector<CbmL1HitId>& vHitIds) { fpvHitIds = &vHitIds; } + void RegisterHitIndexContainer(ca::Vector<CbmL1HitId>& vHitIds) { fpvHitIds = &vHitIds; } /// @brief Registers CA parameters object /// @param pParameters A shared pointer to the parameters object @@ -158,7 +158,7 @@ namespace cbm::ca /// @brief Registers debug hit container /// @param vQaHits Reference to debug hit container - void RegisterQaHitContainer(L1Vector<CbmL1HitDebugInfo>& vQaHits) { fpvQaHits = &vQaHits; } + void RegisterQaHitContainer(ca::Vector<CbmL1HitDebugInfo>& vQaHits) { fpvQaHits = &vQaHits; } private: /// @brief Check class initialization @@ -208,7 +208,7 @@ namespace cbm::ca /// @return true Point is correct and is to be saved to the MCData object /// @return false Point is incorrect and will be ignored template<L1DetectorID DetID> - std::optional<::ca::tools::MCPoint> FillMCPoint(int iExtId, int iEvent, int iFile); + std::optional<tools::MCPoint> FillMCPoint(int iExtId, int iEvent, int iFile); /// @enum EMonitorKey /// @brief Monitor keys @@ -225,7 +225,7 @@ namespace cbm::ca kMissedMatchesTof, ///< Number of missed TOF matches kEND }; - ::ca::Monitor<EMonitorKey> fMonitor {"CA MC Module"}; ///< Monitor + ca::Monitor<EMonitorKey> fMonitor {"CA MC Module"}; ///< Monitor // ------ Flags DetIdArr_t<bool> fvbUseDet = {{false}}; ///< Flag: is detector subsystem used @@ -249,12 +249,12 @@ namespace cbm::ca std::set<std::pair<int, int>> fFileEventIDs; ///< Set of file and event indexes: first - iFile, second - iEvent // ----- Internal MC data - ::ca::tools::MCData* fpMCData = nullptr; ///< MC information (hits and tracks) instance + tools::MCData* fpMCData = nullptr; ///< MC information (hits and tracks) instance // ----- Internal reconstructed data - L1Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to reconstructed track container - L1Vector<CbmL1HitId>* fpvHitIds = nullptr; ///< Pointer to hit index container - L1Vector<CbmL1HitDebugInfo>* fpvQaHits = nullptr; ///< Pointer to QA hit container + ca::Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to reconstructed track container + ca::Vector<CbmL1HitId>* fpvHitIds = nullptr; ///< Pointer to hit index container + ca::Vector<CbmL1HitDebugInfo>* fpvQaHits = nullptr; ///< Pointer to QA hit container /// @brief Pointer to array of first hit indexes in the detector subsystem @@ -262,234 +262,233 @@ namespace cbm::ca /// This array must be initialized in the run initialization function. const std::array<int, constants::size::MaxNdetectors + 1>* fpvFstHitId = nullptr; }; -} // namespace cbm::ca -// ********************************************** -// ** Template function implementation ** -// ********************************************** + // ********************************************** + // ** Template function implementation ** + // ********************************************** + + // ------------------------------------------------------------------------------------------------------------------- + // + template<L1DetectorID DetID> + std::tuple<int, std::vector<int>> MCModule::MatchHitWithMc(int iHitExt) + { + int iPoint = -1; + std::vector<int> vPoints; + const auto* pHitMatch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[DetID]->At(iHitExt)); + if (!pHitMatch) { + LOG(warn) << "Hit match with index " << iHitExt << " is missing for " << kDetName[DetID]; + if constexpr (L1DetectorID::kMvd == DetID) { fMonitor.Increment(EMonitorKey::kMissedMatchesMvd); } + else if constexpr (L1DetectorID::kSts == DetID) { + fMonitor.Increment(EMonitorKey::kMissedMatchesSts); + } + else if constexpr (L1DetectorID::kMuch == DetID) { + fMonitor.Increment(EMonitorKey::kMissedMatchesMuch); + } + else if constexpr (L1DetectorID::kTrd == DetID) { + fMonitor.Increment(EMonitorKey::kMissedMatchesTrd); + } + else if constexpr (L1DetectorID::kTof == DetID) { + fMonitor.Increment(EMonitorKey::kMissedMatchesTof); + } + return std::tuple(iPoint, vPoints); + } -// --------------------------------------------------------------------------------------------------------------------- -// -template<L1DetectorID DetID> -std::tuple<int, vector<int>> cbm::ca::MCModule::MatchHitWithMc(int iHitExt) -{ - int iPoint = -1; - std::vector<int> vPoints; - const auto* pHitMatch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[DetID]->At(iHitExt)); - if (!pHitMatch) { - LOG(warn) << "Hit match with index " << iHitExt << " is missing for " << kDetName[DetID]; - if constexpr (L1DetectorID::kMvd == DetID) { fMonitor.Increment(EMonitorKey::kMissedMatchesMvd); } + for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); ++iLink) { + const auto& link = pHitMatch->GetLink(iLink); + int iPointExt = link.GetIndex(); + int iEvent = link.GetEntry(); + int iFile = link.GetFile(); + if (iPointExt < 0) continue; + int id = fpMCData->FindInternalPointIndex(DetID, iPointExt, iEvent, iFile); + vPoints.push_back(id); + } + + if (pHitMatch->GetNofLinks() > 0) { + const auto& link = pHitMatch->GetMatchedLink(); + if (link.GetIndex() > -1) { + int index = link.GetIndex(); + int event = link.GetEntry(); + int file = link.GetFile(); + iPoint = fpMCData->FindInternalPointIndex(DetID, index, event, file); + } + } + + return std::tuple(iPoint, vPoints); + } + + // ------------------------------------------------------------------------------------------------------------------- + // + template<L1DetectorID DetID> + std::optional<tools::MCPoint> MCModule::FillMCPoint(int iExtId, int iEvent, int iFile) + { + auto oPoint = std::make_optional<tools::MCPoint>(); + + // ----- Get positions, momenta, time and track ID + TVector3 posIn; // Position at entrance to station [cm] + TVector3 posOut; // Position at exist of station [cm] + TVector3 momIn; // 3-momentum at entrance to station [GeV/c] + TVector3 momOut; // 3-momentum at exit of station [GeV/c] + + auto* pBrPoints = fvpBrPoints[DetID]; + + PointTypes_t::at<DetID>* pExtPoint = dynamic_cast<PointTypes_t::at<DetID>*>(pBrPoints->Get(iFile, iEvent, iExtId)); + if (!pExtPoint) { + LOG(warn) << "CbmCaMCModule: " << kDetName[DetID] << " MC point with iExtId = " << iExtId + << ", iEvent = " << iEvent << ", iFile = " << iFile << " does not exist"; + return std::nullopt; + } + if constexpr (L1DetectorID::kMvd == DetID) { + pExtPoint->Position(posIn); + pExtPoint->PositionOut(posOut); + pExtPoint->Momentum(momIn); + pExtPoint->MomentumOut(momOut); + } + // STS else if constexpr (L1DetectorID::kSts == DetID) { - fMonitor.Increment(EMonitorKey::kMissedMatchesSts); + pExtPoint->Position(posIn); + pExtPoint->PositionOut(posOut); + pExtPoint->Momentum(momIn); + pExtPoint->MomentumOut(momOut); } + // MuCh else if constexpr (L1DetectorID::kMuch == DetID) { - fMonitor.Increment(EMonitorKey::kMissedMatchesMuch); + pExtPoint->Position(posIn); + pExtPoint->PositionOut(posOut); + pExtPoint->Momentum(momIn); + pExtPoint->Momentum(momOut); } + // TRD else if constexpr (L1DetectorID::kTrd == DetID) { - fMonitor.Increment(EMonitorKey::kMissedMatchesTrd); + pExtPoint->Position(posIn); + pExtPoint->PositionOut(posOut); + pExtPoint->Momentum(momIn); + pExtPoint->MomentumOut(momOut); } + // TOF else if constexpr (L1DetectorID::kTof == DetID) { - fMonitor.Increment(EMonitorKey::kMissedMatchesTof); + pExtPoint->Position(posIn); + pExtPoint->Position(posOut); + pExtPoint->Momentum(momIn); + pExtPoint->Momentum(momOut); } - return std::tuple(iPoint, vPoints); - } - - for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); ++iLink) { - const auto& link = pHitMatch->GetLink(iLink); - int iPointExt = link.GetIndex(); - int iEvent = link.GetEntry(); - int iFile = link.GetFile(); - if (iPointExt < 0) continue; - int id = fpMCData->FindInternalPointIndex(DetID, iPointExt, iEvent, iFile); - vPoints.push_back(id); - } + double time = pExtPoint->GetTime(); + int iTmcExt = pExtPoint->GetTrackID(); - if (pHitMatch->GetNofLinks() > 0) { - const auto& link = pHitMatch->GetMatchedLink(); - if (link.GetIndex() > -1) { - int index = link.GetIndex(); - int event = link.GetEntry(); - int file = link.GetFile(); - iPoint = fpMCData->FindInternalPointIndex(DetID, index, event, file); + if (iTmcExt < 0) { + LOG(warn) << "CbmCaMCModule: For MC point with iExtId = " << iExtId << ", iEvent = " << iEvent + << ", iFile = " << iFile << " MC track is undefined (has ID = " << iTmcExt << ')'; + return std::nullopt; + } + TVector3 posMid = 0.5 * (posIn + posOut); + TVector3 momMid = 0.5 * (momIn + momOut); + + // // ----- Get station index + int iStLoc = fvpDetInterface[DetID]->GetTrackingStationIndex(pExtPoint); + int stationID = fpParameters->GetStationIndexActive(iStLoc, DetID); + if (stationID == -1) { return std::nullopt; } // Skip points from inactive stations + + // Update point time with event time + time += fpMCEventList->GetEventTime(iEvent, iFile); + + // ----- Reject MC points falling out of the time slice + // STS, MuCh, TRD, TOF + if constexpr (DetID != L1DetectorID::kMvd) { + double startT = fpTimeSlice->GetStartTime(); + double endT = fpTimeSlice->GetEndTime(); + + if ((startT > 0. && time < startT) || (endT > 0. && time > endT)) { + LOG(warn) << "CbmCaMCModule: MC point with iExtId = " << iExtId << ", iEvent = " << iEvent + << ", iFile = " << iFile << " and det id " << int(DetID) << " fell out of the TS duration [" << startT + << ", " << endT << "] with measured time = " << time << " [ns]"; + return std::nullopt; + } } - } - - return std::tuple(iPoint, vPoints); -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<L1DetectorID DetID> -std::optional<::ca::tools::MCPoint> cbm::ca::MCModule::FillMCPoint(int iExtId, int iEvent, int iFile) -{ - auto oPoint = std::make_optional<::ca::tools::MCPoint>(); - - // ----- Get positions, momenta, time and track ID - TVector3 posIn; // Position at entrance to station [cm] - TVector3 posOut; // Position at exist of station [cm] - TVector3 momIn; // 3-momentum at entrance to station [GeV/c] - TVector3 momOut; // 3-momentum at exit of station [GeV/c] - - auto* pBrPoints = fvpBrPoints[DetID]; - - PointTypes_t::at<DetID>* pExtPoint = dynamic_cast<PointTypes_t::at<DetID>*>(pBrPoints->Get(iFile, iEvent, iExtId)); - if (!pExtPoint) { - LOG(warn) << "CbmCaMCModule: " << kDetName[DetID] << " MC point with iExtId = " << iExtId << ", iEvent = " << iEvent - << ", iFile = " << iFile << " does not exist"; - return std::nullopt; - } - if constexpr (L1DetectorID::kMvd == DetID) { - pExtPoint->Position(posIn); - pExtPoint->PositionOut(posOut); - pExtPoint->Momentum(momIn); - pExtPoint->MomentumOut(momOut); - } - // STS - else if constexpr (L1DetectorID::kSts == DetID) { - pExtPoint->Position(posIn); - pExtPoint->PositionOut(posOut); - pExtPoint->Momentum(momIn); - pExtPoint->MomentumOut(momOut); - } - // MuCh - else if constexpr (L1DetectorID::kMuch == DetID) { - pExtPoint->Position(posIn); - pExtPoint->PositionOut(posOut); - pExtPoint->Momentum(momIn); - pExtPoint->Momentum(momOut); - } - // TRD - else if constexpr (L1DetectorID::kTrd == DetID) { - pExtPoint->Position(posIn); - pExtPoint->PositionOut(posOut); - pExtPoint->Momentum(momIn); - pExtPoint->MomentumOut(momOut); - } - // TOF - else if constexpr (L1DetectorID::kTof == DetID) { - pExtPoint->Position(posIn); - pExtPoint->Position(posOut); - pExtPoint->Momentum(momIn); - pExtPoint->Momentum(momOut); - } - double time = pExtPoint->GetTime(); - int iTmcExt = pExtPoint->GetTrackID(); - if (iTmcExt < 0) { - LOG(warn) << "CbmCaMCModule: For MC point with iExtId = " << iExtId << ", iEvent = " << iEvent - << ", iFile = " << iFile << " MC track is undefined (has ID = " << iTmcExt << ')'; - return std::nullopt; - } - TVector3 posMid = 0.5 * (posIn + posOut); - TVector3 momMid = 0.5 * (momIn + momOut); - - // // ----- Get station index - int iStLoc = fvpDetInterface[DetID]->GetTrackingStationIndex(pExtPoint); - int stationID = fpParameters->GetStationIndexActive(iStLoc, DetID); - if (stationID == -1) { return std::nullopt; } // Skip points from inactive stations - - // Update point time with event time - time += fpMCEventList->GetEventTime(iEvent, iFile); - - // ----- Reject MC points falling out of the time slice - // STS, MuCh, TRD, TOF - if constexpr (DetID != L1DetectorID::kMvd) { - double startT = fpTimeSlice->GetStartTime(); - double endT = fpTimeSlice->GetEndTime(); - - if ((startT > 0. && time < startT) || (endT > 0. && time > endT)) { - LOG(warn) << "CbmCaMCModule: MC point with iExtId = " << iExtId << ", iEvent = " << iEvent - << ", iFile = " << iFile << " and det id " << int(DetID) << " fell out of the TS duration [" << startT - << ", " << endT << "] with measured time = " << time << " [ns]"; + // ----- Fill MC point + oPoint->SetExternalId(fpMCData->GetPointGlobExtIndex(DetID, iExtId)); + oPoint->SetEventId(iEvent); + oPoint->SetFileId(iFile); + oPoint->SetTime(time); + oPoint->SetX(posMid.X()); + oPoint->SetY(posMid.Y()); + oPoint->SetZ(posMid.Z()); + oPoint->SetXIn(posIn.X()); + oPoint->SetYIn(posIn.Y()); + oPoint->SetZIn(posIn.Z()); + oPoint->SetXOut(posOut.X()); + oPoint->SetYOut(posOut.Y()); + oPoint->SetZOut(posOut.Z()); + oPoint->SetPx(momMid.X()); + oPoint->SetPy(momMid.Y()); + oPoint->SetPz(momMid.Z()); + oPoint->SetPxIn(momIn.X()); + oPoint->SetPyIn(momIn.Y()); + oPoint->SetPzIn(momIn.Z()); + oPoint->SetPxOut(momOut.X()); + oPoint->SetPyOut(momOut.Y()); + oPoint->SetPzOut(momOut.Z()); + + // Select current number of points as a local id of points + oPoint->SetId(fpMCData->GetNofPoints()); + + // Match MC track and point to each other + int iTmcInt = fpMCData->FindInternalTrackIndex(iTmcExt, iEvent, iFile); + + oPoint->SetTrackId(iTmcInt); + if (iTmcInt > -1) { fpMCData->GetTrack(iTmcInt).AddPointIndex(oPoint->GetId()); } + + oPoint->SetStationId(stationID); + oPoint->SetDetectorId(DetID); + + auto* pExtTrk = L1_DYNAMIC_CAST<CbmMCTrack*>(fpMCTracks->Get(iFile, iEvent, iTmcExt)); + if (!pExtTrk) { + LOG(warn) << "CbmCaMCModule: MC track with iTmcExt = " << iTmcExt << ", iEvent = " << iEvent + << ", iFile = " << iFile << " MC track is undefined (nullptr)"; return std::nullopt; } - } - - // ----- Fill MC point - oPoint->SetExternalId(fpMCData->GetPointGlobExtIndex(DetID, iExtId)); - oPoint->SetEventId(iEvent); - oPoint->SetFileId(iFile); - oPoint->SetTime(time); - oPoint->SetX(posMid.X()); - oPoint->SetY(posMid.Y()); - oPoint->SetZ(posMid.Z()); - oPoint->SetXIn(posIn.X()); - oPoint->SetYIn(posIn.Y()); - oPoint->SetZIn(posIn.Z()); - oPoint->SetXOut(posOut.X()); - oPoint->SetYOut(posOut.Y()); - oPoint->SetZOut(posOut.Z()); - oPoint->SetPx(momMid.X()); - oPoint->SetPy(momMid.Y()); - oPoint->SetPz(momMid.Z()); - oPoint->SetPxIn(momIn.X()); - oPoint->SetPyIn(momIn.Y()); - oPoint->SetPzIn(momIn.Z()); - oPoint->SetPxOut(momOut.X()); - oPoint->SetPyOut(momOut.Y()); - oPoint->SetPzOut(momOut.Z()); - - // Select current number of points as a local id of points - oPoint->SetId(fpMCData->GetNofPoints()); - - // Match MC track and point to each other - int iTmcInt = fpMCData->FindInternalTrackIndex(iTmcExt, iEvent, iFile); - - oPoint->SetTrackId(iTmcInt); - if (iTmcInt > -1) { fpMCData->GetTrack(iTmcInt).AddPointIndex(oPoint->GetId()); } - - oPoint->SetStationId(stationID); - oPoint->SetDetectorId(DetID); - - auto* pExtTrk = L1_DYNAMIC_CAST<CbmMCTrack*>(fpMCTracks->Get(iFile, iEvent, iTmcExt)); - if (!pExtTrk) { - LOG(warn) << "CbmCaMCModule: MC track with iTmcExt = " << iTmcExt << ", iEvent = " << iEvent - << ", iFile = " << iFile << " MC track is undefined (nullptr)"; - return std::nullopt; - } - oPoint->SetPdgCode(pExtTrk->GetPdgCode()); - oPoint->SetMotherId(pExtTrk->GetMotherId()); + oPoint->SetPdgCode(pExtTrk->GetPdgCode()); + oPoint->SetMotherId(pExtTrk->GetMotherId()); - auto* pPdgDB = TDatabasePDG::Instance()->GetParticle(oPoint->GetPdgCode()); - oPoint->SetMass(pPdgDB ? pPdgDB->Mass() : 0.); /// TODO: Set from track - oPoint->SetCharge(pPdgDB ? pPdgDB->Charge() / 3. : 0.); + auto* pPdgDB = TDatabasePDG::Instance()->GetParticle(oPoint->GetPdgCode()); + oPoint->SetMass(pPdgDB ? pPdgDB->Mass() : 0.); /// TODO: Set from track + oPoint->SetCharge(pPdgDB ? pPdgDB->Charge() / 3. : 0.); - return oPoint; -} + return oPoint; + } -// --------------------------------------------------------------------------------------------------------------------- -// -template<L1DetectorID DetID> -void cbm::ca::MCModule::MatchPointsAndHits() -{ - if (!fvbUseDet[DetID]) { return; } - - for (int iH = (*fpvFstHitId)[static_cast<int>(DetID)]; iH < (*fpvFstHitId)[static_cast<int>(DetID) + 1]; ++iH) { - auto& hit = (*fpvQaHits)[iH]; - auto [iBestP, vAllP] = MatchHitWithMc<DetID>(hit.ExtIndex); - if (iBestP >= 0) { hit.SetBestMcPointId(iBestP); } - for (auto iP : vAllP) { - if (iP >= 0) { - hit.AddMcPointId(iP); - fpMCData->GetPoint(iP).AddHitID(iH); + // ------------------------------------------------------------------------------------------------------------------- + // + template<L1DetectorID DetID> + void MCModule::MatchPointsAndHits() + { + if (!fvbUseDet[DetID]) { return; } + + for (int iH = (*fpvFstHitId)[static_cast<int>(DetID)]; iH < (*fpvFstHitId)[static_cast<int>(DetID) + 1]; ++iH) { + auto& hit = (*fpvQaHits)[iH]; + auto [iBestP, vAllP] = MatchHitWithMc<DetID>(hit.ExtIndex); + if (iBestP >= 0) { hit.SetBestMcPointId(iBestP); } + for (auto iP : vAllP) { + if (iP >= 0) { + hit.AddMcPointId(iP); + fpMCData->GetPoint(iP).AddHitID(iH); + } } } } -} -// --------------------------------------------------------------------------------------------------------------------- -// NOTE: template is used, because another template function FillMCPoint is used inside -template<L1DetectorID DetID> -void cbm::ca::MCModule::ReadMCPointsForDetector() -{ - for (const auto& [iFile, iEvent] : fFileEventIDs) { - int nPointsEvent = fvpBrPoints[DetID]->Size(iFile, iEvent); - for (int iP = 0; iP < nPointsEvent; ++iP) { - std::optional<::ca::tools::MCPoint> oPoint = FillMCPoint<DetID>(iP, iEvent, iFile); - if (oPoint) { fpMCData->AddPoint(*oPoint); } - } // iP: end - } // key: end -} - + // ------------------------------------------------------------------------------------------------------------------- + // NOTE: template is used, because another template function FillMCPoint is used inside + template<L1DetectorID DetID> + void MCModule::ReadMCPointsForDetector() + { + for (const auto& [iFile, iEvent] : fFileEventIDs) { + int nPointsEvent = fvpBrPoints[DetID]->Size(iFile, iEvent); + for (int iP = 0; iP < nPointsEvent; ++iP) { + std::optional<tools::MCPoint> oPoint = FillMCPoint<DetID>(iP, iEvent, iFile); + if (oPoint) { fpMCData->AddPoint(*oPoint); } + } // iP: end + } // key: end + } +} // namespace cbm::ca #endif // CbmCaMCModule_h diff --git a/reco/L1/CbmCaTimeSliceReader.cxx b/reco/L1/CbmCaTimeSliceReader.cxx index 19b89ec680578436ace9826ae2396335940e5bdd..e63cc0b62cff134f3e018c97b7c35789a213c09d 100644 --- a/reco/L1/CbmCaTimeSliceReader.cxx +++ b/reco/L1/CbmCaTimeSliceReader.cxx @@ -26,7 +26,7 @@ #include "L1InputData.h" #include "L1Parameters.h" -using ca::tools::HitRecord; +using cbm::ca::tools::HitRecord; using namespace cbm::algo::ca::constants; @@ -308,7 +308,7 @@ void TimeSliceReader::RegisterIODataManager(std::shared_ptr<L1IODataManager>& pI void TimeSliceReader::SortQaHits() { int nStationsActive = fpParameters->GetNstationsActive(); - L1Vector<CbmL1HitDebugInfo> vNewHits(fpvQaHits->size()); + ca::Vector<CbmL1HitDebugInfo> vNewHits(fpvQaHits->size()); std::vector<int> vHitFstIndexes(nStationsActive + 1, 0); std::vector<int> vNofHitsStored(nStationsActive, 0); diff --git a/reco/L1/CbmCaTimeSliceReader.h b/reco/L1/CbmCaTimeSliceReader.h index 70d02fa5f69b52abe12444bd8fb7f009a1db3d1a..01e260936158e26653fe1887a6f079de74df6780 100644 --- a/reco/L1/CbmCaTimeSliceReader.h +++ b/reco/L1/CbmCaTimeSliceReader.h @@ -31,7 +31,7 @@ #include "CaConstants.h" #include "CaToolsHitRecord.h" -#include "L1Vector.h" +#include "CaVector.h" class CbmTimeSlice; @@ -40,8 +40,7 @@ class L1Parameters; namespace cbm::ca { - using ::ca::tools::HitRecord; - + namespace ca = cbm::algo::ca; /// @brief A reader of time slice for CA tracker /// @@ -101,12 +100,12 @@ namespace cbm::ca /// @brief Registers hit debug info container /// @param vQaHits Reference to Qa hit container /// @note If no container is registered, all related routines are omitted - void RegisterQaHitContainer(L1Vector<CbmL1HitDebugInfo>& vQaHits) { fpvQaHits = &vQaHits; } + void RegisterQaHitContainer(ca::Vector<CbmL1HitDebugInfo>& vQaHits) { fpvQaHits = &vQaHits; } /// @brief Registers hit index container /// @param vHitIds Reference to hits indexes container /// @note If no container is registered, all related routines are omitted - void RegisterHitIndexContainer(L1Vector<CbmL1HitId>& vHitIds) { fpvHitIds = &vHitIds; } + void RegisterHitIndexContainer(ca::Vector<CbmL1HitId>& vHitIds) { fpvHitIds = &vHitIds; } /// @brief Registers CA parameters object /// @param pParameters A shared pointer to the parameters object @@ -120,7 +119,7 @@ namespace cbm::ca /// @brief Register the reconstructed tracks container /// @param vTracks Reference to reconstructed tracks container /// @note If no container is registered, all related routines are omitted - void RegisterTracksContainer(L1Vector<CbmL1Track>& vTracks) { fpvTracks = &vTracks; } + void RegisterTracksContainer(ca::Vector<CbmL1Track>& vTracks) { fpvTracks = &vTracks; } /// @brief Sets used detector subsystems /// @param detID Id of detector @@ -158,7 +157,7 @@ namespace cbm::ca /// @param hitRecord Filled hit record /// /// Stores recorded hit information into registered hit containers - void StoreHitRecord(const HitRecord& hitRecord); + void StoreHitRecord(const tools::HitRecord& hitRecord); bool fbReadTracks = true; ///< flag to read reconstructed tracks from reco.root @@ -180,14 +179,14 @@ namespace cbm::ca TClonesArray* fpBrTofTracks = nullptr; ///< Input branch for reconstructed TOF tracks ("TofTrack") // Pointers to output data containers - L1Vector<CbmL1HitId>* fpvHitIds = nullptr; ///< Pointer to array of hit index objects - L1Vector<CbmL1HitDebugInfo>* fpvQaHits = nullptr; ///< Pointer to array of debug hits - L1Vector<CbmL1Track>* fpvTracks = nullptr; ///< Pointer to array of reconstructed tracks + ca::Vector<CbmL1HitId>* fpvHitIds = nullptr; ///< Pointer to array of hit index objects + ca::Vector<CbmL1HitDebugInfo>* fpvQaHits = nullptr; ///< Pointer to array of debug hits + ca::Vector<CbmL1Track>* fpvTracks = nullptr; ///< Pointer to array of reconstructed tracks std::shared_ptr<L1IODataManager> fpIODataManager = nullptr; ///< Pointer to input data manager std::shared_ptr<L1Parameters> fpParameters = nullptr; ///< Pointer to tracking parameters object // Maps of hit indexes: ext -> int - DetIdArr_t<L1Vector<int>> fvvHitExtToIntIndexMap; ///< Hit index map ext -> int + DetIdArr_t<ca::Vector<int>> fvvHitExtToIntIndexMap; ///< Hit index map ext -> int DetIdArr_t<int> fvNofHitsTotal = {{0}}; ///< Total hit number in detector DetIdArr_t<int> fvNofHitsUsed = {{0}}; ///< Number of used hits in detector @@ -223,7 +222,7 @@ int cbm::ca::TimeSliceReader::ReadHitsForDetector() fFirstHitKey = fNofHitKeys; for (int iHext = 0; iHext < nHitsTot; ++iHext) { - HitRecord hitRecord; + tools::HitRecord hitRecord; CbmPixelHit* pPixelHit = nullptr; // Pointer to hit object int iStGeom = -1; // Geometry station number diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx index 7ac966a811f627603b64098e0454b78a533f4c7d..9a0b09c8a2b56e67f9c9b84a00a7d0b17b2215b3 100644 --- a/reco/L1/CbmL1.cxx +++ b/reco/L1/CbmL1.cxx @@ -1,4 +1,4 @@ -/* Copyright (C) 2006-2021 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt +/* Copyright (C) 2006-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt SPDX-License-Identifier: GPL-3.0-only Authors: Ivan Kisel, Sergey Gorbunov [committer], Valentina Akishina, Igor Kulakov, Maksym Zyzak, Sergei Zharko */ @@ -75,6 +75,7 @@ using std::cout; using std::endl; using std::ios; using CaTrack = cbm::algo::ca::Track; +using cbm::ca::tools::MaterialHelper; ClassImp(CbmL1); @@ -768,9 +769,7 @@ void CbmL1::Reconstruct(CbmEvent* event) int trackFirstHit = 0; - for (L1Vector<CaTrack>::iterator it = fpAlgo->fRecoTracks.begin(); it != fpAlgo->fRecoTracks.end(); - trackFirstHit += it->fNofHits, it++) { - + for (auto it = fpAlgo->fRecoTracks.begin(); it != fpAlgo->fRecoTracks.end(); trackFirstHit += it->fNofHits, it++) { CbmL1Track t; for (int i = 0; i < L1TrackPar::kNparTr; i++) { @@ -823,7 +822,6 @@ void CbmL1::Reconstruct(CbmEvent* event) if constexpr (0) { using std::setw; { - LOG(info) << "---------- Reco track sample"; int id = 0; if (fvRecoTracks.size()) { LOG(info) << setw(4) << "No." << ' ' << fvRecoTracks[0].ToString(3, true); } for (const auto& trk : fvRecoTracks) { @@ -964,7 +962,7 @@ void CbmL1::GenerateMaterialMaps() LOG(info) << "Generating material maps..."; auto timerStart = std::chrono::high_resolution_clock::now(); - ca::tools::MaterialHelper matHelper; + cbm::ca::tools::MaterialHelper matHelper; matHelper.SetSafeMaterialInitialization(fDoSafeMaterialInitialization); if (!fMatBudgetParallelProjection) { matHelper.SetDoRadialProjection(fTargetZ); } @@ -1013,9 +1011,7 @@ void CbmL1::IdealTrackFinder() fpAlgo->fRecoTracks.clear(); fpAlgo->fRecoHits.clear(); - for (L1Vector<CbmL1MCTrack>::iterator i = fvMCTracks.begin(); i != fvMCTracks.end(); ++i) { - CbmL1MCTrack& MC = *i; - + for (auto& MC : fvMCTracks) { if (!MC.IsReconstructable()) continue; if (!(MC.ID >= 0)) continue; @@ -1024,7 +1020,7 @@ void CbmL1::IdealTrackFinder() CaTrack algoTr; algoTr.fNofHits = 0; - L1Vector<int> hitIndices("CbmL1::hitIndices", fpAlgo->GetParameters()->GetNstationsActive(), -1); + ca::Vector<int> hitIndices("CbmL1::hitIndices", fpAlgo->GetParameters()->GetNstationsActive(), -1); for (unsigned int iH = 0; iH < MC.Hits.size(); iH++) { const int hitI = MC.Hits[iH]; @@ -1184,7 +1180,7 @@ void CbmL1::WriteSIMDKFData() McTracksOut.open("mctracksout.dat", std::fstream::out | std::fstream::app); } - for (L1Vector<CbmL1Track>::iterator RecTrack = fvRecoTracks.begin(); RecTrack != fvRecoTracks.end(); ++RecTrack) { + for (auto RecTrack = fvRecoTracks.begin(); RecTrack != fvRecoTracks.end(); ++RecTrack) { if (RecTrack->IsGhost()) continue; CbmL1MCTrack* MCTrack = RecTrack->GetMCTrack(); diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h index c981b36076f98d068353f1414f98dcb4faffe171..7f5f70e38f29f511faa677010b3b7ef204e7f970 100644 --- a/reco/L1/CbmL1.h +++ b/reco/L1/CbmL1.h @@ -54,8 +54,8 @@ #include <utility> #include "CaMonitor.h" +#include "CaVector.h" #include "L1Algo/L1Algo.h" -#include "L1Algo/L1Vector.h" #include "L1EventEfficiencies.h" #include "L1IODataManager.h" #include "L1InitManager.h" @@ -73,6 +73,11 @@ class TProfile2D; class TNtuple; class TGeoNode; +namespace +{ + namespace ca = cbm::algo::ca; +} + /// Internal structure to handle link keys struct CbmL1LinkKey { /// Constructor from links @@ -289,12 +294,12 @@ public: // static bool compareZ(const int &a, const int &b ); inline double Get_Z_vMCPoint(int a) const { return fvMCPoints[a].z; } - const L1Vector<CbmL1MCPoint>& GetMcPoints() const { return fvMCPoints; } + const ca::Vector<CbmL1MCPoint>& GetMcPoints() const { return fvMCPoints; } - const L1Vector<CbmL1MCTrack>& GetMcTracks() const { return fvMCTracks; } + const ca::Vector<CbmL1MCTrack>& GetMcTracks() const { return fvMCTracks; } - const L1Vector<vector<int>>& GetHitMcRefs() const { return fvHitPointIndices; } - const L1Vector<int>& GetHitBestMcRefs() const { return fvHitBestPointIndices; } + const ca::Vector<std::vector<int>>& GetHitMcRefs() const { return fvHitPointIndices; } + const ca::Vector<int>& GetHitBestMcRefs() const { return fvHitBestPointIndices; } static double boundedGaus(double sigma); @@ -495,7 +500,7 @@ public: DFSET fvSelectedMcEvents {}; ///< Set of selected MC events with fileID and eventId - L1Vector<CbmL1Track> fvRecoTracks = {"CbmL1::fvRecoTracks"}; ///< Reconstructed tracks container + ca::Vector<CbmL1Track> fvRecoTracks = {"CbmL1::fvRecoTracks"}; ///< Reconstructed tracks container private: static CbmL1* fpInstance; ///< Instance of CbmL1 @@ -519,7 +524,7 @@ private: int fNpointsTrdAll = 0; ///< Number of MC points for TRD int fNpointsTofAll = 0; ///< Number of MC points for TOF - L1Vector<CbmL1MCPoint> fvMCPoints = {"CbmL1::fvMCPoints"}; ///< Container of MC points + ca::Vector<CbmL1MCPoint> fvMCPoints = {"CbmL1::fvMCPoints"}; ///< Container of MC points int fNStations = 0; ///< number of total active detector stations int fNMvdStations = 0; ///< number of active MVD stations @@ -628,9 +633,11 @@ private: TClonesArray* fpMuchDigiMatches = nullptr; ///< Array of MuCh digi matches (NOTE: currently unsused) - L1Vector<CbmL1MCTrack> fvMCTracks = {"CbmL1::fvMCTracks"}; ///< Array of MC tracks - L1Vector<vector<int>> fvHitPointIndices = {"CbmL1::fvHitPointIndices"}; ///< Indices of MC points vs. hit index - L1Vector<int> fvHitBestPointIndices = {"CbmL1::fvHitBestPointIndices"}; ///< Indices of best MC points vs. hit index + ca::Vector<CbmL1MCTrack> fvMCTracks = {"CbmL1::fvMCTracks"}; ///< Array of MC tracks + ca::Vector<std::vector<int>> fvHitPointIndices = { + "CbmL1::fvHitPointIndices"}; ///< Indices of MC points vs. hit index + ca::Vector<int> fvHitBestPointIndices = { + "CbmL1::fvHitBestPointIndices"}; ///< Indices of best MC points vs. hit index int fEventNo = 0; ///< Current number of event/TS int fNofRecoTracks = 0; ///< Total number of reconstructed tracks @@ -638,13 +645,13 @@ private: public: // ** Repacked input data ** - L1Vector<CbmL1HitId> fvExternalHits = {"CbmL1::fvExternalHits"}; ///< Array of hits + ca::Vector<CbmL1HitId> fvExternalHits = {"CbmL1::fvExternalHits"}; ///< Array of hits private: - L1Vector<CbmL1HitDebugInfo> fvHitDebugInfo = { + ca::Vector<CbmL1HitDebugInfo> fvHitDebugInfo = { "CbmL1::fvHitDebugInfo"}; ///< Container of hits with extended information // indices of MCPoints in fvMCPoints, indexed by index of hit in algo->vHits array. According to StsMatch. Used for IdealResponce - // L1Vector<int> vHitMCRef1; + // ca::Vector<int> vHitMCRef1; // CbmMatch HitMatch; private: std::unordered_map<CbmL1LinkKey, int> fmMCPointsLinksMap = {}; /// Internal MC point index vs. link diff --git a/reco/L1/CbmL1Counters.h b/reco/L1/CbmL1Counters.h index 503a01558903073db00c905aae79e4aa26369150..5d12ef004db9acf6353f96fe5ed9e9a54ec1968d 100644 --- a/reco/L1/CbmL1Counters.h +++ b/reco/L1/CbmL1Counters.h @@ -12,7 +12,12 @@ #include <iostream> #include <map> -#include "L1Vector.h" +#include "CaVector.h" + +namespace +{ + namespace cacore = cbm::algo::ca; +} /// counters used for efficiency calculation template<typename T> @@ -113,7 +118,7 @@ private: double Div(double a, double b) { return (b > 0) ? a / b : -1.; }; public: - L1Vector<T> counters {"TL1TracksCatCounters::counters"}; + cacore::Vector<T> counters {"TL1TracksCatCounters::counters"}; }; struct TL1Efficiencies { diff --git a/reco/L1/CbmL1DetectorID.h b/reco/L1/CbmL1DetectorID.h index d8678178a5ebe076d8c4714c941bcd48c853187f..75979dbbcbc24559c35b227cf98e15aa4894774c 100644 --- a/reco/L1/CbmL1DetectorID.h +++ b/reco/L1/CbmL1DetectorID.h @@ -2,7 +2,7 @@ SPDX-License-Identifier: GPL-3.0-only Authors: Sergei Zharko [committer] */ -/// @file CbmL1DetectorID.h +/// @file CbmL1DetectorID.h TODO: Rename to CbmCaDefs.h ! /// @brief Implementation of L1DetectorID enum class for CBM /// @author S.Zharko /// @since 01.12.2022 @@ -12,8 +12,17 @@ #include <string> +#include "CaConstants.h" #include "L1EArray.h" +namespace cbm::ca +{ + namespace phys = cbm::algo::ca::constants::phys; + namespace undef = cbm::algo::ca::constants::undef; + namespace clrs = cbm::algo::ca::constants::clrs; + namespace ca = cbm::algo::ca; ///< CA core namespace (TO BE USED IN CBM-specific interfaces) +} // namespace cbm::ca + /// Enumeration for the detector subsystems used in L1 tracking /// It is important for the subsystems to be specified in the actual order. The order is used /// for the L1Station array filling. diff --git a/reco/L1/CbmL1Hit.h b/reco/L1/CbmL1Hit.h index 8441879aa2ffe90fd9e2777e9972f7e7836506fb..5b3db374ae5b3fb381f8e3ae0eb57be3a0bcf6db 100644 --- a/reco/L1/CbmL1Hit.h +++ b/reco/L1/CbmL1Hit.h @@ -8,8 +8,9 @@ #include <iomanip> #include <sstream> #include <string> +#include <vector> -#include "L1Vector.h" +#include <cmath> // TODO: SZh: Complete the rule of five // TODO: SZh: Make variables private diff --git a/reco/L1/CbmL1MCPoint.h b/reco/L1/CbmL1MCPoint.h index 5055de0a1c0148519d9c82f4aa9bb8017eea6562..555c67de3edced9d8ce3e75f3b6f4c378d877089 100644 --- a/reco/L1/CbmL1MCPoint.h +++ b/reco/L1/CbmL1MCPoint.h @@ -25,7 +25,12 @@ #include <sstream> #include <string> -#include "L1Vector.h" +#include "CaVector.h" + +namespace +{ + using cbm::algo::ca::Vector; +} struct CbmL1MCPoint { @@ -70,7 +75,7 @@ struct CbmL1MCPoint { int pointId = -1; int file = -1; int event = -1; - L1Vector<int> hitIds {"CbmL1MCPoint::hitIds"}; // indices of CbmL1Hits in L1->vStsHits array + Vector<int> hitIds {"CbmL1MCPoint::hitIds"}; // indices of CbmL1Hits in L1->vStsHits array /// Temporary log function for debugging std::string ToString(int verbose = 3, bool printHeader = false) const diff --git a/reco/L1/CbmL1MCTrack.h b/reco/L1/CbmL1MCTrack.h index b7393a6eed15790875e51b577eef8988c7d2176c..991086c04cbe7adcfb715d034a79cda826d32c45 100644 --- a/reco/L1/CbmL1MCTrack.h +++ b/reco/L1/CbmL1MCTrack.h @@ -29,7 +29,12 @@ #include <iostream> #include <string> -#include "L1Vector.h" +#include "CaVector.h" + +namespace +{ + namespace cacore = cbm::algo::ca; +} class CbmL1Track; @@ -55,7 +60,7 @@ public: void AddRecoTrack(CbmL1Track* rTr) { rTracks.push_back_no_warning(rTr); } void AddRecoTrackIndex(int iT) { rTrackIndexes.push_back_no_warning(iT); } - L1Vector<CbmL1Track*>& GetRecoTracks() { return rTracks; } + cacore::Vector<CbmL1Track*>& GetRecoTracks() { return rTracks; } int GetNClones() const { return rTracks.size() - 1; } bool IsReconstructed() const { return rTracks.size(); } @@ -99,8 +104,8 @@ public: int pdg = -1; unsigned int process_ID = (unsigned int) -1; bool isSignal {0}; - L1Vector<int> Points {"CbmL1MCTrack::Points"}; // indices of pints in CbmL1::fvMCPoints - L1Vector<int> Hits {"CbmL1MCTrack::Hits"}; // indices of hits in algo->vHits or L1::vHits + cacore::Vector<int> Points {"CbmL1MCTrack::Points"}; // indices of pints in CbmL1::fvMCPoints + cacore::Vector<int> Hits {"CbmL1MCTrack::Hits"}; // indices of hits in algo->vHits or L1::vHits private: int nMCContStations = 0; // number of consecutive stations with mcPoints @@ -116,14 +121,14 @@ private: bool isAdditional = false; // is not reconstructable, but stil interesting // next members filled and used in Performance - L1Vector<CbmL1Track*> rTracks {"CbmL1MCTrack::rTracks"}; // array of associated recoTracks - L1Vector<CbmL1Track*> tTracks {"CbmL1MCTrack::tTracks"}; // array of recoTracks - // which aren't associated with this mcTrack, - // but use some hits from it. + cacore::Vector<CbmL1Track*> rTracks {"CbmL1MCTrack::rTracks"}; // array of associated recoTracks + cacore::Vector<CbmL1Track*> tTracks {"CbmL1MCTrack::tTracks"}; // array of recoTracks + // which aren't associated with this mcTrack, + // but use some hits from it. // NOTE: SZh 14.12.2022: on the replacement from rTracks and tTracks - L1Vector<int> rTrackIndexes = {"CbmL1MCTrack::rTrackIndexes"}; // array of associated recoTrack indexes - L1Vector<int> tTrackIndexes = {"CbmL1MCTrack::tTrackIndexes"}; // ..... + cacore::Vector<int> rTrackIndexes = {"CbmL1MCTrack::rTrackIndexes"}; // array of associated recoTrack indexes + cacore::Vector<int> tTrackIndexes = {"CbmL1MCTrack::tTrackIndexes"}; // ..... }; diff --git a/reco/L1/CbmL1Performance.cxx b/reco/L1/CbmL1Performance.cxx index 4189bf933e16eaa7a87a0cf7d10854b67f297f4f..7248aec1eb7d0d0f7a335a62e41de8425fb5fb6c 100644 --- a/reco/L1/CbmL1Performance.cxx +++ b/reco/L1/CbmL1Performance.cxx @@ -58,6 +58,7 @@ #include "L1Algo/L1Def.h" #include "L1Algo/L1Fit.h" +using cbm::ca::tools::Debugger; using std::cout; using std::endl; using std::ios; @@ -313,7 +314,7 @@ void CbmL1::EfficienciesPerformance() ntra.fOutDir = fTableDir; // Setup a pointer for output directory L1_NTRA.fOutDir = fTableDir; // Save average efficiencies - ca::tools::Debugger::Instance().AddNtuple("ghost", "it:ih:p:x:y:z:t:dx:dy"); + Debugger::Instance().AddNtuple("ghost", "it:ih:p:x:y:z:t:dx:dy"); static int statNghost = 0; for (vector<CbmL1Track>::iterator rtraIt = fvRecoTracks.begin(); rtraIt != fvRecoTracks.end(); ++rtraIt) { @@ -337,8 +338,7 @@ void CbmL1::EfficienciesPerformance() const L1Hit& h = fpAlgo->fInputData.GetHit(rtraIt->Hits[i]); const CbmL1HitDebugInfo& s = fvHitDebugInfo[rtraIt->Hits[i]]; cout << " x y z t " << s.x << " " << s.y << " " << h.z << " dx " << s.dx << " dy " << s.dy << std::endl; - ca::tools::Debugger::Instance().FillNtuple("ghost", statNghost, i, fabs(1. / tr.qp[0]), s.x, s.y, h.z, h.t, - s.dx, s.dy); + Debugger::Instance().FillNtuple("ghost", statNghost, i, fabs(1. / tr.qp[0]), s.x, s.y, h.z, h.t, s.dx, s.dy); } tr.Print(0); statNghost++; @@ -353,7 +353,7 @@ void CbmL1::EfficienciesPerformance() sta_nfakes[i] = 0; } - for (vector<CbmL1MCTrack>::iterator mtraIt = fvMCTracks.begin(); mtraIt != fvMCTracks.end(); mtraIt++) { + for (auto mtraIt = fvMCTracks.begin(); mtraIt != fvMCTracks.end(); mtraIt++) { CbmL1MCTrack& mtra = *(mtraIt); // if( !( mtra.pdg == -11 && mtra.mother_ID == -1 ) ) continue; // electrons only @@ -365,7 +365,7 @@ void CbmL1::EfficienciesPerformance() // is track killed. At least one hit of it belong to some recoTrack const bool killed = !reco && mtra.IsDisturbed(); // ration length for current mcTrack - L1Vector<CbmL1Track*>& rTracks = mtra.GetRecoTracks(); // for length calculations + auto& rTracks = mtra.GetRecoTracks(); // for length calculations double ratio_length = 0; double ratio_fakes = 0; double mc_length_hits = mtra.NStations(); diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx index b0f83de98b441c6d699db62f6bf2bee386951b8f..89b62a5721bd1025f46ed07512fdc35a351b4044 100644 --- a/reco/L1/CbmL1ReadEvent.cxx +++ b/reco/L1/CbmL1ReadEvent.cxx @@ -54,7 +54,8 @@ using std::cout; using std::endl; using std::find; - +using namespace cbm::algo; +using cbm::ca::tools::Debugger; static bool compareMcPointZ(const int& a, const int& b) { @@ -154,7 +155,7 @@ void CbmL1::ReadEvent(CbmEvent* event) if (fVerbose >= 10) cout << "ReadEvent: clear is done." << endl; - L1Vector<TmpHit> tmpHits("CbmL1ReadEvent::tmpHits"); + ca::Vector<TmpHit> tmpHits("CbmL1ReadEvent::tmpHits"); { // reserve enough space for hits int nHitsTotal = 0; @@ -556,7 +557,7 @@ void CbmL1::ReadEvent(CbmEvent* event) assert(fpTrdHits); - vector<bool> mcUsed(fNpointsTrd, 0); + std::vector<char> mcUsed(fNpointsTrd, 0); Int_t nEntTrd = (event ? event->GetNofData(ECbmDataType::kTrdHit) : fpTrdHits->GetEntriesFast()); @@ -713,7 +714,7 @@ void CbmL1::ReadEvent(CbmEvent* event) fpIODataManager->SetNhitKeys(NStrips); if (0) { - ca::tools::Debugger::Instance().AddNtuple("hits", "ev:mcpoint:mc:mcMother:nmc:sta:x:y:z"); + Debugger::Instance().AddNtuple("hits", "ev:mcpoint:mc:mcMother:nmc:sta:x:y:z"); for (int iHit = 0; iHit < nHits; ++iHit) { TmpHit& h = tmpHits[iHit]; int mc = -1; @@ -722,8 +723,7 @@ void CbmL1::ReadEvent(CbmEvent* event) mc = fvMCPoints[h.iBestMc].ID; mcMother = fvMCPoints[h.iBestMc].mother_ID; } - ca::tools::Debugger::Instance().FillNtuple("hits", nCalls, h.iBestMc, mc, mcMother, h.vMc.size(), h.iStation, h.x, - h.y, h.z); + Debugger::Instance().FillNtuple("hits", nCalls, h.iBestMc, mc, mcMother, h.vMc.size(), h.iStation, h.x, h.y, h.z); } } diff --git a/reco/L1/CbmL1Track.h b/reco/L1/CbmL1Track.h index a605d60fd0adda18029ec4215340a955febfa70d..509ea0a572a0e2c18b6e3a47c8ef302d94c9f1ee 100644 --- a/reco/L1/CbmL1Track.h +++ b/reco/L1/CbmL1Track.h @@ -32,8 +32,13 @@ #include <map> #include <string> #include <vector> -using std::map; -using std::vector; + +#include "CaVector.h" + +namespace +{ + namespace ca = cbm::algo::ca; +} class CbmL1MCTrack; @@ -91,7 +96,7 @@ public: /// Gets container of pointers to MC tracks // TODO: this function is to be replaced with GetMCTrackIndexes() - vector<CbmL1MCTrack*>& GetMCTracks() { return mcTracks; } + std::vector<CbmL1MCTrack*>& GetMCTracks() { return mcTracks; } /// Gets pointer of the associated MC track CbmL1MCTrack* GetMCTrack() { return mcTracks[0]; } @@ -162,21 +167,21 @@ public: std::array<double, L1TrackPar::kNparTr> TLast; ///< Track parameters in the end of the track std::array<double, L1TrackPar::kNparCov> CLast; ///< Track covariance matrix at the end of the track - vector<int> Hits; ///< Indexes of hits of this track + std::vector<int> Hits; ///< Indexes of hits of this track int nStations; ///< Number of stations with hits of this track int index; ///< Index of this track (TODO: it seems to be not initialized) double fTrackTime; ///< Time of the track [ns] ??? - map<int, int> hitMap; // N hits (second) from each mcTrack (first is a MC track ID) belong to current recoTrack + std::map<int, int> hitMap; // N hits (second) from each mcTrack (first is a MC track ID) belong to current recoTrack // FIXME: SZh 14.12.2022: map => unordered_map private: // next members filled and used in Performance - vector<CbmL1MCTrack*> mcTracks; // array of assosiated recoTracks. Should be only one. + std::vector<CbmL1MCTrack*> mcTracks; // array of assosiated recoTracks. Should be only one. - L1Vector<int> fvMcTrackIndexes = {"CbmL1Track::fvMcTrackIndexes"}; // global indexes of MC tracks + ca::Vector<int> fvMcTrackIndexes = {"CbmL1Track::fvMcTrackIndexes"}; // global indexes of MC tracks // NOTE: mcTracks should be replaced with fvMcTrackIndexes double maxPurity; ///< Maximum persent of hits, which belong to one mcTrack. diff --git a/reco/L1/L1Algo/L1Algo.cxx b/reco/L1/L1Algo/L1Algo.cxx index 45120bef1427e7c0740b53bc12d8cfda8c94e766..719d258018854a53cfd93fad49c48eb111266ab1 100644 --- a/reco/L1/L1Algo/L1Algo.cxx +++ b/reco/L1/L1Algo/L1Algo.cxx @@ -18,10 +18,11 @@ L1Algo::L1Algo() } } +using cbm::ca::tools::Debugger; void L1Algo::Init(const TrackingMode mode) { fTrackingMode = mode; } -void L1Algo::Finish() { ca::tools::Debugger::Instance().Write(); } +void L1Algo::Finish() { Debugger::Instance().Write(); } // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/L1Algo/L1Algo.h b/reco/L1/L1Algo/L1Algo.h index c223882d60bdc53b7d3d571b818f4e64dc0b786b..252578e0f329cccd6f56f57632838754d46c48a2 100644 --- a/reco/L1/L1Algo/L1Algo.h +++ b/reco/L1/L1Algo/L1Algo.h @@ -30,6 +30,7 @@ class L1AlgoDraw; #include "CaConstants.h" #include "CaTrack.h" +#include "CaVector.h" #include "L1Branch.h" #include "L1CloneMerger.h" #include "L1Field.h" @@ -42,11 +43,16 @@ class L1AlgoDraw; #include "L1Station.h" #include "L1TrackPar.h" #include "L1Triplet.h" -#include "L1Utils.h" -#include "L1Vector.h" +#include "L1Utils.h" // ? DEPRECATED ? class CbmL1MCTrack; +namespace +{ + using cbm::algo::ca::Track; // TMP + using cbm::algo::ca::Vector; // TMP +} // namespace + // ******************************* // ** Types definition (global) ** // ******************************* @@ -215,7 +221,7 @@ public: #ifdef DRAW L1AlgoDraw* draw {nullptr}; - void DrawRecoTracksTime(const L1Vector<CbmL1Track>& tracks); + void DrawRecoTracksTime(const Vector<CbmL1Track>& tracks); #endif /// TODO: Move to L1 @@ -277,11 +283,11 @@ private: L1Parameters fParameters; ///< Object of L1Algo parameters class L1InputData fInputData; ///< Tracking input data - L1Vector<unsigned char> fvHitKeyFlags { + Vector<unsigned char> fvHitKeyFlags { "L1Algo::fvHitKeyFlags"}; ///< List of key flags: has been this hit or cluster already used public: - L1Vector<L1HitTimeInfo> fHitTimeInfo; + Vector<L1HitTimeInfo> fHitTimeInfo; L1Grid vGrid[constants::size::MaxNstations]; ///< L1Grid vGridTime[constants::size::MaxNstations]; ///< @@ -292,36 +298,35 @@ public: double fCaRecoTime {0.}; // time of the track finder + fitter - L1Vector<cbm::algo::ca::Track> fRecoTracks {"L1Algo::fRecoTracks"}; ///< reconstructed tracks - L1Vector<L1HitIndex_t> fRecoHits {"L1Algo::fRecoHits"}; ///< packed hits of reconstructed tracks + Vector<Track> fRecoTracks {"L1Algo::fRecoTracks"}; ///< reconstructed tracks + Vector<L1HitIndex_t> fRecoHits {"L1Algo::fRecoHits"}; ///< packed hits of reconstructed tracks - L1Vector<cbm::algo::ca::Track> fSliceRecoTracks { - "L1Algo::fSliceRecoTracks"}; ///< reconstructed tracks in sub-timeslice - L1Vector<L1HitIndex_t> fSliceRecoHits {"L1Algo::fSliceRecoHits"}; ///< packed hits of reconstructed tracks + Vector<Track> fSliceRecoTracks {"L1Algo::fSliceRecoTracks"}; ///< reconstructed tracks in sub-timeslice + Vector<L1HitIndex_t> fSliceRecoHits {"L1Algo::fSliceRecoHits"}; ///< packed hits of reconstructed tracks /// Created triplets vs station index - L1Vector<L1Triplet> fTriplets[constants::size::MaxNstations] {{"L1Algo::fTriplets"}}; + Vector<L1Triplet> fTriplets[constants::size::MaxNstations] {{"L1Algo::fTriplets"}}; /// Track candidates created out of adjacent triplets before the final track selection. /// The candidates may share any amount of hits. - L1Vector<L1Branch> fTrackCandidates {"L1Algo::fTrackCandidates"}; + Vector<L1Branch> fTrackCandidates {"L1Algo::fTrackCandidates"}; ///< indices of the sub-slice hits - L1Vector<L1HitIndex_t> fSliceHitIds[constants::size::MaxNstations] {"L1Algo::fSliceHitIds"}; + Vector<L1HitIndex_t> fSliceHitIds[constants::size::MaxNstations] {"L1Algo::fSliceHitIds"}; - L1Vector<L1Hit> fGridHits {"L1Algo::fGridHits"}; ///< hits, ordered with respect to their grid bins - L1Vector<L1Hit> fGridHitsBuf {"L1Algo::fGridHitsBuf"}; ///< hits, ordered with respect to their grid bins + Vector<L1Hit> fGridHits {"L1Algo::fGridHits"}; ///< hits, ordered with respect to their grid bins + Vector<L1Hit> fGridHitsBuf {"L1Algo::fGridHitsBuf"}; ///< hits, ordered with respect to their grid bins - L1Vector<L1HitIndex_t> fGridHitIds {"L1Algo::fGridHitIds"}; ///< indices of grid hits: iGridHit -> iCaHit - L1Vector<L1HitIndex_t> fGridHitIdsBuf {"L1Algo::fGridHitIdsBuf"}; ///< buffer for a new fGridHitIds + Vector<L1HitIndex_t> fGridHitIds {"L1Algo::fGridHitIds"}; ///< indices of grid hits: iGridHit -> iCaHit + Vector<L1HitIndex_t> fGridHitIdsBuf {"L1Algo::fGridHitIdsBuf"}; ///< buffer for a new fGridHitIds - L1Vector<L1HitPoint> fGridPoints {"L1Algo::fGridPoints"}; ///< grid points parallel to fGridHits - L1Vector<L1HitPoint> fGridPointsBuf {"L1Algo::fGridPointsBuf"}; + Vector<L1HitPoint> fGridPoints {"L1Algo::fGridPoints"}; ///< grid points parallel to fGridHits + Vector<L1HitPoint> fGridPointsBuf {"L1Algo::fGridPointsBuf"}; L1HitIndex_t fGridHitStartIndex[constants::size::MaxNstations + 1] {0}; L1HitIndex_t fGridHitStopIndex[constants::size::MaxNstations + 1] {0}; - L1Vector<int> fStripToTrack {"L1Algo::fStripToTrack"}; // strip to track pointers + Vector<int> fStripToTrack {"L1Algo::fStripToTrack"}; // strip to track pointers TrackingMode fTrackingMode {kSts}; @@ -332,8 +337,8 @@ public: int isec {0}; // iteration TODO: to be dispatched (S.Zharko, 21.06.2022) const L1CAIteration* fpCurrentIteration = nullptr; ///< pointer to the current CA track finder iteration - L1Vector<int> fHitFirstTriplet {"L1Algo::fHitFirstTriplet"}; /// link hit -> first triplet { hit, *, *} - L1Vector<int> fHitNtriplets {"L1Algo::fHitNtriplets"}; /// link hit ->n triplets { hit, *, *} + Vector<int> fHitFirstTriplet {"L1Algo::fHitFirstTriplet"}; /// link hit -> first triplet { hit, *, *} + Vector<int> fHitNtriplets {"L1Algo::fHitNtriplets"}; /// link hit ->n triplets { hit, *, *} private: diff --git a/reco/L1/L1Algo/L1BaseStationInfo.h b/reco/L1/L1Algo/L1BaseStationInfo.h index cba776b6b40675f4e4df3d3bcb0a0d67d4b13812..5c0966f39ccc1db5b2d220e441c462ef3a949072 100644 --- a/reco/L1/L1Algo/L1BaseStationInfo.h +++ b/reco/L1/L1Algo/L1BaseStationInfo.h @@ -16,9 +16,9 @@ #define L1BaseStationInfo_h 1 // L1 Core +#include "L1Material.h" #include "L1ObjectInitController.h" #include "L1Station.h" -#include "L1Material.h" // C++ std #include <bitset> #include <functional> diff --git a/reco/L1/L1Algo/L1Branch.h b/reco/L1/L1Algo/L1Branch.h index 07ff6c13845a2b4ad909279d1ac57d3705675de3..477d1ded22cca47263eb0a90c6d974aaea920fc8 100644 --- a/reco/L1/L1Algo/L1Branch.h +++ b/reco/L1/L1Algo/L1Branch.h @@ -8,9 +8,14 @@ #ifndef L1Branch_h #define L1Branch_h +#include "CaVector.h" #include "L1Def.h" #include "L1Hit.h" -#include "L1Vector.h" + +namespace +{ + using namespace cbm::algo::ca; // TMP!! +} /// /// L1Branch class describes a search branch of the CA tracker @@ -25,7 +30,7 @@ struct L1Branch { fscal chi2 {0.}; int fID {0}; bool fIsAlive {0}; - L1Vector<L1HitIndex_t> fHits {"L1Branch::fHits"}; + Vector<L1HitIndex_t> fHits {"L1Branch::fHits"}; // static bool compareCand(const L1Branch *a, const L1Branch *b){ // diff --git a/reco/L1/L1Algo/L1BranchExtender.cxx b/reco/L1/L1Algo/L1BranchExtender.cxx index 42dc1bbb076473e105b78e6fc8d083fdc25473ab..ad235dced702f3a9345d8a5146b9c0453e4fc9ca 100644 --- a/reco/L1/L1Algo/L1BranchExtender.cxx +++ b/reco/L1/L1Algo/L1BranchExtender.cxx @@ -5,6 +5,7 @@ #include <iostream> #include "CaTrack.h" +#include "CaVector.h" #include "L1Algo.h" #include "L1Branch.h" #include "L1HitArea.h" @@ -12,10 +13,10 @@ #include "L1TrackPar.h" // using namespace std; +using cbm::algo::ca::Vector; // TMP!! using std::cout; using std::endl; - /// Fit track /// t - track with hits /// T - track params @@ -33,7 +34,7 @@ void L1Algo::BranchFitterFast(const L1Branch& t, L1TrackPar& Tout, const bool up L1TrackPar& T = fit.Tr(); // get hits of current track - const L1Vector<L1HitIndex_t>& hits = t.fHits; // array of indeses of hits of current track + const Vector<L1HitIndex_t>& hits = t.fHits; // array of indeses of hits of current track const int nHits = t.NHits; const signed short int step = -2 * static_cast<int>(upstream) + 1; // increment for station index @@ -137,7 +138,7 @@ void L1Algo::BranchFitter(const L1Branch& t, L1TrackPar& T, const bool upstream, /// initialize - should be params ititialized. 1 - yes. void L1Algo::FindMoreHits(L1Branch& t, L1TrackPar& Tout, const bool upstream, const fvec qp0) { - L1Vector<L1HitIndex_t> newHits {"L1TrackExtender::newHits"}; + Vector<L1HitIndex_t> newHits {"L1TrackExtender::newHits"}; newHits.reserve(fParameters.GetNstationsActive()); L1Fit fit; diff --git a/reco/L1/L1Algo/L1CaTrackFinder.cxx b/reco/L1/L1Algo/L1CaTrackFinder.cxx index 1ad290493866bb3ac102d99df4d29c54b5ddefe4..6e45e6fd4b20aaa6f6983287097cc024a4a2c083 100644 --- a/reco/L1/L1Algo/L1CaTrackFinder.cxx +++ b/reco/L1/L1Algo/L1CaTrackFinder.cxx @@ -1,6 +1,6 @@ -/* Copyright (C) 2009-2021 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt +/* Copyright (C) 2009-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt SPDX-License-Identifier: GPL-3.0-only - Authors: Valentina Akishina, Ivan Kisel, Sergey Gorbunov [committer], Igor Kulakov, Maksym Zyzak */ + Authors: Valentina Akishina, Ivan Kisel, Sergey Gorbunov [committer], Igor Kulakov, Sergei Zharko, Maksym Zyzak */ /* *===================================================== @@ -180,9 +180,7 @@ void L1Algo::CaTrackFinder() // TODO: only add those hits from the region before tsStartNew that belong to the not stored tracks int trackFirstHit = 0; - for (L1Vector<Track>::iterator it = fSliceRecoTracks.begin(); it != fSliceRecoTracks.end(); - trackFirstHit += it->fNofHits, it++) { - + for (auto it = fSliceRecoTracks.begin(); it != fSliceRecoTracks.end(); trackFirstHit += it->fNofHits, it++) { bool isTrackCompletelyInOverlap = true; for (int ih = 0; ih < it->fNofHits; ih++) { int caHitId = fSliceRecoHits[trackFirstHit + ih]; diff --git a/reco/L1/L1Algo/L1CaTrackFinderSlice.cxx b/reco/L1/L1Algo/L1CaTrackFinderSlice.cxx index 94dcaf6fb735fc43f20110d6beefa8f111798f43..5efbf005dbbfb83895db88aeef47eb03f92c29c4 100644 --- a/reco/L1/L1Algo/L1CaTrackFinderSlice.cxx +++ b/reco/L1/L1Algo/L1CaTrackFinderSlice.cxx @@ -938,7 +938,7 @@ void L1Algo::CAFindTrack(int ista, L1Branch& best_tr, unsigned char& best_L, fsc } #ifdef DRAW -void L1Algo::DrawRecoTracksTime(const L1Vector<CbmL1Track>& tracks) +void L1Algo::DrawRecoTracksTime(const cbm::algo::ca::Vector<CbmL1Track>& tracks) { // TODO: find where the DrawRecoTracksTime is. It is missing in the git repository. //draw->DrawRecoTracksTime(tracks); diff --git a/reco/L1/L1Algo/L1CloneMerger.cxx b/reco/L1/L1Algo/L1CloneMerger.cxx index 2145c9f355196d18b6dd8c04d75176bd28cdeb4b..ef106c9536269f1348e2d4bdd81afad05eb71cf0 100644 --- a/reco/L1/L1Algo/L1CloneMerger.cxx +++ b/reco/L1/L1Algo/L1CloneMerger.cxx @@ -12,11 +12,13 @@ #include <iostream> #include "CaTrack.h" +#include "CaVector.h" #include "L1Algo.h" #include "L1Fit.h" #include "L1Parameters.h" -using Track = cbm::algo::ca::Track; +using cbm::algo::ca::Track; // TMP +using cbm::algo::ca::Vector; // TMP // --------------------------------------------------------------------------------------------------------------------- // @@ -28,16 +30,16 @@ L1CloneMerger::~L1CloneMerger() {} // --------------------------------------------------------------------------------------------------------------------- // -void L1CloneMerger::Exec(L1Vector<Track>& extTracks, L1Vector<L1HitIndex_t>& extRecoHits) +void L1CloneMerger::Exec(Vector<Track>& extTracks, Vector<L1HitIndex_t>& extRecoHits) { - L1Vector<unsigned short>& firstStation = fTrackFirstStation; - L1Vector<unsigned short>& lastStation = fTrackLastStation; - L1Vector<L1HitIndex_t>& firstHit = fTrackFirstHit; - L1Vector<L1HitIndex_t>& lastHit = fTrackLastHit; - L1Vector<unsigned short>& neighbour = fTrackNeighbour; - L1Vector<float>& trackChi2 = fTrackChi2; - L1Vector<char>& isStored = fTrackIsStored; - L1Vector<char>& isDownstreamNeighbour = fTrackIsDownstreamNeighbour; + Vector<unsigned short>& firstStation = fTrackFirstStation; + Vector<unsigned short>& lastStation = fTrackLastStation; + Vector<L1HitIndex_t>& firstHit = fTrackFirstHit; + Vector<L1HitIndex_t>& lastHit = fTrackLastHit; + Vector<unsigned short>& neighbour = fTrackNeighbour; + Vector<float>& trackChi2 = fTrackChi2; + Vector<char>& isStored = fTrackIsStored; + Vector<char>& isDownstreamNeighbour = fTrackIsDownstreamNeighbour; int nTracks = extTracks.size(); diff --git a/reco/L1/L1Algo/L1CloneMerger.h b/reco/L1/L1Algo/L1CloneMerger.h index f4e89d8e060e6cfc51212901d1033ea438c5f474..de0132f8a217f74e0bd29a06a94b754aadbc4d72 100644 --- a/reco/L1/L1Algo/L1CloneMerger.h +++ b/reco/L1/L1Algo/L1CloneMerger.h @@ -11,9 +11,10 @@ #define L1CloneMerger_h 1 #include "CaConstants.h" // TEMPORARY FOR fvec, fscal +#include "CaVector.h" #include "L1Def.h" #include "L1Hit.h" // For L1HitIndex_t -#include "L1Vector.h" + namespace cbm::algo::ca { @@ -22,6 +23,12 @@ namespace cbm::algo::ca class L1Algo; +namespace +{ // TMP + using cbm::algo::ca::Track; + using cbm::algo::ca::Vector; +} // namespace + /// Class implements a clones merger algorithm for the CA (L1) track finder /// class L1CloneMerger { @@ -49,7 +56,7 @@ public: /// Executes track clones merging algorithm and updates input containers /// \param extTracks Reference to the external container of reconstructed tracks /// \param extRecoHits Reference to the external container of reconstructed hit indexes - void Exec(L1Vector<cbm::algo::ca::Track>& extTracks, L1Vector<L1HitIndex_t>&); + void Exec(Vector<Track>& extTracks, Vector<L1HitIndex_t>&); private: // *************** @@ -91,32 +98,32 @@ private: // *************** /// First station of a track - L1Vector<unsigned short> fTrackFirstStation {"L1CloneMerger::fTrackFirstStation"}; + Vector<unsigned short> fTrackFirstStation {"L1CloneMerger::fTrackFirstStation"}; /// Last station of a track - L1Vector<unsigned short> fTrackLastStation {"L1CloneMerger::fTrackLastStation"}; + Vector<unsigned short> fTrackLastStation {"L1CloneMerger::fTrackLastStation"}; /// Index of the first hit of a track - L1Vector<L1HitIndex_t> fTrackFirstHit {"L1CloneMerger::fTrackFirstHit"}; + Vector<L1HitIndex_t> fTrackFirstHit {"L1CloneMerger::fTrackFirstHit"}; /// Index of the last hit of a track - L1Vector<L1HitIndex_t> fTrackLastHit {"L1CloneMerger::fTrackLastHit"}; + Vector<L1HitIndex_t> fTrackLastHit {"L1CloneMerger::fTrackLastHit"}; /// Index (TODO:??) of a track that can be merge with the given track - L1Vector<unsigned short> fTrackNeighbour {"L1CloneMerger::fTrackNeighbour"}; + Vector<unsigned short> fTrackNeighbour {"L1CloneMerger::fTrackNeighbour"}; /// Chi2 value of the track merging procedure - L1Vector<float> fTrackChi2 {"L1CloneMerger::fTrackChi2"}; + Vector<float> fTrackChi2 {"L1CloneMerger::fTrackChi2"}; /// Flag: is the given track already stored to the output - L1Vector<char> fTrackIsStored {"L1CloneMerger::fTrackIsStored"}; + Vector<char> fTrackIsStored {"L1CloneMerger::fTrackIsStored"}; /// Flag: is the track a downstream neighbour of another track - L1Vector<char> fTrackIsDownstreamNeighbour {"L1CloneMerger::fTrackIsDownstreamNeighbour"}; + Vector<char> fTrackIsDownstreamNeighbour {"L1CloneMerger::fTrackIsDownstreamNeighbour"}; - L1Vector<cbm::algo::ca::Track> fTracksNew {"L1CAMergerClones::fTracksNew"}; ///< vector of tracks after the merge + Vector<Track> fTracksNew {"L1CAMergerClones::fTracksNew"}; ///< vector of tracks after the merge - L1Vector<L1HitIndex_t> fRecoHitsNew {"L1CAMergerClones::fRecoHitsNew"}; ///< vector of track hits after the merge + Vector<L1HitIndex_t> fRecoHitsNew {"L1CAMergerClones::fRecoHitsNew"}; ///< vector of track hits after the merge const L1Algo& frAlgo; ///< Reference to the main track finder algorithm class }; diff --git a/reco/L1/L1Algo/L1Field.h b/reco/L1/L1Algo/L1Field.h index dc66cb30a6275244f42dc11a951471fa9af9fc57..4072708740c5a785afd41abc98c6edae3dc4093f 100644 --- a/reco/L1/L1Algo/L1Field.h +++ b/reco/L1/L1Algo/L1Field.h @@ -5,11 +5,12 @@ #ifndef L1Field_h #define L1Field_h 1 +#include <boost/serialization/access.hpp> + #include <string> #include "CaConstants.h" #include "CaSimd.h" -#include <boost/serialization/access.hpp> class L1TrackPar; diff --git a/reco/L1/L1Algo/L1Grid.cxx b/reco/L1/L1Algo/L1Grid.cxx index 8a291e75e952cd6e816a535d52e7e74f9237ab98..f5aa2dab4b98e07baf6290c591c14e0ae039fd80 100644 --- a/reco/L1/L1Algo/L1Grid.cxx +++ b/reco/L1/L1Algo/L1Grid.cxx @@ -14,7 +14,8 @@ #include "L1Algo.h" -using namespace std; +using namespace std; // !! REMOVE +using cbm::algo::ca::Vector; // TMP; /// Copy to memory block [@dest, @dest+@num] num number of times the value of i of type @T with size @typesize. /// uses binary expansion of copied volume for speed up @@ -33,11 +34,9 @@ inline void memset(T* dest, T i, size_t num) } } - -void L1Grid::UpdateIterGrid(unsigned int Nelements, L1Hit* hits, L1Vector<L1HitIndex_t>& indicesBuf, - L1HitIndex_t* indices, L1Vector<L1Hit>& hitsBuf, L1Vector<L1HitPoint>& pointsBuf, - L1HitPoint* points, int& NHitsOnStation, L1Algo& Algo, - const L1Vector<unsigned char>& vSFlag) +void L1Grid::UpdateIterGrid(unsigned int Nelements, L1Hit* hits, Vector<L1HitIndex_t>& indicesBuf, + L1HitIndex_t* indices, Vector<L1Hit>& hitsBuf, Vector<L1HitPoint>& pointsBuf, + L1HitPoint* points, int& NHitsOnStation, L1Algo& Algo, const Vector<unsigned char>& vSFlag) { //L1_SHOW(vSFlag.size()); fFirstHitInBin.reset(fN + 2, 0); diff --git a/reco/L1/L1Algo/L1Grid.h b/reco/L1/L1Algo/L1Grid.h index eb51e2d3a25143c70c1371035138175b797be3d5..c8189811d78710e8eab8bc25a2a5ef9074177d37 100644 --- a/reco/L1/L1Algo/L1Grid.h +++ b/reco/L1/L1Algo/L1Grid.h @@ -15,14 +15,19 @@ #include <assert.h> #include <string.h> +#include "CaVector.h" #include "L1Def.h" #include "L1Hit.h" #include "L1HitPoint.h" -#include "L1Vector.h" class L1Algo; //class L1HitPoint; +namespace +{ + using cbm::algo::ca::Vector; // TMP!! +} + /** * @class L1Grid @@ -92,9 +97,9 @@ public: // }; - void UpdateIterGrid(unsigned int Nelements, L1Hit* hits, L1Vector<L1HitIndex_t>& indicesBuf, L1HitIndex_t* indices, - L1Vector<L1Hit>& hitsBuf, L1Vector<L1HitPoint>& pointsBuf, L1HitPoint* points, - int& NHitsOnStation, L1Algo& Algo, const L1Vector<unsigned char>& vSFlag); + void UpdateIterGrid(unsigned int Nelements, L1Hit* hits, Vector<L1HitIndex_t>& indicesBuf, L1HitIndex_t* indices, + Vector<L1Hit>& hitsBuf, Vector<L1HitPoint>& pointsBuf, L1HitPoint* points, int& NHitsOnStation, + L1Algo& Algo, const Vector<unsigned char>& vSFlag); private: @@ -110,8 +115,8 @@ private: float fStepTInv {0.f}; //* inverse bin size in T int fBinInGrid {0}; - L1Vector<L1HitIndex_t> fFirstHitInBin {"L1Grid::fFirstHitInBin"}; - L1Vector<L1HitIndex_t> fHitsInBin {"L1Grid::fHitsInBin"}; + Vector<L1HitIndex_t> fFirstHitInBin {"L1Grid::fFirstHitInBin"}; + Vector<L1HitIndex_t> fHitsInBin {"L1Grid::fHitsInBin"}; // vector <omp_lock_t> lock; }; diff --git a/reco/L1/L1Algo/L1InitManager.cxx b/reco/L1/L1Algo/L1InitManager.cxx index b68f493bd2b8cd64b6990f1856ac295621865d8b..73ff58facb345113fb8271cb92606409e5fdaead 100644 --- a/reco/L1/L1Algo/L1InitManager.cxx +++ b/reco/L1/L1Algo/L1InitManager.cxx @@ -398,8 +398,8 @@ void L1InitManager::SetActiveDetectorIDs(const L1DetectorIDSet_t& detectorIDs) // void L1InitManager::SetCAIterationsNumberCrosscheck(int nIterations) { - fCAIterationsNumberCrosscheck = nIterations; - L1Vector<L1CAIteration>& iterationsContainer = fParameters.fCAIterations; + fCAIterationsNumberCrosscheck = nIterations; + auto& iterationsContainer = fParameters.fCAIterations; // NOTE: should be called to prevent multiple copies of objects between the memory reallocations iterationsContainer.reserve(nIterations); diff --git a/reco/L1/L1Algo/L1InitManager.h b/reco/L1/L1Algo/L1InitManager.h index 1a0f5e9fa2103ce1b4853441ab399e0f6b5e67dc..efd6382249edc7299c1504440f982f1c63215e24 100644 --- a/reco/L1/L1Algo/L1InitManager.h +++ b/reco/L1/L1Algo/L1InitManager.h @@ -18,6 +18,7 @@ #include <unordered_map> #include "CaConstants.h" +#include "CaSimd.h" #include "L1BaseStationInfo.h" #include "L1CAIteration.h" #include "L1EArray.h" @@ -25,8 +26,7 @@ #include "L1ObjectInitController.h" #include "L1Parameters.h" #include "L1Utils.h" -#include "L1Vector.h" -#include "CaSimd.h" + class L1ConfigRW; class L1Algo; diff --git a/reco/L1/L1Algo/L1InputData.h b/reco/L1/L1Algo/L1InputData.h index 35afbd0a602f84ebafeb483e90513a7695a67fed..d7c6283670c85aeb5ead8c6dbe0b8d20bb29aed4 100644 --- a/reco/L1/L1Algo/L1InputData.h +++ b/reco/L1/L1Algo/L1InputData.h @@ -14,8 +14,13 @@ #include <boost/serialization/array.hpp> #include "CaConstants.h" +#include "CaVector.h" #include "L1Hit.h" -#include "L1Vector.h" + +namespace +{ + using cbm::algo::ca::Vector; +} /// Class L1InputData represents a block of the input data to the L1 tracking algorithm per event or time slice. /// Filling of the L1InputData is carried out with L1IODataManager class @@ -67,7 +72,7 @@ public: const L1Hit& GetHit(L1HitIndex_t index) const { return fHits[index]; } /// Gets reference to hits vector - const L1Vector<L1Hit>& GetHits() const { return fHits; } + const Vector<L1Hit>& GetHits() const { return fHits; } /// Gets number of hits in the hits vector L1HitIndex_t GetNhits() const { return fHits.size(); } @@ -107,11 +112,11 @@ private: // *************************** /// @brief Sample of input hits - L1Vector<L1Hit> fHits {"L1InputData::fHits"}; + Vector<L1Hit> fHits {"L1InputData::fHits"}; /// @brief Index of the first hit in the sorted hits vector for a given data stream - L1Vector<L1HitIndex_t> fStreamStartIndices {"L1InputData::fStreamStartIndices"}; - L1Vector<L1HitIndex_t> fStreamStopIndices {"L1InputData::fStreamStopIndices"}; + Vector<L1HitIndex_t> fStreamStartIndices {"L1InputData::fStreamStartIndices"}; + Vector<L1HitIndex_t> fStreamStopIndices {"L1InputData::fStreamStopIndices"}; /// @brief Number of hit keys used for rejecting fake STS hits int fNhitKeys = -1; diff --git a/reco/L1/L1Algo/L1Parameters.h b/reco/L1/L1Algo/L1Parameters.h index 077a62ce89ff54b14e3136668a00c7ccf1104af3..719b4271d9a2ec675c8a3c448a47010072367403 100644 --- a/reco/L1/L1Algo/L1Parameters.h +++ b/reco/L1/L1Algo/L1Parameters.h @@ -19,18 +19,18 @@ #include <utility> #include "CaConstants.h" +#include "CaVector.h" #include "L1CAIteration.h" #include "L1Material.h" #include "L1SearchWindow.h" #include "L1Station.h" -#include "L1Utils.h" -#include "L1Vector.h" +#include "L1Utils.h" // IS DEPRECATED? class L1InitManager; enum class L1DetectorID; /// Type definitions for used containers -using L1IterationsContainer_t = L1Vector<L1CAIteration>; +using L1IterationsContainer_t = cbm::algo::ca::Vector<L1CAIteration>; using L1StationsContainer_t = std::array<L1Station, constants::size::MaxNstations>; using L1MaterialContainer_t = std::array<L1Material, constants::size::MaxNstations>; diff --git a/reco/L1/L1Algo/L1Station.h b/reco/L1/L1Algo/L1Station.h index 62389f9d51b7a3c62721dda0c8cf842690e1639f..78d0d3b579265cc7c5e3d703b64f84b765930840 100644 --- a/reco/L1/L1Algo/L1Station.h +++ b/reco/L1/L1Algo/L1Station.h @@ -9,8 +9,8 @@ #include <type_traits> #include "CaConstants.h" -#include "L1Field.h" #include "CaSimd.h" +#include "L1Field.h" #include "L1Utils.h" /// Structure L1Station diff --git a/reco/L1/L1Algo/L1TripletConstructor.cxx b/reco/L1/L1Algo/L1TripletConstructor.cxx index 3acf332ab1e144e108a0b5f5bedeffd34ab362b4..97aa57118a09ec6e7850bb7ca0698a7d4770c646 100644 --- a/reco/L1/L1Algo/L1TripletConstructor.cxx +++ b/reco/L1/L1Algo/L1TripletConstructor.cxx @@ -19,6 +19,8 @@ #include "L1HitPoint.h" #include "L1TrackPar.h" +using cbm::ca::tools::Debugger; + L1TripletConstructor::L1TripletConstructor(L1Algo* algo) { fAlgo = algo; } void L1TripletConstructor::InitStations(int istal, int istam, int istar) @@ -197,7 +199,7 @@ void L1TripletConstructor::FitDoublets() // ---- Add the middle hits to parameters estimation ---- - L1Vector<L1HitIndex_t> hitsMtmp("L1TripletConstructor::hitsMtmp", fHitsM_2); + Vector<L1HitIndex_t> hitsMtmp("L1TripletConstructor::hitsMtmp", fHitsM_2); fHitsM_2.clear(); fTracks_2.clear(); @@ -337,7 +339,7 @@ void L1TripletConstructor::FindRightHit() if (iMC < 0 || iMC != fAlgo->GetMcTrackIdForGridHit(indM)) { continue; } } - L1Vector<L1HitIndex_t> collectedHits; + Vector<L1HitIndex_t> collectedHits; CollectHits(T2, fIstaR, fAlgo->fTripletChi2Cut, iMC, collectedHits, fAlgo->fParameters.GetMaxTripletPerDoublets()); if (collectedHits.size() >= fAlgo->fParameters.GetMaxTripletPerDoublets()) { @@ -376,8 +378,8 @@ void L1TripletConstructor::FitTriplets() /// Refit Triplets if (dumpTriplets) { - ca::tools::Debugger::Instance().AddNtuple("triplets", - "ev:iter:i0:x0:y0:z0:i1:x1:y1:z1:i2:x2:y2:z2:mc:sta:p:vx:vy:vz:chi2:ndf"); + Debugger::Instance().AddNtuple("triplets", + "ev:iter:i0:x0:y0:z0:i1:x1:y1:z1:i2:x2:y2:z2:mc:sta:p:vx:vy:vz:chi2:ndf"); } L1Fit fit; @@ -547,9 +549,9 @@ void L1TripletConstructor::FitTriplets() if ((mc1 >= 0) && (mc1 == mc2) && (mc1 == mc3)) { const CbmL1MCTrack& mctr = CbmL1::Instance()->GetMcTracks()[mc1]; - ca::tools::Debugger::Instance().FillNtuple("triplets", mctr.iEvent, fAlgo->isec, ih0, h0.x, h0.y, h0.z, ih1, - h1.x, h1.y, h1.z, ih2, h2.x, h2.y, h2.z, mc1, fIstaL, mctr.p, mctr.x, - mctr.y, mctr.z, (float) T.chi2[0], (float) T.NDF[0]); + Debugger::Instance().FillNtuple("triplets", mctr.iEvent, fAlgo->isec, ih0, h0.x, h0.y, h0.z, ih1, h1.x, h1.y, + h1.z, ih2, h2.x, h2.y, h2.z, mc1, fIstaL, mctr.p, mctr.x, mctr.y, mctr.z, + (float) T.chi2[0], (float) T.NDF[0]); } } } //i3 @@ -610,7 +612,7 @@ void L1TripletConstructor::StoreTriplets() void L1TripletConstructor::CollectHits(const L1TrackPar& Tr, const int iSta, const double chi2Cut, const int iMC, - L1Vector<L1HitIndex_t>& collectedHits, int maxNhits) + Vector<L1HitIndex_t>& collectedHits, int maxNhits) { /// Collect hits on a station diff --git a/reco/L1/L1Algo/L1TripletConstructor.h b/reco/L1/L1Algo/L1TripletConstructor.h index 3e1e8c64680fcf48554ba06f8094d33088d405b6..d3b2d6269f7f742a0382a93b36a3362f7da87790 100644 --- a/reco/L1/L1Algo/L1TripletConstructor.h +++ b/reco/L1/L1Algo/L1TripletConstructor.h @@ -5,6 +5,7 @@ #ifndef L1TripletConstructor_h #define L1TripletConstructor_h +#include "CaVector.h" #include "L1Algo.h" #include "L1Field.h" #include "L1Fit.h" @@ -12,7 +13,11 @@ #include "L1Station.h" #include "L1TrackPar.h" #include "L1Triplet.h" -#include "L1Vector.h" + +namespace +{ // TMP!! + using cbm::algo::ca::Vector; +} /// Construction of triplets for the CA tracker @@ -47,7 +52,7 @@ public: void CreateTripletsForHit(int istal, int istam, int istar, L1HitIndex_t ihl); - const L1Vector<L1Triplet>& GetTriplets() const { return fTriplets; } + const Vector<L1Triplet>& GetTriplets() const { return fTriplets; } /// Find the doublets. Reformat data in the portion of doublets. void FindDoublets(); @@ -67,7 +72,7 @@ public: void StoreTriplets(); void CollectHits(const L1TrackPar& Tr, const int iSta, const double chi2Cut, const int iMC, - L1Vector<L1HitIndex_t>& collectedHits, int maxNhits); + Vector<L1HitIndex_t>& collectedHits, int maxNhits); private: /// left station @@ -100,14 +105,14 @@ private: L1TrackPar fTrL; L1FieldRegion fFldL; - L1Vector<L1HitIndex_t> fHitsM_2 {"L1TripletConstructor::fHitsM_2"}; - L1Vector<L1TrackPar> fTracks_2 {"L1TripletConstructor::fTracks_2"}; + Vector<L1HitIndex_t> fHitsM_2 {"L1TripletConstructor::fHitsM_2"}; + Vector<L1TrackPar> fTracks_2 {"L1TripletConstructor::fTracks_2"}; - L1Vector<L1TrackPar> fTracks_3 {"L1TripletConstructor::fTracks_3"}; - L1Vector<L1HitIndex_t> fHitsM_3 {"L1TripletConstructor::fHitsM_3"}; - L1Vector<L1HitIndex_t> fHitsR_3 {"L1TripletConstructor::fHitsR_3"}; + Vector<L1TrackPar> fTracks_3 {"L1TripletConstructor::fTracks_3"}; + Vector<L1HitIndex_t> fHitsM_3 {"L1TripletConstructor::fHitsM_3"}; + Vector<L1HitIndex_t> fHitsR_3 {"L1TripletConstructor::fHitsR_3"}; - L1Vector<L1Triplet> fTriplets {"L1TripletConstructor::fTriplets"}; + Vector<L1Triplet> fTriplets {"L1TripletConstructor::fTriplets"}; L1Fit fFit; diff --git a/reco/L1/L1Algo/L1Vector.h b/reco/L1/L1Algo/L1Vector.h deleted file mode 100644 index 258b875b69efa79763eacdd7dacdc604658003ff..0000000000000000000000000000000000000000 --- a/reco/L1/L1Algo/L1Vector.h +++ /dev/null @@ -1,278 +0,0 @@ -/* Copyright (C) 2021 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt - SPDX-License-Identifier: GPL-3.0-only - Authors: Sergey Gorbunov [committer] */ - -#ifndef L1Vector_H -#define L1Vector_H - -/// @file L1Vector.h -/// @author Sergey Gorbunov -/// @date 2021-06-16 - - -#include "L1Def.h" -#ifndef FAST_CODE -#include <Logger.h> -#endif -#include <boost/serialization/access.hpp> -#include <boost/serialization/base_object.hpp> -#include <boost/serialization/string.hpp> -#include <boost/serialization/vector.hpp> - -#include <sstream> - -/// L1Vector class is a wrapper around std::vector. -/// It does the following: -/// 1. gives names to vectors for better debugging -/// 2. controls the out-of-range access in debug mode -/// 3. supresses methods that are currently not controlled -/// 4. warns when slow memory operations are called, -/// i.e. when the preallocated capacity is reached and the entire vector should be copied to a new place -/// 5. blocks usage of boolean vectors, as they have a special -/// space-optimized but slow implementation in std. (Use L1Vector<char> instead). -/// - - -template<class T> -class L1Vector : private std::vector<T> { - friend class boost::serialization::access; - -public: - typedef std::vector<T> Tbase; - - template<typename... Tinput> - L1Vector(Tinput... value) : Tbase(value...) - { - } - - template<typename... Tinput> - L1Vector(const char* name, Tinput... value) : Tbase(value...) - , fName(name) - { - } - - /// Constructor to make initializations like L1Vector<int> myVector {"MyVector", {1, 2, 3}} - L1Vector(const std::string& name, std::initializer_list<T> init) : Tbase(init), fName(name) {} - - /// Copy constructor - L1Vector(const L1Vector& v) : Tbase() - { - //LOG(info) << "\033[1;32mCALL L1Vector copy constructor\033[0m"; - *this = v; - } - - /// Move constructor - L1Vector(L1Vector&& v) noexcept : Tbase(std::move(v)), fName(std::move(v.fName)) - { - //LOG(info) << "\033[1;32mCALL L1Vector move constructor\033[0m"; - } - - /// Copy assignment operator - L1Vector& operator=(const L1Vector& v) - { - //LOG(info) << "\033[1;32mCALL L1Vector copy assignment operator\033[0m"; - if (this != &v) { - fName = v.fName; - Tbase::reserve(v.capacity()); // make sure that the capacity is transmitted - Tbase::assign(v.begin(), v.end()); - } - return *this; - } - - /// Move assignment operator - L1Vector& operator=(L1Vector&& v) noexcept - { - //LOG(info) << "\033[1;32mCALL L1Vector move assignment operator\033[0m"; - if (this != &v) { - std::swap(fName, v.fName); - Tbase::swap(v); - } - return *this; - } - - /// Swap operator - void swap(L1Vector& v) noexcept - { - //LOG(info) << "\033[1;32mCALL L1Vector swap operator\033[0m"; - if (this != &v) { - std::swap(fName, v.fName); - Tbase::swap(v); - } - } - - void SetName(const std::string& s) { fName = s; } - - void SetName(const std::basic_ostream<char>& s) - { - // helps to set a composed name in a single line via: - // SetName(std::stringstream()<<"my name "<<..whatever..); - fName = dynamic_cast<const std::stringstream&>(s).str(); - } - - std::string GetName() const - { - std::string s = " L1Vector<"; - s += fName + "> "; - return s; - } - - template<typename... Tinput> - void reset(std::size_t count, Tinput... value) - { - // does the same as Tbase::assign(), but works with the default T constructor too - // (no second parameter) - Tbase::clear(); - Tbase::resize(count, value...); - } - - template<typename... Tinput> - void enlarge(std::size_t count, Tinput... value) - { - if (count < Tbase::size()) { - LOG(fatal) << "L1Vector \"" << fName << "\"::enlarge(" << count - << "): the new size is smaller than the current one " << Tbase::size() << ", something goes wrong." - << std::endl; - assert(count >= Tbase::size()); - } - if ((!Tbase::empty()) && (count > Tbase::capacity())) { - LOG(warning) << "L1Vector \"" << fName << "\"::enlarge(" << count << "): allocated capacity of " - << Tbase::capacity() << " is reached, the vector of size " << Tbase::size() - << " will be copied to the new place." << std::endl; - } - Tbase::resize(count, value...); - } - - void reduce(std::size_t count) - { - if (count > Tbase::size()) { - LOG(fatal) << "L1Vector \"" << fName << "\"::reduce(" << count - << "): the new size is bigger than the current one " << Tbase::size() << ", something goes wrong." - << std::endl; - assert(count < Tbase::size()); - } - Tbase::resize(count); - } - - void reserve(std::size_t count) - { - if (!Tbase::empty()) { - LOG(fatal) << "L1Vector \"" << fName << "\"::reserve(" << count << "): the vector is not empty; " - << " it will be copied to the new place." << std::endl; - assert(Tbase::empty()); - } - Tbase::reserve(count); - } - - - template<typename Tinput> - void push_back(Tinput value) - { -#ifndef FAST_CODE - if (Tbase::size() >= Tbase::capacity()) { - LOG(warning) << "L1Vector \"" << fName << "\"::push_back(): allocated capacity of " << Tbase::capacity() - << " is reached, re-allocate and copy." << std::endl; - } -#endif - Tbase::push_back(value); - } - - template<typename Tinput> - void push_back_no_warning(Tinput value) - { - Tbase::push_back(value); - } - - template<typename... Tinput> - void emplace_back(Tinput&&... value) - { -#ifndef FAST_CODE - if (Tbase::size() >= Tbase::capacity()) { - LOG(warning) << "L1Vector \"" << fName << "\"::emplace_back(): allocated capacity of " << Tbase::capacity() - << " is reached, re-allocate and copy." << std::endl; - } -#endif - Tbase::emplace_back(value...); - } - - T& operator[](std::size_t pos) - { -#ifndef FAST_CODE - if (pos >= Tbase::size()) { - LOG(fatal) << "L1Vector \"" << fName << "\": trying to access element " << pos - << " outside of the vector of the size of " << Tbase::size() << std::endl; - assert(pos < Tbase::size()); - } -#endif - return Tbase::operator[](pos); - } - - const T& operator[](std::size_t pos) const - { -#ifndef FAST_CODE - if (pos >= Tbase::size()) { - LOG(fatal) << "L1Vector \"" << fName << "\": trying to access element " << pos - << " outside of the vector of the size of " << Tbase::size() << std::endl; - assert(pos < Tbase::size()); - } -#endif - return Tbase::operator[](pos); - } - - T& back() - { -#ifndef FAST_CODE - if (Tbase::size() == 0) { - LOG(fatal) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; - assert(Tbase::size() > 0); - } -#endif - return Tbase::back(); - } - - const T& back() const - { -#ifndef FAST_CODE - if (Tbase::size() == 0) { - LOG(fatal) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; - assert(Tbase::size() > 0); - } -#endif - return Tbase::back(); - } - - using Tbase::begin; - using Tbase::capacity; - using Tbase::clear; - using Tbase::end; - using Tbase::insert; //TODO:: make it private - using Tbase::pop_back; - using Tbase::rbegin; - using Tbase::reserve; - using Tbase::size; - using typename Tbase::iterator; - -private: - std::string fName {"no name"}; - using Tbase::assign; // use reset() instead - using Tbase::at; - using Tbase::resize; - - template<class Archive> - void serialize(Archive& ar, const unsigned int /*version*/) - { - ar& boost::serialization::base_object<Tbase>(*this); - ar& fName; - } -}; - -/// -/// std::vector<bool> has a special implementation that is space-optimized -/// and therefore slow and not thread-safe. -/// That is why one should use L1Vector<char> instead. -/// -template<> -class L1Vector<bool> { /// Make sure that L1Vector<bool> is not used -}; - - -#endif diff --git a/reco/L1/L1Algo/utils/CaAlgoRandom.cxx b/reco/L1/L1Algo/utils/CaAlgoRandom.cxx index 218abc50bcf55c1f94163b24d98839360bb867a4..da8d04f31ec1521cd9e21bddaa84a94596056a12 100644 --- a/reco/L1/L1Algo/utils/CaAlgoRandom.cxx +++ b/reco/L1/L1Algo/utils/CaAlgoRandom.cxx @@ -9,7 +9,7 @@ #include "CaAlgoRandom.h" -using ca::algo::Random; +using cbm::algo::ca::Random; // --------------------------------------------------------------------------------------------------------------------- diff --git a/reco/L1/L1Algo/utils/CaAlgoRandom.h b/reco/L1/L1Algo/utils/CaAlgoRandom.h index ba3dd2a7026c89b86dec35991202c666e1a83a63..f203c6d60dc9bbec5f69cc152a1da38434d7c1a5 100644 --- a/reco/L1/L1Algo/utils/CaAlgoRandom.h +++ b/reco/L1/L1Algo/utils/CaAlgoRandom.h @@ -16,7 +16,7 @@ #include <random> #include <type_traits> -namespace ca::algo +namespace cbm::algo::ca { /// @brief A class, providing ROOT-free access to randomly generated variables class Random { @@ -78,38 +78,39 @@ namespace ca::algo mutable GeneratorType_t fGenerator; ///< Random number generator unsigned int fSeed = 1; ///< Random number seed }; -} // namespace ca::algo -// ********************************************************* -// ** Template and inline function implementation ** -// ********************************************************* - -// --------------------------------------------------------------------------------------------------------------------- -// -template<typename T, std::enable_if_t<std::is_floating_point<T>::value, bool>> -T ca::algo::Random::BoundedGaus(const T& mean, const T& sigma, const T& nSigmas) const -{ - LOG_IF(fatal, !(nSigmas > 0 && std::isfinite(nSigmas))) << "ca::algo::Random::BoundedGaus nSigmas = " << nSigmas; - LOG_IF(fatal, !(sigma > 0 && std::isfinite(sigma))) << "ca::algo::Random::BoundedGaus sigma = " << sigma; - - std::normal_distribution rndGaus {mean, sigma}; - double res = 0; - do { - res = rndGaus(fGenerator); - } while (std::fabs(res - mean) > nSigmas * sigma); - return res; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<typename T, std::enable_if_t<std::is_floating_point<T>::value, bool>> -T ca::algo::Random::Uniform(const T& mean, const T& sigma) const -{ - LOG_IF(fatal, !(sigma > 0 && std::isfinite(sigma))) << "ca::algo::Random::Uniform sigma = " << sigma << " is illegal"; - std::uniform_real_distribution rnd {-sigma * std::sqrt(3) + mean, sigma * std::sqrt(3) + mean}; - return rnd(fGenerator); -} - + // ********************************************************* + // ** Template and inline function implementation ** + // ********************************************************* + + // --------------------------------------------------------------------------------------------------------------------- + // + template<typename T, std::enable_if_t<std::is_floating_point<T>::value, bool>> + T Random::BoundedGaus(const T& mean, const T& sigma, const T& nSigmas) const + { + LOG_IF(fatal, !(nSigmas > 0 && std::isfinite(nSigmas))) << "ca::algo::Random::BoundedGaus nSigmas = " << nSigmas; + LOG_IF(fatal, !(sigma > 0 && std::isfinite(sigma))) << "ca::algo::Random::BoundedGaus sigma = " << sigma; + + std::normal_distribution rndGaus {mean, sigma}; + double res = 0; + do { + res = rndGaus(fGenerator); + } while (std::fabs(res - mean) > nSigmas * sigma); + return res; + } + + // --------------------------------------------------------------------------------------------------------------------- + // + template<typename T, std::enable_if_t<std::is_floating_point<T>::value, bool>> + T Random::Uniform(const T& mean, const T& sigma) const + { + LOG_IF(fatal, !(sigma > 0 && std::isfinite(sigma))) + << "ca::algo::Random::Uniform sigma = " << sigma << " is illegal"; + std::uniform_real_distribution rnd {-sigma * std::sqrt(3) + mean, sigma * std::sqrt(3) + mean}; + return rnd(fGenerator); + } + +} // namespace cbm::algo::ca #endif // CaAlgoRandom_h diff --git a/reco/L1/L1Algo/utils/CaMonitor.h b/reco/L1/L1Algo/utils/CaMonitor.h index 70f908db614b018899956acd061a41c20a60d645..f205f9d2517023cd5ecc728689d59c3657a6cd3e 100644 --- a/reco/L1/L1Algo/utils/CaMonitor.h +++ b/reco/L1/L1Algo/utils/CaMonitor.h @@ -16,7 +16,7 @@ #include "L1EArray.h" -namespace ca +namespace cbm::algo::ca { /// @class Monitor /// @brief Monitor class for CA tracking @@ -70,45 +70,45 @@ namespace ca MonitorableArr_t<std::string> faKeyNames {}; ///< Names of keys MonitorableArr_t<int> faKeyCounters {}; ///< Counters of keys }; -} // namespace ca -// ***************************************** -// ** Template function implementations ** -// ***************************************** + // ***************************************** + // ** Template function implementations ** + // ***************************************** -// --------------------------------------------------------------------------------------------------------------------- -// -template<class EMonitorKey> -std::string ca::Monitor<EMonitorKey>::ToString() const -{ - using std::left; - using std::right; - using std::setfill; - using std::setw; - std::stringstream msg; - constexpr size_t width = 24; - size_t fillWidth = (width + 2) * (1 + fvRatioKeys.size()); - msg << "----- Monitor: " << fsName << ' ' << setw(fillWidth - fsName.size() - 16) << setfill('-') << '-' << '\n'; - msg << setfill(' '); - msg << setw(width) << left << "Key" << ' '; - msg << setw(width) << left << "Total" << ' '; - for (auto key : fvRatioKeys) { - msg << setw(width) << left << std::string("per ") + faKeyNames[key] << ' '; - } - msg << '\n'; - for (int iKey = 0; iKey < static_cast<int>(EMonitorKey::kEND); ++iKey) { - msg << setw(width) << left << faKeyNames[iKey] << ' '; - msg << setw(width) << right << faKeyCounters[iKey] << ' '; - for (auto keyDen : fvRatioKeys) { - auto ratio = static_cast<double>(faKeyCounters[iKey]) / faKeyCounters[keyDen]; - msg << setw(width) << right << ratio << ' '; + // --------------------------------------------------------------------------------------------------------------------- + // + template<class EMonitorKey> + std::string Monitor<EMonitorKey>::ToString() const + { + using std::left; + using std::right; + using std::setfill; + using std::setw; + std::stringstream msg; + constexpr size_t width = 24; + size_t fillWidth = (width + 2) * (1 + fvRatioKeys.size()); + msg << "----- Monitor: " << fsName << ' ' << setw(fillWidth - fsName.size() - 16) << setfill('-') << '-' << '\n'; + msg << setfill(' '); + msg << setw(width) << left << "Key" << ' '; + msg << setw(width) << left << "Total" << ' '; + for (auto key : fvRatioKeys) { + msg << setw(width) << left << std::string("per ") + faKeyNames[key] << ' '; } msg << '\n'; + for (int iKey = 0; iKey < static_cast<int>(EMonitorKey::kEND); ++iKey) { + msg << setw(width) << left << faKeyNames[iKey] << ' '; + msg << setw(width) << right << faKeyCounters[iKey] << ' '; + for (auto keyDen : fvRatioKeys) { + auto ratio = static_cast<double>(faKeyCounters[iKey]) / faKeyCounters[keyDen]; + msg << setw(width) << right << ratio << ' '; + } + msg << '\n'; + } + msg << setw(fillWidth) << setfill('-') << '-' << '\n'; + msg << setfill(' '); + return msg.str(); } - msg << setw(fillWidth) << setfill('-') << '-' << '\n'; - msg << setfill(' '); - return msg.str(); -} +} // namespace cbm::algo::ca -#endif // CaToolsMonitor_h +#endif // CaMonitor_h diff --git a/reco/L1/L1Algo/utils/L1AlgoDraw.cxx b/reco/L1/L1Algo/utils/L1AlgoDraw.cxx index 63aebd7ac79566b35e2f84e66b38e7c0e7e28932..2b458eb4d0c6afd9b4da2631a2e5aa9e3a534e6d 100644 --- a/reco/L1/L1Algo/utils/L1AlgoDraw.cxx +++ b/reco/L1/L1Algo/utils/L1AlgoDraw.cxx @@ -26,9 +26,11 @@ // NOTE: using std::cout, std::endl were removed from one of the headers, and it // caused errors. Please, don't use "using" in headers +using cbm::algo::ca::Track; using std::cout; using std::endl; -using Track = cbm::algo::ca::Track; +using std::map; +using std::vector; L1AlgoDraw::L1AlgoDraw() { @@ -237,7 +239,7 @@ void L1AlgoDraw::DrawRecoTracks() // CbmL1 &L1 = *CbmL1::Instance(); int curRecoHit = 0; - L1Vector<L1HitIndex_t>& recoHits = algo->fSliceRecoHits; + cbm::algo::ca::Vector<L1HitIndex_t>& recoHits = algo->fSliceRecoHits; for (vector<Track>::iterator it = algo->fSliceRecoTracks.begin(); it != algo->fSliceRecoTracks.end(); ++it) { Track& T = *it; int nHits = T.fNofHits; diff --git a/reco/L1/L1LinkDef.h b/reco/L1/L1LinkDef.h index b56c0faefb32c0da82e2c6326e9aed6d70780e1c..b915e9a64789329b9ffb0afaba7da8e66f25f85e 100644 --- a/reco/L1/L1LinkDef.h +++ b/reco/L1/L1LinkDef.h @@ -32,13 +32,13 @@ #pragma link C++ enum cbm::ca::ETrackType; #pragma link C++ enum class L1DetectorID; #pragma link C++ class cbm::ca::OutputQa + ; -#pragma link C++ class ca::tools::WindowFinder + ; +#pragma link C++ class cbm::ca::tools::WindowFinder + ; //#pragma link C++ class cbm::ca::IdealHitProducer < L1DetectorID::kMvd > + ; //#pragma link C++ class cbm::ca::IdealHitProducer < L1DetectorID::kSts > + ; //#pragma link C++ class cbm::ca::IdealHitProducer < L1DetectorID::kMuch > + ; //#pragma link C++ class cbm::ca::IdealHitProducer < L1DetectorID::kTrd > + ; //#pragma link C++ class cbm::ca::IdealHitProducer < L1DetectorID::kTof > + ; #pragma link C++ class cbm::ca::IdealHitProducer + ; -#pragma link C++ class ca::tools::MaterialHelper + ; +#pragma link C++ class cbm::ca::tools::MaterialHelper + ; #endif diff --git a/reco/L1/catools/CaToolsDebugger.cxx b/reco/L1/catools/CaToolsDebugger.cxx index 119d685e847c0ba3b94a98dcf07185fa70a8b848..9394c4993e6e8e00d450b6c5441b78598d9d0909 100644 --- a/reco/L1/catools/CaToolsDebugger.cxx +++ b/reco/L1/catools/CaToolsDebugger.cxx @@ -25,7 +25,7 @@ n->Draw("chi2") #include <stdarg.h> -using namespace ca::tools; +using cbm::ca::tools::Debugger; Debugger& Debugger::Instance() { diff --git a/reco/L1/catools/CaToolsDebugger.h b/reco/L1/catools/CaToolsDebugger.h index 8c4aa8fffa77d91106d0831241989b3fc3722c61..fc9f2c8c3168ca0e6262578f016debc4ec7c5e84 100644 --- a/reco/L1/catools/CaToolsDebugger.h +++ b/reco/L1/catools/CaToolsDebugger.h @@ -16,84 +16,80 @@ class TFile; class TNtuple; -namespace ca +namespace cbm::ca::tools { - namespace tools - { - /// Class ca::tools::Debugger helps to debug CA tracking - /// - class Debugger { - public: - // ***************************************** - // ** Constructors and destructor ** - // ***************************************** - - - /// Default constructor - Debugger(const char* fileName = "CAdebug.root"); - - /// Destructor - ~Debugger(); - - /// Copy constructor - Debugger(const Debugger& other) = delete; - - /// Move constructor - Debugger(Debugger&& other) = delete; - - /// Copy assignment operator - Debugger& operator=(const Debugger& other) = delete; - - /// Move assignment operator - Debugger& operator=(Debugger&& other) = delete; - - /// Instance - static Debugger& Instance(); - - /// Write ntuples to the file - void Write(); - - /// Set new ntuple - void AddNtuple(const char* name, const char* varlist); - - /// Add an entry to ntuple - void FillNtuple(const char* name, float v[]); - - /// Add an entry to ntuple - template<typename... Targs> - void FillNtuple(const char* name, Targs... args) - { - constexpr std::size_t n = sizeof...(Targs); - if (n <= 0) return; - float v[n]; - FillFloatArray(v, args...); - FillNtuple(name, v); - } - - /// Get ntuple index - int GetNtupleIndex(const char* name); - - private: - template<typename T, typename... Targs> - void FillFloatArray(float* v, T val, Targs... args) - { - v[0] = (float) val; - if (sizeof...(args) > 0) { FillFloatArray(v + 1, args...); } - } - - template<typename T> - void FillFloatArray(float* v, T last) - { - v[0] = (float) last; - } - - private: - bool fIsWritten {0}; - TFile* fFile {nullptr}; - std::vector<TNtuple*> fNtuples; - }; - - } // namespace tools -} // namespace ca + /// Class ca::tools::Debugger helps to debug CA tracking + /// + class Debugger { + public: + // ***************************************** + // ** Constructors and destructor ** + // ***************************************** + + + /// Default constructor + Debugger(const char* fileName = "CAdebug.root"); + + /// Destructor + ~Debugger(); + + /// Copy constructor + Debugger(const Debugger& other) = delete; + + /// Move constructor + Debugger(Debugger&& other) = delete; + + /// Copy assignment operator + Debugger& operator=(const Debugger& other) = delete; + + /// Move assignment operator + Debugger& operator=(Debugger&& other) = delete; + + /// Instance + static Debugger& Instance(); + + /// Write ntuples to the file + void Write(); + + /// Set new ntuple + void AddNtuple(const char* name, const char* varlist); + + /// Add an entry to ntuple + void FillNtuple(const char* name, float v[]); + + /// Add an entry to ntuple + template<typename... Targs> + void FillNtuple(const char* name, Targs... args) + { + constexpr std::size_t n = sizeof...(Targs); + if (n <= 0) return; + float v[n]; + FillFloatArray(v, args...); + FillNtuple(name, v); + } + + /// Get ntuple index + int GetNtupleIndex(const char* name); + + private: + template<typename T, typename... Targs> + void FillFloatArray(float* v, T val, Targs... args) + { + v[0] = (float) val; + if (sizeof...(args) > 0) { FillFloatArray(v + 1, args...); } + } + + template<typename T> + void FillFloatArray(float* v, T last) + { + v[0] = (float) last; + } + + private: + bool fIsWritten {0}; + TFile* fFile {nullptr}; + std::vector<TNtuple*> fNtuples; + }; +} // namespace cbm::ca::tools #endif // CaToolsDebugger_h diff --git a/reco/L1/catools/CaToolsDef.h b/reco/L1/catools/CaToolsDef.h new file mode 100644 index 0000000000000000000000000000000000000000..d2607ec3b0da0ce268c35ef97e0d382885fde1ad --- /dev/null +++ b/reco/L1/catools/CaToolsDef.h @@ -0,0 +1,23 @@ +/* Copyright (C) 2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Sergei Zharko [committer] */ + +/// @file CaToolsDef.h +/// @brief Definitions for ca::tools namespace +/// @since 16.09.2023 +/// @author S.Zharko <s.zharko@gsi.de> + +#ifndef CaToolsDef_h +#define CaToolsDef_h 1 + +#include "CaConstants.h" + +namespace cbm::ca::tools +{ + namespace undef = cbm::algo::ca::constants::undef; + namespace phys = cbm::algo::ca::constants::phys; + namespace clrs = cbm::algo::ca::constants::clrs; + namespace ca = cbm::algo::ca; +} // namespace cbm::ca::tools + +#endif // CaToolsDef_h diff --git a/reco/L1/catools/CaToolsHitRecord.cxx b/reco/L1/catools/CaToolsHitRecord.cxx index 2286b1273005a65672b19fd79344ee2d63dad3a6..d18b8b41c59444f100a762fe8e47127e34444439 100644 --- a/reco/L1/catools/CaToolsHitRecord.cxx +++ b/reco/L1/catools/CaToolsHitRecord.cxx @@ -15,7 +15,7 @@ #include <cmath> -using ca::tools::HitRecord; +using cbm::ca::tools::HitRecord; // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/catools/CaToolsHitRecord.h b/reco/L1/catools/CaToolsHitRecord.h index 98d99a982a4609b5f8c1c3808492294654a44aba..180e74604540a9b0a1cb5b8b01380a9ce07b7a0c 100644 --- a/reco/L1/catools/CaToolsHitRecord.h +++ b/reco/L1/catools/CaToolsHitRecord.h @@ -7,12 +7,12 @@ /// \since 15.05.2023 /// \author S.Zharko <s.zharko@gsi.de> -#ifndef CbmCaHitRecord_h -#define CbmCaHitRecord_h 1 +#ifndef CbmCaToolsHitRecord_h +#define CbmCaToolsHitRecord_h 1 #include <cstdint> -namespace ca::tools +namespace cbm::ca::tools { /// @brief A helper structure to store hits information from different detectors in a uniform manner /// @@ -47,4 +47,4 @@ namespace ca::tools }; } // namespace ca::tools -#endif // CbmCaHitRecord_h +#endif // CbmCaToolsHitRecord_h diff --git a/reco/L1/catools/CaToolsLinkKey.h b/reco/L1/catools/CaToolsLinkKey.h index 3cd440f9334eafb70834c966b41c6aa6cff23af4..3689c65bdb7c1fb84b06580e53bd38899591de01 100644 --- a/reco/L1/catools/CaToolsLinkKey.h +++ b/reco/L1/catools/CaToolsLinkKey.h @@ -12,9 +12,9 @@ #include <boost/functional/hash.hpp> -#include "CaConstants.h" +#include "CaToolsDef.h" -namespace ca::tools +namespace cbm::ca::tools { struct LinkKey { /// Constructor from index, event and file of a link @@ -29,9 +29,9 @@ namespace ca::tools return lhs.fFile == rhs.fFile && lhs.fEvent == rhs.fEvent && lhs.fIndex == rhs.fIndex; } - int fIndex = cbm::algo::ca::constants::undef::Int; ///< Index of MC point/track in external data structures - int fEvent = cbm::algo::ca::constants::undef::Int; ///< Index of MC event - int fFile = cbm::algo::ca::constants::undef::Int; ///< Index of MC file + int fIndex = undef::Int; ///< Index of MC point/track in external data structures + int fEvent = undef::Int; ///< Index of MC event + int fFile = undef::Int; ///< Index of MC file }; } // namespace ca::tools @@ -39,8 +39,8 @@ namespace std { /// A hash specialization for ca::tools::LinkKey objects template<> - struct hash<ca::tools::LinkKey> { - std::size_t operator()(const ca::tools::LinkKey& key) const + struct hash<cbm::ca::tools::LinkKey> { + std::size_t operator()(const cbm::ca::tools::LinkKey& key) const { std::size_t res = 0; boost::hash_combine(res, key.fFile); diff --git a/reco/L1/catools/CaToolsMCData.cxx b/reco/L1/catools/CaToolsMCData.cxx index 70d56bf99c08e1f4b5879de5fcd4c58a8aa91277..0aa8c4ad5b055574dbb54c0cd1f5014ae463cc51 100644 --- a/reco/L1/catools/CaToolsMCData.cxx +++ b/reco/L1/catools/CaToolsMCData.cxx @@ -1,4 +1,4 @@ -/* Copyright (C) 2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt +/* Copyright (C) 2022-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt SPDX-License-Identifier: GPL-3.0-only Authors: Sergei Zharko [committer] */ @@ -15,7 +15,7 @@ #include <sstream> #include <utility> // for std::move -using namespace ca::tools; +using cbm::ca::tools::MCData; // --------------------------------------------------------------------------------------------------------------------- // @@ -81,7 +81,7 @@ void MCData::Clear() // --------------------------------------------------------------------------------------------------------------------- // -void MCData::InitTrackInfo(const L1Vector<CbmL1HitDebugInfo>& vHits) +void MCData::InitTrackInfo(const ca::Vector<CbmL1HitDebugInfo>& vHits) { for (auto& aTrk : fvTracks) { // Assign hits to tracks diff --git a/reco/L1/catools/CaToolsMCData.h b/reco/L1/catools/CaToolsMCData.h index d09577ab77fcc85f7c574966e69e1f62c5e3cfeb..d562e9954a1d71b9969de87187f6f0e744c74ce0 100644 --- a/reco/L1/catools/CaToolsMCData.h +++ b/reco/L1/catools/CaToolsMCData.h @@ -13,16 +13,16 @@ #include <numeric> #include <string> -#include "CaConstants.h" +#include "CaToolsDef.h" #include "CaToolsLinkKey.h" #include "CaToolsMCPoint.h" #include "CaToolsMCTrack.h" +#include "CaVector.h" #include "L1EArray.h" -#include "L1Vector.h" enum class L1DetectorID; -namespace ca::tools +namespace cbm::ca::tools { /// This class represents a package for tracking-related data class MCData { @@ -156,7 +156,7 @@ namespace ca::tools /// Initialize tracks: defines indexes of hits and points related to the track, calculates max number of points and /// hits on a station, number of consecutive stations containing a hit or point and number of stations and points /// with hits. - void InitTrackInfo(const L1Vector<CbmL1HitDebugInfo>& vHits); + void InitTrackInfo(const ca::Vector<CbmL1HitDebugInfo>& vHits); /// Reserves memory for tracks to avoid extra allocations void ReserveNofTracks(int nTracks) { fvTracks.reserve(nTracks); } @@ -181,8 +181,8 @@ namespace ca::tools // ** Member variables ** // ****************************** - L1Vector<MCPoint> fvPoints = {"ca::tools::MCData::fvPoints"}; ///< Container of points - L1Vector<MCTrack> fvTracks = {"ca::tools::MCData::fvTracks"}; ///< Container of tracks + ca::Vector<MCPoint> fvPoints = {"ca::tools::MCData::fvPoints"}; ///< Container of points + ca::Vector<MCTrack> fvTracks = {"ca::tools::MCData::fvTracks"}; ///< Container of tracks std::array<int, constants::size::MaxNdetectors> fvNofPointsOrig = {0}; ///< Total number of points by detector std::array<int, constants::size::MaxNdetectors> fvNofPointsUsed = {0}; ///< Number of points used vs. detector @@ -190,28 +190,28 @@ namespace ca::tools std::unordered_map<LinkKey, int> fmPointLinkMap = {}; ///< MC point internal index vs. link std::unordered_map<LinkKey, int> fmTrackLinkMap = {}; ///< MC track internal index vs. link }; -} // namespace ca::tools -// ********************************************************* -// ** Template and inline function implementation ** -// ********************************************************* - -// --------------------------------------------------------------------------------------------------------------------- -// -inline void ca::tools::MCData::AddPoint(const MCPoint& point) -{ - fmPointLinkMap[point.GetLinkKey()] = point.GetId(); - fvPoints.push_back(point); - ++fvNofPointsUsed[static_cast<int>(point.GetDetectorId())]; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -inline void ca::tools::MCData::AddTrack(const MCTrack& track) -{ - fmTrackLinkMap[track.GetLinkKey()] = track.GetId(); - fvTracks.push_back(track); -} + // ********************************************************* + // ** Template and inline function implementation ** + // ********************************************************* + + // ------------------------------------------------------------------------------------------------------------------- + // + inline void MCData::AddPoint(const MCPoint& point) + { + fmPointLinkMap[point.GetLinkKey()] = point.GetId(); + fvPoints.push_back(point); + ++fvNofPointsUsed[static_cast<int>(point.GetDetectorId())]; + } + + // ------------------------------------------------------------------------------------------------------------------- + // + inline void MCData::AddTrack(const MCTrack& track) + { + fmTrackLinkMap[track.GetLinkKey()] = track.GetId(); + fvTracks.push_back(track); + } +} // namespace cbm::ca::tools #endif // CaToolsMCData_h diff --git a/reco/L1/catools/CaToolsMCPoint.cxx b/reco/L1/catools/CaToolsMCPoint.cxx index fc3cc474c44bab32a9e3ef50208dee4f8d644691..f810e8361e830e2dba4187de3edb137d3dc77a61 100644 --- a/reco/L1/catools/CaToolsMCPoint.cxx +++ b/reco/L1/catools/CaToolsMCPoint.cxx @@ -13,7 +13,7 @@ #include <iomanip> #include <sstream> -using namespace ca::tools; +using cbm::ca::tools::MCPoint; // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/catools/CaToolsMCPoint.h b/reco/L1/catools/CaToolsMCPoint.h index d0f7deb010bdf260875cd0a987b1a5befb9e8c3e..bda000cee018d1b200a582cb7b275021f7344a28 100644 --- a/reco/L1/catools/CaToolsMCPoint.h +++ b/reco/L1/catools/CaToolsMCPoint.h @@ -12,17 +12,14 @@ #include <string> -#include "CaConstants.h" +#include "CaToolsDef.h" #include "CaToolsLinkKey.h" -#include "L1Vector.h" +#include "CaVector.h" enum class L1DetectorID; -namespace ca::tools +namespace cbm::ca::tools { - namespace phys = cbm::algo::ca::constants::phys; - namespace undef = cbm::algo::ca::constants::undef; - /// @brief Class describes a unified MC-point, used in CA tracking QA analysis /// class MCPoint { @@ -76,21 +73,18 @@ namespace ca::tools int GetId() const { return fId; } /// @brief Gets inverse speed at reference z of station [ns/cm] - int GetInvSpeed() const - { - return std::sqrt(1. + GetMass() * GetMass() / GetP() / GetP()) * constants::phys::SpeedOfLightInv; - } + int GetInvSpeed() const { return std::sqrt(1. + GetMass() * GetMass() / GetP() / GetP()) * phys::SpeedOfLightInv; } /// @brief Gets inverse speed at entrance to station [ns/cm] int GetInvSpeedIn() const { - return std::sqrt(1. + GetMass() * GetMass() / GetPIn() / GetPIn()) * constants::phys::SpeedOfLightInv; + return std::sqrt(1. + GetMass() * GetMass() / GetPIn() / GetPIn()) * phys::SpeedOfLightInv; } /// @brief Gets inverse speed at exit of station [ns/cm] int GetInvSpeedOut() const { - return std::sqrt(1. + GetMass() * GetMass() / GetPOut() / GetPOut()) * constants::phys::SpeedOfLightInv; + return std::sqrt(1. + GetMass() * GetMass() / GetPOut() / GetPOut()) * phys::SpeedOfLightInv; } /// @brief Gets container of matched hit indexes @@ -336,18 +330,18 @@ namespace ca::tools // **************************** private: - std::array<double, 3> fPos = {undef::Double, undef::Double, - undef::Double}; ///< Position at reference z of station [cm] - std::array<double, 3> fPosIn = {undef::Double, undef::Double, - undef::Double}; ///< Position at entrance to station [cm] - std::array<double, 3> fPosOut = {undef::Double, undef::Double, - undef::Double}; ///< Position at exit of station [cm] - std::array<double, 3> fMom = {undef::Double, undef::Double, - undef::Double}; ///< Momentum at reference z of station [cm] - std::array<double, 3> fMomIn = {undef::Double, undef::Double, - undef::Double}; ///< Momentum at entrance to station [cm] - std::array<double, 3> fMomOut = {undef::Double, undef::Double, - undef::Double}; ///< Momentum at exit of station [cm] + /// \brief Position at reference z of station [cm] + std::array<double, 3> fPos = {undef::Double, undef::Double, undef::Double}; + /// \brief Position at entrance to station [cm] + std::array<double, 3> fPosIn = {undef::Double, undef::Double, undef::Double}; + /// \brief Position at exit of station [cm] + std::array<double, 3> fPosOut = {undef::Double, undef::Double, undef::Double}; + /// \brief Momentum at reference z of station [GeV/c] + std::array<double, 3> fMom = {undef::Double, undef::Double, undef::Double}; + /// \brief Momentum at entrance to station [GeV/c] + std::array<double, 3> fMomIn = {undef::Double, undef::Double, undef::Double}; + /// \brief Momentum at exit of station [cm] + std::array<double, 3> fMomOut = {undef::Double, undef::Double, undef::Double}; double fMass = undef::Double; ///< Particle mass [GeV/c2] double fCharge = undef::Double; ///< Particle charge [e] @@ -364,7 +358,7 @@ namespace ca::tools L1DetectorID fDetectorId; ///< Detector ID of MC point // TODO: SZh 17.05.2023: Check, if there are more then one index can be added - L1Vector<int> fvHitIndexes {"ca::tools::MCPoint::fvHitIndexes"}; ///< Indexes of hits, assigned to this point + ca::Vector<int> fvHitIndexes {"ca::tools::MCPoint::fvHitIndexes"}; ///< Indexes of hits, assigned to this point }; } // namespace ca::tools diff --git a/reco/L1/catools/CaToolsMCTrack.cxx b/reco/L1/catools/CaToolsMCTrack.cxx index 4aed3ac5807e61a6112770f8480fffd3a9043263..36afa3917273dd1da032546b5e196121b5695fa2 100644 --- a/reco/L1/catools/CaToolsMCTrack.cxx +++ b/reco/L1/catools/CaToolsMCTrack.cxx @@ -1,4 +1,4 @@ -/* Copyright (C) 2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt +/* Copyright (C) 2022-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt SPDX-License-Identifier: GPL-3.0-only Authors: Sergey Gorbunov, Sergei Zharko [committer] */ @@ -16,7 +16,7 @@ #include "CaToolsMCPoint.h" -using ca::tools::MCTrack; +using cbm::ca::tools::MCTrack; // --------------------------------------------------------------------------------------------------------------------- @@ -32,7 +32,7 @@ void MCTrack::Clear() // --------------------------------------------------------------------------------------------------------------------- // -void MCTrack::InitHitsInfo(const L1Vector<CbmL1HitDebugInfo>& vHits) +void MCTrack::InitHitsInfo(const ca::Vector<CbmL1HitDebugInfo>& vHits) { // NOTE: vHits must be sorted over stations! fMaxNofHitsOnStation = 0; @@ -81,7 +81,7 @@ void MCTrack::InitHitsInfo(const L1Vector<CbmL1HitDebugInfo>& vHits) // --------------------------------------------------------------------------------------------------------------------- // -void MCTrack::InitPointsInfo(const L1Vector<MCPoint>& vPoints) +void MCTrack::InitPointsInfo(const ca::Vector<MCPoint>& vPoints) { fMaxNofPointsOnStation = 0; fMaxNofPointsOnSensor = 0; diff --git a/reco/L1/catools/CaToolsMCTrack.h b/reco/L1/catools/CaToolsMCTrack.h index 1b50bb7f1b6e15a2c1d4b235ae3f2dd5f3061ade..2db0454b23c9babdcda7d808688016226a5eb3d4 100644 --- a/reco/L1/catools/CaToolsMCTrack.h +++ b/reco/L1/catools/CaToolsMCTrack.h @@ -16,16 +16,14 @@ #include <functional> -#include "CaConstants.h" +#include "CaToolsDef.h" #include "CaToolsLinkKey.h" -#include "L1Vector.h" +#include "CaVector.h" class CbmL1HitDebugInfo; -namespace ca::tools +namespace cbm::ca::tools { - using namespace cbm::algo::ca::constants; - class MCPoint; class MCTrack { public: @@ -220,7 +218,7 @@ namespace ca::tools /// #2) Maximal number of hits within one station /// #3) Number of consecutive stations with a hit in MC track /// \param vHits Vector of hits for a given TS - void InitHitsInfo(const L1Vector<CbmL1HitDebugInfo>& vHits); + void InitHitsInfo(const ca::Vector<CbmL1HitDebugInfo>& vHits); /// \brief Initializes information about MC track points arrangement within stations /// Defines: @@ -229,7 +227,7 @@ namespace ca::tools /// #3) Maximal number of points within one sensor (with same z-position) /// #4) Number of consecutive stations with a point in MC track /// \param vPoints Vector of points for a given TS - void InitPointsInfo(const L1Vector<MCPoint>& vPoints); + void InitPointsInfo(const ca::Vector<MCPoint>& vPoints); // ******************* @@ -363,11 +361,11 @@ namespace ca::tools int fMaxNofPointsOnSensor = undef::Int; ///< Max number of MC points with same Z (means on same sensor) int fMaxNofHitsOnStation = undef::Int; ///< Max number of hits on a station - L1Vector<int> fvPointIndexes = {"ca::tools::fvPointIndexes"}; ///< Indexes of MC points in ext.container - L1Vector<int> fvHitIndexes = {"ca::tools::fvHitIndexes"}; ///< Indexes of hits in int.container - L1Vector<int> fvRecoTrackIndexes = {"ca::tools::fvRecoTrackIndexes"}; ///< Indexes of associated reco tracks - L1Vector<int> fvTouchTrackIndexes = {"ca::tools::fvTouchTrackIndexes"}; ///< Pointers to non-associated tracks, - ///< which use hits from this track + ca::Vector<int> fvPointIndexes = {"ca::tools::fvPointIndexes"}; ///< Indexes of MC points in ext.container + ca::Vector<int> fvHitIndexes = {"ca::tools::fvHitIndexes"}; ///< Indexes of hits in int.container + ca::Vector<int> fvRecoTrackIndexes = {"ca::tools::fvRecoTrackIndexes"}; ///< Indexes of associated reco tracks + ca::Vector<int> fvTouchTrackIndexes = {"ca::tools::fvTouchTrackIndexes"}; ///< Pointers to non-associated tracks, + ///< which use hits from this track }; } // namespace ca::tools diff --git a/reco/L1/catools/CaToolsMaterialHelper.cxx b/reco/L1/catools/CaToolsMaterialHelper.cxx index 5c2577efcc402b52068433ca7d9d0ed7f1c84636..c14ac8f87686c5414b2016ad34bab9385904d41f 100644 --- a/reco/L1/catools/CaToolsMaterialHelper.cxx +++ b/reco/L1/catools/CaToolsMaterialHelper.cxx @@ -16,7 +16,8 @@ #include <iostream> #include <thread> -using namespace ca::tools; + +using cbm::ca::tools::MaterialHelper; ClassImp(MaterialHelper); diff --git a/reco/L1/catools/CaToolsMaterialHelper.h b/reco/L1/catools/CaToolsMaterialHelper.h index 7105b0c996112fef57e0cbd81b958b0f56bd469b..99a0d7ee6a5bed307f7093df7720dfd5410598aa 100644 --- a/reco/L1/catools/CaToolsMaterialHelper.h +++ b/reco/L1/catools/CaToolsMaterialHelper.h @@ -19,79 +19,76 @@ class TGeoNavigator; -namespace ca +namespace cbm::ca::tools { - namespace tools - { - /// Class ca::tools::Debugger helps to debug CA tracking + /// Class ca::tools::Debugger helps to debug CA tracking + /// + /// + /// class to create L1Material material maps form the ROOT geometry + /// + class MaterialHelper : public TObject { + public: + MaterialHelper(); + ~MaterialHelper(); + + /// + /// project rays radially from the targetZ througth the XY-bins at a reference z. + /// + void SetDoRadialProjection(double targetZ) + { + fDoRadialProjection = true; + fTargetZ = targetZ; + } + + /// + /// project rays horisontally along the Z axis (default) + /// + void SetDoHorisontalProjection() { fDoRadialProjection = false; } + + /// + /// shoot nRaysDim * nRaysDim rays for each bin in the map + /// + void SetNraysPerDim(int nRaysDim) { fNraysBinPerDim = (nRaysDim > 0) ? nRaysDim : 1; } + + + /// generate a material map, collecting the material between [zMin, zMax] + /// with radial rays from (0,0,targetZ) througth the XY-bins at z == zRef. + /// + /// It creates a map with [nBinsDim x nBinsDim] bins, of a size of [+-xyMax, +-xyMax] + /// shooting fNraysBinPerDim x fNraysBinPerDim through each bin /// + /// The calculated radiation thickness is a projection of the rad.thickness along the ray onto the Z axis. + /// RadThick = sum of (dZ / radLength) over ray trajectory pieces /// - /// class to create L1Material material maps form the ROOT geometry + /// When doRadialProjection==false the rays are shoot horisontally along the Z axis /// - class MaterialHelper : public TObject { - public: - MaterialHelper(); - ~MaterialHelper(); - - /// - /// project rays radially from the targetZ througth the XY-bins at a reference z. - /// - void SetDoRadialProjection(double targetZ) - { - fDoRadialProjection = true; - fTargetZ = targetZ; - } - - /// - /// project rays horisontally along the Z axis (default) - /// - void SetDoHorisontalProjection() { fDoRadialProjection = false; } - - /// - /// shoot nRaysDim * nRaysDim rays for each bin in the map - /// - void SetNraysPerDim(int nRaysDim) { fNraysBinPerDim = (nRaysDim > 0) ? nRaysDim : 1; } - - - /// generate a material map, collecting the material between [zMin, zMax] - /// with radial rays from (0,0,targetZ) througth the XY-bins at z == zRef. - /// - /// It creates a map with [nBinsDim x nBinsDim] bins, of a size of [+-xyMax, +-xyMax] - /// shooting fNraysBinPerDim x fNraysBinPerDim through each bin - /// - /// The calculated radiation thickness is a projection of the rad.thickness along the ray onto the Z axis. - /// RadThick = sum of (dZ / radLength) over ray trajectory pieces - /// - /// When doRadialProjection==false the rays are shoot horisontally along the Z axis - /// - L1Material GenerateMaterialMap(double zRef, double zMin, double zMax, double xyMax, int nBinsDim); - - void SetSafeMaterialInitialization(bool val = true) { fDoSafeInitialization = val; } - - private: - /// initialise the necessary amount of threads in TGeoManager - void InitThreads(); - - /// Clean up TGeoManager: threadIds, create a defailt navigator - void CleanUpThreads(); - - /// Get navigator for current thread, create it if it doesnt exist. - /// Produce an error when anything goes wrong - TGeoNavigator* GetCurrentNavigator(int iThread); - - private: - int fNthreadsOld {0}; // number of threads in TGeoManager before the helper was created - int fNthreads {0}; // number of threads - bool fDoRadialProjection {false}; // project rays horisontally along the Z axis (special mode) - int fNraysBinPerDim {3}; // shoot fNraysBinPerDim * fNraysBinPerDim rays in each map bin - double fTargetZ {0.}; // z of the target for the radial projection - std::vector<TGeoNavigator*> fNavigators {}; // list of created navigators - bool fDoSafeInitialization {false}; // perform slow but safe initialization - // to get around the crashes in TGeoVoxelFinder - ClassDef(MaterialHelper, 0); - }; - - } // namespace tools + L1Material GenerateMaterialMap(double zRef, double zMin, double zMax, double xyMax, int nBinsDim); + + void SetSafeMaterialInitialization(bool val = true) { fDoSafeInitialization = val; } + + private: + /// initialise the necessary amount of threads in TGeoManager + void InitThreads(); + + /// Clean up TGeoManager: threadIds, create a defailt navigator + void CleanUpThreads(); + + /// Get navigator for current thread, create it if it doesnt exist. + /// Produce an error when anything goes wrong + TGeoNavigator* GetCurrentNavigator(int iThread); + + private: + int fNthreadsOld {0}; // number of threads in TGeoManager before the helper was created + int fNthreads {0}; // number of threads + bool fDoRadialProjection {false}; // project rays horisontally along the Z axis (special mode) + int fNraysBinPerDim {3}; // shoot fNraysBinPerDim * fNraysBinPerDim rays in each map bin + double fTargetZ {0.}; // z of the target for the radial projection + std::vector<TGeoNavigator*> fNavigators {}; // list of created navigators + bool fDoSafeInitialization {false}; // perform slow but safe initialization + // to get around the crashes in TGeoVoxelFinder + ClassDef(MaterialHelper, 0); + }; + } // namespace ca diff --git a/reco/L1/catools/CaToolsQa.cxx b/reco/L1/catools/CaToolsQa.cxx deleted file mode 100644 index 04fc6bd33012b046d3c25a3842248e18ca543efa..0000000000000000000000000000000000000000 --- a/reco/L1/catools/CaToolsQa.cxx +++ /dev/null @@ -1,229 +0,0 @@ -/* Copyright (C) 2022-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt - SPDX-License-Identifier: GPL-3.0-only - Authors: Sergei Zharko [committer] */ - -/// \file CaToolsQa.h -/// \brief Tracking performance class (implementation) -/// \since 23.09.2022 -/// \author S.Zharko <s.zharko@gsi.de> - -#include "CaToolsQa.h" - -#include "CbmL1Hit.h" // TODO: SZh: move CbmL1Hit to ca::tools::Hit -#include "CbmL1Track.h" // TODO: SZh: move CbmL1Track to ca::tools::RecoTrack - -#include "FairRunAna.h" - -#include <utility> // for std::move - -#include "CaToolsMCData.h" -#include "L1Parameters.h" - -using namespace ca::tools; - -// --------------------------------------------------------------------------------------------------------------------- -// -void Qa::InitHistograms() -{ - LOG_IF(fatal, !fpOutFile) << "CA QA: Output file was not defined"; - - TDirectory* currentDir = gDirectory; - gDirectory = fpHistoDir; - - // ----- Initialize histograms and profiles -------------------------------------------------------------------------- - // Distributions of reconstructed tracks by different quantities - fph_reco_p = MakeHisto<TH1F>("ca_reco_p", "", kNbinsRecoP, 0., kAxMaxRecoP); - fph_reco_pt = MakeHisto<TH1F>("ca_reco_pt", "", kNbinsRecoPt, 0., kAxMaxRecoPt); - fph_reco_phi = MakeHisto<TH1F>("ca_reco_phi", "", kNbinsRecoPhi, -kAxMaxRecoPhi, kAxMaxRecoPhi); - fph_reco_nhits = MakeHisto<TH1F>("ca_reco_nhits", "", kNbinsNofHits, -0.5, kNbinsNofHits - 0.5); - fph_reco_fsta = MakeHisto<TH1F>("ca_reco_fsta", "", kNbinsNofSta, -0.5, kNbinsNofSta - 0.5); - fph_reco_purity = MakeHisto<TH1F>("ca_reco_purity", "", kNbinsPurity, -0.5, 100.5); - fph_reco_chi2_ndf = MakeHisto<TH1F>("ca_reco_chi2_ndf", "", kNbinsChi2NDF, -0.5, kAxMaxChi2NDF); - fph_reco_prob = MakeHisto<TH1F>("ca_reco_prob", "", kNbinsProb, 0., 1.01); - fph_rest_chi2_ndf = MakeHisto<TH1F>("ca_rest_chi2_ndf", "", kNbinsChi2NDF, -0.5, kAxMaxChi2NDF); - fph_rest_prob = MakeHisto<TH1F>("ca_rest_prob", "", kNbinsProb, 0., 1.01); - - fph_reco_p->SetTitle("Momentum of reconstructed track;(p/q)_{reco} [GeV/ec]"); - fph_reco_pt->SetTitle("Transverse momentum of reconstructed track;(p_{T}/q)_{reco} [GeV/ec]"); - fph_reco_phi->SetTitle("Azimuthal angle of reconstructed track;#phi_{reco} [rad]"); - fph_reco_nhits->SetTitle("Number of hits in a reconstructed track;N_{hits}"); - fph_reco_fsta->SetTitle("First station of reconstructed track;ID_{f.st.}"); - fph_reco_purity->SetTitle("Purity of reconstructed track;purity [%]"); - fph_reco_chi2_ndf->SetTitle("Fit #chi^{2}/NDF of reconstructed track;(#chi^{2}/NDF)_{reco}"); - fph_reco_chi2_ndf->SetTitle("Fit #chi^{2}/NDF of the rest track;(#chi^{2}/NDF)_{reco}"); - fph_reco_prob->SetTitle("Fit probability of reconstructed track;Prob_{reco}"); - fph_rest_prob->SetTitle("Fit probability of the rest tracks;Prob_{reco}"); - - - // Distributions of reconstructed tracks by different quantities - fph_ghost_p = MakeHisto<TH1F>("ca_ghost_p", "", kNbinsRecoP, 0., kAxMaxRecoP); - fph_ghost_pt = MakeHisto<TH1F>("ca_ghost_pt", "", kNbinsRecoPt, 0., kAxMaxRecoPt); - fph_ghost_phi = MakeHisto<TH1F>("ca_ghost_phi", "", kNbinsRecoPhi, -kAxMaxRecoPhi, kAxMaxRecoPhi); - fph_ghost_nhits = MakeHisto<TH1F>("ca_ghost_nhits", "", kNbinsNofHits, -0.5, kNbinsNofHits - 0.5); - fph_ghost_fsta = MakeHisto<TH1F>("ca_ghost_fsta", "", kNbinsNofSta, -0.5, kNbinsNofSta - 0.5); - fph_ghost_purity = MakeHisto<TH1F>("ca_ghost_purity", "", kNbinsPurity, -0.5, 100.5); - fph_ghost_chi2_ndf = MakeHisto<TH1F>("ca_ghost_chi2_ndf", "", kNbinsChi2NDF, 0., kAxMaxChi2NDF); - fph_ghost_prob = MakeHisto<TH1F>("ca_ghost_prob", "", kNbinsProb, 0., 1.01); - fph_ghost_tx = MakeHisto<TH1F>("ca_ghost_tx", "", kNbinsTx, -2., 2.); - fph_ghost_ty = MakeHisto<TH1F>("ca_ghost_ty", "", kNbinsTy, -2., 2.); - fph_ghost_fhitR = MakeHisto<TH1F>("ca_ghost_fhitR", "", kNbinsHitDist, 0., 150.); - fph_ghost_nhits_vs_p = - MakeHisto<TH2F>("ca_ghost_nhits_vs_p", "", kNbinsRecoP, 0., kAxMaxRecoP, kNbinsNofHits, -0.5, kNbinsNofHits - 0.5); - fph_ghost_fsta_vs_p = - MakeHisto<TH2F>("ca_ghost_fsta_vs_p", "", kNbinsRecoP, 0., kAxMaxRecoP, kNbinsNofSta, -0.5, kNbinsNofSta - 0.5); - fph_ghost_lsta_vs_fsta = MakeHisto<TH2F>("ca_ghost_lsta_vs_fsta", "", kNbinsNofSta, -0.5, kNbinsNofSta - 0.5, - kNbinsNofSta, -0.5, kNbinsNofSta - 0.5); - - fph_ghost_p->SetTitle("Momentum of ghost track;(p/q)_{reco} [GeV/ec]"); - fph_ghost_pt->SetTitle("Transverse momentum of ghost track;(p_{T}/q)_{reco} [GeV/ec]"); - fph_ghost_phi->SetTitle("Azimuthal angle of ghost track;#phi_{reco} [rad]"); - fph_ghost_nhits->SetTitle("Number of hits in a ghost track;N_{hits}"); - fph_ghost_fsta->SetTitle("First station of ghost track;ID_{f.st.}"); - fph_ghost_purity->SetTitle("Purity of ghost track;purity [%]"); - fph_ghost_chi2_ndf->SetTitle("Fit #chi^{2}/NDF of ghost track;(#chi^{2}/NDF)_{reco}"); - fph_ghost_prob->SetTitle("Fit probability of ghost track;Prob_{reco}"); - fph_ghost_tx->SetTitle("Slope along x-axis of ghost tracks;t_{x}"); - fph_ghost_ty->SetTitle("Slope along y-axis of ghost tracks;t_{y}"); - fph_ghost_fhitR->SetTitle("Distance of the first hit from z-axis;R [cm]"); - fph_ghost_nhits_vs_p->SetTitle("Number of hits vs. momentum for ghost track;(p/q)_{reco} [GeV/ec];N_{hits}"); - fph_ghost_fsta_vs_p->SetTitle("First station vs. momentum for ghost track;(p/q)_{reco} [GeV/ec];ID_{f.st.}"); - fph_ghost_lsta_vs_fsta->SetTitle("Last station vs. first station for ghost track;ID_{first st.};ID_{last st.}"); - - // ----- Reset event counter - fNofEvents = 0; - - LOG(info) << "CA QA: Initialized histograms: (directory \033[1;31m" << gDirectory->GetName() << "\033[0m)"; - for (auto& entry : fmpHistoList) { - LOG(info) << "\t- " << entry.first; - } - - gDirectory = currentDir; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -void Qa::FillHistograms() -{ - assert(fpParameters); - assert(fpMCData); - assert(fpvRecoTracks); - assert(fpvHits); - - // Update event counter - fNofEvents++; - - - // ************************************************************************ - // ** Fill distributions for reconstructed tracks (including ghost ones) ** - // ************************************************************************ - - //for (const auto& recoTrk : *fpvRecoTracks) { - // // Reject tracks, which do not contain hits - // if (recoTrk.GetNofHits() < 1) { continue; } - - // const auto& hitFst = (*fpvHits)[recoTrk.GetHitIndexes()[0]]; // first hit of track - // const auto& hitSnd = (*fpvHits)[recoTrk.GetHitIndexes()[1]]; // second hit of track - // const auto& hitLst = (*fpvHits)[recoTrk.GetHitIndexes()[recoTrk.GetNofHits() - 1]]; // last hit of track - - // fph_reco_p->Fill(recoTrk.GetP()); - // fph_reco_pt->Fill(recoTrk.GetPt()); - // fph_reco_phi->Fill(recoTrk.GetPhi()); - // fph_reco_nhits->Fill(recoTrk.GetNofHits()); - // //fph_reco_fsta->Fill(hitFst.GetStationId()); - // fph_reco_purity->Fill(100 * recoTrk.GetMaxPurity()); - - // if (recoTrk.GetNDF() > 0) { - // if (recoTrk.IsGhost()) { - // fph_ghost_chi2_ndf->Fill(recoTrk.GetChiSq() / recoTrk.GetNDF()); - // fph_ghost_prob->Fill(recoTrk.GetProb()); - // } - // else { - // int iTmc = recoTrk.GetMCTrackIndexes()[0]; // Index of first MC-track - // const auto& mcTrk = fpMCData->GetTrack(iTmc); - // // NOTE: reconstructed tracks usually have reference to only one MC-track - // if (mcTrk.IsReconstructable()) { - // fph_reco_chi2_ndf->Fill(recoTrk.GetChiSq() / recoTrk.GetNDF()); - // fph_reco_prob->Fill(recoTrk.GetProb()); - // } - // else { - // fph_rest_chi2_ndf->Fill(recoTrk.GetChiSq() / recoTrk.GetNDF()); - // fph_rest_prob->Fill(recoTrk.GetProb()); - // } - // } - // } // recoTrk.GetNDF() > 0: end - - // if (recoTrk.IsGhost()) { - // fph_ghost_purity->Fill(100 * recoTrk.GetMaxPurity()); - // fph_ghost_p->Fill(recoTrk.GetP()); - // fph_ghost_pt->Fill(recoTrk.GetPt()); - // fph_ghost_nhits->Fill(recoTrk.GetNofHits()); - // fph_ghost_fsta->Fill(hitFst.GetStationId()); - // fph_ghost_fhitR->Fill(hitFst.GetDistFromBP()); - - // // z-positions of the first and second hit stations - // double z0 = fpParameters->GetStation(hitFst.GetStationId()).fZ[0]; - // double z1 = fpParameters->GetStation(hitSnd.GetStationId()).fZ[0]; - - // if (fabs(z1 - z0) > 1.e-4) { // test for z1 != z2 - // fph_ghost_tx->Fill((hitSnd.GetX() - hitFst.GetX()) / (z1 - z0)); - // fph_ghost_ty->Fill((hitSnd.GetY() - hitFst.GetY()) / (z1 - z0)); - // } - - // fph_ghost_nhits_vs_p->Fill(recoTrk.GetP(), recoTrk.GetNofHits()); - // fph_ghost_fsta_vs_p->Fill(recoTrk.GetP(), hitFst.GetStationId()); - // if (hitFst.GetStationId() <= hitLst.GetStationId()) { - // fph_ghost_lsta_vs_fsta->Fill(hitFst.GetStationId(), hitLst.GetStationId()); - // } - // } // recoTrk.IsGhost(): end - //} // loop over recoTrk: end - - - // ************************************* - // ** Fill distributions of MC-tracks ** - // ************************************* - - //for (const auto& mcTrk: fpMCData->GetTrackContainer()) { - // - //} -} - -// --------------------------------------------------------------------------------------------------------------------- -// -void Qa::SaveHistograms() -{ - fpOutFile->cd(); - RecursiveWrite(fpHistoDir); - LOG(info) << "CA QA: Histograms have been written to \033[1;31m" << fpOutFile->GetName() << "\033[0m"; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -void Qa::RecursiveWrite(TObject* pObj) -{ - if (!pObj->IsFolder()) { pObj->Write(); } - else { - TDirectory* pObjDir = dynamic_cast<TDirectory*>(pObj); - LOG_IF(fatal, !pObjDir) << "CA QA: Incorrect dynamic cast from the input pointer to RecursiveWrite function"; - - TDirectory* pCurDir = gDirectory; - TDirectory* pSubDir = pCurDir->mkdir(pObjDir->GetName()); - pSubDir->cd(); - TList* pObjList = pObjDir->GetList(); - TIter it(pObjList); - while (TObject* pNext = it()) { - RecursiveWrite(pNext); - } - pCurDir->cd(); - } -} - -// --------------------------------------------------------------------------------------------------------------------- -// -void Qa::SetOutputFilename(const char* filename) -{ - fpOutFile = std::make_unique<TFile>(filename, "RECREATE"); - if (!fpOutFile.get() || fpOutFile->IsZombie()) { - LOG(fatal) << "CA QA: ROOT histogram output file \"" << filename << "\" cannot be created"; - } - fpHistoDir = gROOT->mkdir("histograms"); -} diff --git a/reco/L1/catools/CaToolsQa.h b/reco/L1/catools/CaToolsQa.h deleted file mode 100644 index e6c05fed338fca5e1049cf40d948729a8ea3b41e..0000000000000000000000000000000000000000 --- a/reco/L1/catools/CaToolsQa.h +++ /dev/null @@ -1,248 +0,0 @@ -/* Copyright (C) 2022-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt - SPDX-License-Identifier: GPL-3.0-only - Authors: Sergei Zharko [committer] */ - -/// \file CaToolsQa.h -/// \brief Tracking performance class (header) -/// \since 23.09.2022 -/// \author S.Zharko <s.zharko@gsi.de> - -#ifndef CaToolsQa_h -#define CaToolsQa_h 1 - -#include "Logger.h" - -#include "TFile.h" -#include "TH1F.h" -#include "TH2F.h" -#include "TProfile.h" -#include "TROOT.h" - -#include <string_view> - -#include "L1Vector.h" - -class L1Parameters; -class TDirectory; -class CbmL1HitDebugInfo; -class CbmL1Track; - -namespace ca::tools -{ - class MCData; - - /// Class ca::tools::Qa defines CA tracking internal performance measurement - /// - class Qa { - public: - // ***************************************** - // ** Constructors and destructor ** - // ***************************************** - - /// Default constructor - Qa() = default; - - /// Destructor - ~Qa() = default; - - /// Copy constructor - Qa(const Qa& other) = delete; - - /// Move constructor - Qa(Qa&& other) = delete; - - /// Copy assignment operator - Qa& operator=(const Qa& other) = delete; - - /// Move assignment operator - Qa& operator=(Qa&& other) = delete; - - /// Initializes histograms - /// \note should be called in the beginning of the run - void InitHistograms(); - - /// Fills histograms - /// \note should be called in the end of the event - void FillHistograms(); - - /// Saves histograms to the file - /// \note Should be called in the end of the run - void SaveHistograms(); - - // *********************** - // ** Accessors ** - // *********************** - - /// Sets pointer to hit container from reference - void SetHitContainer(const L1Vector<CbmL1HitDebugInfo>& vHits) { fpvHits = &vHits; } - - /// Sets pointer to hit container from pointer - void SetHitContainer(const L1Vector<CbmL1HitDebugInfo>* pvHits) - { - assert(pvHits); - fpvHits = pvHits; - } - - - /// Sets pointer to MC data object - void SetMCData(const MCData* pMCData) { fpMCData = pMCData; } - - /// Sets output filename - void SetOutputFilename(const char* filename); - - /// Sets pointer to tracking parameters object - void SetParameters(const L1Parameters* pParameters) { fpParameters = pParameters; } - - /// Sets pointer to the vector of reconstructed tracks from reference - void SetRecoTrackContainer(const L1Vector<CbmL1Track>& vRecoTracks) { fpvRecoTracks = &vRecoTracks; } - - /// Sets pointer to the vector of reconstructed tracks from pointer - void SetRecoTrackContainer(const L1Vector<CbmL1Track>* pvRecoTracks) - { - assert(pvRecoTracks); - fpvRecoTracks = pvRecoTracks; - } - - - private: - /// Function, which creates a histogram and initializes it - /// \tparam T Type of the histogram (\note T should not be a pointer) - /// \param name Name of histogram - /// \param args The rest of the arguments, which will be passed to the histogram constructor - template<typename T, typename... Args> - T* MakeHisto(const char* name, Args... args); - - /// Writes recursively all registered objects inside a current directory - /// \param pObj Pointer to object to write - void RecursiveWrite(TObject* pObj); - - const L1Parameters* fpParameters = nullptr; ///< Instance of the tracking core class parameters - const MCData* fpMCData = nullptr; ///< Container for MC information on the event - const L1Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to external reconstructed track container - const L1Vector<CbmL1HitDebugInfo>* fpvHits = nullptr; ///< Pointer to external hit container - - // *********************** - // ** Histograms ** - // *********************** - - TDirectory* fpHistoDir = nullptr; ///< directory to store QA histograms - - std::unique_ptr<TFile> fpOutFile = nullptr; - - // ** Number of bins for different quantities - - static constexpr int kNbinsRecoP = 50; ///< momentum of reconstructed tracks - static constexpr int kNbinsRecoPt = 50; ///< transverse momentum of reconstructed tracks - static constexpr int kNbinsRecoPhi = 68; ///< azimuthal angle of reconstructed tracks - static constexpr int kNbinsNofHits = 10; ///< number of hits of reconstructed tracks - static constexpr int kNbinsPurity = 100; ///< track purity - static constexpr int kNbinsChi2NDF = 50; ///< fit chi-square over NDF - static constexpr int kNbinsProb = 505; ///< fit probability - static constexpr int kNbinsNofSta = 15; ///< number of stations - static constexpr int kNbinsTx = 50; ///< slope along x-axis - static constexpr int kNbinsTy = 50; ///< slope along y-axis - static constexpr int kNbinsHitDist = 50; ///< distance from transverse hit position to z-axis - - // ** Upper bound for quantity range for histograms - static constexpr double kAxMaxRecoP = 5.; ///< momentum of reconstructed tracks [GeV/c] - static constexpr double kAxMaxRecoPt = 5.; ///< transverse momentum of reconstructed tracks [GeV/c] - static constexpr double kAxMaxRecoPhi = 3.2; ///< azimuthal angle of reconstructed tracks [rad] - static constexpr double kAxMaxChi2NDF = 100.; ///< fit chi-square over NDF - - - // ** Pointers to histogram objects ** - - // Reconstructed tracks - TH1F* fph_reco_p = nullptr; ///< Total momentum over charge of reconstructed tracks - TH1F* fph_reco_pt = nullptr; ///< Transverse momentum over charge of reconstructed tracks - TH1F* fph_reco_phi = nullptr; ///< Azimuthal angle of reconstructed tracks - TH1F* fph_reco_nhits = nullptr; ///< Hit number of reconstructed tracks - TH1F* fph_reco_fsta = nullptr; ///< First station index of reconstructed tracks - TH1F* fph_reco_purity = nullptr; ///< Purity of reconstructed tracks (\note purity requires MC information) - TH1F* fph_reco_chi2_ndf = nullptr; ///< Fit chi2 over NDF of reconstructed tracks - TH1F* fph_reco_prob = nullptr; ///< Fit probability of reconstructed tracks - TH1F* fph_rest_chi2_ndf = nullptr; ///< Fit chi2 over NDF of non-reconstructable tracks - TH1F* fph_rest_prob = nullptr; ///< Fit probability of non-reconstructable tracks - - // Ghost tracks - TH1F* fph_ghost_p = nullptr; ///< Total momentum over charge of ghost tracks - TH1F* fph_ghost_pt = nullptr; ///< Transverse momentum over charge of ghost tracks - TH1F* fph_ghost_phi = nullptr; ///< Azimuthal angle of ghost tracks - TH1F* fph_ghost_nhits = nullptr; ///< Hit number of ghost tracks - TH1F* fph_ghost_fsta = nullptr; ///< First station index of ghost tracks - TH1F* fph_ghost_purity = nullptr; ///< Purity of ghost tracks - TH1F* fph_ghost_chi2_ndf = nullptr; ///< Fit chi2 over NDF of ghost tracks - TH1F* fph_ghost_prob = nullptr; ///< Fit probability of ghost tracks - TH1F* fph_ghost_tx = nullptr; ///< Slope along x-axis of ghost tracks - TH1F* fph_ghost_ty = nullptr; ///< Slope along y-axis of ghost tracks - TH1F* fph_ghost_fhitR = nullptr; ///< Distance of the first hit from z-axis for ghost tracks - TH2F* fph_ghost_nhits_vs_p = nullptr; ///< Hit number vs. total momentum over charge of ghost tracks - TH2F* fph_ghost_fsta_vs_p = nullptr; ///< First station index vs. total momentum over charge for ghost tracks - TH2F* fph_ghost_lsta_vs_fsta = nullptr; ///< Last station index vs. first station index of ghost tracks - - // TODO: Add other histograms - - // Registered primary tracks - // Registered secondary tracks - // Reconstructable primary tracks - // Reconstructable secondary tracks - // Reconstructed primary tracks - // Reconstructed secondary tracks - // Reconstructed tracks - // Ghost tracks - - // ** Histograms appearance attributes - static constexpr double kHtextSize = 0.04; - - - /// A map of the histogram pointers. The key is a name of the histogram - std::unordered_map<std::string_view, TH1*> fmpHistoList; ///< List of pointers to histograms - - // ******************************* - // ** Utility variables ** - // ******************************* - - int fNofEvents = 0; ///< Event counter - }; -} // namespace ca::tools - -// ********************************************************** -// ** Inline and template function implementations ** -// ********************************************************** - -// --------------------------------------------------------------------------------------------------------------------- -// -template<typename T, typename... Args> -T* ca::tools::Qa::MakeHisto(const char* name, Args... args) -{ - // Check, if the histogram with a given name was already created. If so, delete it - if (gROOT->FindObject(name)) { - LOG(warn) << "CA QA: A histogram with name \"" << name << "\" was previously created and will be deleted now " - << "to avoid memory leaks"; - T* pHisto = (T*) gROOT->FindObject(name); - delete pHisto; - } - - // Create a new histogram or profile - T* pHisto = new T(name, args...); - - // Register histogram in the current TDirectory - pHisto->SetDirectory(fpHistoDir); - - // Register histogram in the list - fmpHistoList[pHisto->GetName()] = (TH1*) pHisto; - - pHisto->Sumw2(); - - // Appearance settings - std::array<TAxis*, 2> vpAxes = {pHisto->GetXaxis(), pHisto->GetYaxis()}; - for (auto* pAxis : vpAxes) { - pAxis->SetTitleSize(kHtextSize); - pAxis->SetLabelSize(kHtextSize); - } - - return pHisto; -} - - -#endif // CaToolsQa_h diff --git a/reco/L1/catools/CaToolsWFExpression.cxx b/reco/L1/catools/CaToolsWFExpression.cxx index 7f9fbd2aa39a5365b9f8578ef0991598def5ff57..8c04aff3d227419897775a7b3e0b4ab1489b81a2 100644 --- a/reco/L1/catools/CaToolsWFExpression.cxx +++ b/reco/L1/catools/CaToolsWFExpression.cxx @@ -17,7 +17,7 @@ #include <algorithm> -using namespace ca::tools; +using cbm::ca::tools::WFExpression; // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/catools/CaToolsWFExpression.h b/reco/L1/catools/CaToolsWFExpression.h index 99bd7f7e8bb94d89b8d7ffcbb42b4831925763af..98773ccee450032653e7cd39de495a10fba60d49 100644 --- a/reco/L1/catools/CaToolsWFExpression.h +++ b/reco/L1/catools/CaToolsWFExpression.h @@ -12,120 +12,119 @@ #include <tuple> #include <vector> +#include "CaToolsDef.h" + class TTree; class TH2F; class TPad; -namespace ca +namespace cbm::ca::tools { - namespace tools - { - /// A helper class for ca::tools::WindowFinder to handle a particular expression (dx vs. x0 etc.) and all related - /// methods to work with it. - /// \note DISTANCE (expression, axis) -- dx or dy, signed distance between the extrapolated point and real MC-point - /// PARAMETER (expression, axis) -- x0 or y0, a parameter, vs. which the distance is studied - class WFExpression { - public: - static constexpr int kNpars = 1; /// TMP: number of parameters - - /// Default constructor - WFExpression() = delete; - - /// Constructor - /// \param pChain A pointer to a tree with MC triplets - /// \param exprDist Expression for distance - /// \param exprParam Expression for parameter - WFExpression(TTree* pTree, const char* exprDist, const char* exprParam); - - /// Copy constructor - WFExpression(const WFExpression& other) = delete; - - /// Move constructor - WFExpression(WFExpression&& other) = delete; - - /// Destructor - ~WFExpression() = default; - - /// Copy assignment operator - WFExpression& operator=(const WFExpression& other) = delete; - - /// Move assignment operator - WFExpression& operator=(WFExpression&& other) = delete; - - /// Calculates parameters - /// \param pTree Pointer to a tree with MC triplets - /// \return A tuple: - /// - vector of parameters for upper bound - /// - vector of parameters for lower bound - std::tuple<std::vector<float>, std::vector<float>> CalculateParameters(); - - /// Sets cut, including information on the station and track finder iteration - void SetCut(const TCut& cut) { fCut = cut; } - - /// Sets fraction of the events, which can be left out of boundaries - void SetEps(float eps); - - /// Sets number of slices - void SetNslices(int nSlices); - - /// Sets number of bins - /// \param nBinsDist Number of bins along the distance axis - /// \param nBinsParam Number of bins along the parameter axis - void SetNbins(int nBinsDist, int nBinsParam); - - /// Sets base pad pointer - void SetPadBase(TPad* pad); - - /// Sets slices pad pointer - void SetPadSlices(TPad* pad); - - /// Sets title of the histograms - void SetTitle(const char* title) { fsTitle = title; } - - private: - // ***************************** - // ** Private class functions ** - // ***************************** - - /// Process slice - /// \param iBinMax Max bin of the parameter axis to start a projection (date from the bin are included) - /// \return Tuple: lower bound, upper bound, slice center - std::tuple<float, float, float> ProcessSlice(int iBinMin, int iBinMax); - - /// Gets window parameterisations assuming there is no dependence from parameter - void GetConstWindowParams(); - // TODO: use other functions for other window shapes: GetParabWindowParams, GetEllipticWindowParams etc. - - - // ********************* - // ** Class variables ** - // ********************* - TTree* fpTree = nullptr; ///< Tree to be analyzed - TH2F* fpHistBase = nullptr; ///< Base histogram (distance vs. parameter (x0 or y0)) - - TString fsExprDist = ""; ///< Expression along the distance axis - TString fsExprParam = ""; ///< Expression along the parameter axis - - TCut fCut = ""; ///< Cut used to draw and expression - int fNslices = 8; ///< Number of slices along the parameter axis - float fEps = 0.0005; ///< A fraction of triplets, which can be lost - - - std::vector<float> fvUpSBoundaries = std::vector<float>(fNslices); ///< Upper boundaries for diff. slices - std::vector<float> fvLoSBoundaries = std::vector<float>(fNslices); ///< Lower boundaries for diff. slices - std::vector<float> fvSCenters = std::vector<float>(fNslices); ///< Slice centers - std::vector<float> fvUpParams = std::vector<float>(kNpars); ///< Parameters for max - std::vector<float> fvLoParams = std::vector<float>(kNpars); ///< Parameters for min - - // ----- Plotting options - int fNbinsParam = 400; ///< Number of bins along the parameter axis - int fNbinsDist = 400; ///< Number of bins along the distance axis - TString fsTitle = ""; ///< Title of expression - TString fsName = ""; ///< Name of the expression (expr + cut, must be unique!!, TODO: make a check) - TPad* fpPadBase = nullptr; ///< Pointer to a pad for base histogram - TPad* fpPadSlices = nullptr; ///< Pointer to a pad for slices - }; - } // namespace tools + /// A helper class for ca::tools::WindowFinder to handle a particular expression (dx vs. x0 etc.) and all related + /// methods to work with it. + /// \note DISTANCE (expression, axis) -- dx or dy, signed distance between the extrapolated point and real MC-point + /// PARAMETER (expression, axis) -- x0 or y0, a parameter, vs. which the distance is studied + class WFExpression { + public: + static constexpr int kNpars = 1; /// TMP: number of parameters + + /// Default constructor + WFExpression() = delete; + + /// Constructor + /// \param pChain A pointer to a tree with MC triplets + /// \param exprDist Expression for distance + /// \param exprParam Expression for parameter + WFExpression(TTree* pTree, const char* exprDist, const char* exprParam); + + /// Copy constructor + WFExpression(const WFExpression& other) = delete; + + /// Move constructor + WFExpression(WFExpression&& other) = delete; + + /// Destructor + ~WFExpression() = default; + + /// Copy assignment operator + WFExpression& operator=(const WFExpression& other) = delete; + + /// Move assignment operator + WFExpression& operator=(WFExpression&& other) = delete; + + /// Calculates parameters + /// \param pTree Pointer to a tree with MC triplets + /// \return A tuple: + /// - vector of parameters for upper bound + /// - vector of parameters for lower bound + std::tuple<std::vector<float>, std::vector<float>> CalculateParameters(); + + /// Sets cut, including information on the station and track finder iteration + void SetCut(const TCut& cut) { fCut = cut; } + + /// Sets fraction of the events, which can be left out of boundaries + void SetEps(float eps); + + /// Sets number of slices + void SetNslices(int nSlices); + + /// Sets number of bins + /// \param nBinsDist Number of bins along the distance axis + /// \param nBinsParam Number of bins along the parameter axis + void SetNbins(int nBinsDist, int nBinsParam); + + /// Sets base pad pointer + void SetPadBase(TPad* pad); + + /// Sets slices pad pointer + void SetPadSlices(TPad* pad); + + /// Sets title of the histograms + void SetTitle(const char* title) { fsTitle = title; } + + private: + // ***************************** + // ** Private class functions ** + // ***************************** + + /// Process slice + /// \param iBinMax Max bin of the parameter axis to start a projection (date from the bin are included) + /// \return Tuple: lower bound, upper bound, slice center + std::tuple<float, float, float> ProcessSlice(int iBinMin, int iBinMax); + + /// Gets window parameterisations assuming there is no dependence from parameter + void GetConstWindowParams(); + // TODO: use other functions for other window shapes: GetParabWindowParams, GetEllipticWindowParams etc. + + + // ********************* + // ** Class variables ** + // ********************* + TTree* fpTree = nullptr; ///< Tree to be analyzed + TH2F* fpHistBase = nullptr; ///< Base histogram (distance vs. parameter (x0 or y0)) + + TString fsExprDist = ""; ///< Expression along the distance axis + TString fsExprParam = ""; ///< Expression along the parameter axis + + TCut fCut = ""; ///< Cut used to draw and expression + int fNslices = 8; ///< Number of slices along the parameter axis + float fEps = 0.0005; ///< A fraction of triplets, which can be lost + + + std::vector<float> fvUpSBoundaries = std::vector<float>(fNslices); ///< Upper boundaries for diff. slices + std::vector<float> fvLoSBoundaries = std::vector<float>(fNslices); ///< Lower boundaries for diff. slices + std::vector<float> fvSCenters = std::vector<float>(fNslices); ///< Slice centers + std::vector<float> fvUpParams = std::vector<float>(kNpars); ///< Parameters for max + std::vector<float> fvLoParams = std::vector<float>(kNpars); ///< Parameters for min + + // ----- Plotting options + int fNbinsParam = 400; ///< Number of bins along the parameter axis + int fNbinsDist = 400; ///< Number of bins along the distance axis + TString fsTitle = ""; ///< Title of expression + TString fsName = ""; ///< Name of the expression (expr + cut, must be unique!!, TODO: make a check) + TPad* fpPadBase = nullptr; ///< Pointer to a pad for base histogram + TPad* fpPadSlices = nullptr; ///< Pointer to a pad for slices + }; } // namespace ca #endif // CaToolsWFExpression_h diff --git a/reco/L1/catools/CaToolsWindowFinder.cxx b/reco/L1/catools/CaToolsWindowFinder.cxx index b837294204e5f71f556b3ec189ed53f724287a67..f8a5d998f09d120c277b6c33a96b435855e5925f 100644 --- a/reco/L1/catools/CaToolsWindowFinder.cxx +++ b/reco/L1/catools/CaToolsWindowFinder.cxx @@ -29,10 +29,10 @@ #include "L1ConfigRW.h" #include "L1SearchWindow.h" -using namespace ca::tools; +using namespace cbm::ca::tools; using namespace cbm::algo::ca::constants; // for colored logs -ClassImp(ca::tools::WindowFinder); +ClassImp(cbm::ca::tools::WindowFinder); // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/catools/CaToolsWindowFinder.h b/reco/L1/catools/CaToolsWindowFinder.h index 9ab2dcc0fb7b4e18f4c2bd051caa38b44a6e0ce4..aad528d84ab13bdc6e4f55afce62819b64bb9a4f 100644 --- a/reco/L1/catools/CaToolsWindowFinder.h +++ b/reco/L1/catools/CaToolsWindowFinder.h @@ -25,135 +25,132 @@ class L1SearchWindow; class TPad; class TCanvas; -namespace ca +namespace cbm::ca::tools { - namespace tools + /// Enumeration to handle processed expressions + enum EExpr { - /// Enumeration to handle processed expressions - enum EExpr - { - kDxVsX0, - kDxVsY0, - kDyVsX0, - kDyVsY0, - kEND - }; - - /// TODO: ... write an instruction ... - class WindowFinder : public TObject { - public: - // TODO: TEMPORARY CONSTANT EXPRESSIONS (TO BE MOVED TO A SEPARATE HEADER) - static constexpr const char* kTreeName = "t"; ///< Name of the input MC triplets tree - - public: - /// Default constructor - WindowFinder(); - - /// Destructor - virtual ~WindowFinder() = default; - - /// Copy and move are forbidden - WindowFinder(const WindowFinder&) = delete; - WindowFinder(WindowFinder&&) = delete; - WindowFinder& operator=(const WindowFinder&) = delete; - WindowFinder& operator=(WindowFinder&&) = delete; - - /// Adds an input file with a tree object of MC triplets - /// \note Multiple file input is possible - void AddInputFile(const char* filename); - - /// Saves canvases to a set of canvases to pdf - void DumpCanvasesToPdf(const char* filename = "WFLog") const; - - /// Process a tree (TEST) - /// \param opt Define options to process: - /// 'T' - triplets are used instead of doublets - void Process(Option_t* opt = ""); - - /// Reads the iterations from YAML config - void ReadTrackingIterationsFromYAML(const char* filename); - - /// Sets binning of the dx (dy) distribution - /// \param nBinsX Number of bins for the x-axis - /// \param nBinsY Number of bins for the y-axis - void SetBinning(int nBinsX, int nBinsY); - - /// Sets a fraction of triplets (doublets), which can be omitted by the window - void SetEpsilon(float eps); - - /// Sets additional cut (can be useful to reduce distribution contamination by outliers) - /// \param cut A cut object - void SetExtraCut(const TCut& cut) { fExtraCut = cut; } - - /// Sets number of slices along the x-axiso - /// If the number of slices larger then 1, the window bounds will be fitted with a function (...which TODO). - /// Otherwise, the bounds will be represented with constants (independent from x0 or y0) - void SetNslices(int nSlices); - - /// Sets name of output file with search windows array - /// \note Format of output file: - /// \note <number of parameters> <number of search windows stored> <array of serialized L1SearchWindow objects> - /// \param filename Name of the file - void SetOutputName(const char* filename) { fsOutputName = filename; } - - /// Define indexes of stations for which windows are needed - /// \note: A stution - /// TODO: Get number of stations from a tree ... - void SetStationIndexes(const std::vector<int>& vStationIndexes); - - /// Sets components of the target center position - /// \param x Position x-component [cm] - /// \param y Position y-component [cm] - /// \param z Position z-component [cm] - void SetTarget(double x, double y, double z); - - - private: - /// Creates a search window for a selected station and iteration - L1SearchWindow CreateSW(int iStation, const L1CAIteration& caIter); - - /// Returns expression for dx or dy to be drawn in a tree - const char* GetDistExpr(EExpr expr) const; - - /// Gets a cut for doublets/triplets defined by station and a tracking iteration - /// \param iStation Global index of an active station - /// \param caIter CA track finder iteration object - TCut GetTrackSelectionCut(int iStation, const L1CAIteration& caIter) const; - - /// Prints information on the dx and dy expression as well as used cuts on the pad - /// \param pPad A pad to print the information - /// \param iStation Global index of an active station - /// \param caIter CA track finder iteration object - void PrintCaseInformation(TPad* pPad, int iStation, const L1CAIteration& caIter) const; - - - // ********************* - // ** Class variables ** - // ********************* - int fNparams = 1; ///< number of parameters of the searching window - - std::string fsOutputName = "SearchWindows.dat"; ///< Name for output file with estimated search windows - - // ----- Input parameters (iterations and stations) - std::vector<L1CAIteration> fvCaIters = {}; ///< Tracking iterations - std::vector<int> fvStationIndexes = {}; ///< Global indexes of active stations to find the windows - std::array<double, 3> fTargetPos = {0}; ///< Target position {x, y, z} [cm] + kDxVsX0, + kDxVsY0, + kDyVsX0, + kDyVsY0, + kEND + }; + + /// TODO: ... write an instruction ... + class WindowFinder : public TObject { + public: + // TODO: TEMPORARY CONSTANT EXPRESSIONS (TO BE MOVED TO A SEPARATE HEADER) + static constexpr const char* kTreeName = "t"; ///< Name of the input MC triplets tree + + public: + /// Default constructor + WindowFinder(); + + /// Destructor + virtual ~WindowFinder() = default; + + /// Copy and move are forbidden + WindowFinder(const WindowFinder&) = delete; + WindowFinder(WindowFinder&&) = delete; + WindowFinder& operator=(const WindowFinder&) = delete; + WindowFinder& operator=(WindowFinder&&) = delete; + + /// Adds an input file with a tree object of MC triplets + /// \note Multiple file input is possible + void AddInputFile(const char* filename); + + /// Saves canvases to a set of canvases to pdf + void DumpCanvasesToPdf(const char* filename = "WFLog") const; + + /// Process a tree (TEST) + /// \param opt Define options to process: + /// 'T' - triplets are used instead of doublets + void Process(Option_t* opt = ""); + + /// Reads the iterations from YAML config + void ReadTrackingIterationsFromYAML(const char* filename); + + /// Sets binning of the dx (dy) distribution + /// \param nBinsX Number of bins for the x-axis + /// \param nBinsY Number of bins for the y-axis + void SetBinning(int nBinsX, int nBinsY); + + /// Sets a fraction of triplets (doublets), which can be omitted by the window + void SetEpsilon(float eps); + + /// Sets additional cut (can be useful to reduce distribution contamination by outliers) + /// \param cut A cut object + void SetExtraCut(const TCut& cut) { fExtraCut = cut; } + + /// Sets number of slices along the x-axiso + /// If the number of slices larger then 1, the window bounds will be fitted with a function (...which TODO). + /// Otherwise, the bounds will be represented with constants (independent from x0 or y0) + void SetNslices(int nSlices); + + /// Sets name of output file with search windows array + /// \note Format of output file: + /// \note <number of parameters> <number of search windows stored> <array of serialized L1SearchWindow objects> + /// \param filename Name of the file + void SetOutputName(const char* filename) { fsOutputName = filename; } + + /// Define indexes of stations for which windows are needed + /// \note: A stution + /// TODO: Get number of stations from a tree ... + void SetStationIndexes(const std::vector<int>& vStationIndexes); + + /// Sets components of the target center position + /// \param x Position x-component [cm] + /// \param y Position y-component [cm] + /// \param z Position z-component [cm] + void SetTarget(double x, double y, double z); + + + private: + /// Creates a search window for a selected station and iteration + L1SearchWindow CreateSW(int iStation, const L1CAIteration& caIter); + + /// Returns expression for dx or dy to be drawn in a tree + const char* GetDistExpr(EExpr expr) const; + + /// Gets a cut for doublets/triplets defined by station and a tracking iteration + /// \param iStation Global index of an active station + /// \param caIter CA track finder iteration object + TCut GetTrackSelectionCut(int iStation, const L1CAIteration& caIter) const; + + /// Prints information on the dx and dy expression as well as used cuts on the pad + /// \param pPad A pad to print the information + /// \param iStation Global index of an active station + /// \param caIter CA track finder iteration object + void PrintCaseInformation(TPad* pPad, int iStation, const L1CAIteration& caIter) const; + + + // ********************* + // ** Class variables ** + // ********************* + int fNparams = 1; ///< number of parameters of the searching window + + std::string fsOutputName = "SearchWindows.dat"; ///< Name for output file with estimated search windows + + // ----- Input parameters (iterations and stations) + std::vector<L1CAIteration> fvCaIters = {}; ///< Tracking iterations + std::vector<int> fvStationIndexes = {}; ///< Global indexes of active stations to find the windows + std::array<double, 3> fTargetPos = {0}; ///< Target position {x, y, z} [cm] - TCut fExtraCut = TCut(""); ///< Optional cut on the triplets/doublets + TCut fExtraCut = TCut(""); ///< Optional cut on the triplets/doublets - // Window extraction settings - int fNbinsX = 400; ///< Number of bins for the X axis of the distribution - int fNbinsY = 400; ///< Number of bins for the Y axis of the distribution - int fNslices = 1; ///< Number of slices along the X axis - float fEps = 0.001; ///< Fraction of triplets (doublets), which can be omitted + // Window extraction settings + int fNbinsX = 400; ///< Number of bins for the X axis of the distribution + int fNbinsY = 400; ///< Number of bins for the Y axis of the distribution + int fNslices = 1; ///< Number of slices along the X axis + float fEps = 0.001; ///< Fraction of triplets (doublets), which can be omitted - TChain* fpMcTripletsTree = nullptr; ///< Chain of trees containing MC triplets, generated in CbmL1 Performance + TChain* fpMcTripletsTree = nullptr; ///< Chain of trees containing MC triplets, generated in CbmL1 Performance - std::vector<TCanvas*> fvpCanvases; // Vector of pointer to cavnases + std::vector<TCanvas*> fvpCanvases; // Vector of pointer to cavnases - ClassDef(WindowFinder, 0); - }; - } // namespace tools + ClassDef(WindowFinder, 0); + }; } // namespace ca #endif // CaToolsWindowsFinder_h diff --git a/reco/L1/qa/CbmCaOutputQa.cxx b/reco/L1/qa/CbmCaOutputQa.cxx index 200de0a9dadf33ff862911b42162acdb31908bca..4d9d11ff4a65751146833e84a1eece427ba145e9 100644 --- a/reco/L1/qa/CbmCaOutputQa.cxx +++ b/reco/L1/qa/CbmCaOutputQa.cxx @@ -23,10 +23,10 @@ #include "L1InitManager.h" -using ::ca::tools::Debugger; -using ::ca::tools::MCTrack; using cbm::ca::MCModule; using cbm::ca::OutputQa; +using cbm::ca::tools::Debugger; +using cbm::ca::tools::MCTrack; // --------------------------------------------------------------------------------------------------------------------- // diff --git a/reco/L1/qa/CbmCaOutputQa.h b/reco/L1/qa/CbmCaOutputQa.h index 5aab6a9d13dd9eb1d2c72ff1c21adb13aca25334..3e3921cabcb45da44f655059243c980a745f01e0 100644 --- a/reco/L1/qa/CbmCaOutputQa.h +++ b/reco/L1/qa/CbmCaOutputQa.h @@ -23,6 +23,7 @@ #include "CaMonitor.h" #include "CaToolsDebugger.h" +#include "CaVector.h" #include "L1Parameters.h" namespace cbm::ca @@ -132,7 +133,7 @@ namespace cbm::ca /// Creates a debugger and enables its usage inside a QA task void EnableDebugger(const char* filename) { - if (!fpDebugger.get()) { fpDebugger = std::make_shared<::ca::tools::Debugger>(filename); } + if (!fpDebugger.get()) { fpDebugger = std::make_shared<tools::Debugger>(filename); } } /// @brief Reads defined parameters object from file @@ -265,13 +266,13 @@ namespace cbm::ca std::unique_ptr<TimeSliceReader> fpTSReader = nullptr; ///< Reader of the time slice std::shared_ptr<MCModule> fpMCModule = nullptr; ///< MC module std::shared_ptr<L1IODataManager> fpDataManager = nullptr; ///< Data manager - std::shared_ptr<::ca::tools::Debugger> fpDebugger = nullptr; ///< Debugger + std::shared_ptr<tools::Debugger> fpDebugger = nullptr; ///< Debugger std::shared_ptr<L1Parameters> fpParameters = nullptr; ///< Tracking parameters object - L1Vector<CbmL1HitId> fvHitIds {"CbmCaOutputQa::fvHitIds"}; - L1Vector<CbmL1HitDebugInfo> fvHits {"CbmCaOutputQa::fvHits"}; - L1Vector<CbmL1Track> fvRecoTracks {"CbmCaOutputQa::fvRecoTracks"}; - ::ca::tools::MCData fMCData; ///< Input MC data (points and tracks) + ca::Vector<CbmL1HitId> fvHitIds {"CbmCaOutputQa::fvHitIds"}; + ca::Vector<CbmL1HitDebugInfo> fvHits {"CbmCaOutputQa::fvHits"}; + ca::Vector<CbmL1Track> fvRecoTracks {"CbmCaOutputQa::fvRecoTracks"}; + tools::MCData fMCData; ///< Input MC data (points and tracks) /// @enum EMonitorKey /// @brief QA monitor counters @@ -285,7 +286,7 @@ namespace cbm::ca kEND }; - ::ca::Monitor<EMonitorKey> fMonitor {"Output tracking QA"}; + ca::Monitor<EMonitorKey> fMonitor {"Output tracking QA"}; std::set<ETrackType> fmSummaryTableEntries; ///< Which track types should be listed in the summary table diff --git a/reco/L1/qa/CbmCaTrackFitQa.cxx b/reco/L1/qa/CbmCaTrackFitQa.cxx index ba6655a205b8319ec3e2492ae945d924035eb3f8..f80f0c138e1dab702bd84648217b6b307f8b3be2 100644 --- a/reco/L1/qa/CbmCaTrackFitQa.cxx +++ b/reco/L1/qa/CbmCaTrackFitQa.cxx @@ -31,6 +31,7 @@ // ******************************************************* using cbm::ca::TrackFitQa; +using cbm::ca::tools::MCPoint; // --------------------------------------------------------------------------------------------------------------------- // @@ -128,7 +129,7 @@ void TrackFitQa::Init() // --------------------------------------------------------------------------------------------------------------------- // -void TrackFitQa::Fill(const L1TrackPar& trPar, const ::ca::tools::MCPoint& mcPoint, bool bTimeMeasured, double weight) +void TrackFitQa::Fill(const L1TrackPar& trPar, const tools::MCPoint& mcPoint, bool bTimeMeasured, double weight) { // Probably, a bottleneck L1Fit fitter; diff --git a/reco/L1/qa/CbmCaTrackFitQa.h b/reco/L1/qa/CbmCaTrackFitQa.h index 90c1f116997012ce363265b457cafbc57b8cd046..b854bdb30dd49920fa7d1e2e20042dc068b2a0a8 100644 --- a/reco/L1/qa/CbmCaTrackFitQa.h +++ b/reco/L1/qa/CbmCaTrackFitQa.h @@ -18,10 +18,9 @@ #include "L1EArray.h" #include "L1Field.h" #include "L1Fit.h" -#include "L1Vector.h" // Forward declarations -namespace ca::tools +namespace cbm::ca::tools { class MCData; class MCPoint; @@ -99,7 +98,7 @@ namespace cbm::ca /// @note To be called in the loop over reconstructed tracks full sample /// @param iTrkReco Index of reconstructed track /// @param weight Weight - void Fill(const L1TrackPar& trPar, const ::ca::tools::MCPoint& mcPoint, bool bTimeMeasured, double weight = 1); + void Fill(const L1TrackPar& trPar, const tools::MCPoint& mcPoint, bool bTimeMeasured, double weight = 1); /// @brief Sets particle mass, used for fitting a track /// @param mass Particle mass [GeV/c2] diff --git a/reco/L1/qa/CbmCaTrackTypeQa.cxx b/reco/L1/qa/CbmCaTrackTypeQa.cxx index 8241263cab123f368a524434e9d2be2c2ed3a45a..f25c0717d78e160ac9774a94b01b9c4bf90b69fd 100644 --- a/reco/L1/qa/CbmCaTrackTypeQa.cxx +++ b/reco/L1/qa/CbmCaTrackTypeQa.cxx @@ -14,9 +14,9 @@ #include "CaToolsMCData.h" -using ca::tools::MCPoint; -using ca::tools::MCTrack; using cbm::ca::TrackTypeQa; +using cbm::ca::tools::MCPoint; +using cbm::ca::tools::MCTrack; // --------------------------------------------------------------------------------------------------------------------- diff --git a/reco/L1/qa/CbmCaTrackTypeQa.h b/reco/L1/qa/CbmCaTrackTypeQa.h index f64ef78bb82d6340b05d92e66f9975193dc63636..b3b9bbeefae0e3aa25de43cc93ce0b5c0f2aab54 100644 --- a/reco/L1/qa/CbmCaTrackTypeQa.h +++ b/reco/L1/qa/CbmCaTrackTypeQa.h @@ -11,6 +11,7 @@ #define CbmCaTrackTypeQa_h 1 #include "CbmCaTrackFitQa.h" +#include "CbmL1DetectorID.h" #include "CbmL1Hit.h" #include "CbmQaIO.h" @@ -20,7 +21,7 @@ #include "L1Parameters.h" // Forward declarations -namespace ca::tools +namespace cbm::ca::tools { class MCData; class MCTrack; @@ -146,15 +147,15 @@ namespace cbm::ca /// @brief Registers the reconstructed tracks container /// @param vTracks Reference to reconstructed tracks container - void RegisterRecoTracks(L1Vector<CbmL1Track>& vTracks) { fpvRecoTracks = &vTracks; } + void RegisterRecoTracks(ca::Vector<CbmL1Track>& vTracks) { fpvRecoTracks = &vTracks; } /// @brief Register the reconstructed hits container /// @param vHits Reference to reconstructed hits container - void RegisterRecoHits(L1Vector<CbmL1HitDebugInfo>& vHits) { fpvHits = &vHits; } + void RegisterRecoHits(ca::Vector<CbmL1HitDebugInfo>& vHits) { fpvHits = &vHits; } /// @brief Registers the MC data /// @param mcData Reference to MC data object - void RegisterMCData(::ca::tools::MCData& mcData) { fpMCData = &mcData; } + void RegisterMCData(tools::MCData& mcData) { fpMCData = &mcData; } /// @brief Registers CA parameters object /// @param pParameters A shared pointer to the parameters object @@ -175,7 +176,7 @@ namespace cbm::ca /// @param cut Functor <bool(const MCTrack&)> defining cut on MC track /// /// When function returns false, - void SetMCTrackCut(std::function<bool(const ::ca::tools::MCTrack&)> cut) { fMCTrackCut = cut; } + void SetMCTrackCut(std::function<bool(const tools::MCTrack&)> cut) { fMCTrackCut = cut; } /// @brief Sets selection cuts on reconstructed tracks /// @param cut Functor <bool(const CbmL1Track&)> defining cut on reconstructed track @@ -285,15 +286,15 @@ namespace cbm::ca bool fbUseMC = false; ///< Flag: true - MC information is used TString fsTitle = ""; ///< Title of the track category - L1Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to vector of reconstructed tracks - L1Vector<CbmL1HitDebugInfo>* fpvHits = nullptr; ///< Pointer to vector of reconstructed hits - ::ca::tools::MCData* fpMCData = nullptr; ///< Pointer to MC data object + ca::Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to vector of reconstructed tracks + ca::Vector<CbmL1HitDebugInfo>* fpvHits = nullptr; ///< Pointer to vector of reconstructed hits + tools::MCData* fpMCData = nullptr; ///< Pointer to MC data object std::shared_ptr<L1Parameters> fpParameters = nullptr; ///< Pointer to parameters object // ** Cuts on tracks for a given track class ** /// Cut function on MC tracks - std::function<bool(const ::ca::tools::MCTrack&)> fMCTrackCut = [](const ::ca::tools::MCTrack&) { return true; }; + std::function<bool(const tools::MCTrack&)> fMCTrackCut = [](const tools::MCTrack&) { return true; }; /// Cut function on reconstructed tracks std::function<bool(const CbmL1Track&)> fRecoTrackCut = [](const CbmL1Track&) { return true; }; diff --git a/reco/L1/utils/CbmCaIdealHitProducerDet.h b/reco/L1/utils/CbmCaIdealHitProducerDet.h index 4916c6fe631074671f340cb3c15285a678809d1e..9b58a81df48da6a1a00bc5f8ade3a29b8e83472b 100644 --- a/reco/L1/utils/CbmCaIdealHitProducerDet.h +++ b/reco/L1/utils/CbmCaIdealHitProducerDet.h @@ -60,8 +60,6 @@ namespace cbm::ca { - using namespace cbm::algo::ca::constants; - /// @brief Ideal hit producer class /// /// @@ -177,7 +175,7 @@ namespace cbm::ca std::vector<HitParameters> fvStationPars; ///< Parameters, stored for each station - ::ca::algo::Random fRandom {1}; ///< Random generator + ca::Random fRandom {1}; ///< Random generator /// @brief Management flag, which does run the routine if there was a request bool fbRunTheRoutine = true; diff --git a/reco/detectors/rich/CMakeLists.txt b/reco/detectors/rich/CMakeLists.txt index 6e9d3f88911b98b3bbd1f231f75258f40dcdb5df..7fecbf83e010114b46eaceebb273453eff64791c 100644 --- a/reco/detectors/rich/CMakeLists.txt +++ b/reco/detectors/rich/CMakeLists.txt @@ -90,6 +90,7 @@ set(PUBLIC_DEPENDENCIES CbmBase CbmData CbmRichBase + CaCore L1 Littrack FairRoot::Base diff --git a/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.cxx b/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.cxx index 506262a95c91500417037dd75b4312c3316c00ed..a94fc78b009e33868ee6caea04773cbcdee2cb86 100644 --- a/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.cxx +++ b/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.cxx @@ -121,8 +121,8 @@ Int_t CbmL1RichENNRingFinderParallel::DoFind(CbmEvent* event, TClonesArray* HitA sort(Down.begin(), Down.end(), ENNHit::Compare); // save local-out indices correspondece - L1Vector<Int_t> outIndicesUp; // indices in HitArray indexed by index in Up - L1Vector<Int_t> outIndicesDown; // indices in HitArray indexed by index in Down + cbm::algo::ca::Vector<Int_t> outIndicesUp; // indices in HitArray indexed by index in Up + cbm::algo::ca::Vector<Int_t> outIndicesDown; // indices in HitArray indexed by index in Down outIndicesUp.reset(Up.size()); outIndicesDown.reset(Down.size()); for (THitIndex k = 0; k < Up.size(); k++) { @@ -143,8 +143,8 @@ Int_t CbmL1RichENNRingFinderParallel::DoFind(CbmEvent* event, TClonesArray* HitA #endif // PRINT_TIMING // save local-out indices correspondece - L1Vector<ENNHitV> UpV; - L1Vector<ENNHitV> DownV; + cbm::algo::ca::Vector<ENNHitV> UpV; + cbm::algo::ca::Vector<ENNHitV> DownV; UpV.reset((Up.size() + fvec::size() - 1) / fvec::size()); DownV.reset((Down.size() + fvec::size() - 1) / fvec::size()); for (THitIndex k = 0; k < Up.size(); k++) { @@ -234,8 +234,9 @@ Int_t CbmL1RichENNRingFinderParallel::DoFind(CbmEvent* event, TClonesArray* HitA } -void CbmL1RichENNRingFinderParallel::ENNRingFinder(const int NHits, L1Vector<ENNHitV>& HitsV, vector<ENNRing>& Rings, - float HitSize, THitIndex MinRingHits, fvec /*RMin*/, fvec RMax) +void CbmL1RichENNRingFinderParallel::ENNRingFinder(const int NHits, cbm::algo::ca::Vector<ENNHitV>& HitsV, + vector<ENNRing>& Rings, float HitSize, THitIndex MinRingHits, + fvec /*RMin*/, fvec RMax) { #ifdef PRINT_TIMING GetTimer("All").Start(0); @@ -267,8 +268,8 @@ void CbmL1RichENNRingFinderParallel::ENNRingFinder(const int NHits, L1Vector<ENN #ifdef PRINT_TIMING GetTimer("Ring finding").Start(0); #endif // PRINT_TIMING - L1Vector<ENNSearchHitV> SearchArea; - L1Vector<ENNHitV> PickUpArea; + cbm::algo::ca::Vector<ENNSearchHitV> SearchArea; + cbm::algo::ca::Vector<ENNHitV> PickUpArea; const int MaxAreaHits = 200; // TODO take into account NHits SearchArea.reset(MaxAreaHits, ENNSearchHitV()); @@ -276,7 +277,7 @@ void CbmL1RichENNRingFinderParallel::ENNRingFinder(const int NHits, L1Vector<ENN // TODO 1 #if 0 - L1Vector<ENNRingV> rings_tmp; // simd hits for tmp store + cbm::algo::ca::Vector<ENNRingV> rings_tmp; // simd hits for tmp store rings_tmp.reset(100); // TODO use NRings int nRings_tmp = 0; #endif @@ -518,7 +519,7 @@ void CbmL1RichENNRingFinderParallel::ENNRingFinder(const int NHits, L1Vector<ENN ringV.localIHits.push_back(iHit.localIndex); ringV.NHits = 1; ringV.chi2 = 0; - L1Vector<fvec> Shadow; // save hit indices of hits, which's quality will be changed + cbm::algo::ca::Vector<fvec> Shadow; // save hit indices of hits, which's quality will be changed Shadow.reserve(25); Shadow.push_back(iHit.localIndex); for( THitIndex ih = 0; ih < MaxSearchAreaSize; ih++){ diff --git a/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.h b/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.h index 9daca1c082716bf1b0b231a1431abe8b495e5131..977c4b73979234e8b0f69716e15ddb66fae3c34d 100644 --- a/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.h +++ b/reco/detectors/rich/finder/CbmL1RichENNRingFinderParallel.h @@ -24,17 +24,22 @@ #include "CbmRichRingFinder.h" -#include "CaSimd.h" -using namespace cbm::algo::ca; // TODO: remove "using" from headers - #include "TStopwatch.h" #include "TString.h" -#include "L1Vector.h" +#include "CaSimd.h" +#include "CaVector.h" +#include "L1Def.h" class ENNHit; class ENNRing; +namespace +{ + using cbm::algo::ca::fmask; + using cbm::algo::ca::fscal; + using cbm::algo::ca::fvec; +} // namespace class CbmL1RichENNRingFinderParallel : public CbmRichRingFinder { @@ -178,8 +183,8 @@ class CbmL1RichENNRingFinderParallel : public CbmRichRingFinder { }; - void ENNRingFinder(const int NHits, L1Vector<ENNHitV>& HitsV, std::vector<ENNRing>& Rings, float HitSize = 1., - THitIndex MinRingHits = 5, fvec RMin = 2., fvec RMax = 6.); + void ENNRingFinder(const int NHits, cbm::algo::ca::Vector<ENNHitV>& HitsV, std::vector<ENNRing>& Rings, + float HitSize = 1., THitIndex MinRingHits = 5, fvec RMin = 2., fvec RMax = 6.); public: /** Standard constructor **/