From 5368c76eb518ea6814cebb0e785752a231abfd39 Mon Sep 17 00:00:00 2001
From: Dominik Smith <d.smith@gsi.de>
Date: Fri, 4 Feb 2022 16:24:55 +0100
Subject: [PATCH] CbmSeedFinderQa: Added histograms and canvas to output.

---
 reco/eventbuilder/digis/CbmSeedFinderQa.cxx | 106 +++++++++++++++++++-
 reco/eventbuilder/digis/CbmSeedFinderQa.h   |  56 ++++++++---
 2 files changed, 145 insertions(+), 17 deletions(-)

diff --git a/reco/eventbuilder/digis/CbmSeedFinderQa.cxx b/reco/eventbuilder/digis/CbmSeedFinderQa.cxx
index b16ca70a6c..eabe9e48fd 100644
--- a/reco/eventbuilder/digis/CbmSeedFinderQa.cxx
+++ b/reco/eventbuilder/digis/CbmSeedFinderQa.cxx
@@ -4,6 +4,54 @@
 
 #include "CbmSeedFinderQa.h"
 
+#include "CbmQaCanvas.h"
+
+#include "FairRootManager.h"
+#include <Logger.h>
+
+#include "TH1.h"
+#include "TH2.h"
+
+CbmSeedFinderQa::CbmSeedFinderQa() : fOutFolder("SeedFinderQA", "Seed finder QA")
+{
+  // --- Init histogram folder
+  histFolder = fOutFolder.AddFolder("hist", "Histogramms");
+
+  // --- Init histograms
+  fhCorrectDigiRatio = new TH1F("fhCorrectDigiRatio", "Correct digis per seed [pct]", 416, -2, 102);
+  fhCorrectDigiRatioNoNoise =
+    new TH1F("fhCorrectDigiRatioNoNoise", "Correct digis per seed [pct], disregarding noise", 416, -2, 102);
+  fhNoiseDigiRatio = new TH1F("fhNoiseDigiRatio", "Noise digis per seed [pct]", 416, -2, 102);
+  fhFoundDigiRatio = new TH1F("fhFoundDigiRatio", "Found digis per seed [pct]", 416, -2, 102);
+  fhCorrectVsFound = new TH2I("fhCorrectVsFound", "Correct digis  [pct] vs. Found digis [pct]; Correct; Found ", 110,
+                              -5., 105., 110, -5., 105.);
+  fhCorrectVsFoundNoNoise =
+    new TH2I("fhCorrectVsFoundNoNoise", "Correct digis  [pct] vs. Found digis [pct], no noise; Correct; Found ", 110,
+             -5., 105., 110, -5., 105.);
+
+  histFolder->Add(fhCorrectDigiRatio);
+  histFolder->Add(fhCorrectDigiRatioNoNoise);
+  histFolder->Add(fhNoiseDigiRatio);
+  histFolder->Add(fhFoundDigiRatio);
+  histFolder->Add(fhCorrectVsFound);
+  histFolder->Add(fhCorrectVsFoundNoNoise);
+
+  fCanv = new CbmQaCanvas("cSummary", "", 3 * 400, 2 * 400);
+  fCanv->Divide2D(6);
+  fOutFolder.Add(fCanv);
+}
+
+CbmSeedFinderQa::~CbmSeedFinderQa()
+{
+  delete fhCorrectDigiRatio;
+  delete fhCorrectDigiRatioNoNoise;
+  delete fhNoiseDigiRatio;
+  delete fhFoundDigiRatio;
+  delete fhCorrectVsFound;
+  delete fhCorrectVsFoundNoNoise;
+  delete fCanv;
+}
+
 void CbmSeedFinderQa::FillQaInfo(const int32_t WinStart, const int32_t WinEnd, const std::vector<CbmMatch>* vDigiMatch,
                                  const double seedTime)
 {
@@ -64,9 +112,12 @@ void CbmSeedFinderQa::FillQaInfo(const int32_t WinStart, const int32_t WinEnd, c
       if (entry == seedMatch.GetMatchedLink().GetEntry()) { correctDigiCount++; }
     }
   }
-  const double correctDigiRatio = (double) correctDigiCount / (digiCount - noiseDigiCount);
+  const double correctDigiRatio = (double) correctDigiCount / digiCount;
   fvCorrectDigiRatio.push_back(correctDigiRatio);
 
+  const double correctDigiRatioNoNoise = (double) correctDigiCount / (digiCount - noiseDigiCount);
+  fvCorrectDigiRatioNoNoise.push_back(correctDigiRatioNoNoise);
+
   //found digis of matched event in seed window
   for (uint32_t iDigi = 0; iDigi < vDigiMatch->size(); iDigi++) {
     const CbmMatch* digiMatch = &(vDigiMatch->at(iDigi));
@@ -78,6 +129,24 @@ void CbmSeedFinderQa::FillQaInfo(const int32_t WinStart, const int32_t WinEnd, c
   fvFoundDigiRatio.push_back(foundDigiRatio);
 }
 
+void CbmSeedFinderQa::FillHistos()
+{
+  for (uint32_t iEvent = 0; iEvent < fvEventMatches.size(); iEvent++) {
+    const int32_t digiCount                 = fvFullDigiCount.at(iEvent);
+    const int32_t noiseCount                = fvNoiseDigiCount.at(iEvent);
+    const double noiseDigisPercent          = 100. * (double) noiseCount / digiCount;
+    const double correctDigisPercent        = 100. * fvCorrectDigiRatio.at(iEvent);
+    const double correctDigisPercentNoNoise = 100. * fvCorrectDigiRatioNoNoise.at(iEvent);
+    const double foundDigisPercent          = 100. * fvFoundDigiRatio.at(iEvent);
+    fhCorrectDigiRatio->Fill(correctDigisPercent);
+    fhCorrectDigiRatioNoNoise->Fill(correctDigisPercentNoNoise);
+    fhNoiseDigiRatio->Fill(noiseDigisPercent);
+    fhFoundDigiRatio->Fill(foundDigisPercent);
+    fhCorrectVsFound->Fill(correctDigisPercent, foundDigisPercent);
+    fhCorrectVsFoundNoNoise->Fill(correctDigisPercentNoNoise, foundDigisPercent);
+  }
+}
+
 void CbmSeedFinderQa::OutputQa()
 {
   for (uint32_t iEvent = 0; iEvent < fvEventMatches.size(); iEvent++) {
@@ -94,9 +163,42 @@ void CbmSeedFinderQa::OutputQa()
     }
     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: " << fvCorrectDigiRatio.at(iEvent) << std::endl;
     std::cout << "Fraction of correctly matched digis in seed window (disregarding noise): "
-              << fvCorrectDigiRatio.at(iEvent) << std::endl;
+              << fvCorrectDigiRatioNoNoise.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;
   }
+  FillHistos();
+  WriteHistos();
+}
+
+void CbmSeedFinderQa::WriteHistos()
+{
+  //output histograms
+  if (!FairRootManager::Instance() || !FairRootManager::Instance()->GetSink()) {
+    LOG(error) << "No sink found";
+    return;
+  }
+
+  fCanv->cd(1);
+  fhCorrectDigiRatio->DrawCopy("colz", "");
+
+  fCanv->cd(2);
+  fhCorrectDigiRatioNoNoise->DrawCopy("colz", "");
+
+  fCanv->cd(3);
+  fhNoiseDigiRatio->DrawCopy("colz", "");
+
+  fCanv->cd(4);
+  fhFoundDigiRatio->DrawCopy("colz", "");
+
+  fCanv->cd(5);
+  fhCorrectVsFound->DrawCopy("colz", "");
+
+  fCanv->cd(6);
+  fhCorrectVsFoundNoNoise->DrawCopy("colz", "");
+
+  FairSink* sink = FairRootManager::Instance()->GetSink();
+  sink->WriteObject(&fOutFolder, nullptr);
 }
diff --git a/reco/eventbuilder/digis/CbmSeedFinderQa.h b/reco/eventbuilder/digis/CbmSeedFinderQa.h
index 415d1dd9d4..1bb924d983 100644
--- a/reco/eventbuilder/digis/CbmSeedFinderQa.h
+++ b/reco/eventbuilder/digis/CbmSeedFinderQa.h
@@ -19,33 +19,27 @@
 
 #include "CbmMatch.h"
 
+#include <TFolder.h>
+
 #include <cstdint>
 #include <iostream>
 #include <vector>
 
+class TH1F;
+class TH2I;
+class CbmQaCanvas;
+
 class CbmSeedFinderQa {
 public:
   /**
    * @brief Create the CbmSeedFinderQa object
    */
-  CbmSeedFinderQa() {};
+  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() {};
+  ~CbmSeedFinderQa();
 
   /** @brief Matches that link constructed event seeds to MC events. */
   std::vector<CbmMatch> fvEventMatches;
@@ -55,11 +49,43 @@ public:
   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). */
+  /** @brief Ratio of digis from matched MC event. */
   std::vector<double> fvCorrectDigiRatio;
+  /** @brief Ratio of digis from matched MC event (disregarding noise). */
+  std::vector<double> fvCorrectDigiRatioNoNoise;
   /** @brief Ratio of digis of matched events that were included in event seed. */
   std::vector<double> fvFoundDigiRatio;
 
+  /** @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 Finalize histograms and canvases and write to file. */
+  void WriteHistos();
+
 private:
+  TFolder* histFolder = nullptr;  /// subfolder for histograms
+  TFolder fOutFolder;             /// output folder with histos and canvases
+
+  /** Histograms **/
+  TH1F* fhCorrectDigiRatio        = nullptr;  /// correct digis per event
+  TH1F* fhCorrectDigiRatioNoNoise = nullptr;  /// correct digis per event, disregarding noise
+  TH1F* fhNoiseDigiRatio          = nullptr;  /// noise digis per event
+  TH1F* fhFoundDigiRatio          = nullptr;  /// digis found per event
+  TH2I* fhCorrectVsFound          = nullptr;  ///  correct digis per event vs found digis per event
+  TH2I* fhCorrectVsFoundNoNoise   = nullptr;  ///  correct digis per event vs found digis per event, disregarding noise
+
+  CbmQaCanvas* fCanv;  ///summary canvas
+
+  /** @brief Fill output histograms with data from current timeslice. */
+  void FillHistos();
 };
 #endif  //CbmSeedFinderQa_h
-- 
GitLab