From 8a4c0a889adcaa4f1907460be78ef200cf33cd65 Mon Sep 17 00:00:00 2001 From: "s.zharko@gsi.de" <s.zharko@gsi.de> Date: Tue, 27 Jun 2023 13:16:44 +0200 Subject: [PATCH] CA: Matching schemes fixes for MVD, MUCH and TRD --- reco/L1/CbmCaMCModule.cxx | 114 ----------------------- reco/L1/CbmCaMCModule.h | 37 ++++++++ reco/L1/CbmL1.h | 63 +++++++++++++ reco/L1/CbmL1ReadEvent.cxx | 183 ------------------------------------- 4 files changed, 100 insertions(+), 297 deletions(-) diff --git a/reco/L1/CbmCaMCModule.cxx b/reco/L1/CbmCaMCModule.cxx index 549d5d8b07..d70f501990 100644 --- a/reco/L1/CbmCaMCModule.cxx +++ b/reco/L1/CbmCaMCModule.cxx @@ -202,120 +202,6 @@ void MCModule::MatchRecoAndMC() this->InitTrackInfo(); } - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int MCModule::MatchHitWithMc<L1DetectorID::kMvd>(int iHitExt) const -{ - int iPoint = -1; - const auto* hitMatch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[L1DetectorID::kMvd]->At(iHitExt)); - assert(hitMatch); - if (hitMatch->GetNofLinks() > 0 && hitMatch->GetLink(0).GetIndex() < fpMCData->GetNofPointsUsed(L1DetectorID::kMvd)) { - iPoint = hitMatch->GetLink(0).GetIndex(); - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int MCModule::MatchHitWithMc<L1DetectorID::kSts>(int iHitExt) const -{ - int iPoint = -1; - const auto* pHitMatch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[L1DetectorID::kSts]->At(iHitExt)); - assert(pHitMatch); - if (pHitMatch->GetNofLinks() > 0) { - const auto& link = pHitMatch->GetMatchedLink(); - if (link.GetIndex() > 0) { - int index = link.GetIndex(); - int event = link.GetEntry(); - int file = link.GetFile(); - if (index > -1) { - iPoint = fpMCData->FindInternalPointIndex(L1DetectorID::kSts, index, event, file); - assert(iPoint != -1); - } - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int MCModule::MatchHitWithMc<L1DetectorID::kMuch>(int iHitExt) const -{ - int iPoint = -1; - const auto* hitMatchMuch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[L1DetectorID::kMuch]->At(iHitExt)); - if (hitMatchMuch) { - for (int iLink = 0; iLink < hitMatchMuch->GetNofLinks(); ++iLink) { - if (hitMatchMuch->GetLink(iLink).GetIndex() < fpMCData->GetNofPointsUsed(L1DetectorID::kMuch)) { - int iMc = hitMatchMuch->GetLink(iLink).GetIndex(); - if (iMc < 0) { - iPoint = -1; - break; - } - int iFile = hitMatchMuch->GetLink(iLink).GetFile(); - int iEvent = hitMatchMuch->GetLink(iLink).GetEntry(); - iPoint = fpMCData->FindInternalPointIndex(L1DetectorID::kMuch, iMc, iEvent, iFile); - assert(iPoint != -1); - } - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int MCModule::MatchHitWithMc<L1DetectorID::kTrd>(int iHitExt) const -{ - int iPoint = -1; - const auto* hitMatch = dynamic_cast<const CbmMatch*>(fvpBrHitMatches[L1DetectorID::kTrd]->At(iHitExt)); - if (hitMatch) { - int iMC = -1; - if (hitMatch->GetNofLinks() > 0) { - iMC = hitMatch->GetLink(0).GetIndex(); - assert(iMC >= 0 && iMC < fpMCData->GetNofPointsUsed(L1DetectorID::kTrd)); - - int iMc = hitMatch->GetLink(0).GetIndex(); - if (iMc < 0) return iPoint; - int iFile = hitMatch->GetLink(0).GetFile(); - int iEvent = hitMatch->GetLink(0).GetEntry(); - - iPoint = fpMCData->FindInternalPointIndex(L1DetectorID::kTrd, iMc, iEvent, iFile); - assert(iPoint != -1); - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int MCModule::MatchHitWithMc<L1DetectorID::kTof>(int iHitExt) const -{ - // TODO: SZh 23.05.2023: Use TOF interactions class - int iPoint = -1; - const auto* hitMatch = dynamic_cast<const CbmMatch*>(fvpBrHitMatches[L1DetectorID::kTof]->At(iHitExt)); - if (hitMatch) { - for (int iLink = 0; iLink < hitMatch->GetNofLinks(); ++iLink) { - int iMc = hitMatch->GetLink(iLink).GetIndex(); - int iFile = hitMatch->GetLink(iLink).GetFile(); - int iEvent = hitMatch->GetLink(iLink).GetEntry(); - if (iMc < 0) return iPoint; - - assert(iMc >= 0 && iMc < fvpBrPoints[L1DetectorID::kTof]->Size(iFile, iEvent)); - int iPointPrelim = fpMCData->FindInternalPointIndex(L1DetectorID::kTof, iMc, iEvent, iFile); - if (iPointPrelim == -1) { continue; } - iPoint = iPointPrelim; - } - } - return iPoint; -} - - - // --------------------------------------------------------------------------------------------------------------------- // void MCModule::MatchRecoAndMCTracks() diff --git a/reco/L1/CbmCaMCModule.h b/reco/L1/CbmCaMCModule.h index 767e7385c6..b46bf53fc6 100644 --- a/reco/L1/CbmCaMCModule.h +++ b/reco/L1/CbmCaMCModule.h @@ -242,6 +242,43 @@ namespace cbm::ca // ** Template function implementation ** // ********************************************** +// --------------------------------------------------------------------------------------------------------------------- +// +template<L1DetectorID DetID> +int cbm::ca::MCModule::MatchHitWithMc(int iHitExt) const +{ + int iPoint = -1; + const auto* pHitMatch = dynamic_cast<CbmMatch*>(fvpBrHitMatches[DetID]->At(iHitExt)); + assert(pHitMatch); + + if constexpr (L1DetectorID::kTof == DetID) { + for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); ++iLink) { + const auto& link = pHitMatch->GetLink(iLink); + if (link.GetIndex() > -1) { + int iPointExt = link.GetIndex(); + int iEvent = link.GetEntry(); + int iFile = link.GetFile(); + iPoint = fpMCData->FindInternalPointIndex(DetID, iPointExt, iEvent, iFile); + } + else { + return -1; + } + } // iLink + } + else { + if (pHitMatch->GetNofLinks() > 0) { + const auto& link = pHitMatch->GetMatchedLink(); + if (link.GetIndex() > -1) { + int index = link.GetIndex(); + int event = link.GetEntry(); + int file = link.GetFile(); + iPoint = fpMCData->FindInternalPointIndex(DetID, index, event, file); + } + } + } + return iPoint; +} + // --------------------------------------------------------------------------------------------------------------------- // template<L1DetectorID DetID> diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h index 8d663caccc..31efe477b3 100644 --- a/reco/L1/CbmL1.h +++ b/reco/L1/CbmL1.h @@ -660,5 +660,68 @@ private: ClassDef(CbmL1, 0); }; +// --------------------------------------------------------------------------------------------------------------------- +// +template<L1DetectorID DetID> +int CbmL1::MatchHitWithMc(int iHitExt) const +{ + const CbmMatch* pHitMatch = nullptr; + int indexShift = 0; // Number of points in detectors, standing in front + if constexpr (L1DetectorID::kMvd == DetID) { + pHitMatch = dynamic_cast<CbmMatch*>(fpMvdHitMatches->At(iHitExt)); + indexShift = 0; + } + else if constexpr (L1DetectorID::kSts == DetID) { + pHitMatch = dynamic_cast<CbmMatch*>(fpStsHitMatches->At(iHitExt)); + indexShift = fNpointsMvdAll; + } + else if constexpr (L1DetectorID::kMuch == DetID) { + pHitMatch = dynamic_cast<CbmMatch*>(fpMuchHitMatches->At(iHitExt)); + indexShift = fNpointsMvdAll + fNpointsStsAll; + } + else if constexpr (L1DetectorID::kTrd == DetID) { + pHitMatch = dynamic_cast<CbmMatch*>(fpTrdHitMatches->At(iHitExt)); + indexShift = fNpointsMvdAll + fNpointsStsAll + fNpointsMuchAll; + } + else if constexpr (L1DetectorID::kTof == DetID) { + pHitMatch = dynamic_cast<CbmMatch*>(fpTofHitMatches->At(iHitExt)); + indexShift = fNpointsMvdAll + fNpointsStsAll + fNpointsMuchAll + fNpointsTrdAll; + } + + assert(pHitMatch); + + int iPointInt = -1; + if constexpr (L1DetectorID::kTof == DetID) { + for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); ++iLink) { + const auto& link = pHitMatch->GetLink(iLink); + if (link.GetIndex() > -1) { + int iPointExt = link.GetIndex(); + int iEvent = link.GetEntry(); + int iFile = link.GetFile(); + auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iPointExt + indexShift, iEvent, iFile)); + if (itPoint == fmMCPointsLinksMap.cend()) { continue; } + iPointInt = itPoint->second; + } + else { + return -1; + } + } // iLink + } + else { + if (pHitMatch->GetNofLinks() > 0) { + const auto& link = pHitMatch->GetMatchedLink(); + if (link.GetIndex() > -1) { + int iPointExt = link.GetIndex(); + int iEvent = link.GetEntry(); + int iFile = link.GetFile(); + auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iPointExt + indexShift, iEvent, iFile)); + if (itPoint == fmMCPointsLinksMap.cend()) { return -1; } + iPointInt = itPoint->second; + } + } + } + return iPointInt; +} + #endif //_CbmL1_h_ diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx index 7920c9cf9f..6a616ba207 100644 --- a/reco/L1/CbmL1ReadEvent.cxx +++ b/reco/L1/CbmL1ReadEvent.cxx @@ -139,189 +139,6 @@ struct TmpHit { }; -// ***************************************** -// ** Hits and MC-points matching ** -// ***************************************** - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int CbmL1::MatchHitWithMc<L1DetectorID::kMvd>(int iHit) const -{ - int iPoint = -1; - if (fpMvdHitMatches) { - int iHitExt = iHit; - const auto* hitMatch = dynamic_cast<CbmMatch*>(fpMvdHitMatches->At(iHitExt)); - assert(hitMatch); - if (hitMatch->GetNofLinks() > 0 && hitMatch->GetLink(0).GetIndex() < fNpointsMvd) { - iPoint = hitMatch->GetLink(0).GetIndex(); - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int CbmL1::MatchHitWithMc<L1DetectorID::kSts>(int iHit) const -{ - //int iPointCHECK = -1; - //const auto* pHitMatch = dynamic_cast<CbmMatch*>(fpStsHitMatches->At(iHit)); - //assert(pHitMatch); // Note: it's guaranteed, that a match object (empty or full) exists for every hit in the tree - //if (pHitMatch->GetNofLinks() > 0) { - // const auto& link = pHitMatch->GetMatchedLink(); - // if (link.GetIndex() > -1) { - // int index = link.GetIndex(); - // int file = link.GetFile(); - // int event = link.GetEntry(); - // auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(index + fNpointsMvdAll, event, file)); - // assert(itPoint != fmMCPointsLinksMap.cend()); - // iPointCHECK = itPoint->second; - // } - //} - - //float bestWeight = 0.f; - //for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); ++iLink) { - // const CbmLink& link = pHitMatch->GetLink(iLink); - - // int index = link.GetIndex(); - // if (index < 0) { - // iPoint = -1; - // break; - // } - // int file = link.GetFile(); - // int event = link.GetEntry(); - - // auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(index + fNpointsMvdAll, event, file)); - // assert(itPoint != fmMCPointsLinksMap.cend()); - - // if (link.GetWeight() > bestWeight) { - // bestWeight = link.GetWeight(); - // iPoint = itPoint->second; - // } - //} - - int iPoint = -1; - const auto* sh = dynamic_cast<CbmStsHit*>(fpStsHits->At(iHit)); - - // Match MC point - if (fpStsClusterMatches) { - const auto* clusterMatchF = static_cast<const CbmMatch*>(fpStsClusterMatches->At(sh->GetFrontClusterId())); - const auto* clusterMatchB = static_cast<const CbmMatch*>(fpStsClusterMatches->At(sh->GetBackClusterId())); - CbmMatch hitMatch; - - for (int iLinkF = 0; iLinkF < clusterMatchF->GetNofLinks(); ++iLinkF) { - const auto& linkF = clusterMatchF->GetLink(iLinkF); - for (int iLinkB = 0; iLinkB < clusterMatchB->GetNofLinks(); ++iLinkB) { - const auto& linkB = clusterMatchB->GetLink(iLinkB); - if (linkF == linkB) { - hitMatch.AddLink(linkF); - hitMatch.AddLink(linkB); - } - } - } - - - float bestWeight = 0.f; - for (int iLink = 0; iLink < hitMatch.GetNofLinks(); ++iLink) { - const CbmLink& link = hitMatch.GetLink(iLink); - - int iIndex = link.GetIndex(); - if (iIndex < 0) { - iPoint = -1; - break; - } - int iFile = link.GetFile(); - int iEvent = link.GetEntry(); - - auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iIndex + fNpointsMvdAll, iEvent, iFile)); - assert(itPoint != fmMCPointsLinksMap.cend()); - - if (link.GetWeight() > bestWeight) { - bestWeight = link.GetWeight(); - iPoint = itPoint->second; - } - } - } // Match MC point - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int CbmL1::MatchHitWithMc<L1DetectorID::kMuch>(int iHit) const -{ - int iPoint = -1; - const auto* hitMatchMuch = dynamic_cast<CbmMatch*>(fpMuchHitMatches->At(iHit)); - if (hitMatchMuch) { - for (int iLink = 0; iLink < hitMatchMuch->GetNofLinks(); ++iLink) { - if (hitMatchMuch->GetLink(iLink).GetIndex() < fNpointsMuch) { - int iMc = hitMatchMuch->GetLink(iLink).GetIndex(); - if (iMc < 0) { - iPoint = -1; - break; - } - - int iFile = hitMatchMuch->GetLink(iLink).GetFile(); - int iEvent = hitMatchMuch->GetLink(iLink).GetEntry(); - - auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iMc + fNpointsMvdAll + fNpointsStsAll, iEvent, iFile)); - assert(itPoint != fmMCPointsLinksMap.cend()); - iPoint = itPoint->second; - } - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int CbmL1::MatchHitWithMc<L1DetectorID::kTrd>(int iHit) const -{ - int iPoint = -1; - const auto* hitMatch = dynamic_cast<const CbmMatch*>(fpTrdHitMatches->At(iHit)); - if (hitMatch) { - if (hitMatch->GetNofLinks() > 0) { - int iMc = hitMatch->GetLink(0).GetIndex(); - if (iMc < 0) return iPoint; - int iIndex = iMc + fNpointsMvdAll + fNpointsStsAll + fNpointsMuchAll; - int iFile = hitMatch->GetLink(0).GetFile(); - int iEvent = hitMatch->GetLink(0).GetEntry(); - - auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iIndex, iEvent, iFile)); - assert(itPoint != fmMCPointsLinksMap.cend()); - iPoint = itPoint->second; - } - } - return iPoint; -} - -// --------------------------------------------------------------------------------------------------------------------- -// -template<> -int CbmL1::MatchHitWithMc<L1DetectorID::kTof>(int iHit) const -{ - int iPoint = -1; - const auto* hitMatch = dynamic_cast<const CbmMatch*>(fpTofHitMatches->At(iHit)); - if (hitMatch) { - for (int iLink = 0; iLink < hitMatch->GetNofLinks(); ++iLink) { - int iMc = hitMatch->GetLink(iLink).GetIndex(); - int iFile = hitMatch->GetLink(iLink).GetFile(); - int iEvent = hitMatch->GetLink(iLink).GetEntry(); - if (iMc < 0) return iPoint; - assert(iMc >= 0 && iMc < fpTofPoints->Size(iFile, iEvent)); - int iIndex = iMc + fNpointsMvdAll + fNpointsStsAll + fNpointsMuchAll + fNpointsTrdAll; - - auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iIndex, iEvent, iFile)); - if (itPoint == fmMCPointsLinksMap.cend()) { continue; } - iPoint = itPoint->second; - } - } - return iPoint; -} - - // ************************** // ** Event reader ** // ************************** -- GitLab