diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx index 10863478a5c0bff9725438569bdeb70b7c2f0331..aa3de257a24c159cd280b0626094ed797805456e 100644 --- a/reco/L1/CbmL1.cxx +++ b/reco/L1/CbmL1.cxx @@ -1028,7 +1028,7 @@ void CbmL1::Reconstruct(CbmEvent* event) // save reconstructed tracks if (fLegacyEventMode) vRTracksCur.clear(); - int start_hit = 0; + int trackFirstHit = 0; float TsStart_new = TsStart + TsLength - TsOverlap; @@ -1050,68 +1050,31 @@ void CbmL1::Reconstruct(CbmEvent* event) t.chi2 = it->chi2; t.NDF = it->NDF; - //t.T[4] = it->Momentum; t.Hits.clear(); - t.fTrackTime = it->fTrackTime; - - L1Vector<int> HitsLocal("CbmL1::HitsLocal"); - HitsLocal.reserve(it->NHits); - - for (int i = 0; i < it->NHits; i++) { - int start_hit1 = start_hit; - if (fpAlgo->fRecoHits[start_hit1] > fvExternalHits.size() - 1) start_hit++; - else if (!fLegacyEventMode) { - t.Hits.push_back((fpAlgo->GetInputData()->GetHit(fpAlgo->fRecoHits[start_hit]).ID)); - } - else { - t.Hits.push_back(fpAlgo->fRecoHits[start_hit]); - } - HitsLocal.push_back(fpAlgo->fRecoHits[start_hit++]); - } - t.mass = fpAlgo->fDefaultMass; // pion mass t.is_electron = 0; - t.SetId(vRTracksCur.size()); - // CbmL1Track* prtra = &t; - int indd = 0; - bool isInOverlap = 0; + bool isTrackInOverlap = 0; - - for (unsigned int i = 0; i < HitsLocal.size(); i++) { - // if ((*ih) > int(fvExternalHits.size() - 1)) { - // indd = 1; - // break; - // } - - if (fvExternalHits[HitsLocal[i]].t >= (TsStart + TsLength - TsOverlap)) { - isInOverlap = 1; - if (TsStart_new > fvExternalHits[HitsLocal[i]].t) TsStart_new = fvExternalHits[HitsLocal[i]].t; + for (int i = 0; i < it->NHits; i++) { + int caHitId = fpAlgo->fRecoHits[trackFirstHit + i]; + int cbmHitID = fpAlgo->GetInputData()->GetHit(caHitId).ID; + double time = fpAlgo->GetInputData()->GetHit(caHitId).t; + if (!fLegacyEventMode) { t.Hits.push_back(cbmHitID); } + else { + t.Hits.push_back(caHitId); } - - int nMCPoints = fvExternalHits[HitsLocal[i]].mcPointIds.size(); - for (int iP = 0; iP < nMCPoints; iP++) { - int iMC = fvExternalHits[HitsLocal[i]].mcPointIds[iP]; - if (iMC > int(fvMCPoints.size() - 1)) { - // cout << " iMC " << iMC << " fvMCPoints.size() " << fvMCPoints.size() << endl; - indd = 1; - } + if (time >= (TsStart + TsLength - TsOverlap)) { + isTrackInOverlap = 1; + if (TsStart_new > time) { TsStart_new = time; } } } + trackFirstHit += it->NHits; - if (indd == 1) continue; + // Discard tracks from overlap region + if ((!fLegacyEventMode) && (isTrackInOverlap == 1)) { continue; } - if ((!fLegacyEventMode) && (isInOverlap == 1)) { - - continue; ///Discard tracks from overlap region - - /// set strips as unused - //for (unsigned int i = 0; i < HitsLocal.size(); i++) { - // fpAlgo->SetFUnUsed(const_cast<unsigned char&>((*fpAlgo->fStripFlag)[fvExternalHits[HitsLocal[i]].f])); - // fpAlgo->SetFUnUsed(const_cast<unsigned char&>((*fpAlgo->fStripFlag)[fvExternalHits[HitsLocal[i]].b])); - //} - } vRTracksCur.push_back(t); } @@ -1293,8 +1256,8 @@ void CbmL1::IdealTrackFinder() L1Vector<int> hitIndices("CbmL1::hitIndices", fpAlgo->GetParameters()->GetNstationsActive(), -1); for (unsigned int iH = 0; iH < MC.Hits.size(); iH++) { - const int hitI = MC.Hits[iH]; - const CbmL1Hit& hit = fvExternalHits[hitI]; + const int hitI = MC.Hits[iH]; + const CbmL1HitDebugInfo& hit = fvHitDebugInfo[hitI]; const int iStation = fvMCPoints[hit.mcPointIds[0]].iStation; @@ -1319,8 +1282,6 @@ void CbmL1::IdealTrackFinder() fpAlgo->fRecoHits.push_back(hitI); } - - algoTr.Momentum = MC.p; algoTr.TFirst[0] = MC.x; algoTr.TFirst[1] = MC.y; algoTr.TFirst[2] = MC.px / MC.pz; @@ -1470,7 +1431,7 @@ void CbmL1::WriteSIMDKFData() int st[20]; int jHit = 0; for (int iHit = 0; iHit < NHits; iHit++) { - CbmL1HitStore& h = fvHitStore[RecTrack->Hits[iHit]]; + CbmL1HitDebugInfo& h = fvHitDebugInfo[RecTrack->Hits[iHit]]; st[jHit] = h.iStation; if (h.ExtIndex < 0) { CbmMvdHit* MvdH = (CbmMvdHit*) fpMvdHits->At(-h.ExtIndex - 1); diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h index 20995c11b8714546506c5ad27b26cfd5f6687ffe..fbdb116684f61587e05fde3ee400160e4c98c525 100644 --- a/reco/L1/CbmL1.h +++ b/reco/L1/CbmL1.h @@ -71,20 +71,6 @@ class CbmEvent; class TProfile2D; class TNtuple; -/// TODO: SZh 21.09.2022: Replace instances of this class with L1Hit -class CbmL1HitStore { -public: - int ExtIndex; - int iStation; - double x; - double y; - double time; - double dx; - double dy; - double dt; - double dxy; - int Det; -}; /// Internal structure to handle link keys struct CbmL1LinkKey { @@ -626,10 +612,6 @@ private: vector<vector<int>> fTofPointToTrack; ///< - // ** Repacked input data ** - - L1Vector<CbmL1Hit> fvExternalHits = {"CbmL1::fvExternalHits"}; ///< Array of hits with MC information - /// Indexes of hits after hits sorting (used only with fLegacyEventMode = true) L1Vector<int> fvSortedHitsIndexes = {"CbmL1::fvSortedHitsIndexes"}; @@ -639,7 +621,13 @@ private: L1Vector<int> fvHitPointIndexes = {"CbmL1::fvHitPointIndexes"}; ///< Indexes of MC points vs. hit index public: - L1Vector<CbmL1HitStore> fvHitStore = {"CbmL1::fvHitStore"}; ///< Container of hits with extended information + // ** Repacked input data ** + + L1Vector<CbmL1HitId> fvExternalHits = {"CbmL1::fvExternalHits"}; ///< Array of hits + +private: + L1Vector<CbmL1HitDebugInfo> fvHitDebugInfo = { + "CbmL1::fvHitDebugInfo"}; ///< Container of hits with extended information // indices of MCPoints in fvMCPoints, indexed by index of hit in algo->vHits array. According to StsMatch. Used for IdealResponce // L1Vector<int> vHitMCRef1; // CbmMatch HitMatch; diff --git a/reco/L1/CbmL1Hit.h b/reco/L1/CbmL1Hit.h index 64b274f70690b3bd47d10ce3318b0c24403ca420..94f3facaa9bef5bbe01cb464d2d15fbfd37599c0 100644 --- a/reco/L1/CbmL1Hit.h +++ b/reco/L1/CbmL1Hit.h @@ -8,24 +8,34 @@ #include "L1Vector.h" /// -/// hits with hit-mcpoint match information +/// Identificator for cbm hits with their detector and index in cbm arrays /// -class CbmL1Hit { +class CbmL1HitId { public: - CbmL1Hit() = default; + CbmL1HitId() = default; + CbmL1HitId(int det, int index) : detId(det), hitId(index) {}; + + int detId {0}; ///< detector ID (mvd/sts/etc) + int hitId {0}; ///< index of hit in the TClonesArray array +}; - CbmL1Hit(int hitId_, int extIndex_, int Det_) : hitId(hitId_), extIndex(extIndex_), Det(Det_) {}; - int hitId = 0; ///< index of L1Hit in fInputData::fvHits array. Should be equal to index of this in L1->vHits - int extIndex = 0; ///< index of hit in the TClonesArray array - int Det = 0; ///< detector ID (mvd/sts/etc) - float x = 0.f; ///< measured X coordinate - float y = 0.f; ///< measured Y coordinate - float t = 0.f; ///< measured time - int f = 0; ///< front strip index - int b = 0; ///< back strip index - int ID = 0; ///< TODO: check if this ID is redundant - L1Vector<int> mcPointIds {"CbmL1Hit::mcPointIds"}; ///< indices of CbmL1MCPoint in L1->fvMCPoints array +/// +/// a helper class for performance evaluation that contains useful info about cbm hits with hit-mcpoint match information +/// +class CbmL1HitDebugInfo { // TODO: SZh 21.09.2022: Replace instances of this class with L1Hit +public: + int ExtIndex; + int iStation; + double x; + double y; + double time; + double dx; + double dy; + double dt; + double dxy; + int Det; + L1Vector<int> mcPointIds {"CbmL1HitDebugInfo::mcPointIds"}; ///< indices of CbmL1MCPoint in L1->fvMCPoints array }; #endif diff --git a/reco/L1/CbmL1MCTrack.cxx b/reco/L1/CbmL1MCTrack.cxx index 00c6385bae4be13072ac381551fa0575eb8bc5e6..2b75033cb1cd84077c59559bf00003fbdc72b2b7 100644 --- a/reco/L1/CbmL1MCTrack.cxx +++ b/reco/L1/CbmL1MCTrack.cxx @@ -157,7 +157,7 @@ void CbmL1MCTrack::CalculateMaxNStaHits() int lastSta = -1; int cur_maxNStaHits = 0; for (unsigned int iH = 0; iH < Hits.size(); iH++) { - CbmL1HitStore& sh = L1->fvHitStore[Hits[iH]]; + CbmL1HitDebugInfo& sh = L1->fvHitDebugInfo[Hits[iH]]; if (sh.iStation == lastSta) { cur_maxNStaHits++; } else { // new station if (!(sh.iStation > lastSta)) { // tracks going in backward direction are not reconstructable diff --git a/reco/L1/CbmL1Performance.cxx b/reco/L1/CbmL1Performance.cxx index 5a2479c3227d94e65a9e6478dd5d6021ceb5ccfc..c2cf4d0a2ef7e0f07b3757f201b61de0cb891459 100644 --- a/reco/L1/CbmL1Performance.cxx +++ b/reco/L1/CbmL1Performance.cxx @@ -93,9 +93,9 @@ void CbmL1::TrackMatch() map<int, int>& hitmap = prtra->hitMap; // how many hits from each mcTrack belong to current recoTrack for (vector<int>::iterator ih = (prtra->Hits).begin(); ih != (prtra->Hits).end(); ++ih) { - const int nMCPoints = fvExternalHits[*ih].mcPointIds.size(); + const int nMCPoints = fvHitDebugInfo[*ih].mcPointIds.size(); for (int iP = 0; iP < nMCPoints; iP++) { - int iMC = fvExternalHits[*ih].mcPointIds[iP]; + int iMC = fvHitDebugInfo[*ih].mcPointIds[iP]; // cout<<iMC<<" iMC"<<endl; int ID = -1; @@ -293,7 +293,7 @@ void CbmL1::EfficienciesPerformance() for (vector<CbmL1Track>::iterator rtraIt = fvRecoTracks.begin(); rtraIt != fvRecoTracks.end(); ++rtraIt) { ntra.ghosts += rtraIt->IsGhost(); - if (rtraIt->IsGhost()) { // Debug. + if (0 && rtraIt->IsGhost()) { // Debug. cout << " L1: ghost track: nhits " << rtraIt->GetNOfHits() << " p " << 1. / rtraIt->T[5] << " purity " << rtraIt->GetMaxPurity() << " | hits "; for (map<int, int>::iterator posIt = rtraIt->hitMap.begin(); posIt != rtraIt->hitMap.end(); posIt++) { @@ -347,7 +347,7 @@ void CbmL1::EfficienciesPerformance() if (reco) nclones = mtra.GetNClones(); // if (nclones){ // Debug. Look at clones // for (int irt = 0; irt < rTracks.size(); irt++){ - // const int ista = fvHitStore[rTracks[irt]->Hits[0]].iStation; + // const int ista = fvHitDebugInfo[rTracks[irt]->Hits[0]].iStation; // cout << rTracks[irt]->GetNOfHits() << "(" << ista << ") "; // } // cout << mtra.NStations() << endl; @@ -768,7 +768,7 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa // //hit density calculation: h_hit_density[10] // - // for (vector<CbmL1HitStore>::iterator hIt = fvHitStore.begin(); hIt != fvHitStore.end(); ++hIt){ + // for (vector<CbmL1HitDebugInfo>::iterator hIt = fvHitDebugInfo.begin(); hIt != fvHitDebugInfo.end(); ++hIt){ // float x = hIt->x; // float y = hIt->y; // float r = sqrt(x*x+y*y); @@ -788,7 +788,7 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa h_reco_phi->Fill(TMath::ATan2(-prtra->T[3], -prtra->T[2])); // TODO: What is precision? h_reco_nhits->Fill((prtra->Hits).size()); - CbmL1HitStore& mh = fvHitStore[prtra->Hits[0]]; + CbmL1HitDebugInfo& mh = fvHitDebugInfo[prtra->Hits[0]]; h_reco_station->Fill(mh.iStation); } @@ -821,8 +821,8 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa h_ghost_Rmom->Fill(fabs(1.0 / prtra->T[4])); } h_ghost_nhits->Fill((prtra->Hits).size()); - CbmL1HitStore& h1 = fvHitStore[prtra->Hits[0]]; - CbmL1HitStore& h2 = fvHitStore[prtra->Hits[1]]; + CbmL1HitDebugInfo& h1 = fvHitDebugInfo[prtra->Hits[0]]; + CbmL1HitDebugInfo& h2 = fvHitDebugInfo[prtra->Hits[1]]; h_ghost_fstation->Fill(h1.iStation); h_ghost_r->Fill(sqrt(fabs(h1.x * h1.x + h1.y * h1.y))); double z1 = fpAlgo->GetParameters()->GetStation(h1.iStation).fZ[0]; @@ -833,8 +833,8 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa } if (fabs(prtra->T[4]) > 1.e-10) h2_ghost_nhits_vs_mom->Fill(fabs(1.0 / prtra->T[4]), (prtra->Hits).size()); - CbmL1HitStore& hf = fvHitStore[prtra->Hits[0]]; - CbmL1HitStore& hl = fvHitStore[prtra->Hits[(prtra->Hits).size() - 1]]; + CbmL1HitDebugInfo& hf = fvHitDebugInfo[prtra->Hits[0]]; + CbmL1HitDebugInfo& hl = fvHitDebugInfo[prtra->Hits[(prtra->Hits).size() - 1]]; if (fabs(prtra->T[4]) > 1.e-10) h2_ghost_fstation_vs_mom->Fill(fabs(1.0 / prtra->T[4]), hf.iStation + 1); if (hl.iStation >= hf.iStation) h2_ghost_lstation_vs_fstation->Fill(hf.iStation + 1, hl.iStation + 1); } @@ -860,8 +860,8 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa int nSta = mtra.NStations(); - CbmL1HitStore& fh = fvHitStore[*(mtra.Hits.begin())]; - CbmL1HitStore& lh = fvHitStore[*(mtra.Hits.rbegin())]; + CbmL1HitDebugInfo& fh = fvHitDebugInfo[*(mtra.Hits.begin())]; + CbmL1HitDebugInfo& lh = fvHitDebugInfo[*(mtra.Hits.rbegin())]; h_reg_MCmom->Fill(momentum); if (mtra.IsPrimary()) { @@ -926,8 +926,8 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa } - int iph = mtra.Hits[0]; - CbmL1HitStore& ph = fvHitStore[iph]; + int iph = mtra.Hits[0]; + CbmL1HitDebugInfo& ph = fvHitDebugInfo[iph]; h_sec_r->Fill(sqrt(fabs(ph.x * ph.x + ph.y * ph.y))); if (fabs(mtra.pz) > 1.e-8) { @@ -1024,8 +1024,8 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa h_notfound_ty->Fill(pMC.py / pMC.pz); } - // CbmL1HitStore &ph21 = fvHitStore[mtra.Hits[0]]; - // CbmL1HitStore &ph22 = fvHitStore[mtra.Hits[1]]; + // CbmL1HitDebugInfo &ph21 = fvHitDebugInfo[mtra.Hits[0]]; + // CbmL1HitDebugInfo &ph22 = fvHitDebugInfo[mtra.Hits[1]]; // double z21 = fpAlgo->GetParameters()->GetStation(ph21.iStation).fZ[0]; // double z22 = fpAlgo->GetParameters()->GetStation(ph22.iStation).fZ[0]; @@ -1250,7 +1250,7 @@ void CbmL1::TrackFitPerformance() { int nTimeMeasurenments = 0; for (unsigned int ih = 0; ih < it->Hits.size(); ih++) { - int ista = fvHitStore[it->Hits[ih]].iStation; + int ista = fvHitDebugInfo[it->Hits[ih]].iStation; if (fpAlgo->GetParameters()->GetStation(ista).timeInfo) { nTimeMeasurenments++; } } isTimeFitted = (nTimeMeasurenments > 1); @@ -1335,7 +1335,7 @@ void CbmL1::TrackFitPerformance() fld.SetUseOriginalField(); fit.Extrapolate(mc.z, fld); // add material - const int fSta = fvHitStore[it->Hits[0]].iStation; + const int fSta = fvHitDebugInfo[it->Hits[0]].iStation; const int dir = int((mc.z - fpAlgo->GetParameters()->GetStation(fSta).fZ[0]) / fabs(mc.z - fpAlgo->GetParameters()->GetStation(fSta).fZ[0])); // if (abs(mc.z - fpAlgo->GetParameters()->GetStation(fSta).fZ[0]) > 10.) continue; // can't extrapolate on large distance @@ -1390,7 +1390,7 @@ void CbmL1::TrackFitPerformance() fld.SetUseOriginalField(); // add material - const int fSta = fvHitStore[it->Hits[0]].iStation; + const int fSta = fvHitDebugInfo[it->Hits[0]].iStation; const int dir = (mc.z - fpAlgo->GetParameters()->GetStation(fSta).fZ[0]) / abs(mc.z - fpAlgo->GetParameters()->GetStation(fSta).fZ[0]); @@ -1898,11 +1898,9 @@ void CbmL1::InputPerformance() map<unsigned int, unsigned int>::iterator it; if (fpStsHits && fpStsHitMatches) { - for (unsigned int iH = 0; iH < fvExternalHits.size(); iH++) { - const CbmL1Hit& h = fvExternalHits[iH]; + for (int iH = 0; iH < fpStsHits->GetEntriesFast(); iH++) { - if (h.Det != 1) continue; // not sts hit - const CbmStsHit* sh = L1_DYNAMIC_CAST<CbmStsHit*>(fpStsHits->At(h.extIndex)); + const CbmStsHit* sh = L1_DYNAMIC_CAST<CbmStsHit*>(fpStsHits->At(iH)); // int iMCPoint = -1; CbmLink link; @@ -2012,14 +2010,10 @@ void CbmL1::InputPerformance() if (fpMuchPixelHits && fpMuchHitMatches) { - for (unsigned int iH = 0; iH < fvExternalHits.size(); iH++) { - const CbmL1Hit& h = fvExternalHits[iH]; - - if (h.Det != 2) continue; // mvd hit - - const CbmMuchPixelHit* sh = L1_DYNAMIC_CAST<CbmMuchPixelHit*>(fpMuchPixelHits->At(h.extIndex)); - CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpMuchHitMatches->At(h.extIndex)); + for (int iH = 0; iH < fpMuchPixelHits->GetEntriesFast(); iH++) { + const CbmMuchPixelHit* sh = L1_DYNAMIC_CAST<CbmMuchPixelHit*>(fpMuchPixelHits->At(iH)); + const CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpMuchHitMatches->At(iH)); if (hm->GetNofLinks() == 0) continue; Float_t bestWeight = 0.f; @@ -2061,25 +2055,23 @@ void CbmL1::InputPerformance() mcPos.SetZ(hitPos.Z()); { // qa errors - if (sh->GetDx() > 1.e-8) pullXmuch->Fill((h.x - mcPos.X()) / sh->GetDx()); - if (sh->GetDy() > 1.e-8) pullYmuch->Fill((h.y - mcPos.Y()) / sh->GetDy()); - pullTmuch->Fill((h.t - mcTime) / sh->GetTimeError()); + if (sh->GetDx() > 1.e-8) pullXmuch->Fill((sh->GetX() - mcPos.X()) / sh->GetDx()); + if (sh->GetDy() > 1.e-8) pullYmuch->Fill((sh->GetY() - mcPos.Y()) / sh->GetDy()); + pullTmuch->Fill((sh->GetTime() - mcTime) / sh->GetTimeError()); } - resXmuch->Fill((h.x - mcPos.X()) * 10 * 1000); - resYmuch->Fill((h.y - mcPos.Y()) * 10 * 1000); - resTmuch->Fill((h.t - mcTime)); + resXmuch->Fill((sh->GetX() - mcPos.X()) * 10 * 1000); + resYmuch->Fill((sh->GetY() - mcPos.Y()) * 10 * 1000); + resTmuch->Fill((sh->GetTime() - mcTime)); } } // much if (fpTrdHits && fpTrdHitMatches) { - for (unsigned int iH = 0; iH < fvExternalHits.size(); iH++) { - const CbmL1Hit& h = fvExternalHits[iH]; + for (int iH = 0; iH < fpTrdHits->GetEntriesFast(); iH++) { - if (h.Det != 3) continue; // mvd hit - const CbmTrdHit* sh = L1_DYNAMIC_CAST<CbmTrdHit*>(fpTrdHits->At(h.extIndex)); - CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpTrdHitMatches->At(h.extIndex)); + const CbmTrdHit* sh = L1_DYNAMIC_CAST<CbmTrdHit*>(fpTrdHits->At(iH)); + const CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpTrdHitMatches->At(iH)); if (hm->GetNofLinks() == 0) continue; @@ -2124,28 +2116,24 @@ void CbmL1::InputPerformance() mcPos.SetZ(hitPos.Z()); { // qa errors - if (sh->GetDx() > 1.e-8) pullXtrd->Fill((h.x - mcPos.X()) / sh->GetDx()); - if (sh->GetDy() > 1.e-8) pullYtrd->Fill((h.y - mcPos.Y()) / sh->GetDy()); - pullTtrd->Fill((h.t - mcTime) / sh->GetTimeError()); + if (sh->GetDx() > 1.e-8) pullXtrd->Fill((sh->GetX() - mcPos.X()) / sh->GetDx()); + if (sh->GetDy() > 1.e-8) pullYtrd->Fill((sh->GetY() - mcPos.Y()) / sh->GetDy()); + pullTtrd->Fill((sh->GetTime() - mcTime) / sh->GetTimeError()); } - resXtrd->Fill((h.x - mcPos.X()) * 10 * 1000); - resYtrd->Fill((h.y - mcPos.Y()) * 10 * 1000); - resTtrd->Fill((h.t - mcTime)); + resXtrd->Fill((sh->GetX() - mcPos.X()) * 10 * 1000); + resYtrd->Fill((sh->GetY() - mcPos.Y()) * 10 * 1000); + resTtrd->Fill((sh->GetTime() - mcTime)); trdMCPointsZ->Fill(mcPos.Z()); } } // much if (fpTofHits && fpTofHitMatches) { - for (unsigned int iH = 0; iH < fvExternalHits.size(); iH++) { - const CbmL1Hit& h = fvExternalHits[iH]; - - if (h.Det != 4) continue; // mvd hit - - CbmTofHit* sh = L1_DYNAMIC_CAST<CbmTofHit*>(fpTofHits->At(h.extIndex)); - CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpTofHitMatches->At(h.extIndex)); + for (int iH = 0; iH < fpTofHits->GetEntriesFast(); iH++) { + const CbmTofHit* sh = L1_DYNAMIC_CAST<CbmTofHit*>(fpTofHits->At(iH)); + const CbmMatch* hm = L1_DYNAMIC_CAST<CbmMatch*>(fpTofHitMatches->At(iH)); if (hm->GetNofLinks() == 0) continue; Float_t bestWeight = 0.f; @@ -2188,13 +2176,13 @@ void CbmL1::InputPerformance() mcPos.SetZ(hitPos.Z()); { // qa errors - if (sh->GetDx() > 1.e-8) pullXtof->Fill((h.x - mcPos.X()) / sh->GetDx()); - if (sh->GetDy() > 1.e-8) pullYtof->Fill((h.y - mcPos.Y()) / sh->GetDy()); + if (sh->GetDx() > 1.e-8) pullXtof->Fill((sh->GetX() - mcPos.X()) / sh->GetDx()); + if (sh->GetDy() > 1.e-8) pullYtof->Fill((sh->GetY() - mcPos.Y()) / sh->GetDy()); pullTtof->Fill((sh->GetTime() - mcTime) / sh->GetTimeError()); } - resXtof->Fill((h.x - mcPos.X()) * 10 * 1000); - resYtof->Fill((h.y - mcPos.Y()) * 10 * 1000); + resXtof->Fill((sh->GetX() - mcPos.X()) * 10 * 1000); + resYtof->Fill((sh->GetY() - mcPos.Y()) * 10 * 1000); resTtof->Fill((sh->GetTime() - mcTime)); } } // much diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx index 22a86e698267d1e03ae30a22cd016a2628bde890..68e4a66eafab71fa7e5246adf19fb202e83d00f5 100644 --- a/reco/L1/CbmL1ReadEvent.cxx +++ b/reco/L1/CbmL1ReadEvent.cxx @@ -67,28 +67,28 @@ static bool compareMcPointZ(const int& a, const int& b) /// The structure is used to sort hits before writing them into L1 input arrays /// struct TmpHit { - int iStripF; ///< index of front strip - int iStripB; ///< index of back strip - int iStation; ///< index of station - int ExtIndex; ///< index of hit in the external TClonesArray array (NOTE: negative for MVD) - double u; ///< position of hit along axis orthogonal to front strips [cm] - double v; ///< position of hit along axis orthogonal to back strips [cm] - double x; ///< position of hit along x axis [cm] - double y; ///< position of hit along y axis [cm] - double z; ///< position of hit along z axis [cm] - double dx; ///< hit position uncertainty along x axis [cm] - double dy; ///< hit position uncertainty along y axis [cm] - double dxy; ///< hit position covariance along x and y axes [cm2] - double du; ///< hit position uncertainty along axis orthogonal to front strips [cm] - double dv; ///< hit position uncertainty along axis orthogonal to back strips [cm] - double xmc; ///< MC position of hit [cm] - double ymc; ///< MC position of hit [cm] - double p; ///< MC total momentum [GeV/c] - double tx; ///< MC slopes of the mc point - double ty; ///< MC slopes of the mc point - int iMC; ///< index of MCPoint in the fvMCPoints array - double time = 0.; ///< time of the hit [ns] - double dt = 1.e4; ///< time error of the hit [ns] + int iStripF; ///< index of front strip + int iStripB; ///< index of back strip + int iStation; ///< index of station + int ExtIndex; ///< index of hit in the external TClonesArray array (NOTE: negative for MVD) + double u; ///< position of hit along axis orthogonal to front strips [cm] + double v; ///< position of hit along axis orthogonal to back strips [cm] + double x; ///< position of hit along x axis [cm] + double y; ///< position of hit along y axis [cm] + double z; ///< position of hit along z axis [cm] + double dx; ///< hit position uncertainty along x axis [cm] + double dy; ///< hit position uncertainty along y axis [cm] + double dxy; ///< hit position covariance along x and y axes [cm2] + double du; ///< hit position uncertainty along axis orthogonal to front strips [cm] + double dv; ///< hit position uncertainty along axis orthogonal to back strips [cm] + double xmc; ///< MC position of hit [cm] + double ymc; ///< MC position of hit [cm] + double p; ///< MC total momentum [GeV/c] + double tx; ///< MC slopes of the mc point + double ty; ///< MC slopes of the mc point + int iMC; ///< index of MCPoint in the fvMCPoints array + double time = 0.; ///< time of the hit [ns] + double dt = 1.e4; ///< time error of the hit [ns] int Det; int id; ///< index of hit before hits sorting int track; @@ -343,7 +343,7 @@ void CbmL1::ReadEvent(float& TsStart, float& TsLength, float& /*TsOverlap*/, int fvExternalHits.clear(); /* <CbmL1Hit> */ fvRecoTracks.clear(); /* <CbmL1Track> */ // TODO: Move from this function to more suitable place (S.Zharko) fvHitPointIndexes.clear(); /* <int>: indexes of MC-points in fvMCPoints (by index of fpAlgo->vHits) */ - fvHitStore.clear(); /* <CbmL1HitStore> */ + fvHitDebugInfo.clear(); /* <CbmL1HitDebugInfo> */ fmMCPointsLinksMap.clear(); fmMCTracksLinksMap.clear(); @@ -1248,12 +1248,11 @@ void CbmL1::ReadEvent(float& TsStart, float& TsLength, float& /*TsOverlap*/, int if (fVerbose >= 10) { cout << "ReadEvent: strips are read." << endl; } - // ----- Fill and save fvExternalHits, fvHitStore and fvHitPointIndexes vectors as well as fpData->vHits + // ----- Fill and save fvExternalHits, fvHitDebugInfo and fvHitPointIndexes vectors as well as fpData->vHits fvSortedHitsIndexes.reset(maxHitIndex + 1); fvExternalHits.reserve(nHits); - - fvHitStore.reserve(nHits); + fvHitDebugInfo.reserve(nHits); fvHitPointIndexes.reserve(nHits); fIODataManager.ResetInputData(); @@ -1265,7 +1264,7 @@ void CbmL1::ReadEvent(float& TsStart, float& TsLength, float& /*TsOverlap*/, int for (int iHit = 0; iHit < nHits; ++iHit) { TmpHit& th = tmpHits[iHit]; - CbmL1HitStore s; + CbmL1HitDebugInfo s; s.Det = th.Det; s.ExtIndex = th.ExtIndex; s.iStation = th.iStation; @@ -1296,22 +1295,13 @@ void CbmL1::ReadEvent(float& TsStart, float& TsLength, float& /*TsOverlap*/, int // save hit - fvExternalHits.push_back(CbmL1Hit(iHit, th.ExtIndex, th.Det)); - - fvExternalHits[iHit].x = th.x; - fvExternalHits[iHit].y = th.y; - fvExternalHits[iHit].t = th.time; - - fvExternalHits[iHit].ID = th.id; - - fvExternalHits[iHit].f = th.iStripF; - fvExternalHits[iHit].b = th.iStripB; + fvExternalHits.push_back(CbmL1HitId(th.Det, th.ExtIndex)); // TODO: Here one should fill in the fvExternalHits[iHit].mcPointIds fIODataManager.PushBackHit(h); - fvHitStore.push_back(s); + fvHitDebugInfo.push_back(s); fvHitPointIndexes.push_back(th.iMC); } @@ -1571,10 +1561,8 @@ void CbmL1::HitMatch() { const int NHits = fvExternalHits.size(); for (int iH = 0; iH < NHits; iH++) { - CbmL1Hit& hit = fvExternalHits[iH]; - - - int iP = fvHitPointIndexes[iH]; + CbmL1HitDebugInfo& hit = fvHitDebugInfo[iH]; + int iP = fvHitPointIndexes[iH]; if (iP >= 0) { hit.mcPointIds.push_back_no_warning(iP); fvMCPoints[iP].hitIds.push_back_no_warning(iH); diff --git a/reco/L1/L1Algo/L1CATrackFinder.cxx b/reco/L1/L1Algo/L1CATrackFinder.cxx index 0a9f7e8efb27f2d8e5a3cc7a54d9556dbd341db3..2ca00b2142bf79cddcb9d05dd52c621da11ee7d0 100644 --- a/reco/L1/L1Algo/L1CATrackFinder.cxx +++ b/reco/L1/L1Algo/L1CATrackFinder.cxx @@ -2212,8 +2212,6 @@ void L1Algo::CATrackFinder() } t.NHits = tr.NHits; - // t.Momentum = tr.Momentum; - t.fTrackTime = sumTime / t.NHits; fTracks_local[num_thread].push_back(t); if (0) { // SG debug cout << "store track " << iCandidate << " chi2= " << tr.chi2 << endl; diff --git a/reco/L1/L1Algo/L1Track.h b/reco/L1/L1Algo/L1Track.h index e6223f1132fec881708f653cf076fb9c1de9491f..4732851a98358c9ba766690618a750969416b5fc 100644 --- a/reco/L1/L1Algo/L1Track.h +++ b/reco/L1/L1Algo/L1Track.h @@ -31,10 +31,7 @@ class L1Track { public: static constexpr float kNaN {std::numeric_limits<float>::signaling_NaN()}; - unsigned char NHits; ///< Number of hits in track - unsigned char n; ///< TODO: ?? - float Momentum {kNaN}; ///< TODO: ?? - float fTrackTime {kNaN}; ///< Track time + unsigned char NHits {0}; ///< Number of hits in track fscal TFirst[8] {kNaN}; ///< Track parameters on the first station fscal CFirst[28] {kNaN}; ///< Track parameter covariation matrix elements on the first station fscal TLast[8] {kNaN}; ///< Track parameters on the last station @@ -42,28 +39,7 @@ public: fscal Tpv[8] {kNaN}; ///< Track parameters in the primary vertex fscal Cpv[28] {kNaN}; ///< Track parameter covariation matrix elements in the primary vertex fscal chi2 {kNaN}; ///< Track fit chi-square value - short int NDF; ///< Track fit NDF value - - // TODO: shouldn't it be the L1HitIndex_t type instead of int? (S.Zharko) - int FirstHitIndex; ///< Track first hit index in the hits vector - int LastHitIndex; ///< Track last hit index in the hits vector - int index; ///< Index of the track - int ista; ///< TODO: ?? - - /// Provides comparison of two L1Track objects - /// If two tracks have different number of hits, the smallest track has the largest number of hits. - /// If two tracks have the same numbers of hits and ... - static bool compareCand(const L1Track& a, const L1Track& b) - { - if (a.NHits != b.NHits) { return (a.NHits > b.NHits); } - if (a.ista != b.ista) { return (a.ista < b.ista); } - else { - return (a.chi2 < b.chi2); - } - } - - /// Provides comparison for two tracks by the time variance - static bool compare(const L1Track& a, const L1Track& b) { return (a.Cpv[20] <= b.Cpv[20]); } + short int NDF {0}; ///< Track fit NDF value }; diff --git a/reco/L1/L1Algo/utils/L1AlgoEfficiencyPerformance.h b/reco/L1/L1Algo/utils/L1AlgoEfficiencyPerformance.h index 5791f205674a341b4ec32f1fbeb38aa844e040c5..512d6ce217825b0178ccc4f9148844f67dc288d7 100644 --- a/reco/L1/L1Algo/utils/L1AlgoEfficiencyPerformance.h +++ b/reco/L1/L1Algo/utils/L1AlgoEfficiencyPerformance.h @@ -224,9 +224,9 @@ bool L1AlgoEfficiencyPerformance<NHits>::AddOne(L1HitIndex_t* iHits) // obtain mc data for (int iih = 0; iih < NHits; iih++) { - int nMC = fL1->fvExternalHits[iHits[iih]].mcPointIds.size(); + int nMC = fL1->fvHitDebugInfo[iHits[iih]].mcPointIds.size(); for (int iMC = 0; iMC < nMC; iMC++) { - const int iMCP = fL1->fvExternalHits[iHits[iih]].mcPointIds[iMC]; + const int iMCP = fL1->fvHitDebugInfo[iHits[iih]].mcPointIds[iMC]; int mcId = fL1->fvMCPoints[iMCP].ID; mcIds[iih].push_back(mcId); } // for mcPoints @@ -254,7 +254,7 @@ bool L1AlgoEfficiencyPerformance<NHits>::AddOne(L1HitIndex_t* iHits) for (unsigned int i = 0; i < mcsN.size(); i++) { if (mcsN[i] >= 0) { trlet.mcTrackId = mcsN[i]; - trlet.iStation = fL1->fvMCPoints[fL1->fvExternalHits[iHits[0]].mcPointIds[0]].iStation; + trlet.iStation = fL1->fvMCPoints[fL1->fvHitDebugInfo[iHits[0]].mcPointIds[0]].iStation; break; } } @@ -341,7 +341,7 @@ void L1AlgoEfficiencyPerformance<NHits>::CalculateEff() // const int NPointHits = mtra.hitIds[iSta2].size(); // for (int iH = 0; iH < NPointHits; iH++){ // cout << mtra.hitIds[iSta2][iH] << ""; - // cout << "(" << fL1->fvHitStore[mtra.hitIds[iSta2][iH]].x << "\\" << fL1->fvHitStore[mtra.hitIds[iSta2][iH]].y << "= " << fL1->fvHitStore[mtra.hitIds[iSta2][iH]].x/fL1->fvHitStore[mtra.hitIds[iSta2][iH]].y << " ) "; + // cout << "(" << fL1->fvHitDebugInfo[mtra.hitIds[iSta2][iH]].x << "\\" << fL1->fvHitDebugInfo[mtra.hitIds[iSta2][iH]].y << "= " << fL1->fvHitDebugInfo[mtra.hitIds[iSta2][iH]].x/fL1->fvHitDebugInfo[mtra.hitIds[iSta2][iH]].y << " ) "; // } // cout << endl; // } diff --git a/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx b/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx index 6d5f8b2540a3e37e226b1c61f79099887684e0a6..ceef89a2046afadac6ba6b1c6dfc23d5ca789077 100644 --- a/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx +++ b/reco/L1/OffLineInterface/CbmL1GlobalTrackFinder.cxx @@ -109,44 +109,44 @@ Int_t CbmL1GlobalTrackFinder::CopyL1Tracks(CbmEvent* event) bool hasTrdHits = false; bool hasTofHits = false; for (vector<int>::iterator ih = it->Hits.begin(); ih != it->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; + const CbmL1HitId& h = L1->fvExternalHits[*ih]; - if (h.Det == 1 && hasStsHits == false) { + if (h.detId == 1 && hasStsHits == false) { hasStsHits = true; CbmStsTrack* track = new ((*fStsTracks)[stsTrackIndex]) CbmStsTrack(); t->SetStsTrackIndex(stsTrackIndex); if (event) event->AddData(ECbmDataType::kStsTrack, stsTrackIndex); - CbmL1TrackToCbmStsTrack(T, track, h.Det); + CbmL1TrackToCbmStsTrack(T, track, h.detId); stsTrackIndex++; } - if (h.Det == 2 && hasMuchHits == false) { + if (h.detId == 2 && hasMuchHits == false) { hasMuchHits = true; CbmMuchTrack* track = new ((*fMuchTracks)[muchTrackIndex]) CbmMuchTrack(); t->SetMuchTrackIndex(muchTrackIndex); if (event) event->AddData(ECbmDataType::kMuchTrack, muchTrackIndex); - CbmL1TrackToCbmMuchTrack(T, track, h.Det); + CbmL1TrackToCbmMuchTrack(T, track, h.detId); muchTrackIndex++; } - if (h.Det == 3 && hasTrdHits == false) { + if (h.detId == 3 && hasTrdHits == false) { hasTrdHits = true; CbmTrdTrack* track = new ((*fTrdTracks)[trdTrackIndex]) CbmTrdTrack(); t->SetTrdTrackIndex(trdTrackIndex); if (event) event->AddData(ECbmDataType::kTrdTrack, trdTrackIndex); - CbmL1TrackToCbmTrdTrack(T, track, h.Det); + CbmL1TrackToCbmTrdTrack(T, track, h.detId); trdTrackIndex++; } - if (h.Det == 4 && hasTofHits == false) { + if (h.detId == 4 && hasTofHits == false) { hasTofHits = true; CbmTofTrack* track = new ((*fTofTracks)[tofTrackIndex]) CbmTofTrack(); t->SetTofTrackIndex(tofTrackIndex); if (event) event->AddData(ECbmDataType::kTofTrack, tofTrackIndex); - CbmL1TrackToCbmTofTrack(T, track, h.Det); + CbmL1TrackToCbmTofTrack(T, track, h.detId); tofTrackIndex++; - // if (event) event->AddData(ECbmDataType::kTofHit, h.ExtIndex); + // if (event) event->AddData(ECbmDataType::kTofHit, h.hitId); } } //END create detector tracks if needed @@ -166,8 +166,8 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTrack(CbmL1Track l1track, CbmTrack* CbmL1* L1 = CbmL1::Instance(); for (vector<int>::iterator ih = T->Hits.begin(); ih != T->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; - if (h.Det != systemIdT) continue; + CbmL1HitId& h = L1->fvExternalHits[*ih]; + if (h.detId != systemIdT) continue; } ndf -= 5; if (ndf <= 0) ndf = 1; @@ -192,9 +192,9 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmStsTrack(CbmL1Track l1track, CbmStsT CbmL1* L1 = CbmL1::Instance(); for (vector<int>::iterator ih = T->Hits.begin(); ih != T->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; - if (h.Det != systemIdT) continue; - track->AddHit(h.ExtIndex, kSTSHIT); + CbmL1HitId& h = L1->fvExternalHits[*ih]; + if (h.detId != systemIdT) continue; + track->AddHit(h.hitId, kSTSHIT); } ndf -= 5; @@ -228,9 +228,9 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmMuchTrack(CbmL1Track l1track, CbmMuc CbmL1* L1 = CbmL1::Instance(); for (vector<int>::iterator ih = T->Hits.begin(); ih != T->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; - if (h.Det != systemIdT) continue; - track->AddHit(h.ExtIndex, kMUCHPIXELHIT); + CbmL1HitId& h = L1->fvExternalHits[*ih]; + if (h.detId != systemIdT) continue; + track->AddHit(h.hitId, kMUCHPIXELHIT); } ndf -= 5; if (ndf <= 0) ndf = 1; @@ -256,9 +256,9 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTrdTrack(CbmL1Track l1track, CbmTrdT CbmL1* L1 = CbmL1::Instance(); for (vector<int>::iterator ih = T->Hits.begin(); ih != T->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; - if (h.Det != systemIdT) continue; - track->AddHit(h.ExtIndex, kTRDHIT); + CbmL1HitId& h = L1->fvExternalHits[*ih]; + if (h.detId != systemIdT) continue; + track->AddHit(h.hitId, kTRDHIT); } ndf -= 5; if (ndf <= 0) ndf = 1; @@ -285,9 +285,9 @@ void CbmL1GlobalTrackFinder::CbmL1TrackToCbmTofTrack(CbmL1Track l1track, CbmTofT CbmL1* L1 = CbmL1::Instance(); for (vector<int>::iterator ih = T->Hits.begin(); ih != T->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; - if (h.Det != systemIdT) continue; - track->AddHit(h.ExtIndex, kTOFHIT); + CbmL1HitId& h = L1->fvExternalHits[*ih]; + if (h.detId != systemIdT) continue; + track->AddHit(h.hitId, kTOFHIT); } ndf -= 5; if (ndf <= 0) ndf = 1; diff --git a/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx b/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx index cce5fe7aec9d325778e25a85f4bfe66a77f89a85..dfd28eec7de63e290bcc51e85f8e237f5a3eaeb3 100644 --- a/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx +++ b/reco/L1/OffLineInterface/CbmL1StsTrackFinder.cxx @@ -92,20 +92,20 @@ Int_t CbmL1StsTrackFinder::CopyL1Tracks(CbmEvent* event) t->SetLastHitTimeError(T.CLast[20]); for (vector<int>::iterator ih = it->Hits.begin(); ih != it->Hits.end(); ++ih) { - CbmL1HitStore& h = L1->fvHitStore[*ih]; + CbmL1HitId& h = L1->fvExternalHits[*ih]; // double zref = L1->fpAlgo->vStations[h.iStation].z[0]; - if (h.Det > 1) { // not MVD or STS hit + if (h.detId > 1) { // not MVD or STS hit continue; } - if (h.ExtIndex < 0) { + if (h.hitId < 0) { // CbmMvdHit tmp; // tmp.SetZ(zref); - t->AddMvdHit(-h.ExtIndex - 1); //, &tmp ); + t->AddMvdHit(-h.hitId - 1); //, &tmp ); } else { //CbmStsHit tmp; //tmp.SetZ(zref); - t->AddHit(h.ExtIndex, kSTSHIT); //, &tmp ); + t->AddHit(h.hitId, kSTSHIT); //, &tmp ); } } }