Digitization parameters as ConfigurableParam and unit test (#5680)

* DigitizationParameters derived from ConfigurableParams

* unit test for raw 2 digits conversion

* clang

* remove test_raw_conversion_ft0

Author: AllaMaevskaya

Detectors/CTF/test/test_ctf_io_ft0.cxx

@@ -14,6 +14,7 @@
 #include <boost/test/unit_test.hpp>
 #include "DetectorsCommonDataFormats/NameConf.h"
 #include "FT0Reconstruction/CTFCoder.h"
+#include "FT0Simulation/DigitizationParameters.h"
 #include "Framework/Logger.h"
 #include <TFile.h>
 #include <TRandom.h>
@@ -46,7 +47,7 @@ BOOST_AUTO_TEST_CASE(CTFTest)
       uint16_t q = gRandom->Integer(4096);
       uint8_t chain = gRandom->Rndm() > 0.5 ? 0 : 1;
       channels.emplace_back(ich, t, q, chain);
-      if (std::abs(t) < Geometry::mTime_trg_gate) {
+      if (std::abs(t) < DigitizationParameters::mTime_trg_gate) {
         if (ich < 4 * uint8_t(Geometry::NCellsA)) {
           trig.nChanA++;
           ampTotA += q;

Detectors/FIT/FT0/base/include/FT0Base/Geometry.h

@@ -37,15 +37,16 @@ class Geometry
   ///
   TVector3 centerMCP(int imcp) { return mMCP[imcp]; }
 
-  static constexpr int Nchannels = 208;       // number od channels
-  static constexpr int NCellsA = 24;          // number of radiatiors on A side
-  static constexpr int NCellsC = 28;          // number of radiatiors on C side
-  static constexpr float ZdetA = 335;         // number of radiatiors on A side
-  static constexpr float ZdetC = 82;          // number of radiatiors on C side
-  static constexpr float ChannelWidth = 13.02;       // channel width in ps
-  static constexpr float MV_2_Nchannels = 2.2857143; //amplitude channel 7 mV ->16channels
-  static constexpr float MV_2_NchannelsInverse = 0.437499997; //inverse amplitude channel 7 mV ->16channels
-  static constexpr int mTime_trg_gate = 192;                  // #channels
+  static constexpr int Nchannels = 208;                     // number od channels
+  static constexpr int NCellsA = 24;                        // number of radiatiors on A side
+  static constexpr int NCellsC = 28;                        // number of radiatiors on C side
+  static constexpr float ZdetA = 335;                       // number of radiatiors on A side
+  static constexpr float ZdetC = 82;                        // number of radiatiors on C side
+  static constexpr float ChannelWidth = 13.02;              // channel width in ps
+  static constexpr float ChannelWidthInverse = 0.076804916; // channel width in ps inverse
+                                                            // static constexpr float MV_2_Nchannels = 2.2857143;          //amplitude channel 7 mV ->16channels
+                                                            //  static constexpr float MV_2_NchannelsInverse = 0.437499997; //inverse amplitude channel 7 mV ->16channels
+                                                            //  static constexpr int mTime_trg_gate = 192;                  // #channels
 
  private:
   TVector3 mMCP[52];

Detectors/FIT/FT0/reconstruction/CMakeLists.txt

@@ -16,6 +16,7 @@ o2_add_library(FT0Reconstruction
 
   PUBLIC_LINK_LIBRARIES O2::SimulationDataFormat O2::Framework
   O2::FT0Base
+  O2::FT0Simulation
   O2::DataFormatsFT0
   O2::DetectorsRaw
   O2::CommonDataFormat

Detectors/FIT/FT0/reconstruction/include/FT0Reconstruction/CTFCoder.h

@@ -23,6 +23,7 @@
 #include "DataFormatsFT0/Digit.h"
 #include "DataFormatsFT0/ChannelData.h"
 #include "DetectorsCommonDataFormats/DetID.h"
+#include "FT0Simulation/DigitizationParameters.h"
 #include "DetectorsBase/CTFCoderBase.h"
 #include "rANS/rans.h"
 
@@ -60,6 +61,8 @@ class CTFCoder : public o2::ctf::CTFCoderBase
 
   void appendToTree(TTree& tree, CTF& ec);
   void readFromTree(TTree& tree, int entry, std::vector<Digit>& digitVec, std::vector<ChannelData>& channelVec);
+  // DigitizationParameters const &mParameters;
+  //  o2::ft0::Geometry mGeometry;
 
   ClassDefNV(CTFCoder, 1);
 };
@@ -164,7 +167,9 @@ void CTFCoder::decompress(const CompressedDigits& cd, VDIG& digitVec, VCHAN& cha
       const auto& chan = channelVec.emplace_back((chID += cd.idChan[icc]), cd.cfdTime[icc], cd.qtcAmpl[icc], cd.qtcChain[icc]);
       //
       // rebuild digit
-      if (std::abs(chan.CFDTime) < Geometry::mTime_trg_gate) {
+      int isInTrig = DigitizationParameters::mTime_trg_gate;
+      //      if (std::abs(chan.CFDTime) < DigitizationParameters::mTime_trg_gate) {
+      if (std::abs(chan.CFDTime) < isInTrig) {
         if (chan.ChId < 4 * uint8_t(Geometry::NCellsA)) { // A side
           amplA += chan.QTCAmpl;
           timeA += chan.CFDTime;

Detectors/FIT/FT0/reconstruction/include/FT0Reconstruction/CollisionTimeRecoTask.h

@@ -40,6 +40,8 @@ class CollisionTimeRecoTask
                              gsl::span<const o2::ft0::ChannelData> inChData,
                              gsl::span<o2::ft0::ChannelDataFloat> outChData);
   void FinishTask();
+  // DigitizationParameters const &mParameters;
+  //  o2::ft0::Geometry mGeometry;
 
  private:
   ClassDefNV(CollisionTimeRecoTask, 1);

Detectors/FIT/FT0/reconstruction/src/CollisionTimeRecoTask.cxx

@@ -15,6 +15,7 @@
 #include "FairLogger.h" // for LOG
 #include "DataFormatsFT0/RecPoints.h"
 #include "FT0Base/Geometry.h"
+#include "FT0Simulation/DigitizationParameters.h"
 #include <DataFormatsFT0/ChannelData.h>
 #include <DataFormatsFT0/Digit.h>
 #include <cmath>
@@ -48,7 +49,7 @@ o2::ft0::RecPoints CollisionTimeRecoTask::process(o2::ft0::Digit const& bcd,
   for (int ich = 0; ich < nch; ich++) {
     outChData[ich] = o2::ft0::ChannelDataFloat{inChData[ich].ChId,
                                                inChData[ich].CFDTime * Geometry::ChannelWidth,
-                                               (double)inChData[ich].QTCAmpl * Geometry::MV_2_NchannelsInverse,
+                                               (double)inChData[ich].QTCAmpl * DigitizationParameters::mMV_2_NchannelsInverse,
                                                inChData[ich].ChainQTC};
 
     //  only signals with amplitude participate in collision time

Detectors/FIT/FT0/simulation/include/FT0Simulation/DigitizationParameters.h

@@ -11,44 +11,39 @@
 #ifndef ALICEO2_FT0_DIGITIZATION_PARAMETERS
 #define ALICEO2_FT0_DIGITIZATION_PARAMETERS
 #include <FT0Base/Geometry.h>
+#include <CommonUtils/ConfigurableParamHelper.h>
 
 namespace o2::ft0
 {
 
-struct DigitizationParameters {
-  int nCellsA = Geometry::NCellsA;
-  int nCellsC = Geometry::NCellsC;
-  float zDetA = Geometry::ZdetA;
-  float zDetC = Geometry::ZdetC;
-  float bunchWidth = 25; //ns
-  float channelWidth = Geometry::ChannelWidth;
-  float ChannelWidthInverse = 0.076804916;                   // channel width in ps inverse
+struct DigitizationParameters
+  : o2::conf::ConfigurableParamHelper<DigitizationParameters> {
+  float mBunchWidth = 25;                                    //ns
+  float mChannelWidthInverse = 0.076804916;                  // channel width in ps inverse
   float mMCPs = (Geometry::NCellsA + Geometry::NCellsC) * 4; //number of MCPs
   float mCFD_trsh = 3.;                                      // [mV]
   float mAmp_trsh = 100;                                     // [ph.e]
   float mAmpRecordLow = -4;                                  // integrate charge from
-
   float mAmpRecordUp = 15;                                   // to [ns]
-  int mTime_trg_gate = 192;                                  // #channels
-  // int mTime_trg_gate = 2000;                                  //2000/13; (+-2ns)  #channels no limits
-  float mTimeDiffAC = (Geometry::ZdetA - Geometry::ZdetC) * TMath::C();
-  float C_side_cable_cmps = 2.86;   //ns
-  float A_side_cable_cmps = 11.110; //ns
-  int mtrg_central_trh = 600.;      // channels
-  int mtrg_semicentral_trh = 300.;  // channels
+  //  float mTimeDiffAC = (Geometry::ZdetA - Geometry::ZdetC) * TMath::C();
+  float mC_side_cable_cmps = 2.86;   //ns
+  float mA_side_cable_cmps = 11.110; //ns
+  int mtrg_central_trh = 600.;       // channels
+  int mtrg_semicentral_trh = 300.;   // channels
 
-  float mMip_in_V = 7;              //MIP to mV
-  float mPe_in_mip = 0.004;         // invserse Np.e. in MIP 1./250.
-  float mCfdShift = 1.66;           //ns
-  float mCFDShiftPos = 1.47;        //// shift positive part of CFD signal; distance between 0.3 of max amplitude  to max
-  float mCFDdeadTime = 15.6;        // ns
-  float AmpIntegrationTime = 19;    //ns
-  float IntegWindowDelayA = 6;      // ns, A side
-  float IntegWindowDelayC = -1.6;   // ns, C side
-  float charge2amp = 0.22;
+  float mMip_in_V = 7;       //MIP to mV
+  float mPe_in_mip = 0.004;  // invserse Np.e. in MIP 1./250.
+  float mCfdShift = 1.66;    //ns
+  float mCFDShiftPos = 1.47; //// shift positive part of CFD signal; distance between 0.3 of max amplitude  to max
+  float mCFDdeadTime = 15.6; // ns
+  float mCharge2amp = 0.22;
   float mNoiseVar = 0.1;            //noise level
   float mNoisePeriod = 1 / 0.9;     // GHz low frequency noise period;
-  float mV_2_Nchannels = 2.2857143; //7 mV ->16channels
+  static constexpr float mMV_2_Nchannels = 2.2857143;          //amplitude channel 7 mV ->16channels
+  static constexpr float mMV_2_NchannelsInverse = 0.437499997; //inverse amplitude channel 7 mV ->16channels
+  static constexpr int mTime_trg_gate = 192;                   // #channels
+
+  O2ParamDef(DigitizationParameters, "DigitizationParameters");
 };
 } // namespace o2::ft0
 #endif

Detectors/FIT/FT0/simulation/include/FT0Simulation/Digitizer.h

@@ -100,7 +100,7 @@ class Digitizer
     if (x <= 0.0f) {
       return 0.0f;
     }
-    float const y = x / mParameters.bunchWidth * DP::SIGNAL_TABLE_SIZE;
+    float const y = x / mParameters.mBunchWidth * DP::SIGNAL_TABLE_SIZE;
     int const index = std::floor(y);
     if (index + 1 >= DP::SIGNAL_TABLE_SIZE) {
       return mSignalTable.back();
@@ -111,7 +111,7 @@ class Digitizer
 
   inline Vc::float_v signalFormVc(Vc::float_v x) const
   { // table lookup for the signal shape (SIMD version)
-    auto const y = x / mParameters.bunchWidth * DP::SIGNAL_TABLE_SIZE;
+    auto const y = x / mParameters.mBunchWidth * DP::SIGNAL_TABLE_SIZE;
     Vc::float_v::IndexType const index = Vc::floor(y);
     auto const rem = y - index;
     Vc::float_v val(0);
@@ -141,7 +141,8 @@ class Digitizer
   std::deque<BCCache> mCache;
   std::array<GoodInteractionTimeRecord, 208> mDeadTimes;
 
-  DigitizationParameters mParameters;
+  DigitizationParameters const& mParameters;
+  o2::ft0::Geometry mGeometry;
 
   NoiseRandomRingType mRndGaus;
   int mNumNoiseSamples; // number of noise samples in one BC
@@ -156,7 +157,7 @@ class Digitizer
                std::vector<o2::ft0::DetTrigInput>& digitsTrig,
                o2::dataformats::MCTruthContainer<o2::ft0::MCLabel>& labels);
 
-  ClassDefNV(Digitizer, 1);
+  ClassDefNV(Digitizer, 2);
 };
 
 // signal shape function

Detectors/FIT/FT0/simulation/src/Digitizer.cxx

@@ -10,6 +10,7 @@
 
 #include "FT0Simulation/Digitizer.h"
 #include "FT0Simulation/DigitizationConstants.h"
+#include "FT0Simulation/DigitizationParameters.h"
 #include "SimulationDataFormat/MCTruthContainer.h"
 #include <CommonDataFormat/InteractionRecord.h>
 
@@ -22,7 +23,6 @@
 #include <optional>
 
 using namespace o2::ft0;
-//using o2::ft0::Geometry;
 
 ClassImp(Digitizer);
 
@@ -90,22 +90,22 @@ Digitizer::CFDOutput Digitizer::get_time(const std::vector<float>& times, float
   };
   auto const min_time = std::max(deadTime, *std::min_element(std::begin(times),
                                                              std::end(times)));
-  CFDOutput result{std::nullopt, -0.5f * mParameters.bunchWidth};
+  CFDOutput result{std::nullopt, -0.5f * mParameters.mBunchWidth};
   bool is_positive = true;
 
   // reset the chache
   std::fill_n(std::begin(mSignalCache), std::size(mSignalCache), -1.0f);
 
   // we need double precision for time in order to match previous behaviour
-  for (double time = min_time; time < 0.5 * mParameters.bunchWidth; time += DP::SIGNAL_CACHE_DT) {
+  for (double time = min_time; time < 0.5 * mParameters.mBunchWidth; time += DP::SIGNAL_CACHE_DT) {
     float const val = value_at(time);
-    int const index = std::lround((time + 0.5 * mParameters.bunchWidth) / DP::SIGNAL_CACHE_DT);
+    int const index = std::lround((time + 0.5 * mParameters.mBunchWidth) / DP::SIGNAL_CACHE_DT);
     if (index >= 0 && index < mSignalCache.size()) { // save the value for later use
       mSignalCache[index] = val;
     }
     // look up the time-shifted signal value from the past
     float val_prev = 0.0f;
-    int const index_prev = std::lround((time - mParameters.mCFDShiftPos + 0.5f * mParameters.bunchWidth) / DP::SIGNAL_CACHE_DT);
+    int const index_prev = std::lround((time - mParameters.mCFDShiftPos + 0.5f * mParameters.mBunchWidth) / DP::SIGNAL_CACHE_DT);
     val_prev = ((index_prev < 0 || index_prev >= mSignalCache.size() || mSignalCache[index_prev] < 0.0f)
                   ? value_at(time - mParameters.mCFDShiftPos) //  was not computed before
                   : mSignalCache[index_prev]);                //  is available in the cache
@@ -167,8 +167,8 @@ void Digitizer::process(const std::vector<o2::ft0::HitType>* hits,
     }
 
     Int_t hit_ch = hit.GetDetectorID();
-    Bool_t is_A_side = (hit_ch < 4 * mParameters.nCellsA);
-    Float_t time_compensate = is_A_side ? mParameters.A_side_cable_cmps : mParameters.C_side_cable_cmps;
+    Bool_t is_A_side = (hit_ch < 4 * mGeometry.NCellsA);
+    Float_t time_compensate = is_A_side ? mParameters.mA_side_cable_cmps : mParameters.mC_side_cable_cmps;
     Double_t hit_time = hit.GetTime() - time_compensate;
     if (hit_time > 250) {
       continue; //not collect very slow particles
@@ -224,20 +224,20 @@ void Digitizer::storeBC(BCCache& bc,
     if (!cfd.particle) {
       continue;
     }
-    int smeared_time = 1000. * (*cfd.particle - mParameters.mCfdShift) * mParameters.ChannelWidthInverse;
+    int smeared_time = 1000. * (*cfd.particle - mParameters.mCfdShift) * mParameters.mChannelWidthInverse;
     bool is_time_in_signal_gate = (smeared_time > -mParameters.mTime_trg_gate && smeared_time < mParameters.mTime_trg_gate);
-    float charge = measure_amplitude(channel_times) * mParameters.charge2amp;
-    float amp = is_time_in_signal_gate ? mParameters.mV_2_Nchannels * charge : 0;
+    float charge = measure_amplitude(channel_times) * mParameters.mCharge2amp;
+    float amp = is_time_in_signal_gate ? mParameters.mMV_2_Nchannels * charge : 0;
     if (amp > 4095) {
       amp = 4095;
     }
-    LOG(DEBUG) << mEventID << " bc " << firstBCinDeque.bc << " orbit " << firstBCinDeque.orbit << ", ipmt " << ipmt << ", smeared_time " << smeared_time << " nStored " << nStored;
+    LOG(INFO) << mEventID << " bc " << firstBCinDeque.bc << " orbit " << firstBCinDeque.orbit << ", ipmt " << ipmt << ", smeared_time " << smeared_time << " nStored " << nStored;
     digitsCh.emplace_back(ipmt, smeared_time, int(amp), chain);
     nStored++;
 
     // fill triggers
 
-    Bool_t is_A_side = (ipmt <= 4 * mParameters.nCellsA);
+    Bool_t is_A_side = (ipmt <= 4 * mGeometry.NCellsA);
     if (!is_time_in_signal_gate) {
       continue;
     }
@@ -322,7 +322,7 @@ void Digitizer::initParameters()
   auto const sinc = [](double x) { x *= TMath::Pi(); return (std::abs(x) < 1e-12) ? 1.0 : std::sin(x) / x; };
 
   // number of noise samples in one BC
-  mNumNoiseSamples = std::ceil(mParameters.bunchWidth / mParameters.mNoisePeriod);
+  mNumNoiseSamples = std::ceil(mParameters.mBunchWidth / mParameters.mNoisePeriod);
   mNoiseSamples.resize(mNumNoiseSamples);
 
   // set up tables with sinc function values (times noiseVar)
@@ -332,23 +332,24 @@ void Digitizer::initParameters()
     // we make a table of sinc values between -num_noise_samples and 2*num_noise_samples
     mSincTable[i].resize(3 * mNumNoiseSamples);
     for (int j = -mNumNoiseSamples; j < 2 * mNumNoiseSamples; ++j) {
-      mSincTable[i][mNumNoiseSamples + j] = mParameters.mNoiseVar * sinc((time + 0.5f * mParameters.bunchWidth) / mParameters.mNoisePeriod - j);
+      mSincTable[i][mNumNoiseSamples + j] = mParameters.mNoiseVar * sinc((time + 0.5f * mParameters.mBunchWidth) / mParameters.mNoisePeriod - j);
     }
   }
   // set up the lookup table for the signal form
   for (size_t i = 0; i < DP::SIGNAL_TABLE_SIZE; ++i) {
-    float const x = float(i) / float(DP::SIGNAL_TABLE_SIZE) * mParameters.bunchWidth;
+    float const x = float(i) / float(DP::SIGNAL_TABLE_SIZE) * mParameters.mBunchWidth;
     mSignalTable[i] = signalForm_i(x);
   }
 
   // cache for signal time series used by the CFD -BC/2 .. +3BC/2
-  mSignalCache.resize(std::lround(mParameters.bunchWidth / DP::SIGNAL_CACHE_DT));
+  mSignalCache.resize(std::lround(mParameters.mBunchWidth / DP::SIGNAL_CACHE_DT));
 }
 //_______________________________________________________________________
 void Digitizer::init()
 {
   LOG(INFO) << " @@@ Digitizer::init " << std::endl;
   mDeadTimes.fill({InteractionRecord(0), -100.});
+  printParameters();
 }
 //_______________________________________________________________________
 void Digitizer::finish()
@@ -364,5 +365,7 @@ void Digitizer::printParameters() const
             << mParameters.mCfdShift << "CFD distance between 0.3 of max amplitude  to max \n"
             << mParameters.mCFDShiftPos << "MIP -> mV " << mParameters.mMip_in_V << " Pe in MIP \n"
             << mParameters.mPe_in_mip << "noise level " << mParameters.mNoiseVar << " noise frequency \n"
-            << mParameters.mNoisePeriod << std::endl;
+            << mParameters.mNoisePeriod << " mMCPs " << mParameters.mMCPs;
 }
+
+O2ParamImpl(DigitizationParameters);

Steer/DigitizerWorkflow/src/FT0DigitizerSpec.cxx

@@ -44,7 +44,7 @@ class FT0DPLDigitizerTask : public o2::base::BaseDPLDigitizer
   using GRP = o2::parameters::GRPObject;
 
  public:
-  FT0DPLDigitizerTask() : o2::base::BaseDPLDigitizer(), mDigitizer(DigitizationParameters{}) {}
+  FT0DPLDigitizerTask() : o2::base::BaseDPLDigitizer(), mDigitizer(DigitizationParameters::Instance()) {}
   explicit FT0DPLDigitizerTask(o2::ft0::DigitizationParameters const& parameters)
     : o2::base::BaseDPLDigitizer(), mDigitizer(parameters){};
   ~FT0DPLDigitizerTask() override = default;