PAL local changes

macro/sts/plot_module_PN_digi_correlation.C

+
+// U0
+// 0x10000002 0x10010002
+// U1
+// 0x10000012 0x10010012
+// U2
+// 0x10000022 0x10010022
+// pulser
+// 0x99999998 0x99999999
+// 0x11E004E2 0x11E104E2
+
+void plot_module_PN_digi_correlation( uint32_t uModuleAddress = 0x10010022,
+                                      TString sFile = "../run/data/1588_node9_5_0000.digi.root",
+                                      uint32_t uCorrelationThrNs = 50) {
+
+  gROOT->cd();
+  TFile* file = new TFile( sFile, "READ");
+  TTree* tree = (TTree*) (file->Get("cbmsim"));
+
+  std::vector< CbmStsDigi > * vDigisSts = new std::vector< CbmStsDigi >();
+
+  TBranch *pBranchSts  = tree->GetBranch("CbmStsDigi");
+  tree->SetBranchAddress( "StsDigi", &vDigisSts );
+
+  Int_t nentries = 50;//tree->GetEntries();
+  cout << "Entries: " << nentries << endl;
+
+  TH2 *hAutoCoincP = new TH2I("hAutoCoincP","Correlations between P channels;Chan A;Chan B", 1024, 0, 1024, 1024, 0, 1024);
+  TH2 *hAutoCoincN = new TH2I("hAutoCoincN","Correlations between N channels;Chan A;Chan B", 1024, 1024, 2048, 1024, 1024, 2048);
+  TH2 *hCoincPN    = new TH2I("hCoincPN","Correlations between P and N channels;Chan A;Chan B", 1024, 0, 1024, 1024, 1024, 2048);
+
+  ULong64_t prevTime_ch[2048];
+  Int_t prevADC_ch[2048];
+  for (int j=0;j<2048;j++){
+    prevTime_ch[j]=0.0;
+    prevADC_ch[j]=0;
+  }
+
+  for (Int_t iEntry = 0; iEntry < 1; iEntry++) {
+
+    tree->GetEntry(iEntry);
+    UInt_t nDigis  = vDigisSts->size();
+
+    if (iEntry % 10 == 0 ) std::cout << "STS Digis: " << nDigis  << std::endl;
+
+    // STS
+    Int_t iFirstGoodDigi = 0;
+    for (Int_t iDigi = 0; iDigi < nDigis; iDigi++) {
+      CbmStsDigi & pDigi = vDigisSts->at( iDigi );
+      if (iDigi==0) cout << iDigi << " / " <<nDigis << " " << (ULong64_t)pDigi.GetTime() << endl;
+
+      if (pDigi.GetAddress() == uModuleAddress) {
+
+	      Int_t ch = pDigi.GetChannel();
+        ULong64_t timeA = pDigi.GetTime();
+
+        Bool_t bFirst = kFALSE;
+        for (Int_t iDigiB = iFirstGoodDigi; iDigiB < nDigis; iDigiB++) {
+          CbmStsDigi & pDigiB = vDigisSts->at( iDigiB );
+          if (pDigiB.GetAddress() == uModuleAddress) {
+
+            Int_t chB = pDigiB.GetChannel();
+            ULong64_t timeB = pDigiB.GetTime();
+	          ULong64_t diffTime = timeB < timeA ? timeA - timeB : timeB - timeA;
+
+            if (!bFirst) cout << iDigiB << " / " <<nDigis << " " << diffTime << " " << ch << " " << chB << endl;
+
+            if( diffTime < uCorrelationThrNs ) {
+              if( !bFirst ) {
+                bFirst = kTRUE;
+                iFirstGoodDigi = iDigiB;
+              }
+              if( ch < 1024 && chB < 1024 )
+                hAutoCoincP->Fill( ch, chB );
+              else if( 1024 < ch && 1024 < chB )
+                hAutoCoincN->Fill( ch, chB );
+              else if( ch < chB )
+                hCoincPN->Fill( ch, chB );
+              else
+                hCoincPN->Fill( chB, ch );
+            }
+            else if( timeA < timeB ) break;
+          }
+        }
+      }  // feb condition
+    }
+  }
+
+  TCanvas* c = new TCanvas();
+  c->Divide(3);
+
+  c->cd(1);
+  gPad->SetGridx();
+  gPad->SetGridy();
+  gPad->SetLogz();
+  hAutoCoincP->Draw("colz");
+  hAutoCoincP->SetStats(kFALSE);
+
+  c->cd(2);
+  gPad->SetGridx();
+  gPad->SetGridy();
+  gPad->SetLogz();
+  hAutoCoincN->Draw("colz");
+  hAutoCoincN->SetStats(kFALSE);
+
+  c->cd(3);
+  gPad->SetGridx();
+  gPad->SetGridy();
+  gPad->SetLogz();
+  hCoincPN->Draw("colz");
+  hCoincPN->SetStats(kFALSE);
+}

macro/sts/sts_noise_ch_finder.C

@@ -5,13 +5,25 @@ void sts_noise_ch_finder(std::string sInputFile, double_t dShareThrMean = 1, dou
   TFile* pFile = new TFile(sInputFile.data(), "READ");
   gROOT->cd();
 
+  TH1* tempH1                                = nullptr;
   TProfile* tempP1                           = nullptr;
   TProfile2D* tempP2                         = nullptr;
+  TH1* hStsDigisTimeInRun                    = new TH1I();
   std::vector<TProfile*> vFebRawChRatio      = {};
   std::vector<TProfile*> vFebRawChRatioInSys = {};
   std::vector<TProfile*> vFebHitChRatioInSys = {};
   std::vector<TProfile2D*> vFebRawChRatioEvo = {};
+  std::vector<TProfile*> vFebRawChRatioInSysSpillOn = {};
+  std::vector<TProfile*> vFebRawChRatioInSysSpillOff = {};
   for (uint32_t uFeb = 0;; ++uFeb) {
+
+    tempH1 = (TProfile*) (pFile->FindObjectAny("hStsDigisTimeInRun"));
+    if (nullptr != tempH1) {
+      tempH1->Copy(*hStsDigisTimeInRun);
+    }  //  if( NULL != tempH1 )
+    else
+      break;
+
     // If needed to speed up search, folder is "perFeb"
     tempP1 = (TProfile*) (pFile->FindObjectAny(Form("hRawChRatio_%03u", uFeb)));
     if (nullptr != tempP1) {
@@ -46,12 +58,147 @@ void sts_noise_ch_finder(std::string sInputFile, double_t dShareThrMean = 1, dou
     }  //  if( NULL != tempH1 )
     else
       break;
+
+    vFebRawChRatioInSysSpillOn.push_back(new TProfile(Form("FebRawChRatioInSysSpillOn_%03u", uFeb),
+            "In Spill;Channel []; Share of FEB raw msg [Prct]", 1024, -0.5, 1023.5));
+    vFebRawChRatioInSysSpillOff.push_back(new TProfile(Form("FebRawChRatioInSysSpillOff_%03u", uFeb),
+            "Out of Spill;Channel []; Share of FEB raw msg [Prct]", 1024, -0.5, 1023.5));
   }  // for (uint32_t uFeb = 0; ; ++uFeb)
 
+  double_t dMinRate = hStsDigisTimeInRun->GetMinimum();
+  double_t dMaxRate = hStsDigisTimeInRun->GetMaximum();
+  if( dMinRate < 1 ) dMinRate = dMaxRate / 100;
+  uint32_t uNbRateBins = dMaxRate/dMinRate;
+  std::cout << dMinRate << " " << dMaxRate << " " << uNbRateBins << std::endl;
+  TH1* hRateDistrib = new TH1I("hRateDistrib", "hRateDistrib;;", uNbRateBins, 0, dMaxRate );
+  for( uint32_t uBin = 0; uBin < hStsDigisTimeInRun->GetNbinsX(); ++uBin) {
+    double_t dRate = hStsDigisTimeInRun->GetBinContent( 1 + uBin );
+    if( 0 < dRate )
+      hRateDistrib->Fill( dRate );
+  }
+  TF1* gg_rate = new TF1("gg_rate", "gaus(0)+gaus(3)+pol0(6)", 0, dMaxRate);
+  gg_rate->SetParameters( hRateDistrib->GetMaximum(),
+                          hRateDistrib->GetBinCenter(hRateDistrib->GetMaximumBin()),
+                          1 * dMaxRate / uNbRateBins,
+                          hRateDistrib->GetMaximum()/2,
+                          hRateDistrib->GetBinCenter(hRateDistrib->GetMaximumBin()) + 10 * dMaxRate / uNbRateBins,
+                          3 * dMaxRate / uNbRateBins,
+                          1
+                          );
+  gg_rate->SetParLimits(3, 1, hRateDistrib->GetMaximum());
+  gg_rate->SetParLimits(4, hRateDistrib->GetBinCenter(hRateDistrib->GetMaximumBin()), dMaxRate);
+  gg_rate->SetParLimits(5, 1 * dMaxRate / uNbRateBins, dMaxRate);
+  hRateDistrib->Fit("gg_rate", "R");
+  TF1* fitresult = hRateDistrib->GetFunction("gg_rate");
+  std::cout << " parameters from fit, " << std::endl
+            << "Gauss A = " << fitresult->GetParameter(0) << ", " << fitresult->GetParameter(1) << ", "
+            << fitresult->GetParameter(2) << std::endl
+            << "Gauss B = " << fitresult->GetParameter(3) << ", " << fitresult->GetParameter(4) << ", "
+            << fitresult->GetParameter(5) << std::endl;
+  hRateDistrib->Draw();
+
+  double dSpillThreshold = fitresult->GetParameter(1) + 3 * fitresult->GetParameter(2);
+  bool bSpill = false;
+  std::vector<uint32_t> vuSpillStart = {};
+  std::vector<uint32_t> vuSpillStop = {};
+  uint32_t uFirstEmpty = hStsDigisTimeInRun->GetNbinsX();
+  double_t dRate = hStsDigisTimeInRun->GetBinContent( 1 );
+  if ( dRate < dSpillThreshold ) {
+    std::cout << "Starting with Spill OFF" << std::endl;
+    bSpill = false;
+    vuSpillStop.push_back(0);
+  }
+  else {
+    std::cout << "Starting with Spill ON" << std::endl;
+    bSpill = true;
+    vuSpillStart.push_back(0);
+  }
+  for ( uint32_t uBin = 1; uBin < hStsDigisTimeInRun->GetNbinsX(); ++uBin) {
+    dRate = hStsDigisTimeInRun->GetBinContent( 1 + uBin );
+    if ( 0 < dRate ) {
+      if ( bSpill && dRate < dSpillThreshold ) {
+        bSpill = false;
+        vuSpillStop.push_back(uBin);
+        std::cout << "Stop at " << uBin << std::endl;
+      }
+      else if ( !bSpill && dSpillThreshold <= dRate ) {
+        bSpill = true;
+        vuSpillStart.push_back(uBin);
+        std::cout << "Start at " << uBin << std::endl;
+      }
+    }
+    else {
+      std::cout << "Empty at " << uBin << std::endl;
+      uFirstEmpty = uBin;
+      break;
+    }
+  }
+
   uint32_t uNbBadChannelsCandidates = 0;
   double_t dShareRawBadChCandInSts  = 0.0;
   double_t dShareHitBadChCandInSts  = 0.0;
   for (uint32_t uFeb = 0; uFeb < vFebRawChRatio.size(); ++uFeb) {
+    if ( 0 < vuSpillStart.size() && 0 < vuSpillStop.size() ) {
+      uint32_t uOn = 0;
+      uint32_t uOff = 0;
+      auto itStart = vuSpillStart.begin();
+      auto itStop  = vuSpillStop.begin();
+      if ( *itStart < *itStop ) {
+        while ( itStop != vuSpillStop.end() ) {
+          TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_on_%03d", uOn), *itStart, *itStop - 1);
+          vFebRawChRatioInSysSpillOn[uFeb]->Add(temp);
+          delete temp;
+          ++uOn;
+          ++itStart;
+          if( itStart != vuSpillStart.end() ) {
+            TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_off_%03d", uOff), *itStop, *itStart - 1);
+            vFebRawChRatioInSysSpillOff[uFeb]->Add(temp);
+            delete temp;
+          }
+          else {
+            TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_off_%03d", uOff), *itStop, uFirstEmpty - 1);
+            vFebRawChRatioInSysSpillOff[uFeb]->Add(temp);
+            delete temp;
+          }
+          ++uOff;
+          ++itStop;
+        }
+        if( itStart != vuSpillStart.end() ) {
+          TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_on_%03d", uOn), *itStart, uFirstEmpty - 1);
+          vFebRawChRatioInSysSpillOn[uFeb]->Add(temp);
+          delete temp;
+          ++uOn;
+        }
+      }
+      else {
+        while ( itStart != vuSpillStart.end() ) {
+          TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_off_%03d", uOff), *itStop, *itStart - 1);
+          vFebRawChRatioInSysSpillOff[uFeb]->Add(temp);
+          delete temp;
+          ++uOff;
+          ++itStop;
+          if( itStop != vuSpillStop.end() ) {
+            TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_on_%03d", uOn), *itStart, *itStop - 1);
+            vFebRawChRatioInSysSpillOn[uFeb]->Add(temp);
+            delete temp;
+          }
+          else {
+            TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_on_%03d", uOn), *itStart, uFirstEmpty - 1);
+            vFebRawChRatioInSysSpillOn[uFeb]->Add(temp);
+            delete temp;
+          }
+          ++uOn;
+          ++itStart;
+        }
+        if( itStop != vuSpillStop.end() ) {
+          TProfile * temp = vFebRawChRatioEvo[uFeb]->ProfileY( Form("proj_off_%03d", uOff), *itStop, uFirstEmpty - 1);
+          vFebRawChRatioInSysSpillOff[uFeb]->Add(temp);
+          delete temp;
+          ++uOff;
+        }
+      }
+    }
+
     for (uint32_t uCh = 0; uCh < 1024; ++uCh) {
       if (dShareThrMean <= vFebRawChRatioInSys[uFeb]->GetBinContent(uCh + 1)
           || dFebShareThrMean <= vFebRawChRatio[uFeb]->GetBinContent(uCh + 1)) {
@@ -59,6 +206,9 @@ void sts_noise_ch_finder(std::string sInputFile, double_t dShareThrMean = 1, dou
                           uFeb, uCh, vFebRawChRatio[uFeb]->GetBinContent(uCh + 1),
                           vFebRawChRatioInSys[uFeb]->GetBinContent(uCh + 1),
                           vFebHitChRatioInSys[uFeb]->GetBinContent(uCh + 1))
+                  << Form(" In-Spill %6.2f %% in FEB, Out of Spill %6.2f",
+                          vFebRawChRatioInSysSpillOn[uFeb]->GetBinContent(uCh + 1),
+                          vFebRawChRatioInSysSpillOff[uFeb]->GetBinContent(uCh + 1))
                   << std::endl;
         uNbBadChannelsCandidates++;
         dShareRawBadChCandInSts += vFebRawChRatioInSys[uFeb]->GetBinContent(uCh + 1);

reco/detectors/sts/unpack/CbmStsUnpackAlgo.cxx

@@ -512,7 +512,10 @@ void CbmStsUnpackAlgo::processHitInfo(const stsxyter::Message& mess)
                            fviFebSide[uFebIdx]);
       }
 
-      if (fMonitor) fMonitor->CountDigi(uFebIdx, uChanInFeb);
+      if (fMonitor) {
+        fMonitor->CountDigi(uFebIdx, uChanInFeb);
+        fMonitor->FillDigisTimeInRun( ulTimeInNs + fTsStartTime );
+      }
 
       /// Catch the pulser digis and either save them to their own output or drop them
       if (fvbFebPulser[uFebIdx]) {
@@ -779,9 +782,6 @@ bool CbmStsUnpackAlgo::unpack(const fles::Timeslice* ts, std::uint16_t icomp, UI
   //Output debugging info
   if (fMonitor) {
     if (fMonitor->GetDebugMode()) { fMonitor->PrintDebugInfo(fMsStartTime, fNrProcessedTs, msDescriptor.flags, uSize); }
-    for (auto itHit = fOutputVec.begin(); itHit != fOutputVec.end(); ++itHit) {
-      fMonitor->FillDigisTimeInRun(itHit->GetTime());
-    }
     fMonitor->FillVectorSize(ts->index(), fOutputVec.size());
     //fMonitor->DrawCanvases();
   }

reco/detectors/sts/unpack/CbmStsUnpackMonitor.cxx

@@ -69,8 +69,9 @@ Bool_t CbmStsUnpackMonitor::CreateHistograms(CbmMcbm2018StsPar* pUnpackPar)
   }
 
   /// Create general unpacking histograms
-  fhDigisTimeInRun = new TH1I("hStsDigisTimeInRun", "Digis Nb vs Time in Run; Time in run [s]; Digis Nb []", 10, 0, 1);
-  fhDigisTimeInRun->SetCanExtend(TH1::kAllAxes);
+  fhDigisTimeInRun = new TH1I("hStsDigisTimeInRun", "Digis Nb vs Time in Run; Time in run [s]; Digis Nb []",
+                              600, -0.5, 599.5);
+//  fhDigisTimeInRun->SetCanExtend(TH1::kAllAxes);
   AddHistoToVector(fhDigisTimeInRun, "");
 
   fhVectorSize = new TH1I("fhVectorSize", "Size of the vector VS TS index; TS index; Size [bytes]", 10, 0, 10);
@@ -913,6 +914,8 @@ void CbmStsUnpackMonitor::Finish()
 {
   DrawCanvases();
 
+  LOG(info) << "Digis/Raws ratio: " << ((100.0 * fuNbDigisTsSts) / fuNbRawTsSts );
+
   /// Obtain vector of pointers on each histo (+ optionally desired folder)
   std::vector<std::pair<TNamed*, std::string>> vHistos = GetHistoVector();
 

reco/detectors/sts/unpack/CbmStsUnpackMonitor.h

@@ -6,6 +6,7 @@
 #define CbmStsUnpackMonitor_H
 
 #include "CbmMcbm2018StsPar.h"
+#include <Logger.h>
 
 #include "Rtypes.h"
 #include "TH1.h"
@@ -79,7 +80,11 @@ public:
   void FillVectorCapacity(ULong64_t TsIdx, UInt_t Capacity) { fhVectorCapacity->Fill(TsIdx, Capacity); }
   void FillMsCntEvo(ULong64_t MsIdx) { fhMsCntEvo->Fill(MsIdx * 1e-9); }
   void FillMsErrorsEvo(ULong64_t MsIdx, uint16_t MsErrorType) { fhMsErrorsEvo->Fill(1e-9 * MsIdx, MsErrorType); }
-  void FillDigisTimeInRun(Double_t Time) { fhDigisTimeInRun->Fill(Time * 1e-9); }
+  void FillDigisTimeInRun(Double_t Time) {
+    if( 0 < dFirstTsStartTime) {
+      fhDigisTimeInRun->Fill(Time * 1e-9 - dFirstTsStartTime);
+    }
+  }
   void FillStsAllFebsHitRateEvo(Double_t dTimeSinceStartSec, UInt_t uFebIdx)
   {
     fhStsAllFebsHitRateEvo->Fill(dTimeSinceStartSec, uFebIdx);