Beginning time of TOFPulse is now aligned to the beginning of clock c…

…ycle. Units are also added to TOFHitDigiConfig.

src/algorithms/digi/TOFHitDigiConfig.h

@@ -8,32 +8,29 @@
 namespace eicrecon {
 
 struct TOFHitDigiConfig {
-  // single hit energy deposition threshold
-  double threshold{1.0 * dd4hep::keV};
-  double tRes = 0.1; /// TODO 8 of what units??? Same TODO in juggler. Probably [ns]
-                     // digitization settings
-  double resolutionTDC{1};
-  double resolutionADC{1};
-
   // Parameters of AC-LGAD signal generation
   double gain           = 113.755;
-  double risetime       = 0.45; //in ns
-  double sigma_analog   = 0.293951;
-  double sigma_sharingx = 0.1;
-  double sigma_sharingy = 0.5;
+  double risetime       = 0.45 * dd4hep::ns; 
+  double sigma_analog   = 0.293951 * dd4hep::ns;
+  double sigma_sharingx = 0.1 * dd4hep::cm;
+  double sigma_sharingy = 0.5 * dd4hep::cm;
   double Vm             = -1e-4 * dd4hep::GeV; // Vm = voltage maximum. When EDep = 1e-4 GeV, voltage corresponds to ADC = adc_max
-  double t_thres        = 0.1 * Vm;
+  double t_thres        = 0.1 * Vm; // TDC value = time when pulse exceed t_thres
   double ignore_thres   = 0.01 * Vm; // If EDep below this value, digitization for the cell will be ignored. Speed up calculation
-                                     //
-  double tMin     = 0.1;
-  double tMax     = 25;// 25 ns is the period of 40MHz EIC clock
-  int total_time  = ceil(tMax - tMin);
+
   int adc_bit     = 8;
   int tdc_bit     = 10;
 
+  // total number of TDC/ADC values 
+  // Since digitizer starts at zero, max ADC value = adc_range - 1
+  // Similar for TDC
   int adc_range = pow(2, adc_bit);
   int tdc_range = pow(2, tdc_bit);
 
+  // period of the sensor clock. Time internal to sensor will all be digitized to integer multiple of tInterval
+  double tInterval = 25 * dd4hep::ns/(tdc_range - 1);  
+  double tMin      = 0.;
+  double tMax      = 25 * dd4hep::ns; // 25 ns is the period of 40MHz EIC clock
   std::string readout = "TOFBarrelHits";
 };
 

src/algorithms/digi/TOFPulseGeneration.cc

@@ -23,14 +23,20 @@ double TOFPulseGeneration::_Landau(double amp, double x, double mean, double std
   return amp*TMath::Landau(x, mean, std, kTRUE);
 }
 
+double TOFPulseGeneration::_DigitizeTime(double time) const {
+  // digitization always round down to tbe previous bin
+  return static_cast<int>(time/m_cfg.tInterval + 1e-3)*m_cfg.tInterval;
+}
+
 void TOFPulseGeneration::process(const TOFPulseGeneration::Input& input,
                                  const TOFPulseGeneration::Output& output) const {
   const auto [simhits] = input;
   auto [rawADCs] = output;
 
+  double interval = m_cfg.tInterval; 
   double Vm = m_cfg.Vm;
-  double tMin = m_cfg.tMin;
-  double tMax = m_cfg.tMax;
+  double tMin = _DigitizeTime(m_cfg.tMin);
+  double tMax = _DigitizeTime(m_cfg.tMax);
   int adc_range = m_cfg.adc_range;
   int tdc_range = m_cfg.tdc_range;
   int nBins = m_cfg.tdc_range;
@@ -39,14 +45,13 @@ void TOFPulseGeneration::process(const TOFPulseGeneration::Input& input,
   // NOTE: we take the cellID of the most energetic hit in this group so it is a real cellID from an
   // MC hit
   std::unordered_map<dd4hep::rec::CellID, std::vector<double>> adc_sum;
-  double interval = (tMax - tMin) / (nBins - 1);
 
   for (const auto& hit : *simhits) {
     auto cellID          = hit.getCellID();
     double sum_charge = 0.0;
     double mpv_analog = 0.0;
 
-    double  time       = hit.getTime();
+    double  time       = hit.getTime() * dd4hep::ns;
     double  charge     = hit.getEDep();
     // reduce computation power by not simulating low-charge hits
     if(charge < m_cfg.ignore_thres) continue;
@@ -63,11 +68,11 @@ void TOFPulseGeneration::process(const TOFPulseGeneration::Input& input,
     const double x_when_landau_min = -0.22278;
     double landau_min = this -> _Landau(-m_cfg.gain, x_when_landau_min, 0, 1)/m_cfg.sigma_analog;
     double scalingFactor = 1. / Vm / landau_min * adc_range;
-
-    for (int j = 0; j < nBins; ++j) {
-      double x = tMin + j * interval;
-      double y = charge * this -> _Landau(-m_cfg.gain, x, mpv_analog, m_cfg.sigma_analog) * scalingFactor;
-      ADCs[j] += y;;
+ 
+    {
+        int j;
+        for (double t = tMin, j = 0; t < tMax; ++j, t += interval) 
+          ADCs[j] += charge * this -> _Landau(-m_cfg.gain, t, mpv_analog, m_cfg.sigma_analog) * scalingFactor;
     }
 
   }
@@ -76,7 +81,7 @@ void TOFPulseGeneration::process(const TOFPulseGeneration::Input& input,
  for(const auto &[cellID, ADCs] : adc_sum) {
     auto time_series = rawADCs -> create();
     time_series.setCellID(cellID);
-    time_series.setTime(1.); // placeholder. Don't know what to assign when there are two or more hits
+    time_series.setTime(tMin); 
     time_series.setCharge(1.); // placeholder. Don't know what to assign when there are two or more hits
     time_series.setInterval(interval);
 

src/algorithms/digi/TOFPulseGeneration.h

@@ -36,6 +36,7 @@ class TOFPulseGeneration : public TOFPulseGenerationAlgorithm,
 
 protected:
   double _Landau(double amp, double x, double mean, double std) const;
+  double _DigitizeTime(double time) const;
 };
 
 } // namespace eicrecon

src/tests/algorithms_test/digi_TOFPulseDigitization.cc

@@ -62,6 +62,7 @@ TEST_CASE("the BTOF charge sharing algorithm runs", "[TOFPulseDigitization]") {
   algo.init();
 
   SECTION("TDC vs analytic solution scan") {
+    logger->info("Begin TDC vs analytic solution scan");
 
     // test pulse with gaussian shape
     for (double tdc_frac = 0.4; tdc_frac < 1; tdc_frac += 0.1) {
@@ -109,6 +110,7 @@ TEST_CASE("the BTOF charge sharing algorithm runs", "[TOFPulseDigitization]") {
   }
 
   SECTION("ADC scan") {
+    logger->info("Begin ADC scan");
 
     // test pulse with gaussian shape
     for (double adc_frac = 0.4; adc_frac < 1; adc_frac += 0.1) {