diff --git a/reco/eventbuilder/CMakeLists.txt b/reco/eventbuilder/CMakeLists.txt index 6e49500b79f6e259c5e97b59ea4a2ed7c1e8707b..1008551332c6b6ea9fb269713b9f1998d8735126 100644 --- a/reco/eventbuilder/CMakeLists.txt +++ b/reco/eventbuilder/CMakeLists.txt @@ -66,6 +66,7 @@ digis/CbmBuildEventsIdealNew.cxx digis/CbmBuildEventsQa.cxx digis/CbmBuildEventsSimple.cxx digis/CbmSeedFinderSlidingWindow.cxx +digis/CbmSeedFinderQa.cxx #digis/CbmEvBuildSource.cxx tracks/CbmBuildEventsFromTracksIdeal.cxx diff --git a/reco/eventbuilder/digis/CbmSeedFinderQa.cxx b/reco/eventbuilder/digis/CbmSeedFinderQa.cxx new file mode 100644 index 0000000000000000000000000000000000000000..b16ca70a6c7219c3246ee287664cfd919f5747ea --- /dev/null +++ b/reco/eventbuilder/digis/CbmSeedFinderQa.cxx @@ -0,0 +1,102 @@ +/* Copyright (C) 2007-2022 Facility for Antiproton and Ion Research in Europe, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Dominik Smith [committer] */ + +#include "CbmSeedFinderQa.h" + +void CbmSeedFinderQa::FillQaInfo(const int32_t WinStart, const int32_t WinEnd, const std::vector<CbmMatch>* vDigiMatch, + const double seedTime) +{ + fvSeedTimesFull.push_back(seedTime); + + int32_t digiCount = 0; + int32_t noiseDigiCount = 0; + int32_t correctDigiCount = 0; + int32_t matchedEventDigiCount = 0; + CbmMatch seedMatch; + + for (int32_t iDigi = WinStart; iDigi <= WinEnd; iDigi++) { + const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); + digiCount++; + if (digiMatch->GetNofLinks() == 0) { + //skip digis with no links to avoid T0 pollution + noiseDigiCount++; + continue; + } + if (digiMatch->GetMatchedLink().GetEntry() == -1) { + noiseDigiCount++; + continue; //disregard noise digis + } + + // Update seed match with digi links + for (int32_t iLink = 0; iLink < digiMatch->GetNofLinks(); iLink++) { + int32_t entry = digiMatch->GetLink(iLink).GetEntry(); + if (entry == -1) { continue; } //ignore noise links + int32_t file = digiMatch->GetLink(iLink).GetFile(); + double weight = digiMatch->GetLink(iLink).GetWeight(); + // LOG(info) << "Adding link (number, weight, entry, file): " << iLink << " " << weight << " " + // << entry << " " << file; + seedMatch.AddLink(weight, 0, entry, file); + } + } + fvFullDigiCount.push_back(digiCount); + fvNoiseDigiCount.push_back(noiseDigiCount); + + if (seedMatch.GetNofLinks() == 0) //seed is only noise digis + { + seedMatch.AddLink(1.0, 0, -1, 0); + fvEventMatches.push_back(seedMatch); + fvCorrectDigiRatio.push_back(-1.0); + fvFoundDigiRatio.push_back(-1.0); + return; + } + else { + fvEventMatches.push_back(seedMatch); + } + + //correct digis in seed window + for (int32_t iDigi = WinStart; iDigi <= WinEnd; iDigi++) { + const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); + if (digiMatch->GetNofLinks() == 0) { continue; } //skip digis with no links to avoid T0 pollution + const int32_t entry = digiMatch->GetMatchedLink().GetEntry(); + if (entry != -1) // disregarding noise + { + if (entry == seedMatch.GetMatchedLink().GetEntry()) { correctDigiCount++; } + } + } + const double correctDigiRatio = (double) correctDigiCount / (digiCount - noiseDigiCount); + fvCorrectDigiRatio.push_back(correctDigiRatio); + + //found digis of matched event in seed window + for (uint32_t iDigi = 0; iDigi < vDigiMatch->size(); iDigi++) { + const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); + if (digiMatch->GetNofLinks() == 0) { continue; } //skip digis with no links to avoid T0 pollution + const int matchedEvent = digiMatch->GetMatchedLink().GetEntry(); + if (matchedEvent == seedMatch.GetMatchedLink().GetEntry()) { matchedEventDigiCount++; } + } + const double foundDigiRatio = (double) correctDigiCount / matchedEventDigiCount; + fvFoundDigiRatio.push_back(foundDigiRatio); +} + +void CbmSeedFinderQa::OutputQa() +{ + for (uint32_t iEvent = 0; iEvent < fvEventMatches.size(); iEvent++) { + const CbmMatch* match = &(fvEventMatches.at(iEvent)); + const int32_t mcEventNr = match->GetMatchedLink().GetEntry(); + + std::cout << "QA for seed # " << iEvent << std::endl; + std::cout << "Seed time: " << fvSeedTimesFull.at(iEvent) << std::endl; + std::cout << "Links to MC events: " << match->GetNofLinks() << ", matched MC event number " << mcEventNr + << std::endl; + if (mcEventNr == -1) { + std::cout << "Warning: Seed was constructed from noise digis only (MC event = -1)!" << std::endl; + std::cout << " Please increase your noise threshold!" << std::endl; + } + std::cout << "Total digis in seed window: " << fvFullDigiCount.at(iEvent); + std::cout << ", Noise digis in seed window: " << fvNoiseDigiCount.at(iEvent) << std::endl; + std::cout << "Fraction of correctly matched digis in seed window (disregarding noise): " + << fvCorrectDigiRatio.at(iEvent) << std::endl; + std::cout << "Fraction of digis of matched event found in seed window: " << fvFoundDigiRatio.at(iEvent); + std::cout << " (only from this timeslice)" << std::endl; + } +} diff --git a/reco/eventbuilder/digis/CbmSeedFinderQa.h b/reco/eventbuilder/digis/CbmSeedFinderQa.h new file mode 100644 index 0000000000000000000000000000000000000000..415d1dd9d490678d813b5808c4d178426bf4d521 --- /dev/null +++ b/reco/eventbuilder/digis/CbmSeedFinderQa.h @@ -0,0 +1,65 @@ +/* Copyright (C) 2007-2022 Facility for Antiproton and Ion Research in Europe, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: Dominik Smith [committer] */ + +/** + * @file CbmSeedFinderQa.h + * @author Dominik Smith (d.smith@gsi.de) + * @brief Class for sliding window seed finder + * @version 0.1 + * @date 2022-02-04 + * + * @copyright Copyright (c) 2022 + * + * Qa output for the seed finder. + */ + +#ifndef CbmSeedFinderQa_h +#define CbmSeedFinderQa_h + +#include "CbmMatch.h" + +#include <cstdint> +#include <iostream> +#include <vector> + +class CbmSeedFinderQa { +public: + /** + * @brief Create the CbmSeedFinderQa object + */ + CbmSeedFinderQa() {}; + CbmSeedFinderQa(const CbmSeedFinderQa&) = delete; + CbmSeedFinderQa operator=(const CbmSeedFinderQa&) = delete; + + /** @brief Gather QA Information. + * @params WinStart Starting position of seed window. + * @params WinStart End position of seed window. + * @params vDigiMatch Input vector of digi matches (should match input data to MC events). + * @params seedTime Current seed. + **/ + void FillQaInfo(const int32_t WinStart, const int32_t WinEnd, const std::vector<CbmMatch>* vDigiMatch, + const double seedTime); + + /** @brief Output QA Information. */ + void OutputQa(); + + /** @brief Destructor **/ + ~CbmSeedFinderQa() {}; + + /** @brief Matches that link constructed event seeds to MC events. */ + std::vector<CbmMatch> fvEventMatches; + /** @brief Full vector of all event seeds that is not cleared at the end of a timeslice. */ + std::vector<double> fvSeedTimesFull; + /** @brief Counts how many digis contributed to a seed. */ + std::vector<int32_t> fvFullDigiCount; + /** @brief Counts how many noise digis contributed to a seed. */ + std::vector<int32_t> fvNoiseDigiCount; + /** @brief Ratio of digis from matched MC event (disregarding noise). */ + std::vector<double> fvCorrectDigiRatio; + /** @brief Ratio of digis of matched events that were included in event seed. */ + std::vector<double> fvFoundDigiRatio; + +private: +}; +#endif //CbmSeedFinderQa_h diff --git a/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.cxx b/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.cxx index 1bc8c3e29fb5aaf117652fed2ebd181d6d1e8e76..8dd3e7dbf760069875f0ea52c868f1c812e21713 100644 --- a/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.cxx +++ b/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.cxx @@ -9,18 +9,14 @@ #include "CbmMuchDigi.h" #include "CbmPsdDigi.h" #include "CbmRichDigi.h" +#include "CbmSeedFinderQa.h" #include "CbmStsDigi.h" #include "CbmTofDigi.h" #include "CbmTrdDigi.h" CbmSeedFinderSlidingWindow::~CbmSeedFinderSlidingWindow() { - if (fvEventMatches != nullptr) { delete fvEventMatches; } - if (fvSeedTimesFull != nullptr) { delete fvSeedTimesFull; } - if (fvFullDigiCount != nullptr) { delete fvFullDigiCount; } - if (fvNoiseDigiCount != nullptr) { delete fvNoiseDigiCount; } - if (fvCorrectDigiRatio != nullptr) { delete fvCorrectDigiRatio; } - if (fvFoundDigiRatio != nullptr) { delete fvFoundDigiRatio; } + if (fQa != nullptr) { delete fQa; } } template<> @@ -65,12 +61,12 @@ void CbmSeedFinderSlidingWindow::FillSeedTimes(const std::vector<inType>* vIn, c if (nDigisWin >= fminDigis) { // Reached required number of digis - fvSeedTimes->push_back((currentT + winStartT) - / 2.); //one possibility. perhaps better to place seed at end of window. + const double seedTime = + (currentT + winStartT) / 2.; //one possibility. perhaps better to place seed at end of window. + fvSeedTimes->push_back(seedTime); - if (vDigiMatch && fvEventMatches) { // QA mode - FillEventMatch(winStartN, i, vDigiMatch); - fvSeedTimesFull->push_back((currentT + winStartT) / 2.); + if (vDigiMatch && fQa) { // QA mode + fQa->FillQaInfo(winStartN, i, vDigiMatch, seedTime); } nDigisWin = 0; @@ -88,7 +84,7 @@ void CbmSeedFinderSlidingWindow::FillSeedTimes(const std::vector<inType>* vIn, c } } - if (vDigiMatch && fvEventMatches) { // QA mode + if (fQa && vDigiMatch) { // QA mode std::cout << "CbmSeedFinderSlidingWindow::FillSeedTimes(): Found " << GetNofSeeds() << " seeds for this timeslice." << std::endl; } @@ -120,128 +116,17 @@ double CbmSeedFinderSlidingWindow::GetTime(const std::vector<double>* vIn, int32 return vIn->at(i); } -void CbmSeedFinderSlidingWindow::FillEventMatch(int32_t WinStart, int32_t WinEnd, - const std::vector<CbmMatch>* vDigiMatch) -{ - int32_t digiCount = 0; - int32_t noiseDigiCount = 0; - int32_t correctDigiCount = 0; - int32_t matchedEventDigiCount = 0; - CbmMatch seedMatch; - - for (int32_t iDigi = WinStart; iDigi <= WinEnd; iDigi++) { - const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); - digiCount++; - if (digiMatch->GetNofLinks() == 0) { - //skip digis with no links to avoid T0 pollution - noiseDigiCount++; - continue; - } - if (digiMatch->GetMatchedLink().GetEntry() == -1) { - noiseDigiCount++; - continue; //disregard noise digis - } - - // Update seed match with digi links - for (int32_t iLink = 0; iLink < digiMatch->GetNofLinks(); iLink++) { - int32_t entry = digiMatch->GetLink(iLink).GetEntry(); - if (entry == -1) { continue; } //ignore noise links - int32_t file = digiMatch->GetLink(iLink).GetFile(); - double weight = digiMatch->GetLink(iLink).GetWeight(); - // LOG(info) << "Adding link (number, weight, entry, file): " << iLink << " " << weight << " " - // << entry << " " << file; - seedMatch.AddLink(weight, 0, entry, file); - } - } - fvFullDigiCount->push_back(digiCount); - fvNoiseDigiCount->push_back(noiseDigiCount); - - if (seedMatch.GetNofLinks() == 0) //seed is only noise digis - { - seedMatch.AddLink(1.0, 0, -1, 0); - fvEventMatches->push_back(seedMatch); - fvCorrectDigiRatio->push_back(-1.0); - fvFoundDigiRatio->push_back(-1.0); - return; - } - else { - fvEventMatches->push_back(seedMatch); - } - - //correct digis in seed window - for (int32_t iDigi = WinStart; iDigi <= WinEnd; iDigi++) { - const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); - if (digiMatch->GetNofLinks() == 0) { continue; } //skip digis with no links to avoid T0 pollution - const int32_t entry = digiMatch->GetMatchedLink().GetEntry(); - if (entry != -1) // disregarding noise - { - if (entry == seedMatch.GetMatchedLink().GetEntry()) { correctDigiCount++; } - } - } - const double correctDigiRatio = (double) correctDigiCount / (digiCount - noiseDigiCount); - fvCorrectDigiRatio->push_back(correctDigiRatio); - - //found digis of matched event in seed window - for (uint32_t iDigi = 0; iDigi < vDigiMatch->size(); iDigi++) { - const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi)); - if (digiMatch->GetNofLinks() == 0) { continue; } //skip digis with no links to avoid T0 pollution - const int matchedEvent = digiMatch->GetMatchedLink().GetEntry(); - if (matchedEvent == seedMatch.GetMatchedLink().GetEntry()) { matchedEventDigiCount++; } - } - const double foundDigiRatio = (double) correctDigiCount / matchedEventDigiCount; - fvFoundDigiRatio->push_back(foundDigiRatio); -} - void CbmSeedFinderSlidingWindow::SetQa(bool doQA) { if (doQA == true) { - if (fvEventMatches == nullptr) { fvEventMatches = new std::vector<CbmMatch>; }; - if (fvSeedTimesFull == nullptr) { fvSeedTimesFull = new std::vector<double>; }; - if (fvFullDigiCount == nullptr) { fvFullDigiCount = new std::vector<int32_t>; }; - if (fvNoiseDigiCount == nullptr) { fvNoiseDigiCount = new std::vector<int32_t>; }; - if (fvCorrectDigiRatio == nullptr) { fvCorrectDigiRatio = new std::vector<double>; }; - if (fvFoundDigiRatio == nullptr) { fvFoundDigiRatio = new std::vector<double>; }; - fvEventMatches->clear(); - fvFullDigiCount->clear(); - fvNoiseDigiCount->clear(); - fvSeedTimesFull->clear(); - fvCorrectDigiRatio->clear(); - fvFoundDigiRatio->clear(); + if (fQa == nullptr) { fQa = new CbmSeedFinderQa(); } } else { - if (fvEventMatches != nullptr) { delete fvEventMatches; } - if (fvSeedTimesFull != nullptr) { delete fvSeedTimesFull; } - if (fvFullDigiCount != nullptr) { delete fvFullDigiCount; } - if (fvNoiseDigiCount != nullptr) { delete fvNoiseDigiCount; } - if (fvCorrectDigiRatio != nullptr) { delete fvCorrectDigiRatio; } - if (fvFoundDigiRatio != nullptr) { delete fvFoundDigiRatio; } + if (fQa != nullptr) { delete fQa; } } } void CbmSeedFinderSlidingWindow::OutputQa() { - if (fvEventMatches == nullptr) { - std::cout << "CbmSeedFinderSlidingWindow::OutputQa(): Error, QA mode not active (no event matches present)." - << std::endl; - exit(1); - } - for (uint32_t iEvent = 0; iEvent < fvEventMatches->size(); iEvent++) { - const CbmMatch* match = &(fvEventMatches->at(iEvent)); - const int32_t mcEventNr = match->GetMatchedLink().GetEntry(); - - std::cout << "QA for seed # " << iEvent << std::endl; - std::cout << "Seed time: " << fvSeedTimesFull->at(iEvent) << std::endl; - std::cout << "Links to MC events: " << match->GetNofLinks() << ", matched MC event number " << mcEventNr - << std::endl; - if (mcEventNr == -1) { - std::cout << "Warning: Seed was constructed from noise digis only (MC event = -1)!" << std::endl; - std::cout << " Please increase your noise threshold!" << std::endl; - } - std::cout << "Total digis in seed window: " << fvFullDigiCount->at(iEvent); - std::cout << ", Noise digis in seed window: " << fvNoiseDigiCount->at(iEvent) << std::endl; - std::cout << "Fraction of correctly matched digis in seed window (disregarding noise): " - << fvCorrectDigiRatio->at(iEvent) << std::endl; - std::cout << "Fraction of digis of matched event found in seed window: " << fvFoundDigiRatio->at(iEvent); - std::cout << " (only from this timeslice)" << std::endl; - } + if (fQa) { fQa->OutputQa(); } } diff --git a/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.h b/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.h index c3991887fac5775955b8aab74c446196914bfbf7..6bcc7484cfd577a70a66a6943a0076f9932bdec4 100644 --- a/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.h +++ b/reco/eventbuilder/digis/CbmSeedFinderSlidingWindow.h @@ -21,12 +21,11 @@ #ifndef CbmSeedFinderSlidingWindow_h #define CbmSeedFinderSlidingWindow_h +#include <cstdint> #include <vector> -#include <stdint.h> -#include <stdio.h> - class CbmMatch; +class CbmSeedFinderQa; class CbmSeedFinderSlidingWindow { public: @@ -66,22 +65,12 @@ public: size_t GetNofSeeds() { return fvSeedTimes->size(); } private: + /** @brief Processes QA info. */ + CbmSeedFinderQa* fQa = nullptr; + /** @brief Output of the algorithm. Stores seed times for current time slice. */ std::vector<double>* fvSeedTimes = nullptr; - /** @brief For QA. Matches that link constructed event seeds to MC events. */ - std::vector<CbmMatch>* fvEventMatches = nullptr; - /** @brief For QA. Full vector of all event seeds that is not cleared at the end of a timeslice. */ - std::vector<double>* fvSeedTimesFull = nullptr; - /** @brief For QA. Counts how many digis contributed to a seed. */ - std::vector<int32_t>* fvFullDigiCount = nullptr; - /** @brief For QA. Counts how many noise digis contributed to a seed. */ - std::vector<int32_t>* fvNoiseDigiCount = nullptr; - /** @brief For QA. Ratio of digis from matched MC event (disregarding noise). */ - std::vector<double>* fvCorrectDigiRatio = nullptr; - /** @brief For QA. Ratio of digis of matched events that were included in event seed. */ - std::vector<double>* fvFoundDigiRatio = nullptr; - /** @brief Minimum number of digis which must be found in the seed window. */ int32_t fminDigis = 0; /** @brief Size of sliding window. */ @@ -92,12 +81,5 @@ private: /** @brief Fetches time at position i of either a digi vector or vector of times. */ template<class inType> double GetTime(const std::vector<inType>* vIn, int32_t i); - - /** @brief For QA. Fills fvEventMatches (matches that link constructed event seeds to MC events). - * @params WinStart Starting position of seed window. - * @params WinStart End position of seed window. - * @params vDigiMatch Input vector of digi matches (should match input data to MC events). - */ - void FillEventMatch(int32_t WinStart, int32_t WinEnd, const std::vector<CbmMatch>* vDigiMatch); }; #endif //CbmSeedFinderSlidingWindow_h