Add tools to study old/new Bmon

- tool to check data from TSA w/ legacy unp
- tool to compare digis found in the old BMon VS those found in the new BMon
  (from TSA with either legacy or online unpacker)
- tool to study BMon bunches size, extension and counts

macro/beamtime/mcbm2024/bmon_checks/check_bmon_digis_bunches.C

+/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+   SPDX-License-Identifier: GPL-3.0-only
+   Authors: Pierre-Alain Loizeau [committer] */
+
+/** @brief Macro for check of BMON digis bunching in direct vector from legacy or online unpacking
+ ** @param input          Name of input file
+ ** @param nTimeSlices    Number of time-slices to process
+ **/
+void check_bmon_digis_bunches(TString inputFileName, uint32_t uRunId, size_t numTimeslices = 0,
+                              double_t dBunchCutOff = 40.0, double_t dBunchRange = 500.0, double_t dBunchStep = 5.0,
+                              bool bNbDigisZ = true)
+{
+  gROOT->cd();
+  TFile* file = new TFile(inputFileName, "READ");
+  TTree* tree = (TTree*) (file->Get("cbmsim"));
+
+  std::vector<CbmBmonDigi>* vDigisBmon = new std::vector<CbmBmonDigi>();
+  tree->SetBranchAddress("BmonDigi", &vDigisBmon);
+
+  uint32_t nentries = tree->GetEntries();
+  cout << "Entries: " << nentries << endl;
+  nentries = (numTimeslices && numTimeslices < nentries ? numTimeslices : nentries);
+
+  gROOT->cd();
+
+  int32_t iNbBinsRange = std::round(dBunchRange / dBunchStep) + 1;
+
+  TString sName  = "histRangeMulOld";
+  TString sTitle = Form("Digis multiplicity vs range old BMON, run %u; Range to first digi [ns]; Digis mul []", uRunId);
+  TH2* fHistRangeMulOld = new TH2I(sName, sTitle, iNbBinsRange, 0.0, dBunchRange + dBunchStep, 300, 1.0, 300.0);
+
+  sName  = "histRangeMulNew";
+  sTitle = Form("Digis multiplicity vs range new BMON, run %u; Range to first digi [ns]; Digis mul []", uRunId);
+  TH2* fHistRangeMulNew = new TH2I(sName, sTitle, iNbBinsRange, 0.0, dBunchRange + dBunchStep, 300, 0.0, 300.0);
+
+  sName  = "histBunchMulOld";
+  sTitle = Form("Bunch multiplicity vs duration old BMON, cutoff %.2f ns, run %u; Bunch range [ns]; Bunch mul []; %s",
+                dBunchCutOff, uRunId, bNbDigisZ ? " Digis in bunches []" : " Bunches []");
+  TH2* fHistBunchMulOld = new TH2I(sName, sTitle, iNbBinsRange, 0.0, dBunchRange + dBunchStep, 300, 1.0, 300.0);
+
+  sName  = "histBunchMulNew";
+  sTitle = Form("Bunch multiplicity vs duration old BMON, cutoff %.2f ns, run %u; Bunch range [ns]; Bunch mul []; %s",
+                dBunchCutOff, uRunId, bNbDigisZ ? " Digis in bunches []" : " Bunches []");
+  TH2* fHistBunchMulNew = new TH2I(sName, sTitle, iNbBinsRange, 0.0, dBunchRange + dBunchStep, 300, 0.0, 300.0);
+
+  for (Int_t iEntry = 0; iEntry < nentries; iEntry++) {
+    //if (iEntry % 10 == 0 ) std::cout << "Entry " << iEntry << " / " << nentries << std::endl;
+
+    tree->GetEntry(iEntry);
+    uint32_t nDigisBmon = vDigisBmon->size();
+
+    std::vector<CbmBmonDigi> vDigisOld;
+    std::vector<CbmBmonDigi> vDigisNew;
+    vDigisOld.reserve(nDigisBmon);
+    vDigisNew.reserve(nDigisBmon);
+    for (auto& digi : *vDigisBmon) {
+      if (1 == CbmTofAddress::GetChannelSide(digi.GetAddress())) {
+        if (CbmTofAddress::GetChannelId(digi.GetAddress()) < 4) {
+          vDigisNew.push_back(digi);
+        }
+        else {
+          LOG(fatal) << "Bad sCVD channel: " << CbmTofAddress::GetChannelId(digi.GetAddress());
+        }
+      }
+      else {
+        vDigisOld.push_back(digi);
+      }
+    }
+    std::cout << "TS " << iEntry << ", BMON  Digis: " << nDigisBmon << " => " << vDigisOld.size() << " Old + "
+              << vDigisNew.size() << " sCVD" << std::endl;
+
+    for (auto itOld = vDigisOld.begin(); itOld != vDigisOld.end(); ++itOld) {
+      auto nextDigiOld = itOld;
+      ++nextDigiOld;
+      uint32_t uRangeMul = 1;
+      while (nextDigiOld != vDigisOld.end()) {
+        double_t dDt = (*nextDigiOld).GetTime() - (*itOld).GetTime();
+        if (dDt < dBunchRange) {
+          uRangeMul++;
+          fHistRangeMulOld->Fill(dDt, uRangeMul);
+          ++nextDigiOld;
+        }
+        else {
+          break;
+        }
+      }
+      if (1 == uRangeMul) {
+        fHistRangeMulOld->Fill(dBunchRange, uRangeMul);
+      }
+    }
+    for (auto itNew = vDigisNew.begin(); itNew != vDigisNew.end(); ++itNew) {
+      auto nextDigiNew = itNew;
+      ++nextDigiNew;
+      uint32_t uRangeMul = 1;
+      while (nextDigiNew != vDigisNew.end()) {
+        double_t dDt = (*nextDigiNew).GetTime() - (*itNew).GetTime();
+        if (dDt < dBunchRange) {
+          uRangeMul++;
+          fHistRangeMulNew->Fill(dDt, uRangeMul);
+          ++nextDigiNew;
+        }
+        else {
+          break;
+        }
+      }
+      if (1 == uRangeMul) {
+        fHistRangeMulNew->Fill(dBunchRange, uRangeMul);
+      }
+    }
+
+    for (auto itOld = vDigisOld.begin(); itOld != vDigisOld.end();) {
+      auto currDigiOld = itOld;
+      auto nextDigiOld = itOld;
+      ++nextDigiOld;
+      uint32_t uBunchMul = 1;
+      while (nextDigiOld != vDigisOld.end()) {
+        double_t dDt      = (*nextDigiOld).GetTime() - (*itOld).GetTime();
+        double_t dDtDigis = (*nextDigiOld).GetTime() - (*currDigiOld).GetTime();
+        if (dDt < dBunchRange && dDtDigis < dBunchCutOff) {
+          uBunchMul++;
+          currDigiOld = nextDigiOld;
+          ++nextDigiOld;
+        }
+        else {
+          /// End of bunch, fill intervall between first and last digi in bunch + multiplicity
+          double_t dDtBunch = (*currDigiOld).GetTime() - (*itOld).GetTime();
+          // For bNbDigisZ = true, Z = Nb Digis in bunches, else Z = Nb bunches
+          fHistBunchMulOld->Fill(dDtBunch, uBunchMul, bNbDigisZ ? uBunchMul : 1);
+          break;
+        }
+      }
+      /// Go to first digi out of bunch
+      itOld = nextDigiOld;
+      if (1 == uBunchMul) {
+        fHistBunchMulOld->Fill(0.0, 1);
+      }
+    }
+    for (auto itNew = vDigisNew.begin(); itNew != vDigisNew.end();) {
+      auto currDigiNew = itNew;
+      auto nextDigiNew = itNew;
+      ++nextDigiNew;
+      uint32_t uBunchMul = 1;
+      while (nextDigiNew != vDigisNew.end()) {
+        double_t dDt      = (*nextDigiNew).GetTime() - (*itNew).GetTime();
+        double_t dDtDigis = (*nextDigiNew).GetTime() - (*currDigiNew).GetTime();
+        if (dDt < dBunchRange && dDtDigis < dBunchCutOff) {
+          uBunchMul++;
+          currDigiNew = nextDigiNew;
+          ++nextDigiNew;
+        }
+        else {
+          /// End of bunch, fill intervall between first and last digi in bunch + multiplicity
+          double_t dDtBunch = (*currDigiNew).GetTime() - (*itNew).GetTime();
+          // For bNbDigisZ = true, Z = Nb Digis in bunches, else Z = Nb bunches
+          fHistBunchMulNew->Fill(dDtBunch, uBunchMul, bNbDigisZ ? uBunchMul : 1);
+          break;
+        }
+      }
+      /// Go to first digi out of bunch
+      itNew = nextDigiNew;
+      if (1 == uBunchMul) {
+        fHistBunchMulNew->Fill(0.0, 1);
+      }
+    }
+  }
+
+  TCanvas* fCanvBunches = new TCanvas("canvBunches", "BMON Digis Bunch multiplicity at various short ranges");
+  fCanvBunches->Divide(2, 2);
+
+  fCanvBunches->cd(1);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistRangeMulOld->Draw("Colz");
+
+  fCanvBunches->cd(2);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistRangeMulNew->Draw("Colz");
+
+  fCanvBunches->cd(3);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistBunchMulOld->Draw("Colz");
+
+  fCanvBunches->cd(4);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistBunchMulNew->Draw("Colz");
+}

macro/beamtime/mcbm2024/bmon_checks/check_bmon_legacy_digis_scvd.C

+/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
+   SPDX-License-Identifier: GPL-3.0-only
+   Authors: Pierre-Alain Loizeau [committer] */
+
+/** @brief Macro for check of sCVD BMON digis in direct vector from legacy unpacking
+ ** @param input          Name of input file
+ ** @param nTimeSlices    Number of time-slices to process
+ **/
+void check_bmon_legacy_digis_scvd(TString inputFileName, size_t numTimeslices = 0, double_t dPileUpThrNs = 20.0)
+{
+  gROOT->cd();
+  TFile* file = new TFile(inputFileName, "READ");
+  TTree* tree = (TTree*) (file->Get("cbmsim"));
+
+  std::vector<CbmBmonDigi>* vDigisBmon = new std::vector<CbmBmonDigi>();
+  tree->SetBranchAddress("BmonDigi", &vDigisBmon);
+
+  uint32_t nentries = tree->GetEntries();
+  cout << "Entries: " << nentries << endl;
+  nentries = (numTimeslices && numTimeslices < nentries ? numTimeslices : nentries);
+
+  gROOT->cd();
+  TH1* fHistMapBmonOld     = new TH1I("histMapBmonOld", "Channel map, old BMON; Strip []", 16, 0., 16.);
+  TH2* fHistMapBmonScvd    = new TH2I("histMapBmonScvd", "Pad map, sCVD BMON; Pad X []; Pad Y []",  //
+                                   2, 0., 2., 2, 0., 2.);
+  TH2* fHistMapEvoBmonOld  = new TH2I("histMapEvoBmonOld", "Pad map, old BMON; TS []; Strip []",  //
+                                     100, 0., 1000., 16, 0., 16.);
+  TH2* fHistMapEvoBmonScvd = new TH2I("histMapEvoBmonScvd", "Pad map, sCVD BMON; TS []; Channel []",  //
+                                      100, 0., 1000., 4, 0., 4.);
+
+  double_t dDtMin   = -500.;
+  double_t dDtMax   = 500.;
+  int32_t iDtNbBins = std::round(dDtMax - dDtMin);
+  TH1* fHistDtBmon =
+    new TH1I("histDtBmon", "Time difference old vs sCVD BMON; dt [ns]", 10 * iDtNbBins, dDtMin, dDtMax);
+  TH2* fHistDtEvoBmon = new TH2I("histDtEvoBmon", "Evolution Time difference old vs sCVD BMON ; TS []; dt [ns]",  //
+                                 100, 0., 1000., iDtNbBins, dDtMin, dDtMax);
+  TH2* fHistDtDxBmon  = new TH2I("histDtDxBmon", "X correlation vs Time diff, old vs sCVD BMON; dt [ns]; dX []",  //
+                                1000, -500., 500., 33, -16.5, 16.5);
+  TH2* fHistDxCorBmon = new TH2I("histDxCorrBmon", "Pad map, old vs sCVD BMON; Strip []; Pad []",  //
+                                 16, 0., 16., 4, 0., 4.);
+
+  TH2* fHistDtChanOld  = new TH2I("histDtChanOld", "Time difference old vs sCVD BMON VS Chan in old ; Ch []; dt [ns]",
+                                 16, 0., 16., iDtNbBins, dDtMin, dDtMax);
+  TH2* fHistDtChanScvd = new TH2I("histDtChanScvd", "Time difference old vs sCVD BMON VS Chan in Scvd ; Ch []; dt [ns]",
+                                  4, 0., 4., iDtNbBins, dDtMin, dDtMax);
+
+  double_t dSelfDtMin   = 0.;
+  double_t dSelfDtMax   = 10000.;
+  int32_t iSelfDtNbBins = std::round(dSelfDtMax - dSelfDtMin);
+  TH1* fHistSelfDtOld =
+    new TH1I("histSelfDtOld", "Time difference old BMON digis; self dt [ns]", iSelfDtNbBins, dSelfDtMin, dSelfDtMax);
+  TH1* fHistSelfDtNew =
+    new TH1I("histSelfDtNew", "Time difference new BMON digis; self dt [ns]", iSelfDtNbBins, dSelfDtMin, dSelfDtMax);
+  TH1* fHistOldInPileup =
+    new TH1I("histOldInPileup", Form("Fractions of old BMON digis in pileup (dt < %.1fns); Pileup?", dPileUpThrNs), 2,
+             -0.5, 1.5);
+  TH1* fHistNewInPileup =
+    new TH1I("histNewInPileup", Form("Fractions of new BMON digis in pileup (dt < %.1fns); Pileup?", dPileUpThrNs), 2,
+             -0.5, 1.5);
+
+  uint8_t ucScvdX[4] = {1, 1, 0, 0};
+  uint8_t ucScvdY[4] = {1, 0, 0, 1};
+
+  for (Int_t iEntry = 0; iEntry < nentries; iEntry++) {
+    //if (iEntry % 10 == 0 ) std::cout << "Entry " << iEntry << " / " << nentries << std::endl;
+
+    tree->GetEntry(iEntry);
+    uint32_t nDigisBmon = vDigisBmon->size();
+
+    std::vector<CbmBmonDigi> vDigisOld;
+    std::vector<CbmBmonDigi> vDigisScvd;
+    vDigisOld.reserve(nDigisBmon);
+    vDigisScvd.reserve(nDigisBmon);
+    for (auto& digi : *vDigisBmon) {
+      if (1 == CbmTofAddress::GetChannelSide(digi.GetAddress())) {
+        if (CbmTofAddress::GetChannelId(digi.GetAddress()) < 4) {
+          fHistMapBmonScvd->Fill(ucScvdX[CbmTofAddress::GetChannelId(digi.GetAddress())],
+                                 ucScvdY[CbmTofAddress::GetChannelId(digi.GetAddress())]);
+          fHistMapEvoBmonScvd->Fill(iEntry, CbmTofAddress::GetChannelId(digi.GetAddress()));
+          vDigisScvd.push_back(digi);
+        }
+        else {
+          LOG(fatal) << "Bad sCVD channel: " << CbmTofAddress::GetChannelId(digi.GetAddress());
+        }
+      }
+      else {
+        fHistMapBmonOld->Fill(CbmTofAddress::GetChannelId(digi.GetAddress()));
+        fHistMapEvoBmonOld->Fill(iEntry, CbmTofAddress::GetChannelId(digi.GetAddress()));
+        vDigisOld.push_back(digi);
+      }
+    }
+    std::cout << "TS " << iEntry << ", BMON  Digis: " << nDigisBmon << " => " << vDigisOld.size() << " Old + "
+              << vDigisScvd.size() << " sCVD" << std::endl;
+
+    auto itScvdStart = vDigisScvd.begin();
+    for (auto& digiOld : vDigisOld) {
+      for (auto itScvd = itScvdStart; itScvd != vDigisScvd.end(); ++itScvd) {
+        double_t dDt = (*itScvd).GetTime() - digiOld.GetTime();
+        if (dDt < dDtMin) {
+          itScvdStart = itScvd;
+          continue;
+        }
+        if (dDtMax < dDt) {
+          break;
+        }
+        fHistDtBmon->Fill(dDt);
+        fHistDtEvoBmon->Fill(iEntry, dDt);
+        fHistDtDxBmon->Fill(dDt, 16 * ucScvdX[CbmTofAddress::GetChannelId((*itScvd).GetAddress())]
+                                   - CbmTofAddress::GetChannelId(digiOld.GetAddress()));
+        if (-20 < dDt && dDt < 20) {
+          fHistDxCorBmon->Fill(CbmTofAddress::GetChannelId(digiOld.GetAddress()),
+                               CbmTofAddress::GetChannelId((*itScvd).GetAddress()));
+          //                     ucScvdX[CbmTofAddress::GetChannelId((*itScvd).GetAddress())]);
+        }
+        fHistDtChanOld->Fill(CbmTofAddress::GetChannelId(digiOld.GetAddress()), dDt);
+        fHistDtChanScvd->Fill(CbmTofAddress::GetChannelId((*itScvd).GetAddress()), dDt);
+      }
+    }
+    auto prevDigiOld = vDigisOld.begin();
+    for (auto itOld = prevDigiOld; itOld != vDigisOld.end(); ++itOld) {
+      if (itOld != prevDigiOld) {
+        double_t dDt = (*itOld).GetTime() - (*prevDigiOld).GetTime();
+        fHistSelfDtOld->Fill(dDt);
+        if (dDt < dPileUpThrNs) {
+          fHistOldInPileup->Fill(1);
+        }
+        else {
+          fHistOldInPileup->Fill(0);
+        }
+      }
+      prevDigiOld = itOld;
+    }
+    auto prevDigiScvd = vDigisScvd.begin();
+    for (auto itScvd = prevDigiScvd; itScvd != vDigisScvd.end(); ++itScvd) {
+      if (itScvd != prevDigiScvd) {
+        double_t dDt = (*itScvd).GetTime() - (*prevDigiScvd).GetTime();
+        fHistSelfDtNew->Fill(dDt);
+        if (dDt < dPileUpThrNs) {
+          fHistNewInPileup->Fill(1);
+        }
+        else {
+          fHistNewInPileup->Fill(0);
+        }
+      }
+      prevDigiScvd = itScvd;
+    }
+  }
+
+  TCanvas* fCanvMap = new TCanvas("canvMap", "Channel counts mapping for old and sCVD BMON");
+  fCanvMap->Divide(2, 2);
+
+  fCanvMap->cd(1);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistMapBmonOld->Draw("hist");
+
+  fCanvMap->cd(2);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistMapBmonScvd->Draw("ColzText");
+
+  fCanvMap->cd(3);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistMapEvoBmonOld->Draw("colz");
+
+  fCanvMap->cd(4);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistMapEvoBmonScvd->Draw("colz");
+
+  TCanvas* fCanvCorr = new TCanvas("canvCorr", "Correlations (T, X) between old and sCVD BMON");
+  fCanvCorr->Divide(3, 2);
+
+  fCanvCorr->cd(1);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDtBmon->Draw("hist");
+
+  fCanvCorr->cd(2);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDtEvoBmon->Draw("colz");
+
+  fCanvCorr->cd(3);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDtChanOld->Draw("colz");
+
+  fCanvCorr->cd(4);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDtDxBmon->Draw("colz");
+
+  fCanvCorr->cd(5);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDxCorBmon->Draw("colz");
+
+  fCanvCorr->cd(6);
+  gPad->SetLogz();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistDtChanScvd->Draw("colz");
+
+  TCanvas* fCanvSelfCorr = new TCanvas("canvSelfCorr", "Self-Correlations in old and sCVD BMON");
+  fCanvSelfCorr->Divide(2, 2);
+
+  fCanvSelfCorr->cd(1);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistSelfDtOld->Draw("hist");
+
+  fCanvSelfCorr->cd(2);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistSelfDtNew->Draw("hist");
+
+  fCanvSelfCorr->cd(3);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistOldInPileup->Draw("hist");
+  gPad->Update();
+  dynamic_cast<TPaveStats*>(fHistOldInPileup->FindObject("stats"))->SetOptStat(110);
+
+  fCanvSelfCorr->cd(4);
+  gPad->SetLogy();
+  gPad->SetGridx();
+  gPad->SetGridy();
+  fHistNewInPileup->Draw("hist");
+  gPad->Update();
+  dynamic_cast<TPaveStats*>(fHistNewInPileup->FindObject("stats"))->SetOptStat(110);
+}