From 35bf9d244443cba1f213cb399e02e3e0c435deda Mon Sep 17 00:00:00 2001 From: sgorbuno <se.gorbunov@gsi.de> Date: Tue, 13 Jul 2021 16:29:01 +0000 Subject: [PATCH] update of L1Vectors --- reco/L1/CMakeLists.txt | 10 +- reco/L1/CbmL1.cxx | 366 +++++------------- reco/L1/CbmL1.h | 274 ++++++------- reco/L1/CbmL1CATrdTrackFinderSA.cxx | 3 +- reco/L1/CbmL1ReadEvent.cxx | 108 +++--- reco/L1/CbmL1TofMerger.cxx | 3 +- reco/L1/CbmL1TrackMerger.cxx | 3 +- reco/L1/CbmL1TrdTrackFinderSts.cxx | 3 +- reco/L1/L1Algo/L1Algo.cxx | 182 ++------- reco/L1/L1Algo/L1Algo.h | 198 +++++----- reco/L1/L1Algo/L1Branch.h | 3 +- reco/L1/L1Algo/L1CAMergeClones.cxx | 23 +- reco/L1/L1Algo/L1CATrackFinder.cxx | 209 +++++----- reco/L1/{CbmL1Def.h => L1Algo/L1Def.h} | 20 +- reco/L1/L1Algo/L1Extrapolation.h | 3 +- reco/L1/L1Algo/L1Field.h | 4 +- reco/L1/L1Algo/L1Filtration.h | 3 +- reco/L1/L1Algo/L1Fit.h | 3 +- reco/L1/L1Algo/L1FitMaterial.h | 3 +- reco/L1/L1Algo/L1Grid.cxx | 16 +- reco/L1/L1Algo/L1Grid.h | 90 ++--- reco/L1/L1Algo/L1HitArea.h | 3 +- reco/L1/L1Algo/L1MaterialInfo.h | 4 +- reco/L1/L1Algo/L1Timer.h | 4 +- reco/L1/L1Algo/L1TrackExtender.cxx | 8 +- reco/L1/L1Algo/L1TrackPar.h | 2 +- reco/L1/L1Algo/L1TrackParFit.h | 3 +- reco/L1/L1Algo/L1Triplet.h | 2 +- reco/L1/L1Algo/L1UMeasurementInfo.h | 2 +- reco/L1/L1Algo/L1Vector.h | 83 ++-- reco/L1/L1Algo/L1XYMeasurementInfo.h | 2 +- reco/L1/L1Algo/utils/L1AlgoDraw.h | 3 +- reco/L1/L1AlgoInputData.h | 31 +- .../CbmL1RichENNRingFinderParallel.h | 3 +- reco/L1/OffLineInterface/CbmL1RichRingQa.cxx | 3 +- reco/L1/ParticleFinder/CbmL1PFFitter.h | 4 +- 36 files changed, 718 insertions(+), 966 deletions(-) rename reco/L1/{CbmL1Def.h => L1Algo/L1Def.h} (89%) diff --git a/reco/L1/CMakeLists.txt b/reco/L1/CMakeLists.txt index 5276883e12..4ff4ab6bb7 100644 --- a/reco/L1/CMakeLists.txt +++ b/reco/L1/CMakeLists.txt @@ -8,6 +8,14 @@ Set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -Wshadow -Weffc++ -Wno-unu #Set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Weverything -Wno-padded -Wno-global-constructors") #Set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wold-style-cast") +IF (SSE_FOUND) + ADD_DEFINITIONS(-DHAVE_SSE) + Message(STATUS "L1 will be compiled with SSE support") +ELSE (SSE_FOUND) + Message(FATAL_ERROR "L1 can not be compiled without SSE support") +ENDIF (SSE_FOUND) + + Set(INCLUDE_DIRECTORIES ${CBMROOT_SOURCE_DIR}/reco/base ${CBMROOT_SOURCE_DIR}/reco/detectors/sts @@ -162,7 +170,7 @@ OffLineInterface/CbmL1GlobalFindTracksEvents.h #OffLineInterface / CbmL1SttHit.h #OffLineInterface / CbmL1SttTrackFinder.h #OffLineInterface / CbmL1SttTrack.h -#L1Algo / L1Algo.h +L1Algo/L1Def.h ) diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx index 4efb8e8ba6..3be7a99a57 100644 --- a/reco/L1/CbmL1.cxx +++ b/reco/L1/CbmL1.cxx @@ -60,7 +60,6 @@ #include <iomanip> #include <iostream> #include <list> -#include <vector> #include "L1Algo/L1Algo.h" #include "L1Algo/L1Branch.h" @@ -72,7 +71,6 @@ using std::cout; using std::endl; using std::ios; -using std::vector; #include "KFTopoPerformance.h" @@ -85,109 +83,7 @@ ClassImp(CbmL1) CbmL1* CbmL1::fInstance = 0; -CbmL1::CbmL1() - : FairTask("L1") - , algo(0) - , // for access to L1 Algorithm from L1::Instance - fDigiFile() - , fUseHitErrors(0) - , fmCBMmode(0) - , fGlobalMode(0) - , vRTracks() - , // reconstructed tracks - vFileEvent() - , vHitStore() - , vMCPoints() - , nMvdPoints(0) - , vMCPoints_in_Time_Slice() - , NStation(0) - , NMvdStations(0) - , NStsStations(0) - , NMuchStations(0) - , NTrdStations(0) - , NTOFStation(0) - , // number of detector stations (all\sts\mvd) - fPerformance(0) - , fSTAPDataMode(0) - , fSTAPDataDir("") - , - - fTrackingLevel(2) - , // really doesn't used - fMomentumCutOff(0.1) - , // really doesn't used - fGhostSuppression(1) - , // really doesn't used - fUseMVD(0) - , // really doesn't used - fUseMUCH(0) - , fUseTRD(0) - , fUseTOF(0) - , - - PrimVtx() - , fTimeSlice(nullptr) - , fEventList(nullptr) - , listStsDigi() - , fStsPoints(0) - , fMCTracks(0) - , fMvdPoints(0) - , - //listMCTracks (0), - //listStsDigi(0), - listStsPts(0) - , listStsDigiMatch(0) - , listStsClusters(0) - , listStsHits(0) - , listStsHitMatch(0) - , listStsClusterMatch(0) - , - - listMvdPts(0) - , listMvdHits(0) - , listMvdDigiMatches(0) - , listMvdHitMatches(0) - , - - nMuchPoints(0) - , fMuchPoints(nullptr) - , listMuchHitMatches(nullptr) - , fDigiMatchesMuch(nullptr) - , fClustersMuch(nullptr) - , fMuchPixelHits(nullptr) - , fDigisMuch(nullptr) - , fTrdDigiPar(nullptr) - , fTrdModuleInfo(nullptr) - , fTrdPoints(nullptr) - , listTrdHits(nullptr) - , fTrdHitMatches(nullptr) - , fTofPoints(nullptr) - , fTofHitDigiMatches(nullptr) - , fTofHits(nullptr) - , fDigiPar(nullptr) - , fTofDigiBdfPar(nullptr) - , - - fPerfFile(nullptr) - , fHistoDir(nullptr) - , vStsHits() - , SortedIndex(0) - , StsIndex(0) - , vMCTracks() - , vHitMCRef() - , dFEI2vMCPoints() - , dFEI2vMCTracks() - , fData(nullptr) - , fFindParticlesMode() - , fStsMatBudgetFileName("") - , fMvdMatBudgetFileName("") - , fMuchMatBudgetFileName("") - , fTrdMatBudgetFileName("") - , fTofMatBudgetFileName("") - , fExtrapolateToTheEndOfSTS(false) - , fTimesliceMode(0) - , fTopoPerformance(nullptr) - , fEventEfficiency() +CbmL1::CbmL1() : FairTask("L1") { if (!fInstance) fInstance = this; } @@ -195,113 +91,17 @@ CbmL1::CbmL1() CbmL1::CbmL1(const char* name, Int_t iVerbose, Int_t _fPerformance, int fSTAPDataMode_, TString fSTAPDataDir_, int findParticleMode_) : FairTask(name, iVerbose) - , algo(0) - , // for access to L1 Algorithm from L1::Instance - fDigiFile() - , fUseHitErrors(0) - , fmCBMmode(0) - , fGlobalMode(0) - , vRTracks() - , // reconstructed tracks - vFileEvent() - , vHitStore() - , vMCPoints() - , nMvdPoints(0) - , vMCPoints_in_Time_Slice() - , NStation(0) - , NMvdStations(0) - , NStsStations(0) - , NMuchStations(0) - , NTrdStations(0) - , NTOFStation(0) - , // number of detector stations (all\sts\mvd) - fPerformance(_fPerformance) + , fPerformance(_fPerformance) , fSTAPDataMode(fSTAPDataMode_) , fSTAPDataDir(fSTAPDataDir_) - , fTrackingLevel(2) - , // really doesn't used - fMomentumCutOff(0.1) - , // really doesn't used - fGhostSuppression(1) - , // really doesn't used - fUseMVD(0) - , // really doesn't used - fUseMUCH(0) - , fUseTRD(0) - , fUseTOF(0) - , - - - PrimVtx() - , fTimeSlice(nullptr) - , fEventList(nullptr) - , listStsDigi(0) - , fStsPoints(0) - , fMCTracks(0) - , fMvdPoints(NULL) - , - //listMCTracks (0), - - listStsPts(0) - , listStsDigiMatch(0) - , listStsClusters(0) - , listStsHits(0) - , listStsHitMatch(0) - , listStsClusterMatch(0) - , - - listMvdPts(0) - , listMvdHits(0) - , listMvdDigiMatches(0) - , listMvdHitMatches(0) - , - - nMuchPoints(0) - , fMuchPoints(nullptr) - , listMuchHitMatches(nullptr) - , fDigiMatchesMuch(nullptr) - , fClustersMuch(nullptr) - , fMuchPixelHits(nullptr) - , fDigisMuch(nullptr) - , fTrdDigiPar(nullptr) - , fTrdModuleInfo(nullptr) - , fTrdPoints(nullptr) - , listTrdHits(nullptr) - , fTrdHitMatches(nullptr) - , fTofPoints(nullptr) - , fTofHitDigiMatches(nullptr) - , fTofHits(nullptr) - , fDigiPar(nullptr) - , fTofDigiBdfPar(nullptr) - , - - fPerfFile(nullptr) - , fHistoDir(nullptr) - , vStsHits() - , SortedIndex(0) - , StsIndex(0) - , vMCTracks() - , vHitMCRef() - , dFEI2vMCPoints() - , dFEI2vMCTracks() - , fData(nullptr) , fFindParticlesMode(findParticleMode_) - , fStsMatBudgetFileName("") - , fMvdMatBudgetFileName("") - , fMuchMatBudgetFileName("") - , fTrdMatBudgetFileName("") - , fTofMatBudgetFileName("") - , fExtrapolateToTheEndOfSTS(false) - , fTimesliceMode(0) - , fTopoPerformance(nullptr) - , fEventEfficiency() { if (!fInstance) fInstance = this; } CbmL1::~CbmL1() { - if (fInstance == this) fInstance = 0; + if (fInstance == this) fInstance = nullptr; } @@ -556,8 +356,8 @@ InitStatus CbmL1::Init() algo = &algo_static; - vector<fscal> geo; - geo.clear(); + L1Vector<fscal> geo("geo"); + geo.reserve(1000); for (int i = 0; i < 3; i++) { Double_t point[3] = {0, 0, 2.5 * i}; @@ -577,7 +377,7 @@ InitStatus CbmL1::Init() TFile* oldFile = gFile; TDirectory* oldDir = gDirectory; - TFile* file = new TFile(fDigiFile); + TFile* file = new TFile(fMuchDigiFile, "READ"); TObjArray* stations = (TObjArray*) file->Get("stations"); fGeoScheme->Init(stations, 0); for (int iStation = 0; iStation < fGeoScheme->GetNStations(); iStation++) { @@ -610,38 +410,34 @@ InitStatus CbmL1::Init() // count ToF parameters - int maxTofStation = 0; - int St = 0; - - vector<int> NHits; - vector<float> Z_pos; - - NHits.resize(10, 0); - Z_pos.resize(10, 0); + // if (fUseTOF) { + // int maxTofStation = 0; + // L1Vector<int> NHits("NHits", 10, 0); + // L1Vector<float> Z_pos("Z_pos", 10, 0.f); - if (fUseTOF) { + // if (fTofHits) { + // for (int j = 0; j < fTofHits->GetEntriesFast(); j++) { - if (fTofHits) { - for (int j = 0; j < fTofHits->GetEntriesFast(); j++) { + // CbmTofHit* mh = L1_DYNAMIC_CAST<CbmTofHit*>(fTofHits->At(j)); + // int St = fTofDigiBdfPar->GetNbTrackingStations(); - CbmTofHit* mh = L1_DYNAMIC_CAST<CbmTofHit*>(fTofHits->At(j)); - St = fTofDigiBdfPar->GetNbTrackingStations(); + // if (maxTofStation < St) maxTofStation = St; - if (maxTofStation < St) maxTofStation = St; - - TVector3 pos; - mh->Position(pos); - Z_pos[St] += pos.Z(); - NHits[St]++; - } - } + // TVector3 pos; + // mh->Position(pos); + // Z_pos[St] += pos.Z(); + // NHits[St]++; + // } + // } - if (fUseTOF) NTOFStation = 3; //fTofDigiBdfPar->GetNbTrackingStations(); + // for (int i = 0; i < (maxTofStation + 1); i++) { + // Z_pos[i] = Z_pos[i] / NHits[i]; + // } + // } - for (int i = 0; i < (maxTofStation + 1); i++) - Z_pos[i] = Z_pos[i] / NHits[i]; - } + // TODO: Read N TOF stations from geometry + if (fUseTOF) NTOFStation = 3; //fTofDigiBdfPar->GetNbTrackingStations(); CbmStsSetup* stsSetup = CbmStsSetup::Instance(); if (!stsSetup->IsInit()) { stsSetup->Init(nullptr); } @@ -951,8 +747,7 @@ InitStatus CbmL1::Init() algo->Init(geo, fUseHitErrors, fmCBMmode); geo.clear(); - - algo->fRadThick.resize(algo->NStations); + algo->fRadThick.reset(algo->NStations); // Read STS MVD TRD MuCh ToF Radiation Thickness table TString stationName = "Radiation Thickness [%], Station"; @@ -1001,7 +796,6 @@ InitStatus CbmL1::Init() LOG(warn) << "No MVD material budget file is found. Homogenious budget " "will be used"; for (int iSta = 0; iSta < algo->NMvdStations; iSta++) { - cout << iSta << endl; algo->fRadThick[iSta].SetBins(1, 100); // mvd should be in +-100 cm square algo->fRadThick[iSta].table.resize(1); @@ -1242,10 +1036,12 @@ void CbmL1::Exec(Option_t* /*option*/) {} void CbmL1::Reconstruct(CbmEvent* event) { + static int nevent = 0; vFileEvent.clear(); - std::vector<std::pair<double, int>> SortStsHits; + L1Vector<std::pair<double, int>> SortStsHits("CbmL1::SortStsHits"); + SortStsHits.reserve(listStsHits->GetEntriesFast()); float start_t = 10000000000; @@ -1267,14 +1063,13 @@ void CbmL1::Reconstruct(CbmEvent* event) std::sort(SortStsHits.begin(), SortStsHits.end()); - StsIndex.resize(0); - + StsIndex.clear(); + StsIndex.reserve(SortStsHits.size()); for (unsigned int i = 0; i < SortStsHits.size(); i++) { int j = SortStsHits[i].second; StsIndex.push_back(j); }; - if (fTimesliceMode && fPerformance) { int nofEvents = fEventList->GetNofEvents(); @@ -1291,15 +1086,27 @@ void CbmL1::Reconstruct(CbmEvent* event) vFileEvent.insert(DFSET::value_type(iFile, iEvent)); } - if (fVerbose > 1) { cout << endl << "CbmL1::Exec event " << ++nevent << " ..." << endl << endl; } #ifdef _OPENMP omp_set_num_threads(1); #endif // repack data - vector<CbmL1Track> vRTracksCur; // reconstructed tracks + L1Vector<CbmL1Track> vRTracksCur("CbmL1::vRTracksCur"); // reconstructed tracks + { + int nHits = 0; + int nSta = 1; + if (listMvdHits) { + nHits += listMvdHits->GetEntriesFast(); + nSta += NMvdStations; + } + if (listStsHits) { + nHits += listStsHits->GetEntriesFast(); + nSta += NStsStations; + } + vRTracksCur.reserve(10 + (2 * nHits) / nSta); + } while (newTS) { @@ -1322,14 +1129,13 @@ void CbmL1::Reconstruct(CbmEvent* event) algo->SetData(fData->GetStsHits(), fData->GetNStsStrips(), fData->GetStsZPos(), fData->GetSFlag(), fData->GetStsHitsStartIndex(), fData->GetStsHitsStopIndex()); } - else + else { ReadEvent(fData, TsStart, TsLength, TsOverlap, FstHitinTs, newTS, event); + } if (0) { // correct hits on MC // dbg TRandom3 random; - vector<int> strips, zP; - strips.clear(); - zP.clear(); + L1Vector<int> strips("CbmL1::strips"), zP("CbmL1::zP"); for (unsigned int iH = 0; iH < (*algo->vStsHits).size(); ++iH) { L1Hit& h = const_cast<L1Hit&>((*algo->vStsHits)[iH]); #ifdef USE_EVENT_NUMBER @@ -1359,8 +1165,8 @@ void CbmL1::Reconstruct(CbmEvent* event) } strips.push_back(h.b); if (std::find(zP.begin(), zP.end(), h.iz) != zP.end()) { // TODO why do we need it??gives prob=0 - h.iz = (*algo->vStsZPos).size(); - const_cast<std::vector<float>*>(algo->vStsZPos)->push_back(0); + h.iz = algo->vStsZPos->size(); + algo->vStsZPos->push_back(0.f); } zP.push_back(h.iz); @@ -1385,7 +1191,7 @@ void CbmL1::Reconstruct(CbmEvent* event) HitMatch(); // calculate the max number of Hits\mcPoints on continuous(consecutive) stations - for (vector<CbmL1MCTrack>::iterator it = vMCTracks.begin(); it != vMCTracks.end(); ++it) + for (L1Vector<CbmL1MCTrack>::iterator it = vMCTracks.begin(); it != vMCTracks.end(); ++it) it->Init(); } @@ -1406,8 +1212,8 @@ void CbmL1::Reconstruct(CbmEvent* event) for (unsigned int iH = 0; iH < (*algo->vStsHits).size(); ++iH) { #ifdef USE_EVENT_NUMBER - L1Hit& h = const_cast<L1Hit&>((*algo->vStsHits)[iH]); - h.n = -1; + L1Hit& h = const_cast<L1Hit&>((*algo->vStsHits)[iH]); + h.n = -1; #endif if (vStsHits[iH].mcPointIds.size() == 0) continue; #ifdef USE_EVENT_NUMBER @@ -1416,18 +1222,18 @@ void CbmL1::Reconstruct(CbmEvent* event) #endif } - for (vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { + for (L1Vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { CbmL1MCTrack& MC = *i; if (!MC.IsReconstructable()) continue; if (!(MC.ID >= 0)) continue; if (MC.StsHits.size() < 4) continue; - vector<int> hitIndices(algo->NStations, -1); + L1Vector<int> hitIndices("hitIndices", algo->NStations, -1); for (unsigned int iH = 0; iH < MC.StsHits.size(); iH++) { - const int hitI = MC.StsHits[iH]; - CbmL1Hit& hit = const_cast<CbmL1Hit&>(vStsHits[hitI]); + const int hitI = MC.StsHits[iH]; + CbmL1Hit& hit = const_cast<CbmL1Hit&>(vStsHits[hitI]); hit.event = MC.iEvent; @@ -1497,16 +1303,18 @@ void CbmL1::Reconstruct(CbmEvent* event) t.StsHits.clear(); t.fTrackTime = it->fTrackTime; - vector<int> StsHitsLocal; + L1Vector<int> StsHitsLocal("CbmL1::StsHitsLocal"); + StsHitsLocal.reserve(it->NHits); for (int i = 0; i < it->NHits; i++) { int start_hit1 = start_hit; if (algo->fRecoHits[start_hit1] > vStsHits.size() - 1) start_hit++; - else if (fTimesliceMode) + else if (fTimesliceMode) { t.StsHits.push_back(((*algo->vStsHits)[algo->fRecoHits[start_hit]]).ID); - else + } + else { t.StsHits.push_back(algo->fRecoHits[start_hit]); - + } StsHitsLocal.push_back(algo->fRecoHits[start_hit++]); } @@ -1585,7 +1393,7 @@ void CbmL1::Reconstruct(CbmEvent* event) HitMatch(); // calculate the max number of Hits\mcPoints on continuous(consecutive) stations - for (vector<CbmL1MCTrack>::iterator it = vMCTracks.begin(); it != vMCTracks.end(); ++it) + for (L1Vector<CbmL1MCTrack>::iterator it = vMCTracks.begin(); it != vMCTracks.end(); ++it) it->Init(); } // @@ -1599,6 +1407,8 @@ void CbmL1::Reconstruct(CbmEvent* event) // ReadSTAPPerfData(); // }; + vRTracks.clear(); + vRTracks.reserve(vRTracksCur.size()); for (unsigned int iTrack = 0; iTrack < vRTracksCur.size(); iTrack++) { for (unsigned int iHit = 0; iHit < vRTracksCur[iTrack].StsHits.size(); iHit++) @@ -1613,8 +1423,8 @@ void CbmL1::Reconstruct(CbmEvent* event) for (unsigned int iH = 0; iH < (*algo->vStsHits).size(); ++iH) { #ifdef USE_EVENT_NUMBER - L1Hit& h = const_cast<L1Hit&>((*algo->vStsHits)[iH]); - h.n = -1; + L1Hit& h = const_cast<L1Hit&>((*algo->vStsHits)[iH]); + h.n = -1; #endif if (vStsHits[iH].mcPointIds.size() == 0) continue; #ifdef USE_EVENT_NUMBER @@ -1623,18 +1433,18 @@ void CbmL1::Reconstruct(CbmEvent* event) #endif } - for (vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { + for (L1Vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { CbmL1MCTrack& MC = *i; if (!MC.IsReconstructable()) continue; if (!(MC.ID >= 0)) continue; if (MC.StsHits.size() < 4) continue; - vector<int> hitIndices(algo->NStations, -1); + L1Vector<int> hitIndices("CbmL1::hitIndices", algo->NStations, -1); for (unsigned int iH = 0; iH < MC.StsHits.size(); iH++) { - const int hitI = MC.StsHits[iH]; - CbmL1Hit& hit = const_cast<CbmL1Hit&>(vStsHits[hitI]); + const int hitI = MC.StsHits[iH]; + CbmL1Hit& hit = const_cast<CbmL1Hit&>(vStsHits[hitI]); hit.event = MC.iEvent; } @@ -1712,7 +1522,7 @@ void CbmL1::IdealTrackFinder() algo->fTracks.clear(); algo->fRecoHits.clear(); - for (vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { + for (L1Vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) { CbmL1MCTrack& MC = *i; if (!MC.IsReconstructable()) continue; @@ -1723,11 +1533,11 @@ void CbmL1::IdealTrackFinder() L1Track algoTr; algoTr.NHits = 0; - vector<int> hitIndices(algo->NStations, -1); + L1Vector<int> hitIndices("CbmL1::hitIndices", algo->NStations, -1); for (unsigned int iH = 0; iH < MC.StsHits.size(); iH++) { - const int hitI = MC.StsHits[iH]; - const CbmL1Hit& hit = vStsHits[hitI]; + const int hitI = MC.StsHits[iH]; + const CbmL1Hit& hit = vStsHits[hitI]; const int iStation = vMCPoints[hit.mcPointIds[0]].iStation; @@ -1768,7 +1578,7 @@ void CbmL1::IdealTrackFinder() /// ----- STandAlone Package service-functions ----------------------------- -void CbmL1::WriteSTAPGeoData(const vector<float>& geo_) +void CbmL1::WriteSTAPGeoData(const L1Vector<float>& geo_) { // write geo in file TString fgeo_name = fSTAPDataDir + "geo_algo.txt"; @@ -2037,9 +1847,9 @@ istream& CbmL1::eatwhite(istream& is) // skip spaces return is; } -//void CbmL1::ReadSTAPGeoData(vector<float> geo_, int &size) -//void CbmL1::ReadSTAPGeoData(vector<fscal> geo_, int &size) -void CbmL1::ReadSTAPGeoData(vector<fscal>& geo_, int& size) +//void CbmL1::ReadSTAPGeoData(L1Vector<float> geo_, int &size) +//void CbmL1::ReadSTAPGeoData(L1Vector<fscal> geo_, int &size) +void CbmL1::ReadSTAPGeoData(L1Vector<fscal>& geo_, int& size) { TString fgeo_name = fSTAPDataDir + "geo_algo.txt"; ifstream fgeo(fgeo_name); @@ -2065,9 +1875,9 @@ void CbmL1::ReadSTAPAlgoData() if (1) { if (nEvent == 1) fadata.open(fadata_name, fstream::in); - if (algo->vStsHits) const_cast<std::vector<L1Hit>*>(algo->vStsHits)->clear(); + if (algo->vStsHits) algo->vStsHits->clear(); algo->NStsStrips = 0; - if (algo->vStsZPos) const_cast<std::vector<float>*>(algo->vStsZPos)->clear(); + if (algo->vStsZPos) algo->vStsZPos->clear(); if (algo->fStripFlag) algo->fStripFlag->clear(); // check correct position in file @@ -2091,7 +1901,7 @@ void CbmL1::ReadSTAPAlgoData() for (int i = 0; i < n; i++) { fscal element; fadata >> element; - const_cast<std::vector<float>*>(algo->vStsZPos)->push_back(element); + algo->vStsZPos->push_back(element); } if (fVerbose >= 4) { cout << "vStsZPos[" << n << "]" @@ -2120,7 +1930,7 @@ void CbmL1::ReadSTAPAlgoData() element.f = static_cast<THitI>(element_f); element.b = static_cast<THitI>(element_b); element.iz = static_cast<TZPosI>(element_iz); - const_cast<std::vector<L1Hit>*>(algo->vStsHits)->push_back(element); + algo->vStsHits->push_back(element); } if (fVerbose >= 4) { cout << "vStsHits[" << n << "]" @@ -2264,6 +2074,7 @@ void CbmL1::ReadSTAPPerfData() } // vHitMCRef fpdata >> n; + vHitMCRef.reserve(n); for (int i = 0; i < n; i++) { int element; fpdata >> element; @@ -2275,6 +2086,7 @@ void CbmL1::ReadSTAPPerfData() } // vHitStore fpdata >> n; + vHitStore.reserve(n); for (int i = 0; i < n; i++) { CbmL1HitStore element; fpdata >> element.ExtIndex; @@ -2507,7 +2319,7 @@ void CbmL1::WriteSIMDKFData() McTracksOut.open("mctracksout.dat", fstream::out | fstream::app); } - for (vector<CbmL1Track>::iterator RecTrack = vRTracks.begin(); RecTrack != vRTracks.end(); ++RecTrack) { + for (L1Vector<CbmL1Track>::iterator RecTrack = vRTracks.begin(); RecTrack != vRTracks.end(); ++RecTrack) { if (RecTrack->IsGhost()) continue; CbmL1MCTrack* MCTrack = RecTrack->GetMCTrack(); diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h index c8c818cade..4b67502e2a 100644 --- a/reco/L1/CbmL1.h +++ b/reco/L1/CbmL1.h @@ -29,6 +29,8 @@ #include "CbmL1Track.h" #include "CbmL1Vtx.h" +#include "L1Algo/L1Vector.h" + //#include "L1Algo/L1Algo.h" #include "CbmEvent.h" #include "CbmL1Hit.h" @@ -73,12 +75,9 @@ #include <map> #include <set> #include <utility> -#include <vector> #include "L1EventEfficiencies.h" -using std::vector; - class L1AlgoInputData; class L1Algo; class L1Event; @@ -118,16 +117,17 @@ private: CbmL1 operator=(const CbmL1&); public: - L1Algo* algo; // for access to L1 Algorithm from L1::Instance - - TString fDigiFile; // Digitization file - bool fUseHitErrors; - bool fmCBMmode; - bool fGlobalMode; - vector<CbmL1Track> vRTracks; // reconstructed tracks typedef std::set<std::pair<int, int>> DFSET; - DFSET vFileEvent; - // set<pair<int, int> > vFileEvent; + + friend class L1AlgoDraw; + friend class L1AlgoPulls; + + template<int NHits> + friend class L1AlgoEfficiencyPerformance; + + friend class CbmL1MCTrack; + friend class CbmL1PFFitter; + static CbmL1* Instance() { return fInstance; } void SetParContainers(); @@ -154,7 +154,7 @@ public: void SetTofMaterialBudgetFileName(TString fileName) { fTofMatBudgetFileName = fileName; } void SetExtrapolateToTheEndOfSTS(bool b) { fExtrapolateToTheEndOfSTS = b; } void SetDataMode(int TimesliceMode) { fTimesliceMode = TimesliceMode; } - void SetMuchPar(TString fileName) { fDigiFile = fileName; } + void SetMuchPar(TString fileName) { fMuchDigiFile = fileName; } void SetUseHitErrors(bool value) { fUseHitErrors = value; } void SetmCBMmode(bool value) { fmCBMmode = value; } @@ -165,16 +165,8 @@ public: // void SetGhostSuppression( Bool_t b ){ fGhostSuppression= b; } // void SetDetectorEfficiency( Double_t eff ){ fDetectorEfficiency = eff; } - vector<CbmL1HitStore> vHitStore; // diff hit information - void Reconstruct(CbmEvent* event = NULL); - friend class L1AlgoDraw; - friend class L1AlgoPulls; - template<int NHits> - friend class L1AlgoEfficiencyPerformance; - friend class CbmL1MCTrack; - friend class CbmL1PFFitter; // bool ReadMCDataFromFile(const char work_dir[100], const int maxNEvent, const int iVerbose); // vector<CbmL1MCPoint> vMCPoints; @@ -182,9 +174,14 @@ public: inline double Get_Z_vMCPoint(int a) const { return vMCPoints[a].z; } private: - vector<CbmL1MCPoint> vMCPoints; - int nMvdPoints; - vector<int> vMCPoints_in_Time_Slice; + typedef std::map<Double_t, Int_t> DFEI2I; + + struct TH1FParameters { + TString name, title; + int nbins; + float xMin, xMax; + }; + void IdealTrackFinder(); // just copy all reconstructable MCTracks into RecoTracks. /// Read information about hits, mcPoints and mcTracks into L1 classes @@ -211,11 +208,11 @@ private: void HistoPerformance(); // fill some histograms and calculate efficiencies /// STandAlone Package service-functions - void WriteSTAPGeoData(const vector<float>& geo); // create geo_algo.dat - void WriteSTAPAlgoData(); // create data_algo.dat - void WriteSTAPPerfData(); // create data_perfo.dat - //void ReadSTAPGeoData(vector<float> geo, int &size); - void ReadSTAPGeoData(vector<float>& geo, int& size); + void WriteSTAPGeoData(const L1Vector<float>& geo); // create geo_algo.dat + void WriteSTAPAlgoData(); // create data_algo.dat + void WriteSTAPPerfData(); // create data_perfo.dat + //void ReadSTAPGeoData(L1Vector<float> geo, int &size); + void ReadSTAPGeoData(L1Vector<float>& geo, int& size); void ReadSTAPAlgoData(); void ReadSTAPPerfData(); /// SIMD KF Banchmark service-functions @@ -226,135 +223,144 @@ private: static std::istream& eatwhite(std::istream& is); // skip spaces static void writedir2current(TObject* obj); // help procedure - int NStation, NMvdStations, NStsStations, NMuchStations, NTrdStations, - NTOFStation; // number of detector stations (all\sts\mvd) - Int_t fPerformance; // 0 - w\o perf. 1 - L1-Efficiency definition. 2 - QA-Eff.definition - int - fSTAPDataMode; // way to work with file for standalone package. 0 (off) , 1 (write), 2 (read data and work only with it), 3 (debug - write and read) - TString fSTAPDataDir; - Int_t fTrackingLevel; // really doesn't used - Double_t fMomentumCutOff; // really doesn't used - Bool_t fGhostSuppression; // really doesn't used - Bool_t fUseMVD, fUseMUCH, fUseTRD, fUseTOF; // - // Double_t fDetectorEfficiency; // really doesn't used + inline Double_t dFEI(int file, int event, int idx) { return (1000 * idx) + file + (0.0001 * event); } - CbmL1Vtx PrimVtx; - // L1FieldSlice *targetFieldSlice _fvecalignment; +public: + L1Algo* algo {nullptr}; // for access to L1 Algorithm from L1::Instance - /// Input data - CbmTimeSlice* fTimeSlice; - CbmMCEventList* fEventList; //! MC event list (all) - //vector<CbmStsDigi> *listStsDigi; - vector<CbmStsDigi> listStsDigi; - CbmMCDataArray* fStsPoints; - CbmMCDataArray* fMCTracks; - CbmMCDataArray* fMvdPoints; - - - //TClonesArray *listMCTracks ; - TClonesArray* listStsPts; // Sts MC points - //TClonesArray *listStsDigi; - TClonesArray* listStsDigiMatch; - TClonesArray* listStsClusters; - TClonesArray* listStsHits; - TClonesArray* listStsHitMatch; - TClonesArray* listStsClusterMatch; - - // Pointers to data arrays - // CbmEvBasedArray* fMCTracks; - // STS - // CbmEvBasedArray* fStsPts; // CbmStsPoint array - // MVD - // CbmEvBasedArray* fMvdPts; //CbmMvdPoint array - - TClonesArray* listMvdPts; // Mvd MC points - TClonesArray* listMvdHits; - TClonesArray* listMvdDigiMatches; - TClonesArray* listMvdHitMatches; + TString fMuchDigiFile {}; // Much digitization file name + bool fUseHitErrors {false}; + bool fmCBMmode {false}; + bool fGlobalMode {false}; + L1Vector<CbmL1Track> vRTracks {"CbmL1::vRTracks"}; // reconstructed tracks + DFSET vFileEvent {}; - //MuCh - int nMuchPoints; - CbmMCDataArray* fMuchPoints; - TClonesArray* listMuchHitMatches; // Output CbmMatch array - TClonesArray* fDigiMatchesMuch; - TClonesArray* fClustersMuch; + L1Vector<CbmL1HitStore> vHitStore {"CbmL1::vHitStore"}; // diff hit information - TClonesArray* fMuchPixelHits; // CbmMuchPixelHit array - TClonesArray* fDigisMuch; +private: + static CbmL1* fInstance; - //TRD + L1AlgoInputData* fData {nullptr}; - CbmTrdParSetDigi* fTrdDigiPar; //! - CbmTrdParModDigi* fTrdModuleInfo; //! + L1Vector<CbmL1MCPoint> vMCPoints {"CbmL1::vMCPoints"}; + int nMvdPoints {0}; + L1Vector<int> vMCPoints_in_Time_Slice {"CbmL1::vMCPoints_in_Time_Slice"}; - CbmMCDataArray* fTrdPoints; - TClonesArray* listTrdHits; - TClonesArray* fTrdHitMatches; + int NStation {0}; // number of all detector stations + int NMvdStations {0}; // number of mvd stations + int NStsStations {0}; // number of sts stations + int NMuchStations {0}; // number of much stations + int NTrdStations {0}; // number of trd stations + int NTOFStation {0}; // number of tof stations - //ToF - CbmMCDataArray* fTofPoints; - TClonesArray* fTofHitDigiMatches; // CbmMatches array - TClonesArray* fTofHits; // CbmMatches array - CbmTofDigiPar* fDigiPar; - CbmTofDigiBdfPar* fTofDigiBdfPar; + Int_t fPerformance {0}; // 0 - w\o perf. 1 - L1-Efficiency definition. 2 - QA-Eff.definition + int fSTAPDataMode {0}; // way to work with file for standalone package. + // 0 (off) , 1 (write), 2 (read data and work only with it), 3 (debug - write and read) - struct TH1FParameters { - TString name, title; - int nbins; - float xMin, xMax; - }; + TString fSTAPDataDir {}; - TFile* fPerfFile; - TDirectory* fHistoDir; + Int_t fTrackingLevel {2}; // currently not used + Double_t fMomentumCutOff {0.1}; // currently not used + Bool_t fGhostSuppression {true}; // currently not used - static const int fNTimeHistos = 22; - TH1F* fTimeHisto[fNTimeHistos]; + Bool_t fUseMVD {false}; // if Mvd data should be processed + Bool_t fUseMUCH {false}; // if Much data should be processed + Bool_t fUseTRD {false}; // if Trd data should be processed + Bool_t fUseTOF {false}; // if Tof data should be processed - static const int fNGhostHistos = 9; - TH1F* fGhostHisto[fNGhostHistos]; + CbmL1Vtx PrimVtx {}; - //CbmMCEventHeader* fEvent; - /// Used data = Repacked input data - vector<CbmL1Hit> vStsHits; // hits with hit-mcpoint match information - // vector<CbmL1MCPoint> vMCPoints; - vector<int> SortedIndex; - vector<int> StsIndex; - vector<CbmL1MCTrack> vMCTracks; - vector<int> - vHitMCRef; // indices of MCPoints in vMCPoints, indexed by index of hit in algo->vStsHits array. According to StsMatch. Used for IdealResponce - // vector<int> vHitMCRef1; - // CbmMatch stsHitMatch; - //std::map<Double_t, Int_t> FEI2vMCPoints; - typedef std::map<Double_t, Int_t> DFEI2I; - DFEI2I dFEI2vMCPoints; - DFEI2I dFEI2vMCTracks; - inline Double_t dFEI(int file, int event, int idx) { return (1000 * idx) + file + (0.0001 * event); } - // DFEI2I::iterator map_it; - L1AlgoInputData* fData; + /// Input data + CbmTimeSlice* fTimeSlice {nullptr}; + CbmMCEventList* fEventList {nullptr}; //! MC event list (all) + + L1Vector<CbmStsDigi> listStsDigi {"CbmL1::listStsDigi"}; + CbmMCDataArray* fStsPoints {nullptr}; + CbmMCDataArray* fMvdPoints {nullptr}; + CbmMCDataArray* fMCTracks {nullptr}; + + TClonesArray* listStsPts {nullptr}; // Sts MC points + TClonesArray* listStsDigiMatch {nullptr}; + TClonesArray* listStsClusters {nullptr}; + TClonesArray* listStsHits {nullptr}; + TClonesArray* listStsHitMatch {nullptr}; + TClonesArray* listStsClusterMatch {nullptr}; + + TClonesArray* listMvdPts {nullptr}; // Mvd MC points + TClonesArray* listMvdHits {nullptr}; + TClonesArray* listMvdDigiMatches {nullptr}; + TClonesArray* listMvdHitMatches {nullptr}; - static CbmL1* fInstance; + //MuCh + int nMuchPoints {0}; + CbmMCDataArray* fMuchPoints {nullptr}; + TClonesArray* listMuchHitMatches {nullptr}; // Output CbmMatch array + TClonesArray* fDigiMatchesMuch {nullptr}; + TClonesArray* fClustersMuch {nullptr}; + TClonesArray* fMuchPixelHits {nullptr}; // CbmMuchPixelHit array + TClonesArray* fDigisMuch {nullptr}; -private: - int fFindParticlesMode; + //TRD - TString fStsMatBudgetFileName; - TString fMvdMatBudgetFileName; - TString fMuchMatBudgetFileName; - TString fTrdMatBudgetFileName; - TString fTofMatBudgetFileName; - bool fExtrapolateToTheEndOfSTS; - int fTimesliceMode; + CbmTrdParSetDigi* fTrdDigiPar {nullptr}; //! + CbmTrdParModDigi* fTrdModuleInfo {nullptr}; //! - KFTopoPerformance* fTopoPerformance; - L1EventEfficiencies fEventEfficiency; // average efficiencies + CbmMCDataArray* fTrdPoints {nullptr}; + TClonesArray* listTrdHits {nullptr}; + TClonesArray* fTrdHitMatches {nullptr}; + + //ToF + CbmMCDataArray* fTofPoints {nullptr}; + TClonesArray* fTofHitDigiMatches {nullptr}; // CbmMatches array + TClonesArray* fTofHits {nullptr}; // CbmMatches array + CbmTofDigiPar* fDigiPar {nullptr}; + CbmTofDigiBdfPar* fTofDigiBdfPar {nullptr}; + + TFile* fPerfFile {nullptr}; + TDirectory* fHistoDir {nullptr}; + + static const int fNTimeHistos = 22; + TH1F* fTimeHisto[fNTimeHistos] {nullptr}; + + static const int fNGhostHistos = 9; + TH1F* fGhostHisto[fNGhostHistos] {nullptr}; + + /// Used data = Repacked input data + L1Vector<CbmL1Hit> vStsHits {"CbmL1::vStsHits"}; // hits with hit-mcpoint match information + L1Vector<int> SortedIndex {"CbmL1::SortedIndex"}; + L1Vector<int> StsIndex {"CbmL1::StsIndex"}; + L1Vector<CbmL1MCTrack> vMCTracks {"CbmL1::vMCTracks"}; + L1Vector<int> vHitMCRef { + "CbmL1::vHitMCRef"}; // indices of MCPoints in vMCPoints, indexed by index of hit in algo->vStsHits array. According to StsMatch. Used for IdealResponce + // L1Vector<int> vHitMCRef1; + // CbmMatch stsHitMatch; + + DFEI2I dFEI2vMCPoints {}; + DFEI2I dFEI2vMCTracks {}; + + int fFindParticlesMode {0}; // 0 - don't run FindParticles + // 1 - run, all MC particle is reco-able + // 2 - run, MC particle is reco-able if created from reco-able tracks + // 3 - run, MC particle is reco-able if created from reconstructed tracks + + TString fStsMatBudgetFileName {}; + TString fMvdMatBudgetFileName {}; + TString fMuchMatBudgetFileName {}; + TString fTrdMatBudgetFileName {}; + TString fTofMatBudgetFileName {}; + bool fExtrapolateToTheEndOfSTS {false}; + int fTimesliceMode {0}; + + KFTopoPerformance* fTopoPerformance {nullptr}; + L1EventEfficiencies fEventEfficiency {}; // average efficiencies CbmStsParSetSensor* fStsParSetSensor {nullptr}; CbmStsParSetSensorCond* fStsParSetSensorCond {nullptr}; CbmStsParSetModule* fStsParSetModule {nullptr}; - ClassDef(CbmL1, 1); + ClassDef(CbmL1, 0); }; #endif //_CbmL1_h_ diff --git a/reco/L1/CbmL1CATrdTrackFinderSA.cxx b/reco/L1/CbmL1CATrdTrackFinderSA.cxx index 46ff95b8e7..1906f891d8 100644 --- a/reco/L1/CbmL1CATrdTrackFinderSA.cxx +++ b/reco/L1/CbmL1CATrdTrackFinderSA.cxx @@ -14,7 +14,6 @@ #include "CbmKFTrack.h" #include "CbmKFTrdHit.h" #include "CbmL1.h" -#include "CbmL1Def.h" #include "CbmL1TrdTracklet.h" #include "CbmL1TrdTracklet4.h" #include "CbmMCTrack.h" @@ -45,6 +44,8 @@ #include <cmath> +#include "L1Def.h" + using std::cout; using std::endl; using std::fabs; diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx index ebf5250e8f..4a2bca850c 100644 --- a/reco/L1/CbmL1ReadEvent.cxx +++ b/reco/L1/CbmL1ReadEvent.cxx @@ -43,12 +43,10 @@ #include "TRandom.h" #include <iostream> -#include <vector> using std::cout; using std::endl; using std::find; -using std::vector; //#define MVDIDEALHITS //#define STSIDEALHITS @@ -99,9 +97,28 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, vHitStore.clear(); dFEI2vMCPoints.clear(); dFEI2vMCTracks.clear(); + if (fVerbose >= 10) cout << "ReadEvent: clear is done." << endl; - vector<TmpHit> tmpHits; + L1Vector<TmpHit> tmpHits("CbmL1ReadEvent::tmpHits"); + + { // reserve enough space for hits + int nHitsTotal = 0; + if (listMvdHits) nHitsTotal += listMvdHits->GetEntriesFast(); + Int_t nEntSts = 0; + if (listStsHits) { + if (!fTimesliceMode && event) { nEntSts = event->GetNofData(ECbmDataType::kStsHit); } + else { + nEntSts = listStsHits->GetEntriesFast(); + } + if (nEntSts < 0) nEntSts = 0; // GetNofData() can return -1; + } + nHitsTotal += nEntSts; + if (fMuchPixelHits) nHitsTotal += fMuchPixelHits->GetEntriesFast(); + if (listTrdHits) nHitsTotal += listTrdHits->GetEntriesFast(); + if (fTofHits) nHitsTotal += fTofHits->GetEntriesFast(); + tmpHits.reserve(nHitsTotal); + } // -- produce Sts hits from space points -- @@ -126,10 +143,14 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, nMuchPoints = 0; int nTofPoints = 0; - vector<CbmLink*> ToFPointsMatch; + L1Vector<CbmLink*> ToFPointsMatch("CbmL1ReadEvent::ToFPointsMatch"); if (fPerformance) { Fill_vMCTracks(); + vMCPoints.clear(); + vMCPoints.reserve(5 * vMCTracks.size() * NStation); + vMCPoints_in_Time_Slice.clear(); + vMCPoints_in_Time_Slice.reserve(vMCPoints.capacity()); for (DFSET::iterator set_it = vFileEvent.begin(); set_it != vFileEvent.end(); ++set_it) { Int_t iFile = set_it->first; @@ -261,7 +282,7 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } } - ToFPointsMatch.resize(0); + ToFPointsMatch.clear(); if (fTofPoints) { for (int j = 0; j < fTofHits->GetEntriesFast(); j++) { @@ -369,7 +390,6 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } //fPerformance - int NStrips = 0; // get MVD hits @@ -476,7 +496,6 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } // if listMvdHits if (fVerbose >= 10) cout << "ReadEvent: mvd hits are gotten." << endl; - Int_t nEntSts = 0; if (listStsHits) { @@ -600,7 +619,6 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } // for j } // if listStsHits - if (fMuchPixelHits) { Int_t nEnt = fMuchPixelHits->GetEntriesFast(); @@ -802,6 +820,7 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } // for j } // if listTrdHits + if (fTofHits) { int firstDetStrip = NStrips; @@ -912,9 +931,7 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, // save strips in L1Algo fData_->NStsStrips = NStrips; - fData_->fStripFlag.clear(); - fData_->fStripFlag.reserve(NStrips); - fData_->fStripFlag.assign(NStrips, 0); + fData_->fStripFlag.reset(NStrips, 0); for (int ih = 0; ih < nHits; ih++) { TmpHit& th = tmpHits[ih]; char flag = th.iStation * 4; @@ -926,15 +943,21 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, // -- save hits -- int nEffHits = 0; - int maxHitIndex = 0; - if (fTofHits) maxHitIndex = fTofHits->GetEntriesFast() + nMvdHits + nStsHits + nMuchHits + nTrdHits; - else - maxHitIndex = nMvdHits + nStsHits + nMuchHits + nTrdHits; + int maxHitIndex = nMvdHits + nStsHits + nMuchHits + nTrdHits; + if (fTofHits) maxHitIndex += fTofHits->GetEntriesFast(); + + SortedIndex.reset(max(nEntSts, maxHitIndex)); + + L1Vector<float> vStsZPos_temp( + "CbmL1ReadEvent::vStsZPos_temp"); // temp array for unsorted z positions of detectors segments + vStsZPos_temp.reserve(100 * NStation); - SortedIndex.resize(max(nEntSts, maxHitIndex)); + vStsHits.reserve(nHits); + fData_->vStsHits.reserve(nHits); + vHitStore.reserve(nHits); + vHitMCRef.reserve(nHits); - vector<float> vStsZPos_temp; // temp array for unsorted z positions of detectors segments for (int i = 0; i < nHits; i++) { TmpHit& th = tmpHits[i]; @@ -1076,19 +1099,19 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, if (fVerbose >= 10) cout << "ReadEvent: mvd and sts are saved." << endl; - // sort z-pos if (vStsZPos_temp.size() != 0) { - vector<float> vStsZPos_temp2; - vStsZPos_temp2.clear(); + L1Vector<float> vStsZPos_temp2("CbmL1ReadEvent::vStsZPos_temp2"); + vStsZPos_temp2.reserve(vStsZPos_temp.size()); vStsZPos_temp2.push_back(vStsZPos_temp[0]); - vector<int> newToOldIndex; - newToOldIndex.clear(); + + L1Vector<int> newToOldIndex("CbmL1ReadEvent::newToOldIndex"); + newToOldIndex.reserve(vStsZPos_temp.size()); newToOldIndex.push_back(0); for (unsigned int k = 1; k < vStsZPos_temp.size(); k++) { - vector<float>::iterator itpos = vStsZPos_temp2.begin() + 1; - vector<int>::iterator iti = newToOldIndex.begin() + 1; + L1Vector<float>::iterator itpos = vStsZPos_temp2.begin() + 1; + L1Vector<int>::iterator iti = newToOldIndex.begin() + 1; for (; itpos < vStsZPos_temp2.end(); itpos++, iti++) { if (vStsZPos_temp[k] < *itpos) { vStsZPos_temp2.insert(itpos, vStsZPos_temp[k]); @@ -1102,22 +1125,22 @@ void CbmL1::ReadEvent(L1AlgoInputData* fData_, float& TsStart, float& TsLength, } } // k - if (fVerbose >= 10) cout << "ReadEvent: z-pos are sorted." << endl; - for (unsigned int k = 0; k < vStsZPos_temp2.size(); k++) + fData_->vStsZPos.reserve(vStsZPos_temp2.size()); + for (unsigned int k = 0; k < vStsZPos_temp2.size(); k++) { fData_->vStsZPos.push_back(vStsZPos_temp2[k]); + } int size_nto_tmp = newToOldIndex.size(); - vector<int> oldToNewIndex; - oldToNewIndex.clear(); - oldToNewIndex.resize(size_nto_tmp); - for (int k = 0; k < size_nto_tmp; k++) + L1Vector<int> oldToNewIndex("CbmL1ReadEvent::oldToNewIndex", size_nto_tmp); + for (int k = 0; k < size_nto_tmp; k++) { oldToNewIndex[newToOldIndex[k]] = k; - + } int size_hs_tmp = vHitStore.size(); - for (int k = 0; k < size_hs_tmp; k++) + for (int k = 0; k < size_hs_tmp; k++) { fData_->vStsHits[k].iz = oldToNewIndex[fData_->vStsHits[k].iz]; + } } if (fVerbose >= 10) cout << "ReadEvent: z-pos are saved." << endl; @@ -1145,15 +1168,22 @@ void CbmL1::Fill_vMCTracks() int nvtracks = 0, nvtrackscurr = 0; vMCTracks.clear(); + { + Int_t nMCTracks = 0; + for (DFSET::iterator set_it = vFileEvent.begin(); set_it != vFileEvent.end(); ++set_it) { + Int_t iFile = set_it->first; + Int_t iEvent = set_it->second; + nMCTracks += fMCTracks->Size(iFile, iEvent); + } + vMCTracks.reserve(nMCTracks); + } for (DFSET::iterator set_it = vFileEvent.begin(); set_it != vFileEvent.end(); ++set_it) { Int_t iFile = set_it->first; Int_t iEvent = set_it->second; - Int_t nMCTrack = fMCTracks->Size(iFile, iEvent); - for (Int_t iMCTrack = 0; iMCTrack < nMCTrack; iMCTrack++) { CbmMCTrack* MCTrack = L1_DYNAMIC_CAST<CbmMCTrack*>(fMCTracks->Get(iFile, iEvent, iMCTrack)); if (!MCTrack) continue; @@ -1374,9 +1404,6 @@ void CbmL1::HitMatch() int iP = -1; - vector<int> iEvent1; - - if (listStsClusterMatch) { const CbmMatch* frontClusterMatch = @@ -1423,9 +1450,6 @@ void CbmL1::HitMatch() Int_t iEvent = stsHitMatch.GetLink(iLink).GetEntry(); Int_t iIndex = stsHitMatch.GetLink(iLink).GetIndex(); - iEvent1.push_back(iEvent); - - if (!fTimesliceMode) { iFile = vFileEvent.begin()->first; iEvent = vFileEvent.begin()->second; @@ -1450,10 +1474,6 @@ void CbmL1::HitMatch() if (iP >= 0) { - for (unsigned int i = 0; i < iEvent1.size(); i++) { - hit.event = iEvent1[i]; - } - hit.event = vMCPoints[iP].event; hit.mcPointIds.push_back(iP); vMCPoints[iP].hitIds.push_back(iH); diff --git a/reco/L1/CbmL1TofMerger.cxx b/reco/L1/CbmL1TofMerger.cxx index a8960e1543..1c56302daf 100644 --- a/reco/L1/CbmL1TofMerger.cxx +++ b/reco/L1/CbmL1TofMerger.cxx @@ -10,7 +10,6 @@ #include "CbmGlobalTrack.h" #include "CbmKFTrack.h" -#include "CbmL1Def.h" #include "CbmTofHit.h" #include "CbmTrdTrack.h" @@ -23,6 +22,8 @@ #include <map> #include <utility> +#include "L1Def.h" + using std::cout; using std::endl; using std::make_pair; diff --git a/reco/L1/CbmL1TrackMerger.cxx b/reco/L1/CbmL1TrackMerger.cxx index 4fc07b4298..23d1c818c5 100644 --- a/reco/L1/CbmL1TrackMerger.cxx +++ b/reco/L1/CbmL1TrackMerger.cxx @@ -10,7 +10,6 @@ #include "CbmGlobalTrack.h" #include "CbmKFTrack.h" -#include "CbmL1Def.h" #include "CbmStsTrack.h" #include "CbmTrackMatch.h" #include "CbmTrdTrack.h" @@ -24,6 +23,8 @@ #include <iostream> #include <map> +#include "L1Def.h" + using std::cout; using std::endl; using std::map; diff --git a/reco/L1/CbmL1TrdTrackFinderSts.cxx b/reco/L1/CbmL1TrdTrackFinderSts.cxx index 9311c1a095..2c40c4b6f0 100644 --- a/reco/L1/CbmL1TrdTrackFinderSts.cxx +++ b/reco/L1/CbmL1TrdTrackFinderSts.cxx @@ -11,7 +11,6 @@ #include "CbmKFTrack.h" #include "CbmKFTrdHit.h" -#include "CbmL1Def.h" #include "CbmStsTrack.h" #include "CbmTrackMatch.h" #include "CbmTrdHit.h" @@ -37,6 +36,8 @@ #include <map> #include <vector> +#include "L1Def.h" + using std::cout; using std::endl; using std::map; diff --git a/reco/L1/L1Algo/L1Algo.cxx b/reco/L1/L1Algo/L1Algo.cxx index cbbc855c5f..81233e32dd 100644 --- a/reco/L1/L1Algo/L1Algo.cxx +++ b/reco/L1/L1Algo/L1Algo.cxx @@ -7,102 +7,40 @@ #include "L1Grid.h" #include "L1HitPoint.h" -L1Algo::L1Algo(int nThreads, int ExpectedHits) - : NStations(0) // number of all detector stations - , NMvdStations(0) // number of mvd stations - , NStsStations(0) - , NFStations(0) - , fRadThick() - , NStsStrips(0) // strips positions created from hits - , vStsZPos(0) - , // all possible z-positions of hits - vStsHits(nullptr) // hits as a combination of front-, backstrips and z-position - , fStripFlag(nullptr) // information of hits station & using hits in tracks(), - , CATime(0) // time of trackfinding - , StsHitsStartIndex(nullptr) - , StsHitsStopIndex(nullptr) - , NHitsIsecAll(0) - , vStsDontUsedHits_A(ExpectedHits) - , vStsDontUsedHits_B(ExpectedHits) - , vStsDontUsedHits_Buf(ExpectedHits) - , vStsDontUsedHitsxy_A(ExpectedHits) - , vStsDontUsedHitsxy_buf(ExpectedHits) - , vStsDontUsedHitsxy_B(ExpectedHits) - , RealIHit_v(ExpectedHits) - , RealIHit_v_buf(ExpectedHits) - , RealIHit_v_buf2(ExpectedHits) - , fNThreads(nThreads) - , fUseHitErrors(0) - , fmCBMmode(0) - , fGlobal(0) - , isec(0) - , vStsHitsUnused() - , RealIHitP() - , RealIHitPBuf() - , vStsHitPointsUnused() - , RealIHit(nullptr) - , FIRSTCASTATION() - , threadNumberToCpuMap() - , TRACK_CHI2_CUT(10.) - , TRIPLET_CHI2_CUT(5.) - , DOUBLET_CHI2_CUT(5.) - , TIME_CUT1(0.) - , TIME_CUT2(0.) - , MaxDZ(0.) - , - Pick_gather(0) - , PickNeighbour(0) - , // (PickNeighbour < dp/dp_error) => triplets are neighbours - MaxInvMom(0) - , // max considered q/p for tracks - MaxSlope(0) - , targX(0) - , targY(0) - , targZ(0) - , // target coor - targB() - , // field in the target point - TargetXYInfo() - , // target constraint [cm] - vtxFieldRegion() - , // really doesn't used - vtxFieldValue() - , // field at teh vertex position. - //vTripletsP(), // container for triplets got in finding - fTrackingLevel(0) - , fGhostSuppression(0) - , // really doesn't used - fMomentumCutOff(0) // really doesn't used +L1Algo::L1Algo(int nThreads) { + SetNThreads(nThreads); + for (int i = 0; i < MaxNStations; i++) { + vGridTime[i].AllocateMemory(fNThreads); + } +} - fDupletPortionSize.reserve(100000); - fTracks.reserve(40000); - fRecoHits.reserve(400000); - - fStripToTrack.reserve(ExpectedHits * 4); - fStripToTrackB.reserve(ExpectedHits * 4); +void L1Algo::SetNThreads(int n) +{ + if (n > kMaxNthreads) { + LOG(FATAL) << "L1Algo: n threads " << n << " is greater than the maximum " << kMaxNthreads << std::endl; + } + fNThreads = n; for (int i = 0; i < fNThreads; i++) { fTracks_local[i].SetName(std::stringstream() << "L1Algo::fTracks_local[" << i << "]"); - fTracks_local[i].clear(); - fTracks_local[i].reserve(100000); fRecoHits_local[i].SetName(std::stringstream() << "L1Algo::fRecoHits_local[" << i << "]"); - fRecoHits_local[i].clear(); - fRecoHits_local[i].reserve(400000); - - TripForHit[0].resize(ExpectedHits); - TripForHit[1].resize(ExpectedHits); fTrackCandidates[i].SetName(std::stringstream() << "L1Algo::fTrackCandidates[" << i << "]"); - fTrackCandidates[i].clear(); - fTrackCandidates[i].reserve(10000); + fT_3[i].reserve(MaxPortionTriplets / fvecLen); + + fhitsl_3[i].SetName(std::stringstream() << "L1Algo::fhitsl_3[" << i << "]"); + fhitsm_3[i].SetName(std::stringstream() << "L1Algo::fhitsm_3[" << i << "]"); + fhitsr_3[i].SetName(std::stringstream() << "L1Algo::fhitsr_3[" << i << "]"); + fhitsl_3[i].reserve(MaxPortionTriplets); fhitsm_3[i].reserve(MaxPortionTriplets); fhitsr_3[i].reserve(MaxPortionTriplets); + fu_front3[i].reserve(MaxPortionTriplets / fvecLen); fu_back3[i].reserve(MaxPortionTriplets / fvecLen); fz_pos3[i].reserve(MaxPortionTriplets / fvecLen); @@ -115,39 +53,12 @@ L1Algo::L1Algo(int nThreads, int ExpectedHits) for (int j = 0; j < MaxNStations; j++) { fTriplets[j][i].SetName(std::stringstream() << "L1Algo::fTriplets[" << i << "][" << j << "]"); - fTriplets[j][i].reserve(ExpectedHits); - fTriplets[j][i].clear(); } } - - for (int i = 0; i < MaxNStations; i++) { - vGridTime[i].AllocateMemory(fNThreads); - } - -#ifdef _OPENMP - fHitToBestTrackF.reserve(ExpectedHits * 2); - fHitToBestTrackB.reserve(ExpectedHits * 2); -#endif - - NHitsIsecAll = ExpectedHits; - const int kExpectedTracks = ExpectedHits / 8; - - fMergerTrackFirstStation.reserve(kExpectedTracks); - fMergerTrackLastStation.reserve(kExpectedTracks); - fMergerTrackFirstHit.reserve(kExpectedTracks); - fMergerTrackLastHit.reserve(kExpectedTracks); - fMergerTrackNeighbour.reserve(kExpectedTracks); - fMergerTrackChi2.reserve(kExpectedTracks); - fMergerTrackIsStored.reserve(kExpectedTracks); - fMergerTrackIsDownstreamNeighbour.reserve(kExpectedTracks); - fMergerTracksNew.reserve(kExpectedTracks); - fMergerRecoHitsNew.reserve(ExpectedHits); - - // IsNext.resize(kExpectedTracks); } -void L1Algo::Init(const vector<fscal>& geo, const bool UseHitErrors, const bool mCBMmode) +void L1Algo::Init(const L1Vector<fscal>& geo, const bool UseHitErrors, const bool mCBMmode) { for (int iProc = 0; iProc < 4; iProc++) { @@ -312,7 +223,7 @@ void L1Algo::Init(const vector<fscal>& geo, const bool UseHitErrors, const bool } -void L1Algo::SetData(vector<L1Hit>& StsHits_, int nStsStrips_, const vector<fscal>& StsZPos_, +void L1Algo::SetData(L1Vector<L1Hit>& StsHits_, int nStsStrips_, L1Vector<fscal>& StsZPos_, L1Vector<unsigned char>& SFlag_, const THitI* StsHitsStartIndex_, const THitI* StsHitsStopIndex_) { @@ -328,31 +239,36 @@ void L1Algo::SetData(vector<L1Hit>& StsHits_, int nStsStrips_, const vector<fsca int nHits = vStsHits->size(); - vStsDontUsedHits_A.resize(nHits); - vStsDontUsedHits_B.resize(nHits); - vStsDontUsedHits_Buf.resize(nHits); - vStsDontUsedHitsxy_A.resize(nHits); - vStsDontUsedHitsxy_buf.resize(nHits); - vStsDontUsedHitsxy_B.resize(nHits); - RealIHit_v.resize(nHits); - RealIHit_v_buf.resize(nHits); - RealIHit_v_buf2.resize(nHits); + NHitsIsecAll = nHits; + + vStsDontUsedHits_A.reset(nHits); + vStsDontUsedHits_B.reset(nHits); + vStsDontUsedHits_Buf.reset(nHits); + vStsDontUsedHitsxy_A.reset(nHits); + vStsDontUsedHitsxy_buf.reset(nHits); + vStsDontUsedHitsxy_B.reset(nHits); + RealIHit_v.reset(nHits); + RealIHit_v_buf.reset(nHits); + RealIHit_v_buf2.reset(nHits); #ifdef _OPENMP - fHitToBestTrackF.resize(NStsStrips); - fHitToBestTrackB.resize(NStsStrips); + fHitToBestTrackF.reset(NStsStrips); + fHitToBestTrackB.reset(NStsStrips); for (unsigned int j = 0; j < fHitToBestTrackB.size(); j++) { omp_init_lock(&fHitToBestTrackB[j]); omp_init_lock(&fHitToBestTrackF[j]); } #endif + fStripToTrack.clear(); fStripToTrack.reserve(NStsStrips); + + fStripToTrackB.clear(); fStripToTrackB.reserve(NStsStrips); - TripForHit[0].resize(nHits); - TripForHit[1].resize(nHits); - NHitsIsecAll = nHits; + TripForHit[0].reset(nHits); + TripForHit[1].reset(nHits); + fDupletPortionSize.clear(); fDupletPortionSize.reserve(2 * nHits); @@ -368,24 +284,6 @@ void L1Algo::SetData(vector<L1Hit>& StsHits_, int nStsStrips_, const vector<fsca fTriplets[j][i].reserve(2 * nHits); } } - - /* - vStsHits.resize(StsHits_.size()); - vStsStrips.resize(StsStrips_.size()); - vStsStripsB.resize(StsStripsB_.size()); - vStsZPos.resize(StsZPos_.size()); - fStripFlag.resize(SFlag_.size()); - fStripFlagB.resize(SFlagB_.size()); - - for(Tindex i=0; i< static_cast<Tindex>(StsHits_.size()); ++i ) vStsHits[i] = StsHits_[i]; - for(Tindex i=0; i< static_cast<Tindex>(StsStrips_.size()); ++i ) vStsStrips[i] = StsStrips_[i]; - for(Tindex i=0; i< static_cast<Tindex>(StsStripsB_.size()); ++i ) vStsStripsB[i] = StsStripsB_[i]; - for(Tindex i=0; i< static_cast<Tindex>(StsZPos_.size()); ++i ) vStsZPos[i] = StsZPos_[i]; - for(Tindex i=0; i< static_cast<Tindex>(SFlag_.size()); ++i ) fStripFlag[i] = SFlag_[i]; - for(Tindex i=0; i< static_cast<Tindex>(SFlagB_.size()); ++i ) fStripFlagB[i] = SFlagB_[i]; - - for(Tindex i=0; i<MaxNStations+1; ++i) StsHitsStartIndex[i] = StsHitsStartIndex_[i]; - for(Tindex i=0; i<MaxNStations+1; ++i) StsHitsStopIndex[i] = StsHitsStopIndex_[i];*/ } diff --git a/reco/L1/L1Algo/L1Algo.h b/reco/L1/L1Algo/L1Algo.h index a85bb53507..16833a2976 100644 --- a/reco/L1/L1Algo/L1Algo.h +++ b/reco/L1/L1Algo/L1Algo.h @@ -40,7 +40,6 @@ class L1AlgoDraw; #include <iomanip> #include <iostream> #include <map> -#include <vector> #include "L1Branch.h" #include "L1Field.h" @@ -61,7 +60,6 @@ class L1AlgoDraw; #endif using std::map; -using std::vector; #ifdef PULLS #define TRIP_PERFORMANCE @@ -80,8 +78,7 @@ typedef int Tindex; class L1Algo { public: - // L1Algo(int nThreads=7): - L1Algo(int nThreads = 1, int ExpectedHits = 200000); + L1Algo(int nThreads = 1); L1Algo(const L1Algo&) = delete; L1Algo operator=(const L1Algo&) = delete; @@ -97,18 +94,18 @@ public: /// get default particle mass squared float GetDefaultParticleMass2() const { return fDefaultMass * fDefaultMass; } - static const int nTh = 1; - static const int nSta = 25; + static const int kMaxNthreads = 1; + static const int nSta = 25; float fDefaultMass = 0.10565800; // muon mass - L1Vector<L1Triplet> fTriplets[nSta][nTh]; // created triplets at station + thread + L1Vector<L1Triplet> fTriplets[nSta][kMaxNthreads] {{"L1Algo::fTriplets"}}; // created triplets at station + thread // Track candidates created out of adjacent triplets before the final track selection. // The candidates may share any amount of hits. - L1Vector<L1Branch> fTrackCandidates[nTh]; + L1Vector<L1Branch> fTrackCandidates[kMaxNthreads] {"L1Algo::fTrackCandidates"}; - Tindex fDupletPortionStopIndex[nSta]; // end of the duplet portions for the station + Tindex fDupletPortionStopIndex[nSta] {0}; // end of the duplet portions for the station L1Vector<Tindex> fDupletPortionSize {"L1Algo::fDupletPortionSize"}; // N duplets in a portion // @@ -132,13 +129,13 @@ public: #ifdef DRAW L1AlgoDraw* draw {nullptr}; - void DrawRecoTracksTime(const vector<CbmL1Track>& tracks); + void DrawRecoTracksTime(const L1Vector<CbmL1Track>& tracks); #endif - void Init(const vector<fscal>& geo, const bool UseHitErrors, const bool mCBMmode); + void Init(const L1Vector<fscal>& geo, const bool UseHitErrors, const bool mCBMmode); - void SetData(vector<L1Hit>& StsHits_, int nStsStrips_, const vector<fscal>& StsZPos_, L1Vector<unsigned char>& SFlag_, + void SetData(L1Vector<L1Hit>& StsHits_, int nStsStrips_, L1Vector<fscal>& StsZPos_, L1Vector<unsigned char>& SFlag_, const THitI* StsHitsStartIndex_, const THitI* StsHitsStopIndex_); void PrintHits(); @@ -155,51 +152,53 @@ public: /// ----- Input data ----- // filled in CbmL1::ReadEvent(); - void SetNThreads(int n = 1) { fNThreads = n; } + void SetNThreads(int n); enum { MaxNStations = 25 }; - int NStations, // number of all detector stations - NMvdStations, // number of mvd stations - NStsStations, NFStations; + int NStations {0}; // number of all detector stations + int NMvdStations {0}; // number of mvd stations + int NStsStations {0}; // number of sts stations + int NFStations {0}; // ? + L1Station vStations[MaxNStations] _fvecalignment; // station info - vector<L1Material> fRadThick; // material for each station + L1Vector<L1Material> fRadThick {"fRadThick"}; // material for each station - int NStsStrips; // number of strips - const vector<fscal>* vStsZPos; // all possible z-positions of hits - vector<L1Hit>* vStsHits; // hits as a combination of front-, backstrips and z-position - L1Grid vGrid[MaxNStations]; // hits as a combination of front-, backstrips and z-position + int NStsStrips {0}; // number of strips + L1Vector<fscal>* vStsZPos {nullptr}; // all possible z-positions of hits + L1Vector<L1Hit>* vStsHits {nullptr}; // hits as a combination of front-, backstrips and z-position + L1Grid vGrid[MaxNStations]; // hits as a combination of front-, backstrips and z-position L1Grid vGridTime[MaxNStations]; - L1Vector<unsigned char>* fStripFlag; // information of hits station & using hits in tracks; + L1Vector<unsigned char>* fStripFlag {nullptr}; // information of hits station & using hits in tracks; - double CATime; // time of trackfinding + double CATime {0.}; // time of trackfinding L1Vector<L1Track> fTracks {"L1Algo::fTracks"}; // reconstructed tracks L1Vector<THitI> fRecoHits {"L1Algo::fRecoHits"}; // packed hits of reconstructed tracks - const THitI *StsHitsStartIndex, - *StsHitsStopIndex; // station-bounders in vStsHits array + const THitI* StsHitsStartIndex {nullptr}; // station-bounders in vStsHits array + const THitI* StsHitsStopIndex {nullptr}; // station-bounders in vStsHits array // L1Branch* pointer; - unsigned int NHitsIsecAll; + unsigned int NHitsIsecAll {0}; - vector<L1Hit> vStsDontUsedHits_A; - vector<L1Hit> vStsDontUsedHits_B; - vector<L1Hit> vStsDontUsedHits_Buf; - vector<L1HitPoint> vStsDontUsedHitsxy_A; - vector<L1HitPoint> vStsDontUsedHitsxy_buf; - vector<L1HitPoint> vStsDontUsedHitsxy_B; - L1Vector<L1Track> fTracks_local[nTh]; - L1Vector<THitI> fRecoHits_local[nTh]; + L1Vector<L1Hit> vStsDontUsedHits_A {"L1Algo::vStsDontUsedHits_A"}; + L1Vector<L1Hit> vStsDontUsedHits_B {"L1Algo::vStsDontUsedHits_B"}; + L1Vector<L1Hit> vStsDontUsedHits_Buf {"L1Algo::vStsDontUsedHits_Buf"}; + L1Vector<L1HitPoint> vStsDontUsedHitsxy_A {"L1Algo::vStsDontUsedHitsxy_A"}; + L1Vector<L1HitPoint> vStsDontUsedHitsxy_buf {"L1Algo::vStsDontUsedHitsxy_buf"}; + L1Vector<L1HitPoint> vStsDontUsedHitsxy_B {"L1Algo::vStsDontUsedHitsxy_B"}; + L1Vector<L1Track> fTracks_local[kMaxNthreads] {"L1Algo::fTracks_local"}; + L1Vector<THitI> fRecoHits_local[kMaxNthreads] {"L1Algo::fRecoHits_local"}; - vector<THitI> RealIHit_v; - vector<THitI> RealIHit_v_buf; - vector<THitI> RealIHit_v_buf2; + L1Vector<THitI> RealIHit_v {"L1Algo::RealIHit_v"}; + L1Vector<THitI> RealIHit_v_buf {"L1Algo::RealIHit_v_buf"}; + L1Vector<THitI> RealIHit_v_buf2 {"L1Algo::RealIHit_v_buf2"}; #ifdef _OPENMP L1Vector<omp_lock_t> fHitToBestTrackF {"L1Algo::fHitToBestTrackF"}; @@ -209,13 +208,13 @@ public: L1Vector<int> fStripToTrack {"L1Algo::fStripToTrack"}; // front strip to track pointers L1Vector<int> fStripToTrackB {"L1Algo::fStripToTrackB"}; // back strip to track pointers - int fNThreads; - bool fUseHitErrors; - bool fmCBMmode; - bool fGlobal; + int fNThreads {0}; + bool fUseHitErrors {0}; + bool fmCBMmode {0}; + bool fGlobal {0}; - fvec EventTime[nTh][nTh]; - fvec Err[nTh][nTh]; + fvec EventTime[kMaxNthreads][kMaxNthreads] {{0}}; + fvec Err[kMaxNthreads][kMaxNthreads] {{0}}; /// standard sizes of the arrays @@ -240,29 +239,41 @@ public: /// --- data used during finding iterations - int isec; // iteration - vector<L1Hit>* vStsHitsUnused; - vector<THitI>* RealIHitP; - vector<THitI>* RealIHitPBuf; - vector<L1HitPoint>* vStsHitPointsUnused; - THitI* RealIHit; // index in vStsHits indexed by index in vStsHitsUnused - THitI StsHitsUnusedStartIndex[MaxNStations + 1], StsHitsUnusedStopIndex[MaxNStations + 1]; - THitI StsHitsUnusedStartIndexEnd[MaxNStations + 1], StsHitsUnusedStopIndexEnd[MaxNStations + 1]; + int isec {0}; // iteration + L1Vector<L1Hit>* vStsHitsUnused {nullptr}; + L1Vector<THitI>* RealIHitP {nullptr}; + L1Vector<THitI>* RealIHitPBuf {nullptr}; + L1Vector<L1HitPoint>* vStsHitPointsUnused {nullptr}; + THitI* RealIHit {nullptr}; // index in vStsHits indexed by index in vStsHitsUnused + + THitI StsHitsUnusedStartIndex[MaxNStations + 1] {0}; + THitI StsHitsUnusedStopIndex[MaxNStations + 1] {0}; + THitI StsHitsUnusedStartIndexEnd[MaxNStations + 1] {0}; + THitI StsHitsUnusedStopIndexEnd[MaxNStations + 1] {0}; - vector<int> TripForHit[2]; + L1Vector<int> TripForHit[2] {"L1Algo::TripForHit"}; // fvec u_front[Portion/fvecLen], u_back[Portion/fvecLen]; // fvec zPos[Portion/fvecLen]; // fvec fHitTime[Portion/fvecLen]; - nsL1::vector<L1TrackPar>::TSimd fT_3[nTh]; + nsL1::vector<L1TrackPar>::TSimd fT_3[kMaxNthreads]; - vector<THitI> fhitsl_3[nTh], fhitsm_3[nTh], fhitsr_3[nTh]; + L1Vector<THitI> fhitsl_3[kMaxNthreads] {"L1Algo::fhitsl_3"}; + L1Vector<THitI> fhitsm_3[kMaxNthreads] {"L1Algo::fhitsm_3"}; + L1Vector<THitI> fhitsr_3[kMaxNthreads] {"L1Algo::fhitsr_3"}; - nsL1::vector<fvec>::TSimd fu_front3[nTh], fu_back3[nTh], fz_pos3[nTh], fTimeR[nTh], fTimeER[nTh], dx[nTh], dy[nTh], - du[nTh], dv[nTh]; + nsL1::vector<fvec>::TSimd fu_front3[kMaxNthreads]; + nsL1::vector<fvec>::TSimd fu_back3[kMaxNthreads]; + nsL1::vector<fvec>::TSimd fz_pos3[kMaxNthreads]; + nsL1::vector<fvec>::TSimd fTimeR[kMaxNthreads]; + nsL1::vector<fvec>::TSimd fTimeER[kMaxNthreads]; + nsL1::vector<fvec>::TSimd dx[kMaxNthreads]; + nsL1::vector<fvec>::TSimd dy[kMaxNthreads]; + nsL1::vector<fvec>::TSimd du[kMaxNthreads]; + nsL1::vector<fvec>::TSimd dv[kMaxNthreads]; // Tindex NHits_l[MaxNStations]; @@ -395,12 +406,12 @@ private: Tindex n1, L1Station& stam, L1HitPoint* vStsHits_m, L1TrackPar* T_1, THitI* hitsl_1, // output - Tindex& n2, vector<THitI>& i1_2, + Tindex& n2, L1Vector<THitI>& i1_2, #ifdef DOUB_PERFORMANCE - vector<THitI>& hitsl_2, + L1Vector<THitI>& hitsl_2, #endif // DOUB_PERFORMANCE - vector<THitI>& hitsm_2, fvec* Event, vector<bool>& lmDuplets); + L1Vector<THitI>& hitsm_2, fvec* Event, L1Vector<char>& lmDuplets); /// Add the middle hits to parameters estimation. Propagate to right station. /// Find the triplets (right hit). Reformat data in the portion of triplets. @@ -409,12 +420,12 @@ private: int istam, int istar, L1HitPoint* vStsHits_m, L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, - Tindex n2, vector<THitI>& hitsm_2, vector<THitI>& i1_2, + Tindex n2, L1Vector<THitI>& hitsm_2, L1Vector<THitI>& i1_2, - const vector<bool>& mrDuplets, + const L1Vector<char>& mrDuplets, // output - Tindex& n3, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, vector<THitI>& hitsm_3, - vector<THitI>& hitsr_3, nsL1::vector<fvec>::TSimd& u_front_3, nsL1::vector<fvec>::TSimd& u_back_3, + Tindex& n3, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, L1Vector<THitI>& hitsm_3, + L1Vector<THitI>& hitsr_3, nsL1::vector<fvec>::TSimd& u_front_3, nsL1::vector<fvec>::TSimd& u_back_3, nsL1::vector<fvec>::TSimd& z_Pos_3, // nsL1::vector<fvec>::TSimd& dx_, // nsL1::vector<fvec>::TSimd& dy_, @@ -434,15 +445,15 @@ private: /// Refit Triplets. void f32( // input - Tindex n3, int istal, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, vector<THitI>& hitsm_3, - vector<THitI>& hitsr_3, int nIterations = 0); + Tindex n3, int istal, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, L1Vector<THitI>& hitsm_3, + L1Vector<THitI>& hitsr_3, int nIterations = 0); /// Select triplets. Save them into vTriplets. void f4( // input - Tindex n3, int istal, int istam, int istar, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, - vector<THitI>& hitsm_3, vector<THitI>& hitsr_3, + Tindex n3, int istal, int istam, int istar, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, + L1Vector<THitI>& hitsm_3, L1Vector<THitI>& hitsr_3, // output - Tindex& nstaltriplets, vector<THitI>* hitsn_3 = 0, vector<THitI>* hitsr_5 = 0 + Tindex& nstaltriplets, L1Vector<THitI>* hitsn_3 = 0, L1Vector<THitI>* hitsr_5 = 0 // #ifdef XXX // ,unsigned int &stat_n_trip @@ -463,18 +474,18 @@ private: // output L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, - vector<bool>& lmDuplets, + L1Vector<char>& lmDuplets, - Tindex& n_2, vector<THitI>& i1_2, vector<THitI>& hitsm_2); + Tindex& n_2, L1Vector<THitI>& i1_2, L1Vector<THitI>& hitsm_2); /// Find triplets on station void TripletsStaPort( // input int istal, int istam, int istar, Tindex& nstaltriplets, L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, - Tindex& n_2, vector<THitI>& i1_2, vector<THitI>& hitsm_2, + Tindex& n_2, L1Vector<THitI>& i1_2, L1Vector<THitI>& hitsm_2, - const vector<bool>& mrDuplets + const L1Vector<char>& mrDuplets // output @@ -517,7 +528,7 @@ private: /// ----- Different parameters of CATrackFinder ----- - Tindex FIRSTCASTATION; //first station used in CA + Tindex FIRSTCASTATION {0}; //first station used in CA // fNFindIterations - set number of interation for trackfinding // itetation of finding: @@ -564,35 +575,39 @@ private: }; #endif // FIND_GAPED_TRACKS - map<int, int> threadNumberToCpuMap; + map<int, int> threadNumberToCpuMap {}; - float TRACK_CHI2_CUT; - float TRIPLET_CHI2_CUT; // = 5.0; // cut for selecting triplets before collecting tracks.per one DoF - float DOUBLET_CHI2_CUT; - float TIME_CUT1, TIME_CUT2; + float TRACK_CHI2_CUT {10.f}; + float TRIPLET_CHI2_CUT {5.f}; // cut for selecting triplets before collecting tracks.per one DoF + float DOUBLET_CHI2_CUT {5.f}; + float TIME_CUT1 {0.f}; + float TIME_CUT2 {0.f}; - fvec - MaxDZ; // correction in order to take into account overlaping and iff z. if sort by y then it is max diff between same station's modules (~0.4cm) + fvec MaxDZ { + 0.f}; // correction in order to take into account overlaping and iff z. if sort by y then it is max diff between same station's modules (~0.4cm) /// parameters which are different for different iterations. Set in the begin of CAL1TrackFinder - float Pick_gather; // same for attaching additional hits to track - float PickNeighbour; // (PickNeighbour < dp/dp_error) => triplets are neighbours - fvec MaxInvMom; // max considered q/p for tracks - fvec MaxSlope; // max slope (tx\ty) in prim vertex - fvec targX, targY, targZ; // target coor - L1FieldValue targB _fvecalignment; // field in the target point - L1XYMeasurementInfo TargetXYInfo _fvecalignment; // target constraint [cm] + float Pick_gather {0.f}; // same for attaching additional hits to track + float PickNeighbour {0.f}; // (PickNeighbour < dp/dp_error) => triplets are neighbours + fvec MaxInvMom {0.f}; // max considered q/p for tracks + fvec MaxSlope {0.f}; // max slope (tx\ty) in prim vertex + fvec targX {0.f}; // target coor TODO: set defaults to a crasy value to be sure that the target is initialised later + fvec targY {0.f}; + fvec targZ {0.f}; + L1FieldValue targB _fvecalignment {}; // field in the target point + L1XYMeasurementInfo TargetXYInfo _fvecalignment {}; // target constraint [cm] - L1FieldRegion vtxFieldRegion _fvecalignment; // really doesn't used - L1FieldValue vtxFieldValue _fvecalignment; // field at teh vertex position. + L1FieldRegion vtxFieldRegion _fvecalignment {}; // really doesn't used + L1FieldValue vtxFieldValue _fvecalignment {}; // field at teh vertex position. // int TripNumThread; - int fTrackingLevel, fGhostSuppression; // really doesn't used - float fMomentumCutOff; // really doesn't used + int fTrackingLevel {0}; // currently not used + int fGhostSuppression {0}; // currently not used + float fMomentumCutOff {0}; // currently not used /// ----- Debug features ----- #ifdef PULLS @@ -608,7 +623,6 @@ private: #ifdef DRAW friend class L1AlgoDraw; #endif - } _fvecalignment; #endif diff --git a/reco/L1/L1Algo/L1Branch.h b/reco/L1/L1Algo/L1Branch.h index f37ff197ad..f6967b3062 100644 --- a/reco/L1/L1Algo/L1Branch.h +++ b/reco/L1/L1Algo/L1Branch.h @@ -8,8 +8,7 @@ #ifndef L1Branch_h #define L1Branch_h -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1Hit.h" #include "L1Vector.h" diff --git a/reco/L1/L1Algo/L1CAMergeClones.cxx b/reco/L1/L1Algo/L1CAMergeClones.cxx index 97b82acc86..5b42bae5f4 100644 --- a/reco/L1/L1Algo/L1CAMergeClones.cxx +++ b/reco/L1/L1Algo/L1CAMergeClones.cxx @@ -25,7 +25,6 @@ // using namespace std; using std::cout; using std::endl; -using std::vector; void L1Algo::InvertCholetsky(fvec a[15]) { @@ -225,20 +224,20 @@ void L1Algo::CAMergeClones() fMergerRecoHitsNew.clear(); fMergerRecoHitsNew.reserve(fRecoHits.size()); - firstStation.resize(nTracks); - lastStation.resize(nTracks); - firstHit.resize(nTracks); - lastHit.resize(nTracks); - isStored.resize(nTracks); - trackChi2.resize(nTracks); - neighbour.resize(nTracks); - isDownstreamNeighbour.resize(nTracks); + firstStation.reset(nTracks); + lastStation.reset(nTracks); + firstHit.reset(nTracks); + lastHit.reset(nTracks); + isStored.reset(nTracks); + trackChi2.reset(nTracks); + neighbour.reset(nTracks); + isDownstreamNeighbour.reset(nTracks); THitI start_hit = 0; + #ifdef OMP #pragma omp parallel for #endif - for (int iTr = 0; iTr < nTracks; iTr++) { firstHit[iTr] = start_hit; firstStation[iTr] = (*fStripFlag)[(*vStsHits)[fRecoHits[start_hit]].f] / 4; @@ -417,13 +416,13 @@ void L1Algo::CAMergeClones() } } - fTracks.resize(fMergerTracksNew.size()); + fTracks.reset(fMergerTracksNew.size()); for (unsigned int iTr = 0; iTr < fMergerTracksNew.size(); iTr++) { fTracks[iTr] = fMergerTracksNew[iTr]; } assert(fRecoHits.size() == fMergerRecoHitsNew.size()); - fRecoHits.resize(fMergerRecoHitsNew.size()); + fRecoHits.reset(fMergerRecoHitsNew.size()); for (THitI iHi = 0; iHi < fMergerRecoHitsNew.size(); iHi++) { fRecoHits[iHi] = fMergerRecoHitsNew[iHi]; } diff --git a/reco/L1/L1Algo/L1CATrackFinder.cxx b/reco/L1/L1Algo/L1CATrackFinder.cxx index 661beb849c..5108201d93 100644 --- a/reco/L1/L1Algo/L1CATrackFinder.cxx +++ b/reco/L1/L1Algo/L1CATrackFinder.cxx @@ -69,7 +69,6 @@ // using namespace std; using std::cout; using std::endl; -using std::vector; /// Prepare the portion of data of left hits of a triplet: all hits except the last and the second last station will be procesesed in the algorythm, @@ -386,19 +385,15 @@ inline void L1Algo::f11( /// input 1st stage of singlet search /// Find the doublets. Reformat data in the portion of doublets. -inline void L1Algo::f20( // input - Tindex n1, L1Station& stam, - - L1HitPoint* vStsHits_m, L1TrackPar* T_1, THitI* hitsl_1, - - // output - Tindex& n2, vector<THitI>& i1_2, +inline void L1Algo::f20( + /// input + Tindex n1, L1Station& stam, L1HitPoint* vStsHits_m, L1TrackPar* T_1, THitI* hitsl_1, + /// output + Tindex& n2, L1Vector<THitI>& i1_2, #ifdef DOUB_PERFORMANCE - vector<THitI>& hitsl_2, + L1Vector<THitI>& hitsl_2, #endif // DOUB_PERFORMANCE - vector<THitI>& hitsm_2, - // comment unused parameters, FU, 18.01.21 - fvec* /*Event*/, vector<bool>& lmDuplets) + L1Vector<THitI>& hitsm_2, fvec* /*Event*/, L1Vector<char>& lmDuplets) { n2 = 0; // number of doublet for (Tindex i1 = 0; i1 < n1; ++i1) // for each singlet @@ -501,7 +496,6 @@ inline void L1Algo::f20( // input L1ExtrapolateC10Line(T1, zm, C10); fvec chi2 = T1.chi2; - L1UMeasurementInfo info = stam.frontInfo; if (fUseHitErrors) info.sigma2 = hitm.dU() * hitm.dU(); @@ -550,11 +544,11 @@ inline void L1Algo::f20( // input /// Find the triplets(right hit). Reformat data in the portion of triplets. inline void L1Algo::f30( // input L1HitPoint* vStsHits_r, L1Station& stam, L1Station& star, int istam, int istar, L1HitPoint* vStsHits_m, - L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, Tindex n2, vector<THitI>& hitsm_2, vector<THitI>& i1_2, - const vector<bool>& /*mrDuplets*/, + L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, Tindex n2, L1Vector<THitI>& hitsm_2, L1Vector<THitI>& i1_2, + const L1Vector<char>& /*mrDuplets*/, // output - Tindex& n3, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, vector<THitI>& hitsm_3, - vector<THitI>& hitsr_3, nsL1::vector<fvec>::TSimd& u_front_3, nsL1::vector<fvec>::TSimd& u_back_3, + Tindex& n3, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, L1Vector<THitI>& hitsm_3, + L1Vector<THitI>& hitsr_3, nsL1::vector<fvec>::TSimd& u_front_3, nsL1::vector<fvec>::TSimd& u_back_3, nsL1::vector<fvec>::TSimd& z_Pos_3, // nsL1::vector<fvec>::TSimd& dx_, // nsL1::vector<fvec>::TSimd& dy_, @@ -922,8 +916,8 @@ inline void L1Algo::f31( // input /// Refit Triplets. inline void L1Algo::f32( // input // TODO not updated after gaps introduction - Tindex n3, int istal, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, vector<THitI>& hitsm_3, - vector<THitI>& hitsr_3, int nIterations) + Tindex n3, int istal, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, L1Vector<THitI>& hitsm_3, + L1Vector<THitI>& hitsr_3, int nIterations) { L1Fit fit; fit.SetParticleMass(fDefaultMass); @@ -1083,10 +1077,10 @@ inline void L1Algo::f32( // input // TODO not updated after gaps introduction inline void L1Algo::f4( // input - Tindex n3, int istal, int istam, int istar, nsL1::vector<L1TrackPar>::TSimd& T_3, vector<THitI>& hitsl_3, - vector<THitI>& hitsm_3, vector<THitI>& hitsr_3, + Tindex n3, int istal, int istam, int istar, nsL1::vector<L1TrackPar>::TSimd& T_3, L1Vector<THitI>& hitsl_3, + L1Vector<THitI>& hitsm_3, L1Vector<THitI>& hitsr_3, // output - Tindex& nstaltriplets, vector<THitI>* /*hitsn_3*/, vector<THitI>* /*hitsr_5*/ + Tindex& nstaltriplets, L1Vector<THitI>* /*hitsn_3*/, L1Vector<THitI>* /*hitsr_5*/ ) { @@ -1238,8 +1232,8 @@ inline void L1Algo::f5( // input L1_ASSERT(ihitm < StsHitsUnusedStopIndex[istam], ihitm << " < " << StsHitsUnusedStopIndex[istam]); L1_ASSERT(ihitr < StsHitsUnusedStopIndex[istar], ihitr << " < " << StsHitsUnusedStopIndex[istar]); - vector<unsigned int> neighCands; // save neighbour candidates - neighCands.reserve(8); // ~average is 1-2 for central, up to 5 + L1Vector<unsigned int> neighCands("L1CATrackFinder::neighCands"); // save neighbour candidates + neighCands.reserve(8); // ~average is 1-2 for central, up to 5 unsigned int nNeighbours = TripForHit[1][ihitm] - TripForHit[0][ihitm]; @@ -1286,16 +1280,31 @@ inline void L1Algo::f5( // input /// ------------------- doublets on station ---------------------- -inline void L1Algo:: - DupletsStaPort( /// creates duplets: input: @istal - start station number, @istam - last station number, @ip - index of portion, @&n_g - number of elements in portion, @*portionStopIndex - int istal, int istam, Tindex ip, L1Vector<Tindex>& n_g, Tindex* portionStopIndex_, - /// output: - L1TrackPar* - T_1, /// @*T_1 - singlets perameters, @*fld_1 - field aproximation, @*hitsl_1- left hits of triplets, @&lmDuplets - existance of a doublet starting from the left hit, - L1FieldRegion* - fld_1, /// @&n_2 - number of douplets,@&i1_2 - index of 1st hit in portion indexed by doublet index, @&hitsm_2 - index of middle hit in hits array indexed by doublet index - THitI* hitsl_1, vector<bool>& lmDuplets, Tindex& n_2, vector<THitI>& i1_2, vector<THitI>& hitsm_2) +inline void L1Algo::DupletsStaPort( + /// input: + int istal, int istam, Tindex ip, L1Vector<Tindex>& n_g, Tindex* portionStopIndex_, + /// output: + L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, L1Vector<char>& lmDuplets, Tindex& n_2, L1Vector<THitI>& i1_2, + L1Vector<THitI>& hitsm_2 + /// +) { + /// creates duplets. + /// input: + /// @istal - start station number + /// @istam - last station number + /// @ip - index of portion + /// @&n_g - number of elements in portion + /// @*portionStopIndex + /// output: + /// @*T_1 - singlets perameters + /// @*fld_1 - field aproximation + /// @*hitsl_1- left hits of triplets + /// @&lmDuplets - existance of a doublet starting from the left hit + /// @&n_2 - number of douplets + /// @&i1_2 - index of 1st hit in portion indexed by doublet index + /// @&hitsm_2 - index of middle hit in hits array indexed by doublet index + if (istam < NStations) { L1Station& stam = vStations[istam]; @@ -1336,7 +1345,7 @@ inline void L1Algo:: /// Find the doublets. Reformat data in the portion of doublets. #ifdef DOUB_PERFORMANCE - vector<THitI> hitsl_2; + L1Vector<THitI> hitsl_2("L1CATrackFinder::hitsl_2"); #endif // DOUB_PERFORMANCE f20( // input @@ -1378,9 +1387,9 @@ L1Algo::TripletsStaPort( /// creates triplets: input: @istal - start station nu Tindex& nstaltriplets, L1TrackPar* T_1, L1FieldRegion* fld_1, THitI* hitsl_1, - Tindex& n_2, vector<THitI>& i1_2, vector<THitI>& hitsm_2, + Tindex& n_2, L1Vector<THitI>& i1_2, L1Vector<THitI>& hitsm_2, - const vector<bool>& mrDuplets + const L1Vector<char>& mrDuplets /// output: @*vTriplets_part - array of triplets, @*TripStartIndexH, @*TripStopIndexH - start/stop index of a triplet in the array @@ -1408,9 +1417,9 @@ L1Algo::TripletsStaPort( /// creates triplets: input: @istal - start station nu #endif nsL1::vector<L1TrackPar>::TSimd& T_3 = fT_3[Thread]; - vector<THitI>& hitsl_3 = fhitsl_3[Thread]; - vector<THitI>& hitsm_3 = fhitsm_3[Thread]; - vector<THitI>& hitsr_3 = fhitsr_3[Thread]; + L1Vector<THitI>& hitsl_3 = fhitsl_3[Thread]; + L1Vector<THitI>& hitsm_3 = fhitsm_3[Thread]; + L1Vector<THitI>& hitsr_3 = fhitsr_3[Thread]; nsL1::vector<fvec>::TSimd& u_front3 = fu_front3[Thread]; nsL1::vector<fvec>::TSimd& u_back3 = fu_back3[Thread]; nsL1::vector<fvec>::TSimd& z_pos3 = fz_pos3[Thread]; @@ -1616,18 +1625,21 @@ void L1Algo::CATrackFinder() static Tindex stat_nTrCandidates[fNFindIterations] = {0}; #endif - RealIHitP = &RealIHit_v; - RealIHitPBuf = &RealIHit_v_buf; - vStsHitsUnused = &vStsDontUsedHits_B; /// array of hits used on current iteration - vector<L1Hit>* vStsHitsUnused_buf = &vStsDontUsedHits_A; /// buffer for copy + RealIHitP = &RealIHit_v; + RealIHitPBuf = &RealIHit_v_buf; + vStsHitsUnused = &vStsDontUsedHits_B; /// array of hits used on current iteration + L1Vector<L1Hit>* vStsHitsUnused_buf = &vStsDontUsedHits_A; /// buffer for copy vStsHitPointsUnused = &vStsDontUsedHitsxy_B; /// array of info for hits used on current iteration - vector<L1HitPoint>* vStsHitPointsUnused_buf = &vStsDontUsedHitsxy_A; + L1Vector<L1HitPoint>* vStsHitPointsUnused_buf = &vStsDontUsedHitsxy_A; NHitsIsecAll = 0; fTracks.clear(); fRecoHits.clear(); + fRecoHits.reserve(2 * vStsHits->size()); + fTracks.reserve(2 * vStsHits->size() / NStations); + int nDontUsedHits = 0; // #pragma omp parallel for reduction(+:nDontUsedHits) @@ -1684,11 +1696,15 @@ void L1Algo::CATrackFinder() for (int iS = 0; iS < NStations; ++iS) { vGridTime[iS].BuildBins(-1, 1, -0.6, 0.6, starttime, lasttime, xStep, yStep, (lasttime - starttime + 1)); - - vGridTime[iS].StoreHits(StsHitsUnusedStopIndex[iS] - StsHitsUnusedStartIndex[iS], - &((*vStsHits)[StsHitsUnusedStartIndex[iS]]), iS, *this, StsHitsUnusedStartIndex[iS], - &(vStsDontUsedHits_Buf[StsHitsUnusedStartIndex[iS]]), - &((*vStsHits)[StsHitsUnusedStartIndex[iS]]), &(RealIHit_v[StsHitsUnusedStartIndex[iS]])); + int start = StsHitsUnusedStartIndex[iS]; + int nhits = StsHitsUnusedStopIndex[iS] - start; + if (nhits > 0) { + vGridTime[iS].StoreHits(nhits, &((*vStsHits)[start]), iS, *this, start, &(vStsDontUsedHits_Buf[start]), + &((*vStsHits)[start]), &(RealIHit_v[start])); + } + else { // to avoid out-of-range crash in array[start] + vGridTime[iS].StoreHits(nhits, nullptr, iS, *this, start, nullptr, nullptr, nullptr); + } } @@ -1730,7 +1746,7 @@ void L1Algo::CATrackFinder() // n_g1.assign(n_g1.size(), Portion); - for (int n = 0; n < nTh; n++) { + for (int n = 0; n < fNThreads; n++) { for (int j = 0; j < NStations; j++) { fTriplets[j][n].clear(); } @@ -1741,17 +1757,17 @@ void L1Algo::CATrackFinder() #endif if (isec != 0) { - vector<THitI>* RealIHitPTmp = RealIHitP; - RealIHitP = RealIHitPBuf; - RealIHitPBuf = RealIHitPTmp; + L1Vector<THitI>* RealIHitPTmp = RealIHitP; + RealIHitP = RealIHitPBuf; + RealIHitPBuf = RealIHitPTmp; - vector<L1Hit>* vStsHitsUnused_temp = vStsHitsUnused; - vStsHitsUnused = vStsHitsUnused_buf; - vStsHitsUnused_buf = vStsHitsUnused_temp; + L1Vector<L1Hit>* vStsHitsUnused_temp = vStsHitsUnused; + vStsHitsUnused = vStsHitsUnused_buf; + vStsHitsUnused_buf = vStsHitsUnused_temp; - vector<L1HitPoint>* vStsHitsUnused_temp2 = vStsHitPointsUnused; - vStsHitPointsUnused = vStsHitPointsUnused_buf; - vStsHitPointsUnused_buf = vStsHitsUnused_temp2; + L1Vector<L1HitPoint>* vStsHitsUnused_temp2 = vStsHitPointsUnused; + vStsHitPointsUnused = vStsHitPointsUnused_buf; + vStsHitPointsUnused_buf = vStsHitsUnused_temp2; } for (int ist = 0; ist < NStations; ++ist) { @@ -1925,19 +1941,28 @@ void L1Algo::CATrackFinder() L1FieldRegion fldG_1[vPortion]; THitI hitslG_1[Portion]; - vector<THitI> hitsm_2; /// middle hits indexed by number of doublets in portion(i2) - vector<THitI> i1_2; /// index in portion of singlets(i1) indexed by index in portion of doublets(i2) + L1Vector<THitI> hitsm_2("L1CATrackFinder::hitsm_2"); /// middle hits indexed by number of doublets in portion(i2) + L1Vector<THitI> i1_2( + "L1CATrackFinder::i1_2"); /// index in portion of singlets(i1) indexed by index in portion of doublets(i2) - vector<THitI> hitsmG_2; /// middle hits indexed by number of doublets in portion(i2) - vector<THitI> i1G_2; /// index in portion of singlets(i1) indexed by index in portion of doublets(i2) - vector<bool> lmDuplets[MaxNStations]; // is exist a doublet started from indexed by left hit - vector<bool> lmDupletsG[MaxNStations]; // is exist a doublet started from indexed by left hit - hitsm_2.reserve(3500); - i1_2.reserve(3500); + L1Vector<THitI> hitsmG_2("L1CATrackFinder::hitsmG_2"); /// middle hits indexed by number of doublets in portion(i2) + L1Vector<THitI> i1G_2( + "L1CATrackFinder::i1G_2"); /// index in portion of singlets(i1) indexed by index in portion of doublets(i2) + L1Vector<char> lmDuplets[MaxNStations] { + "L1CATrackFinder::lmDuplets"}; // is exist a doublet started from indexed by left hit + L1Vector<char> lmDupletsG[MaxNStations] { + "L1CATrackFinder::lmDupletsG"}; // is exist a doublet started from indexed by left hit - hitsmG_2.reserve(800); - i1G_2.reserve(800); + for (int i = 0; i < NStations; i++) { + lmDuplets[i].SetName(std::stringstream() << "L1CATrackFinder::lmDuplets[" << i << "]"); + lmDupletsG[i].SetName(std::stringstream() << "L1CATrackFinder::lmDupletsG[" << i << "]"); + } + + hitsm_2.reserve(9000); // TODO: make reasonable cuts on n combinations, put them to the header + i1_2.reserve(9000); // TODO: why that large numbers are needed even for mbias??? something goes wrong sometimes.. + hitsmG_2.reserve(9000); + i1G_2.reserve(9000); for (int istal = NStations - 2; istal >= FIRSTCASTATION; istal--) // //start downstream chambers { @@ -1949,8 +1974,8 @@ void L1Algo::CATrackFinder() for (Tindex ip = fDupletPortionStopIndex[istal + 1]; ip < fDupletPortionStopIndex[istal]; ++ip) { Tindex n_2; /// number of doublets in portion int NHitsSta = StsHitsStopIndex[istal] - StsHitsStartIndex[istal]; - lmDuplets[istal].resize(NHitsSta); - lmDupletsG[istal].resize(NHitsSta); + lmDuplets[istal].reset(NHitsSta); + lmDupletsG[istal].reset(NHitsSta); hitsm_2.clear(); i1_2.clear(); @@ -2075,12 +2100,8 @@ void L1Algo::CATrackFinder() fTrackCandidates[i].clear(); } - fStripToTrack.clear(); - fStripToTrack.resize(NStsStrips, -1); - fStripToTrackB.clear(); - fStripToTrackB.resize(NStsStrips, -1); - fStripToTrack.assign(NStsStrips, -1); - fStripToTrackB.assign(NStsStrips, -1); + fStripToTrack.reset(NStsStrips, -1); + fStripToTrackB.reset(NStsStrips, -1); for (int istaF = FIRSTCASTATION; istaF <= NStations - 3 - ilev; ++istaF) { @@ -2330,26 +2351,26 @@ void L1Algo::CATrackFinder() #endif int nTracks = fTracks.size(); - vector<int> offset_tracks(nTh, nTracks); - vector<int> offset_hits(nTh, NHitsIsecAll); + L1Vector<int> offset_tracks("L1CATrackFinder::offset_tracks", fNThreads, nTracks); + L1Vector<int> offset_hits("L1CATrackFinder::offset_hits", fNThreads, NHitsIsecAll); assert(NHitsIsecAll == fRecoHits.size()); //SG!! nTracks += fTracks_local[0].size(); NHitsIsecAll += fRecoHits_local[0].size(); - for (int i = 1; i < nTh; ++i) { + for (int i = 1; i < fNThreads; ++i) { offset_tracks[i] = offset_tracks[i - 1] + fTracks_local[i - 1].size(); offset_hits[i] = offset_hits[i - 1] + fRecoHits_local[i - 1].size(); nTracks += fTracks_local[i].size(); NHitsIsecAll += fRecoHits_local[i].size(); } - fTracks.resize(nTracks); - fRecoHits.resize(NHitsIsecAll); + fTracks.enlarge(nTracks); + fRecoHits.enlarge(NHitsIsecAll); #ifdef _OPENMP #pragma omp parallel for #endif - for (int i = 0; i < nTh; ++i) { + for (int i = 0; i < fNThreads; ++i) { for (Tindex iC = 0; iC < (Tindex) fTracks_local[i].size(); ++iC) { fTracks[offset_tracks[i] + iC] = fTracks_local[i][iC]; } @@ -2370,12 +2391,22 @@ void L1Algo::CATrackFinder() int NHitsOnStation = 0; for (int ista = 0; ista < NStations; ++ista) { - int Nelements = StsHitsUnusedStopIndex[ista] - StsHitsUnusedStartIndex[ista]; + int start = StsHitsUnusedStartIndex[ista]; + int Nelements = StsHitsUnusedStopIndex[ista] - start; + L1Hit* staHits = nullptr; // to avoid out-of-range error in ..[start] + THitI* staHitIndices = nullptr; + L1HitPoint* staHitPoints = nullptr; + if (Nelements > 0) { + staHits = &((*vStsHitsUnused)[start]); + staHitIndices = &((*RealIHitP)[start]); + staHitPoints = &((*vStsHitPointsUnused)[start]); + } + int NHitsOnStationTmp = NHitsOnStation; - vGridTime[ista].UpdateIterGrid( - Nelements, &((*vStsHitsUnused)[StsHitsUnusedStartIndex[ista]]), RealIHitPBuf, - &((*RealIHitP)[StsHitsUnusedStartIndex[ista]]), vStsHitsUnused_buf, vStsHitPointsUnused_buf, - &((*vStsHitPointsUnused)[StsHitsUnusedStartIndex[ista]]), NHitsOnStation, ista, *this, fStripFlag); + + vGridTime[ista].UpdateIterGrid(Nelements, staHits, RealIHitPBuf, staHitIndices, vStsHitsUnused_buf, + vStsHitPointsUnused_buf, staHitPoints, NHitsOnStation, ista, *this, fStripFlag); + StsHitsUnusedStartIndex[ista] = NHitsOnStationTmp; StsHitsUnusedStopIndex[ista] = NHitsOnStation; } @@ -2690,7 +2721,7 @@ inline void L1Algo::CAFindTrack(int ista, L1Branch& best_tr, unsigned char& best } #ifdef DRAW -void L1Algo::DrawRecoTracksTime(const vector<CbmL1Track>& tracks) +void L1Algo::DrawRecoTracksTime(const L1Vector<CbmL1Track>& tracks) { // TODO: find where the DrawRecoTracksTime is. It is missing in the git repository. //draw->DrawRecoTracksTime(tracks); diff --git a/reco/L1/CbmL1Def.h b/reco/L1/L1Algo/L1Def.h similarity index 89% rename from reco/L1/CbmL1Def.h rename to reco/L1/L1Algo/L1Def.h index 603f1e04f3..24c24288db 100644 --- a/reco/L1/CbmL1Def.h +++ b/reco/L1/L1Algo/L1Def.h @@ -2,26 +2,24 @@ SPDX-License-Identifier: GPL-3.0-only Authors: Maksym Zyzak, Igor Kulakov [committer], Sergey Gorbunov */ -#ifndef CbmL1Def_h -#define CbmL1Def_h +#ifndef L1Def_h +#define L1Def_h // #define FAST_CODE // FAST_CODE = more unsafe - #include "TStopwatch.h" #include <iostream> -#include <vector> #include <assert.h> -#ifdef HAVE_SSE + +//#if defined(HAVE_SSE) || SSE_FOUND #include "vectors/P4_F32vec4.h" -#else -#include "vectors/PSEUDO_F32vec4.h" -#error NoSseFound -#endif // HAVE_SSE +//#else +//#include "vectors/PSEUDO_F32vec4.h" +//#error NoSseFound +//#endif // HAVE_SSE -//#include "vectors/PSEUDO_F64vec1.h" template<typename T> T finite(T x) @@ -67,4 +65,4 @@ using namespace std; typedef int index_type; -#endif // CbmL1Def_h +#endif diff --git a/reco/L1/L1Algo/L1Extrapolation.h b/reco/L1/L1Algo/L1Extrapolation.h index 87030095a8..4958db054d 100644 --- a/reco/L1/L1Algo/L1Extrapolation.h +++ b/reco/L1/L1Algo/L1Extrapolation.h @@ -5,8 +5,7 @@ #ifndef L1Extrapolation_h #define L1Extrapolation_h -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1Field.h" #include "L1TrackPar.h" diff --git a/reco/L1/L1Algo/L1Field.h b/reco/L1/L1Algo/L1Field.h index 674426eb8a..5ee064ffea 100644 --- a/reco/L1/L1Algo/L1Field.h +++ b/reco/L1/L1Algo/L1Field.h @@ -5,9 +5,9 @@ #ifndef L1Field_h #define L1Field_h 1 -#include "CbmL1Def.h" - #include <iostream> + +#include "L1Def.h" using std::cout; using std::endl; using std::ostream; diff --git a/reco/L1/L1Algo/L1Filtration.h b/reco/L1/L1Algo/L1Filtration.h index e3ecad8dc3..08e77162c2 100644 --- a/reco/L1/L1Algo/L1Filtration.h +++ b/reco/L1/L1Algo/L1Filtration.h @@ -5,8 +5,7 @@ #ifndef L1Filtration_h #define L1Filtration_h -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1TrackPar.h" #include "L1UMeasurementInfo.h" #include "L1XYMeasurementInfo.h" diff --git a/reco/L1/L1Algo/L1Fit.h b/reco/L1/L1Algo/L1Fit.h index 020d8d27c2..a40bbc2096 100644 --- a/reco/L1/L1Algo/L1Fit.h +++ b/reco/L1/L1Algo/L1Fit.h @@ -10,8 +10,7 @@ #ifndef L1Fit_H #define L1Fit_H -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1MaterialInfo.h" #include "L1TrackPar.h" diff --git a/reco/L1/L1Algo/L1FitMaterial.h b/reco/L1/L1Algo/L1FitMaterial.h index 9679407ba4..24a7696131 100644 --- a/reco/L1/L1Algo/L1FitMaterial.h +++ b/reco/L1/L1Algo/L1FitMaterial.h @@ -5,8 +5,7 @@ #ifndef L1FitMaterial_h #define L1FitMaterial_h -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1MaterialInfo.h" #include "L1TrackPar.h" diff --git a/reco/L1/L1Algo/L1Grid.cxx b/reco/L1/L1Algo/L1Grid.cxx index 0227b0fa1e..3de8de7d23 100644 --- a/reco/L1/L1Algo/L1Grid.cxx +++ b/reco/L1/L1Algo/L1Grid.cxx @@ -4,15 +4,13 @@ #include "L1Grid.h" -#include "CbmL1Def.h" - #include <algorithm> #include <cstdio> #include <assert.h> #include <string.h> -#include "L1Hit.h" +#include "L1Def.h" #ifdef _OPENMP #include "omp.h" #endif @@ -39,12 +37,12 @@ inline void memset(T* dest, T i, size_t num) } -void L1Grid::UpdateIterGrid(unsigned int Nelements, L1Hit* hits, vector<THitI>* indicesBuf, THitI* indices, - vector<L1Hit>* hits2, vector<L1HitPoint>* pointsBuf, L1HitPoint* points, +void L1Grid::UpdateIterGrid(unsigned int Nelements, L1Hit* hits, L1Vector<THitI>* indicesBuf, THitI* indices, + L1Vector<L1Hit>* hits2, L1Vector<L1HitPoint>* pointsBuf, L1HitPoint* points, int& NHitsOnStation, char iS, L1Algo& Algo, const L1Vector<unsigned char>* vSFlag) { - fFirstHitInBin.assign(fN + 2, 0); + fFirstHitInBin.reset(fN + 2, 0); fscal xs = 0; fscal ys = 0; @@ -129,8 +127,8 @@ void L1Grid::AllocateMemory(int NThreads) int binsGrid = 600000; - fFirstHitInBin.resize(binsGrid, 0); - fHitsInBin.resize(binsGrid, 0); + fFirstHitInBin.reset(binsGrid, 0); + fHitsInBin.reset(binsGrid, 0); // for( int i = 0; i < fNThreads; i++ ) // { @@ -174,7 +172,7 @@ void L1Grid::StoreHits(THitI nhits, const L1Hit* hits, char iS, L1Algo& Algo, TH fscal xs = 0; fscal ys = 0; - fFirstHitInBin.assign(fN + 2, 0); + fFirstHitInBin.reset(fN + 2, 0); #ifdef _OPENMP #pragma omp parallel for firstprivate(xs, ys) diff --git a/reco/L1/L1Algo/L1Grid.h b/reco/L1/L1Algo/L1Grid.h index f4df6800a1..fd72235659 100644 --- a/reco/L1/L1Algo/L1Grid.h +++ b/reco/L1/L1Algo/L1Grid.h @@ -10,13 +10,13 @@ #ifndef L1GRID_H #define L1GRID_H -#include "CbmL1Def.h" - #include <algorithm> #include <cstdio> #include <assert.h> #include <string.h> + +#include "L1Def.h" #ifdef _OPENMP #include "omp.h" #endif @@ -40,50 +40,11 @@ class L1Algo; */ class L1Grid { public: - L1Grid() - : fN(0) - , fNy(0) - , fNz(0) - , fNt(0) - , fYMinOverStep(0.) - , fZMinOverStep(0.) - , fTMinOverStep(0.) - , fStepYInv(0.) - , fStepZInv(0.) - , fStepTInv(0.) - , fBinInGrid(0) - , fFirstHitInBin() - , fHitsInBin() - , fNThreads(0) - { - } + L1Grid() = default; - L1Grid(const L1Grid& grid) - : fN(grid.N()) - , fNy(grid.Ny()) - , fNz(grid.Nz()) - , fNt(grid.Nt()) - , fYMinOverStep(0.) - , fZMinOverStep(0.) - , fTMinOverStep(0.) - , fStepYInv(0.) - , fStepZInv(0.) - , fStepTInv(0.) - , fBinInGrid(0) - , fFirstHitInBin() - , fHitsInBin() - , fNThreads(0) - { - } - // ~L1Grid(){ //if ( fFirstHitInBin ) delete[] fFirstHitInBin; - // - // - // - // for( int i = 0; i < fNThreads; i++ ) { - // // if (fFirstHitInBinArray[i]) delete[] fFirstHitInBinArray[i]; - // if (fFirstHitInBin[i]) delete[] fFirstHitInBin[i]; - // } - // } + L1Grid(const L1Grid& grid) : fN(grid.N()), fNy(grid.Ny()), fNz(grid.Nz()), fNt(grid.Nt()) {} + + ~L1Grid() = default; void StoreHits(THitI nhits, const L1Hit* hits, char iS, L1Algo& Algo, THitI n, L1Hit* hitsBuf1, const L1Hit* hits1, THitI* indices1); @@ -142,32 +103,27 @@ public: // }; - void UpdateIterGrid(unsigned int Nelements, L1Hit* hits, vector<THitI>* indicesBuf, THitI* indices, - vector<L1Hit>* hits2, vector<L1HitPoint>* pointsBuf, L1HitPoint* points, int& NHitsOnStation, + void UpdateIterGrid(unsigned int Nelements, L1Hit* hits, L1Vector<THitI>* indicesBuf, THitI* indices, + L1Vector<L1Hit>* hits2, L1Vector<L1HitPoint>* pointsBuf, L1HitPoint* points, int& NHitsOnStation, char iS, L1Algo& Algo, const L1Vector<unsigned char>* vSFlag); private: - unsigned int fN; //* total N bins - unsigned short fNy; //* N bins in Y - unsigned short fNz; //* N bins in Z - unsigned short fNt; //* N bins in Z - float fYMinOverStep; //* minimal Y value * fStepYInv - float fZMinOverStep; //* minimal Z value * fStepZInv - float fTMinOverStep; //* minimal Z value * fStepZInv - float fStepYInv; //* inverse bin size in Y - float fStepZInv; //* inverse bin size in Z - float fStepTInv; //* inverse bin size in Z - int fBinInGrid; - - vector<THitI> fFirstHitInBin; - vector<THitI> fHitsInBin; - // vector <THitI*> fFirstHitInBinArray; - // vector <THitI> fOffsets; - // vector <THitI> fNumberHitsInBin; - - - unsigned short fNThreads; + unsigned int fN {0}; //* total N bins + unsigned short fNy {0}; //* N bins in Y + unsigned short fNz {0}; //* N bins in Z + unsigned short fNt {0}; //* N bins in Z + float fYMinOverStep {0.f}; //* minimal Y value * fStepYInv + float fZMinOverStep {0.f}; //* minimal Z value * fStepZInv + float fTMinOverStep {0.f}; //* minimal Z value * fStepZInv + float fStepYInv {0.f}; //* inverse bin size in Y + float fStepZInv {0.f}; //* inverse bin size in Z + float fStepTInv {0.f}; //* inverse bin size in Z + int fBinInGrid {0}; + unsigned short fNThreads {0}; + + L1Vector<THitI> fFirstHitInBin {"L1Grid::fFirstHitInBin"}; + L1Vector<THitI> fHitsInBin {"L1Grid::fHitsInBin"}; // vector <omp_lock_t> lock; }; diff --git a/reco/L1/L1Algo/L1HitArea.h b/reco/L1/L1Algo/L1HitArea.h index 4e53fb5732..6d3fda118c 100644 --- a/reco/L1/L1Algo/L1HitArea.h +++ b/reco/L1/L1Algo/L1HitArea.h @@ -5,8 +5,7 @@ #ifndef L1HitArea_H #define L1HitArea_H -#include "CbmL1Def.h" - +#include "L1Def.h" #include "L1Grid.h" class L1Row; diff --git a/reco/L1/L1Algo/L1MaterialInfo.h b/reco/L1/L1Algo/L1MaterialInfo.h index 465cf06abe..d0da3ebbca 100644 --- a/reco/L1/L1Algo/L1MaterialInfo.h +++ b/reco/L1/L1Algo/L1MaterialInfo.h @@ -5,10 +5,10 @@ #ifndef L1MaterialInfo_h #define L1MaterialInfo_h -#include "../CbmL1Def.h" - #include <vector> +#include "L1Def.h" + class L1MaterialInfo { public: L1MaterialInfo() : thick(0), RL(0), RadThick(0), logRadThick(0) {}; diff --git a/reco/L1/L1Algo/L1Timer.h b/reco/L1/L1Algo/L1Timer.h index fb392b9094..3688207ab4 100644 --- a/reco/L1/L1Algo/L1Timer.h +++ b/reco/L1/L1Algo/L1Timer.h @@ -5,13 +5,13 @@ #ifndef _L1Timer_H #define _L1Timer_H -#include "../CbmL1Def.h" - #include <iomanip> #include <iostream> #include <string> #include <vector> +#include "L1Def.h" + using std::cout; using std::endl; using std::ios; diff --git a/reco/L1/L1Algo/L1TrackExtender.cxx b/reco/L1/L1Algo/L1TrackExtender.cxx index 2e917ff1f0..afb5d12740 100644 --- a/reco/L1/L1Algo/L1TrackExtender.cxx +++ b/reco/L1/L1Algo/L1TrackExtender.cxx @@ -17,7 +17,6 @@ // using namespace std; using std::cout; using std::endl; -using std::vector; /// Fit track @@ -197,9 +196,8 @@ void L1Algo::BranchFitter(const L1Branch& t, L1TrackPar& T, const bool dir, cons void L1Algo::FindMoreHits(L1Branch& t, L1TrackPar& T, const bool dir, const fvec qp0) // TODO take into account pipe { - std::vector<THitI> newHits; - newHits.clear(); - newHits.reserve(5); + L1Vector<THitI> newHits {"L1TrackExtender::newHits"}; + newHits.reserve(NStations); L1Fit fit; fit.SetParticleMass(GetDefaultParticleMass()); @@ -357,7 +355,7 @@ void L1Algo::FindMoreHits(L1Branch& t, L1TrackPar& T, const bool dir, // save hits const unsigned int NOldHits = t.NHits; const unsigned int NNewHits = newHits.size(); - t.fStsHits.resize(NNewHits + NOldHits); + t.fStsHits.enlarge(NNewHits + NOldHits); t.NHits = (NNewHits + NOldHits); if (dir) { // backward diff --git a/reco/L1/L1Algo/L1TrackPar.h b/reco/L1/L1Algo/L1TrackPar.h index 03b713c31e..b23458b472 100644 --- a/reco/L1/L1Algo/L1TrackPar.h +++ b/reco/L1/L1Algo/L1TrackPar.h @@ -5,7 +5,7 @@ #ifndef L1TrackPar_h #define L1TrackPar_h 1 -#include "../CbmL1Def.h" +#include "L1Def.h" class L1TrackPar { diff --git a/reco/L1/L1Algo/L1TrackParFit.h b/reco/L1/L1Algo/L1TrackParFit.h index 9dceeb8911..9e327be3b1 100644 --- a/reco/L1/L1Algo/L1TrackParFit.h +++ b/reco/L1/L1Algo/L1TrackParFit.h @@ -5,8 +5,7 @@ #ifndef L1TrackParFit_h #define L1TrackParFit_h -#include "../CbmL1Def.h" - +#include "L1Def.h" #include "L1Field.h" #include "L1MaterialInfo.h" #include "L1TrackPar.h" diff --git a/reco/L1/L1Algo/L1Triplet.h b/reco/L1/L1Algo/L1Triplet.h index 342ef3129c..5c31f8e714 100644 --- a/reco/L1/L1Algo/L1Triplet.h +++ b/reco/L1/L1Algo/L1Triplet.h @@ -10,7 +10,7 @@ // @author Valentina Akishina // @date 2021-05-18 -#include "CbmL1Def.h" +#include "L1Def.h" /// L1Triplet class represents a short 3-hit track segment called a "triplet". /// diff --git a/reco/L1/L1Algo/L1UMeasurementInfo.h b/reco/L1/L1Algo/L1UMeasurementInfo.h index a57555bfdc..bcdf46ee0d 100644 --- a/reco/L1/L1Algo/L1UMeasurementInfo.h +++ b/reco/L1/L1Algo/L1UMeasurementInfo.h @@ -5,7 +5,7 @@ #ifndef L1UMeasurementInfo_h #define L1UMeasurementInfo_h 1 -#include "../CbmL1Def.h" +#include "L1Def.h" class L1UMeasurementInfo { diff --git a/reco/L1/L1Algo/L1Vector.h b/reco/L1/L1Algo/L1Vector.h index 99b6bd4884..80f2cf1609 100644 --- a/reco/L1/L1Algo/L1Vector.h +++ b/reco/L1/L1Algo/L1Vector.h @@ -10,7 +10,7 @@ /// @date 2021-06-16 -#include "CbmL1Def.h" +#include "L1Def.h" #ifndef FAST_CODE #include <FairLogger.h> #endif @@ -33,7 +33,16 @@ class L1Vector : private std::vector<T> { public: typedef std::vector<T> Tbase; - L1Vector(const char* name = "no name") : Tbase(), fName(name) {}; + template<typename... Tinput> + L1Vector(Tinput... value) : Tbase(value...) + { + } + + template<typename... Tinput> + L1Vector(const char* name, Tinput... value) : Tbase(value...) + , fName(name) + { + } L1Vector(const L1Vector& v) : Tbase(), fName(v.fName) { @@ -41,6 +50,7 @@ public: Tbase::assign(v.begin(), v.end()); } + void SetName(const std::string& s) { fName = s; } void SetName(const std::basic_ostream<char>& s) @@ -58,17 +68,42 @@ public: } template<typename... Tinput> - void resize(std::size_t count, Tinput... value) + void reset(std::size_t count, Tinput... value) { -#ifndef FAST_CODE - if (!Tbase::empty() && (count > Tbase::capacity())) { - LOG(WARNING) << "L1Vector \"" << fName << "\"::resize(" << count << "): allocated capacity of " - << Tbase::capacity() << " is reached, re-allocate and copy." << std::endl; + // 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; } -#endif Tbase::resize(count, value...); } + 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) { @@ -97,9 +132,9 @@ public: { #ifndef FAST_CODE if (pos >= Tbase::size()) { - LOG(ERROR) << "L1Vector \"" << fName << "\": trying to access element " << pos + LOG(FATAL) << "L1Vector \"" << fName << "\": trying to access element " << pos << " outside of the vector of the size of " << Tbase::size() << std::endl; - exit(0); + assert(pos < Tbase::size()); } #endif return Tbase::operator[](pos); @@ -109,9 +144,9 @@ public: { #ifndef FAST_CODE if (pos >= Tbase::size()) { - LOG(ERROR) << "L1Vector \"" << fName << "\": trying to access element " << pos + LOG(FATAL) << "L1Vector \"" << fName << "\": trying to access element " << pos << " outside of the vector of the size of " << Tbase::size() << std::endl; - exit(0); + assert(pos < Tbase::size()); } #endif return Tbase::operator[](pos); @@ -121,8 +156,8 @@ public: { #ifndef FAST_CODE if (Tbase::size() == 0) { - LOG(ERROR) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; - exit(0); + LOG(FATAL) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; + assert(Tbase::size() > 0); } #endif return Tbase::back(); @@ -132,36 +167,28 @@ public: { #ifndef FAST_CODE if (Tbase::size() == 0) { - LOG(ERROR) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; - exit(0); + LOG(FATAL) << "L1Vector \"" << fName << "\": trying to access element of an empty vector" << std::endl; + assert(Tbase::size() > 0); } #endif return Tbase::back(); } - void assign(std::size_t count, const T& value) - { -#ifndef FAST_CODE - if (!Tbase::empty() && (count > Tbase::capacity())) { - LOG(WARNING) << "L1Vector \"" << fName << "\"::assign(" << count << "): allocated capacity of " - << Tbase::capacity() << " is reached, re-allocate and copy." << std::endl; - } -#endif - Tbase::assign(count, value); - } - 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::reserve; using Tbase::size; using typename Tbase::iterator; private: - std::string fName; + std::string fName {"no name"}; + using Tbase::assign; // use reset() instead using Tbase::at; + using Tbase::resize; }; /// diff --git a/reco/L1/L1Algo/L1XYMeasurementInfo.h b/reco/L1/L1Algo/L1XYMeasurementInfo.h index 69d9413e8d..c452408fb5 100644 --- a/reco/L1/L1Algo/L1XYMeasurementInfo.h +++ b/reco/L1/L1Algo/L1XYMeasurementInfo.h @@ -5,7 +5,7 @@ #ifndef L1XYMeasurementInfo_h #define L1XYMeasurementInfo_h 1 -#include "../CbmL1Def.h" +#include "L1Def.h" class L1XYMeasurementInfo { diff --git a/reco/L1/L1Algo/utils/L1AlgoDraw.h b/reco/L1/L1Algo/utils/L1AlgoDraw.h index 129c9426ed..ea9d7e03f6 100644 --- a/reco/L1/L1Algo/utils/L1AlgoDraw.h +++ b/reco/L1/L1Algo/utils/L1AlgoDraw.h @@ -5,11 +5,10 @@ #ifndef L1AlgoDraw_h #define L1AlgoDraw_h 1 -#include "CbmL1Def.h" - #include <map> #include <vector> +#include "L1Def.h" #include "L1Hit.h" #include "L1Station.h" #include "L1Triplet.h" diff --git a/reco/L1/L1AlgoInputData.h b/reco/L1/L1AlgoInputData.h index 7152a424fa..2fbfd1d0f1 100644 --- a/reco/L1/L1AlgoInputData.h +++ b/reco/L1/L1AlgoInputData.h @@ -5,35 +5,24 @@ #ifndef _L1AlgoInputData_h #define _L1AlgoInputData_h -#include "CbmL1Def.h" - #include <fstream> #include <iostream> -#include <vector> +#include "L1Def.h" #include "L1Hit.h" #include "L1Vector.h" using std::istream; -using std::vector; class L1AlgoInputData { public: - L1AlgoInputData() : vStsHits(), NStsStrips(0), vStsZPos() - // MaxNStations(12) - - { - for (int i = 0; i < MaxNStations + 1; ++i) - StsHitsStartIndex[i] = 0; - for (int i = 0; i < MaxNStations + 1; ++i) - StsHitsStopIndex[i] = 0; - } - ~L1AlgoInputData() {}; + L1AlgoInputData() = default; + ~L1AlgoInputData() = default; - vector<L1Hit>& GetStsHits() { return vStsHits; } + L1Vector<L1Hit>& GetStsHits() { return vStsHits; } int GetNStsStrips() const { return NStsStrips; } - const vector<fscal>& GetStsZPos() const { return vStsZPos; } + L1Vector<fscal>& GetStsZPos() { return vStsZPos; } L1Vector<unsigned char>& GetSFlag() { return fStripFlag; } const THitI* GetStsHitsStartIndex() const { return StsHitsStartIndex; } const THitI* GetStsHitsStopIndex() const { return StsHitsStopIndex; } @@ -80,15 +69,15 @@ public: { MaxNStations = 25 }; - vector<L1Hit> vStsHits; // hits as a combination of front-, backstrips and z-position - int NStsStrips; // Number of strips in sts - vector<fscal> vStsZPos; // all possible z-positions of hits + L1Vector<L1Hit> vStsHits {"L1AlgoInputData::vStsHits"}; // hits as a combination of front-, backstrips and z-position + int NStsStrips {0}; // Number of strips in sts + L1Vector<fscal> vStsZPos {"L1AlgoInputData::vStsZPos"}; // all possible z-positions of hits L1Vector<unsigned char> fStripFlag { "L1AlgoInputData::fStripFlag"}; // information of hits station & using hits in tracks; - THitI StsHitsStartIndex[MaxNStations + 1], - StsHitsStopIndex[MaxNStations + 1]; // station-bounders in vStsHits array + THitI StsHitsStartIndex[MaxNStations + 1] {0}; // station-bounders in vStsHits array + THitI StsHitsStopIndex[MaxNStations + 1] {0}; // station-bounders in vStsHits array } _fvecalignment; diff --git a/reco/L1/OffLineInterface/CbmL1RichENNRingFinderParallel.h b/reco/L1/OffLineInterface/CbmL1RichENNRingFinderParallel.h index 6ebb5ab0c4..731d8cac98 100644 --- a/reco/L1/OffLineInterface/CbmL1RichENNRingFinderParallel.h +++ b/reco/L1/OffLineInterface/CbmL1RichENNRingFinderParallel.h @@ -22,8 +22,9 @@ #define _CBM_L1_RICH_ENN_RING_FINDER2_H_ -#include "CbmL1Def.h" #include "CbmRichRingFinder.h" + +#include "L1Def.h" //#include "../vectors/PSEUDO_F32vec4.h" // for check //#include "../vectors/PSEUDO_F32vec1.h" // for check diff --git a/reco/L1/OffLineInterface/CbmL1RichRingQa.cxx b/reco/L1/OffLineInterface/CbmL1RichRingQa.cxx index 7226939b3e..d697f5aefb 100644 --- a/reco/L1/OffLineInterface/CbmL1RichRingQa.cxx +++ b/reco/L1/OffLineInterface/CbmL1RichRingQa.cxx @@ -4,7 +4,6 @@ #include "CbmL1RichRingQa.h" -#include "CbmL1Def.h" #include "CbmMCTrack.h" #include "CbmRichHit.h" #include "CbmRichPoint.h" @@ -39,6 +38,8 @@ #include <ctype.h> #include <stdio.h> +#include "L1Def.h" + using namespace std; using std::cout; diff --git a/reco/L1/ParticleFinder/CbmL1PFFitter.h b/reco/L1/ParticleFinder/CbmL1PFFitter.h index 8c7c40157a..17904ac338 100644 --- a/reco/L1/ParticleFinder/CbmL1PFFitter.h +++ b/reco/L1/ParticleFinder/CbmL1PFFitter.h @@ -20,10 +20,10 @@ #ifndef _CbmL1PFFitter_h_ #define _CbmL1PFFitter_h_ -#include "CbmL1Def.h" - #include <vector> +#include "L1Def.h" + class CbmL1Track; class CbmStsTrack; class L1TrackPar; -- GitLab