ForwardDetectors/ZDC/ZdcAnalysis/Root/ZDCDataAnalyzer.cxx

 /*
-  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
 */
 
 #include <ZdcAnalysis/ZDCDataAnalyzer.h>
@@ -8,7 +8,6 @@
 #include <sstream>
 #include <utility>
 
-#include "CxxUtils/trapping_fp.h"
 #include "CxxUtils/checker_macros.h"
 ATLAS_NO_CHECK_FILE_THREAD_SAFETY;
 
@@ -26,12 +25,6 @@ ZDCDataAnalyzer::ZDCDataAnalyzer(ZDCMsg::MessageFunctionPtr msgFunc_p, int nSamp
   m_haveECalib(false),
   m_haveT0Calib(false),
   m_currentLB(-1),
-
-  // Default "calibrations"
-  m_currentECalibCoeff ({{{{1, 1, 1, 1}}, {{1, 1, 1, 1}}}}),
-  m_currentT0OffsetsHG ({{{{0, 0, 0, 0}}, {{0, 0, 0, 0}}}}),
-  m_currentT0OffsetsLG ({{{{0, 0, 0, 0}}, {{0, 0, 0, 0}}}}),
-
   m_moduleMask(0),
   m_moduleSum({{0, 0}}),
   m_moduleSumErrSq({{0, 0}}),
@@ -67,6 +60,12 @@ ZDCDataAnalyzer::ZDCDataAnalyzer(ZDCMsg::MessageFunctionPtr msgFunc_p, int nSamp
   m_pedestals[0] = {{100, 100, 100, 100}};
   m_pedestals[1] = {{100, 100, 100, 100}};
 
+  // Default "calibrations"
+  //
+  m_currentECalibCoeff = {{{{1, 1, 1, 1}}, {{1, 1, 1, 1}}}};
+
+  m_currentT0OffsetsHG = {{{{0, 0, 0, 0}}, {{0, 0, 0, 0}}}};
+  m_currentT0OffsetsLG = {{{{0, 0, 0, 0}}, {{0, 0, 0, 0}}}};
 
   // Construct the per-module pulse analyzers
   //
@@ -75,7 +74,7 @@ ZDCDataAnalyzer::ZDCDataAnalyzer(ZDCMsg::MessageFunctionPtr msgFunc_p, int nSamp
       std::ostringstream moduleTag;
       moduleTag << "_s" << side << "_m" << module;
 
-      m_moduleAnalyzers[side][module].reset (new ZDCPulseAnalyzer(m_msgFunc_p, moduleTag.str(), m_nSample, m_deltaTSample, m_preSampleIdx,
+      m_moduleAnalyzers[side][module].reset (new ZDCPulseAnalyzer(m_msgFunc_p, moduleTag.str().c_str(), m_nSample, m_deltaTSample, m_preSampleIdx,
                                              m_pedestals[side][module], m_HGGains[side][module], m_fitFunction,
                                              peak2ndDerivMinSamples[side][module],
                                              peak2ndDerivMinThresholdsHG[side][module],
@@ -328,9 +327,29 @@ void ZDCDataAnalyzer::SetNonlinCorrParams(float refADC, float refScale,
   }
 }
 
+void ZDCDataAnalyzer::enableFADCCorrections(bool correctPerSample,
+					    std::array<std::array<std::unique_ptr<const TH1>, 4>, 2>& corrHistHG,
+					    std::array<std::array<std::unique_ptr<const TH1>, 4>, 2>& corrHistLG)
+{
+  for (size_t side : {0, 1}) {
+    for (size_t module : {0, 1, 2, 3}) {
+      m_moduleAnalyzers[side][module]->enableFADCCorrections(correctPerSample, corrHistHG[side][module], corrHistLG[side][module]);
+    }
+  }
+}
+
+void ZDCDataAnalyzer::disableFADCCorrections()
+{
+  for (size_t side : {0, 1}) {
+    for (size_t module : {0, 1, 2, 3}) {
+      m_moduleAnalyzers[side][module]->disableFADCCorrections();
+    }
+  }
+}
+
 void ZDCDataAnalyzer::enableTimeSigCut(bool AND, float sigCut, std::string TF1String,
-			const std::array<std::array<std::vector<double>, 4>, 2>& parsHGArr, 
-			const std::array<std::array<std::vector<double>, 4>, 2>& parsLGArr)
+				       const std::array<std::array<std::vector<double>, 4>, 2>& parsHGArr, 
+				       const std::array<std::array<std::vector<double>, 4>, 2>& parsLGArr)
 {
   for (size_t side : {0, 1}) {
     for (size_t module : {0, 1, 2, 3}) {
@@ -569,12 +588,9 @@ bool ZDCDataAnalyzer::FinishEvent()
       // subtract the fraction of LGOverflow events if we have fraction available (<0 means unavailable)
       if (pulseAna_p->LGOverflow() && m_moduleAmpFractionLG[side][module] > 0) {tempFraction -= m_moduleAmpFractionLG[side][module];}
     }
-
-    {
-      CXXUTILS_TRAPPING_FP;
-      if (m_moduleSum[side] > 0) m_moduleSumBkgdFrac[side] = sumAmpTimesBkgdFrac/m_moduleSum[side];
-      else m_moduleSumBkgdFrac[side] = 0;
-    }
+    
+    if (m_moduleSum[side] > 0) m_moduleSumBkgdFrac[side] = sumAmpTimesBkgdFrac/m_moduleSum[side];
+    else m_moduleSumBkgdFrac[side] = 0;
     
     if (m_calibModuleSum[side] > 1e-6) {
       m_averageTime[side] /= m_calibModuleSum[side];

ForwardDetectors/ZDC/ZdcAnalysis/Root/ZDCPulseAnalyzer.cxx

 /*
-  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
 */
 
 #include "ZdcAnalysis/ZDCPulseAnalyzer.h"
@@ -7,16 +7,18 @@
 #include "TFitResult.h"
 #include "TFitResultPtr.h"
 #include "TVirtualFitter.h"
+#include "TFitter.h"
 #include "TList.h"
 #include "TMinuit.h"
 
 #include <algorithm>
 #include <sstream>
 #include <cmath>
+#include <memory>
 #include <numeric>
 #include <iomanip>
-#include <stdexcept>
 
+#include "CxxUtils/checker_macros.h"
 ATLAS_NO_CHECK_FILE_THREAD_SAFETY;
 
 extern int gErrorIgnoreLevel;
@@ -87,7 +89,7 @@ void ZDCPulseAnalyzer::CombinedPulsesFCN(int& /*numParam*/, double*, double& f,
 ZDCPulseAnalyzer::ZDCPulseAnalyzer(ZDCMsg::MessageFunctionPtr msgFunc_p, const std::string& tag, int Nsample, float deltaTSample, size_t preSampleIdx, int pedestal,
                                    float gainHG, const std::string& fitFunction, int peak2ndDerivMinSample,
                                    float peak2ndDerivMinThreshHG, float peak2ndDerivMinThreshLG) :
-  m_msgFunc_p(std::move(msgFunc_p)),
+  m_msgFunc_p(msgFunc_p),
   m_tag(tag), m_Nsample(Nsample),
   m_preSampleIdx(preSampleIdx),
   m_deltaTSample(deltaTSample),
@@ -208,6 +210,7 @@ void ZDCPulseAnalyzer::SetDefaults()
   m_useDelayed = false;
   m_enablePreExcl = false;
   m_enablePostExcl = false;
+  m_haveUserFilter = false;
   
   m_timingCorrMode = NoTimingCorr;
   m_haveNonlinCorr = false;
@@ -436,6 +439,31 @@ void ZDCPulseAnalyzer::enableTimeSigCut(bool AND, float sigCut, std::string TF1S
 
 }
 
+void ZDCPulseAnalyzer::enableFADCCorrections(bool correctPerSample, std::unique_ptr<const TH1>& correHistHG, std::unique_ptr<const TH1>& correHistLG)
+{
+  m_haveFADCCorrections = true;
+  m_FADCCorrPerSample = correctPerSample;
+  
+  // check for appropriate limits
+  //
+  if (std::abs(correHistHG->GetXaxis()->GetXmin()) > 1e-3 ||
+      std::abs(correHistHG->GetXaxis()->GetXmax() - 4096) > 1e-3) {
+    (*m_msgFunc_p)(ZDCMsg::Error, ("ZDCPulseAnalyzer::enableFADCCorrections:: invalid high gain correction histogram range: xmin, xmax = " +
+				   std::to_string(correHistHG->GetXaxis()->GetXmin()) + ", " + std::to_string(correHistHG->GetXaxis()->GetXmax())) );
+  }
+  else {
+    m_FADCCorrHG = std::move(correHistHG);
+  }
+  
+  if (std::abs(correHistLG->GetXaxis()->GetXmin()) > 1e-3 ||
+      std::abs(correHistLG->GetXaxis()->GetXmin() - 4096) > 1e-3) {
+    (*m_msgFunc_p)(ZDCMsg::Error, ("ZDCPulseAnalyzer::enableFADCCorrections:: invalid low gain correction histogram range: xmin, xmax = " +
+				   std::to_string(correHistLG->GetXaxis()->GetXmin()) + ", " + std::to_string(correHistLG->GetXaxis()->GetXmax())) );
+  }
+  else {
+    m_FADCCorrLG = std::move(correHistLG);
+  }
+}
 
 std::vector<float>  ZDCPulseAnalyzer::GetFitPulls(bool refitLG) const
 {
@@ -474,6 +502,46 @@ std::vector<float>  ZDCPulseAnalyzer::GetFitPulls(bool refitLG) const
   }
 }
 
+double ZDCPulseAnalyzer::getAmplitudeCorrection(bool highGain)
+{
+  double amplCorrFactor  = 1;
+  
+  // If we have FADC correction and we aren't applying it per-sample, do so here 
+  //
+  if (m_haveFADCCorrections && !m_FADCCorrPerSample) {
+    double fadcCorr = highGain ? m_FADCCorrHG->Interpolate(m_fitAmplitude) : m_FADCCorrLG->Interpolate(m_fitAmplitude);
+    amplCorrFactor *= fadcCorr;
+  }
+
+  double amplCorr = m_fitAmplitude * amplCorrFactor;
+  
+  // If we have a non-linear correction, apply it here   
+  //   We apply it as an inverse correction - i.e. we divide by a correction
+  //     term tha is a sum of coefficients times the ADC minus a reference
+  //     to a power. The lowest power is 1, the highest is deteremined by
+  //     the number of provided coefficients 
+
+  if (m_haveNonlinCorr) {
+    float invNLCorr = 1.0;
+    float nlPolyArg = (amplCorr - m_nonLinCorrRefADC) / m_nonLinCorrRefScale;
+	
+    if (m_useLowGain) {
+      for (size_t power = 1; power <= m_nonLinCorrParamsHG.size(); power++) {
+	invNLCorr += m_nonLinCorrParamsHG[power - 1]*pow(nlPolyArg, power);
+      }
+    }
+    else {
+      for (size_t power = 1; power <= m_nonLinCorrParamsLG.size(); power++) {
+	invNLCorr += m_nonLinCorrParamsLG[power - 1]*pow(nlPolyArg, power);
+      }
+    }
+
+    amplCorrFactor /= invNLCorr;  
+  }
+  
+  return amplCorrFactor;
+}
+
 void ZDCPulseAnalyzer::SetupFitFunctions()
 {
   float prePulseTMin = 0;
@@ -850,34 +918,21 @@ bool ZDCPulseAnalyzer::DoAnalysis(bool repass)
       //
       // --BAC
       //
-      
-       // If we have a non-linear correction, apply it here   
-       //   We apply it as an inverse correction - i.e. we divide by a correction
-       //     term tha is a sum of coefficients times the ADC minus a reference
-       //     to a power. The lowest power is 1, the highest is deteremined by
-       //     the number of provided coefficients 
-
-      double invNLCorr = 1.0;
-      if (m_haveNonlinCorr) {
-        float ampCorrFact = (m_fitAmplitude - m_nonLinCorrRefADC) / m_nonLinCorrRefScale;
-	
-	for (size_t power = 1; power <= m_nonLinCorrParamsLG.size(); power++) {
-	  invNLCorr += m_nonLinCorrParamsLG[power - 1]*pow(ampCorrFact, power);
-	}
-      }
+
+      double amplCorrFactor = getAmplitudeCorrection(false);
       
       //
       // Multiply amplitude by gain factor
       //
       m_ampNoNonLin   = m_fitAmplitude * m_gainFactorLG;
-      m_amplitude     = m_fitAmplitude / invNLCorr * m_gainFactorLG;
-      m_ampError      = m_fitAmpError / invNLCorr * m_gainFactorLG;
+      m_amplitude     = m_fitAmplitude * amplCorrFactor * m_gainFactorLG;
+      m_ampError      = m_fitAmpError * amplCorrFactor * m_gainFactorLG;
       m_preSampleAmp  = m_preSample    * m_gainFactorLG;
       m_preAmplitude  = m_fitPreAmp    * m_gainFactorLG;
       m_postAmplitude = m_fitPostAmp   * m_gainFactorLG;
       m_expAmplitude  = m_fitExpAmp    * m_gainFactorLG;
 
-      // BAC: also scale up the 2nd derivative so low and high gain can be treated on the same footing
+      // BAC: also scale up the 2nd derivative by the gain factor so low and high gain can be treated on the same footing
       //
       m_minDeriv2nd *= m_gainFactorLG;
     }
@@ -916,17 +971,10 @@ bool ZDCPulseAnalyzer::DoAnalysis(bool repass)
       //      to a power. The lowest power is 1, the highest is deteremined by
       //      the number of provided coefficients 
       //
-      if (m_haveNonlinCorr) {
-        float ampCorrFact = (m_fitAmplitude - m_nonLinCorrRefADC) / m_nonLinCorrRefScale;
-	
-	float invNLCorr = 1.0;
-	for (size_t power = 1; power <= m_nonLinCorrParamsHG.size(); power++) {
-	  invNLCorr += m_nonLinCorrParamsHG[power - 1]*pow(ampCorrFact, power);
-	}
+      double amplCorrFactor = getAmplitudeCorrection(true);
 
-        m_amplitude /= invNLCorr;
-        m_ampError /= invNLCorr;
-      }
+      m_amplitude *= amplCorrFactor;
+      m_ampError *= amplCorrFactor;
     }
 
     // If LG refit has been requested, do it now
@@ -935,16 +983,8 @@ bool ZDCPulseAnalyzer::DoAnalysis(bool repass)
       prepareLGRefit(m_ADCSamplesLGSub, m_ADCSSampSigLG, m_useSampleLG);
       DoFit(true);
 
-      if (m_haveNonlinCorr) {
-        float ampCorrFact = (m_refitLGAmpl - m_nonLinCorrRefADC) / m_nonLinCorrRefScale;
-	
-	float invNLCorr = 1.0;
-	for (size_t power = 1; power <= m_nonLinCorrParamsHG.size(); power++) {
-	  invNLCorr += m_nonLinCorrParamsHG[power - 1]*pow(ampCorrFact, power);
-	}
-
-	m_refitLGAmplCorr = m_refitLGAmpl/invNLCorr;
-      }
+      double amplCorrFactor = getAmplitudeCorrection(false);
+      m_refitLGAmplCorr = m_refitLGAmpl*amplCorrFactor;
     }
     
     return result;
@@ -1143,10 +1183,15 @@ bool ZDCPulseAnalyzer::AnalyzeData(size_t nSamples, size_t preSampleIdx,
 	if (sampleSig > m_preExpSig) m_preExpSig = sampleSig;
       }
     }
-    
+
+    // Now we search for maxima before the main pulse
+    //
     int loopLimit = (m_havePulse ? m_minDeriv2ndIndex - 2 : m_peak2ndDerivMinSample - 2);
     int loopStart = m_minSampleEvt == 0 ? 1 : m_minSampleEvt;
     
+    float maxPrepulseSig = 0;
+    unsigned int maxPrepulseSample = 0;
+    
     for (int isample = loopStart; isample <= loopLimit; isample++) {
       if (!useSample[isample]) continue;
       
@@ -1154,30 +1199,42 @@ bool ZDCPulseAnalyzer::AnalyzeData(size_t nSamples, size_t preSampleIdx,
       // If any of the second derivatives prior to the peak are significantly negative, we have a an extra pulse
       //   prior to the main one -- as opposed to just an expnential tail
       //
-      m_prePulseSig = -m_samplesDeriv2nd[isample]/(std::sqrt(6.0)*noiseSig);
+      float prePulseSig = -m_samplesDeriv2nd[isample]/(std::sqrt(6.0)*noiseSig);
+
+      // std::cout << m_tag << ", for sample " << isample << ", 2nd derivative = " << m_samplesDeriv2nd[isample]
+      // 		<< " 0.05 * m_minDeriv2nd) = " << 0.05 * m_minDeriv2nd
+      // 		<< ", prePulseSig = " << prePulseSig << std::endl;
       
-      if ((m_prePulseSig > 6 && m_samplesDeriv2nd[isample] < 0.05 * m_minDeriv2nd) ||
+      if ((prePulseSig > 6 && m_samplesDeriv2nd[isample] < 0.05 * m_minDeriv2nd) ||
 	  m_samplesDeriv2nd[isample]  < 0.5*m_minDeriv2nd)
       {
-	if (m_preExpTail) {
-	  //
-	  // We have a prepulse. If we already indicated an negative exponential,
-	  //   if the prepulse has greater significance, we override the negative exponential
-	  //
-	  if (m_prePulseSig >  m_preExpSig) {
-	    m_prePulse = true;
-	    m_preExpTail = false;
-	  }
+	 m_prePulse = true;
+	if (prePulseSig > maxPrepulseSig) {
+	  maxPrepulseSig = prePulseSig;
+	  maxPrepulseSample = isample;
 	}
-	else {
-	  m_prePulse = true;
+      }
+    }
+
+    if (m_prePulse) {
+      //      std::cout << m_tag << ": prepulse sigma = " << m_prePulseSig << ", pre exp sigma = " << m_preExpSig << std::endl;
+      m_prePulseSig = maxPrepulseSig;
+      
+      if (m_preExpTail) {
+	//
+	// We have a prepulse. If we already indicated an negative exponential,
+	//   if the prepulse has greater significance, we override the negative exponential
+	//
+	if (m_prePulseSig >  m_preExpSig) {
+	  m_preExpTail = false;
 	}
-	
-	if (m_prePulse && m_samplesSub[isample] > m_initialPrePulseAmp) {
-	  m_initialPrePulseAmp = m_samplesSub[isample];
-	  m_initialPrePulseT0 = m_deltaTSample * (isample);
+	else {
+	  m_prePulse = false;
 	}
       }
+
+      m_initialPrePulseAmp = m_samplesSub[maxPrepulseSample];
+      m_initialPrePulseT0 = m_deltaTSample * (maxPrepulseSample);
     }
 
     //    if (m_preExpTail) m_prePulse = true;    
@@ -1884,7 +1941,7 @@ void ZDCPulseAnalyzer::dump() const
 void ZDCPulseAnalyzer::dumpTF1(const TF1* func) const
 {
   std::string message = "Dump of TF1: " + std::string(func->GetName());
-  bool continueDump = (*m_msgFunc_p)(ZDCMsg::Verbose, std::move(message));
+  bool continueDump = (*m_msgFunc_p)(ZDCMsg::Verbose, message);
   if (!continueDump) return;
   
   unsigned int npar = func->GetNpar();
@@ -2072,7 +2129,7 @@ std::vector<float> ZDCPulseAnalyzer::Calculate2ndDerivative(const std::vector <f
 //
 float ZDCPulseAnalyzer::obtainDelayedBaselineCorr(const std::vector<float>& samples)
 {
-  const unsigned int nsamples = samples.size();
+  unsigned int nsamples = samples.size();
   
   std::vector<float> derivVec = CalculateDerivative(samples, 2);
   std::vector<float> deriv2ndVec = Calculate2ndDerivative(samples, 2);
@@ -2108,10 +2165,6 @@ float ZDCPulseAnalyzer::obtainDelayedBaselineCorr(const std::vector<float>& samp
   // Because of the way the above analysis is done, we can always
   // Go back one even and one odd sample and forward one odd sample.
   //
-  //if minIndex is < 2 or >samples.size() the result is undefined; prevent this:
-  if (minIndex<2 or (minIndex+1) >=nsamples){
-    throw std::out_of_range("minIndex out of range in ZDCPulseAnalyzer::obtainDelayedBaselineCorr");
-  }
   float sample0 = samples[minIndex - 2];
   float sample1 = samples[minIndex - 1];
   float sample2 = samples[minIndex];

ForwardDetectors/ZDC/ZdcAnalysis/Root/ZdcAnalysisTool.cxx

@@ -51,6 +51,9 @@ ZdcAnalysisTool::ZdcAnalysisTool(const std::string& name)
     declareProperty("ForceCalibRun", m_forceCalibRun = -1); // last run of Pb+Pb 2015
     declareProperty("ForceCalibLB", m_forceCalibLB = 814); // last LB of Pb+Pb 2015
 
+    declareProperty("DoFADCCorr", m_doFADCCorr = false); 
+    declareProperty("DoFADCCorrPerSample", m_doFADCCorrPerSample = false); 
+    
     // The following parameters are primarily used for the "default" configuration, but also may be
     //   use to modify/tailor other configurations
     //
@@ -2439,6 +2442,91 @@ void ZdcAnalysisTool::setTimeCalibrations(unsigned int runNumber)
     }
 }
 
+void ZdcAnalysisTool::setFADCCorrections(unsigned int runNumber)
+{
+  std::string filename;
+  
+  if (m_LHCRun==3) {
+    std::string runString;
+    
+    if (runNumber == 0) runString = "ZdcFADCCorr_" + m_configuration + "_default.root";
+    else runString = ("ZdcFADCCorr_Run"+TString::Itoa(runNumber,10)+".root").Data();
+    
+    filename = PathResolverFindCalibFile("ZdcAnalysis/" + runString );
+  }
+  else {
+    ATH_MSG_WARNING("setFADCCorrections: FADC corrections not implemented for Run 2");
+    return;
+  }
+  
+  ATH_MSG_INFO("Opening FADC corrections file " << filename);
+  std::unique_ptr<TFile> fFADCCorr(TFile::Open(filename.c_str(), "READ"));
+  
+  if (!fFADCCorr->IsOpen()) {
+    ATH_MSG_INFO ("setFADCCorrections: failed to open file: " << filename);
+    throw std::runtime_error ("ZdcAnalysisTool failed to open FADCCorrections file " + filename);
+  }
+  
+  // Attempt to read histograms with corrections from file
+  //
+  bool readSuccess = true;
+  std::array<std::array<std::unique_ptr<const TH1>, 4>, 2> histogramsHG;
+  std::array<std::array<std::unique_ptr<const TH1>, 4>, 2> histogramsLG;
+  
+  for (size_t side : {0, 1}) {
+    for (int module : {0, 1, 2, 3}) {
+      std::string histNameHG = "ZDC_FADCCorr_s" + std::to_string(side) + "_m_HG" + std::to_string(module);
+      std::string histNameLG = "ZDC_FADCCorr_s" + std::to_string(side) + "_m_LG" + std::to_string(module);
+
+      ATH_MSG_DEBUG("setFADCCorrections: Searching for histograms HG and LG: " << histNameHG << ", " << histNameLG);
+      
+      TH1* histHG_ptr = static_cast<TH1*>(fFADCCorr->GetObjectChecked(histNameHG.c_str(), "TH1"));
+      TH1* histLG_ptr = static_cast<TH1*>(fFADCCorr->GetObjectChecked(histNameLG.c_str(), "TH1"));
+
+      if (!histHG_ptr || !histLG_ptr) {
+	std::string errMsg = "setFADCCorrections: unable to read FADC correction histogram(s) ";
+	if (!histHG_ptr) errMsg += histNameHG + " ";
+	if (!histLG_ptr) errMsg += histNameLG;
+
+	ATH_MSG_ERROR(errMsg);
+	readSuccess = false;
+	break;
+      }
+      else {
+	//
+	//  Check for valid range
+	//
+	if (std::abs(histHG_ptr->GetXaxis()->GetXmin()) > 1e-3 || std::abs(histHG_ptr->GetXaxis()->GetXmax() - 4096) > 1e-3) {
+	  ATH_MSG_ERROR("setFADCCorrections: invalid axis range for HG FADC corrections in histogram with name " << histNameHG);
+	  readSuccess = false;
+	  break;
+	}
+	if (std::abs(histLG_ptr->GetXaxis()->GetXmin()) > 1e-3 || std::abs(histLG_ptr->GetXaxis()->GetXmax() - 4096) > 1e-3) {
+	  ATH_MSG_ERROR("setFADCCorrections: invalid axis range for HG FADC corrections in histogram with name " << histNameLG);
+	  readSuccess = false;
+	  break;
+	}
+	
+	histogramsHG[side][module].reset(histHG_ptr);
+	histogramsLG[side][module].reset(histLG_ptr);
+      }
+    }
+  }
+    
+  fFADCCorr->Close();
+
+  if (readSuccess) {
+    m_zdcDataAnalyzer->enableFADCCorrections(m_doFADCCorrPerSample, histogramsHG, histogramsLG);
+  }
+  else {
+    ATH_MSG_ERROR("setFADCCorrections: due to at least one error, FADC corrections are not implemented");
+    m_doFADCCorr = false;
+  }
+  
+  return;
+}
+
+
 StatusCode ZdcAnalysisTool::reprocessZdc()
 {
     if (!m_init)

ForwardDetectors/ZDC/ZdcAnalysis/ZdcAnalysis/ZDCDataAnalyzer.h

@@ -179,6 +179,12 @@ public:
 			const std::array<std::array<std::vector<double>, 4>, 2>& parsHGArr, 
 			const std::array<std::array<std::vector<double>, 4>, 2>& parsLGArr);
 
+  void enableFADCCorrections(bool correctPerSample,
+			     std::array<std::array<std::unique_ptr<const TH1>, 4>, 2>& correHistHG,
+			     std::array<std::array<std::unique_ptr<const TH1>, 4>, 2>& correHistLG);
+
+  void disableFADCCorrections();
+
   void LoadEnergyCalibrations(std::array<std::array<std::unique_ptr<TSpline>, 4>, 2>& calibSplines)
   {
     (*m_msgFunc_p)(ZDCMsg::Verbose, "Loading energy calibrations");

ForwardDetectors/ZDC/ZdcAnalysis/ZdcAnalysis/ZDCPulseAnalyzer.h

 /*
-  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef ZDCANALYSIS_ZDCPulseAnalyzer_h
@@ -9,6 +9,7 @@
 #include "ZdcAnalysis/ZDCMsg.h"
 #include "TGraphErrors.h"
 #include "TFitter.h"
+#include "TList.h"
 #include "TF1.h"
 #include "TH1.h"
 
@@ -162,6 +163,12 @@ private:
   unsigned int m_postExclLGADCThresh{0};
   unsigned int m_maxSamplesPostExcl{0};
 
+  // Enable of a user-provided filter on the FADC samples
+  //
+  bool m_haveUserFilter{false};
+  void (*m_userFilterHG)(std::vector<float>& FADCSamples, std::vector<bool> useSamples){};
+  void (*m_userFilterLG)(std::vector<float>& FADCSamples, std::vector<bool> useSamples){};
+  
   //
   unsigned int m_timingCorrMode{NoTimingCorr};
   float m_timingCorrRefADC{500};
@@ -175,6 +182,11 @@ private:
   std::vector<float> m_nonLinCorrParamsHG;
   std::vector<float> m_nonLinCorrParamsLG;
 
+  bool m_haveFADCCorrections{false};
+  bool m_FADCCorrPerSample{false};
+  std::unique_ptr<const TH1> m_FADCCorrHG{};
+  std::unique_ptr<const TH1> m_FADCCorrLG{};
+  
   // Histogram used to perform the fits and function wrappers
   //
   mutable std::unique_ptr<TH1> m_fitHist;
@@ -357,6 +369,8 @@ private:
                   );
 
 
+  double getAmplitudeCorrection(bool highGain);
+    
   static std::vector<float> Calculate2ndDerivative(const std::vector <float>& inputData, unsigned int step);
   static std::vector<float> CalculateDerivative(const std::vector <float>& inputData, unsigned int step);
   static float obtainDelayedBaselineCorr(const std::vector<float>& samples);
@@ -518,6 +532,12 @@ public:
     m_haveNonlinCorr = true;
   }
 
+  // Provide a historam that provides per-ADC channel correction factors for integral and differential
+  //   non-linearities
+  //
+  void enableFADCCorrections(bool correctPerSample, std::unique_ptr<const TH1>& correHistHG, std::unique_ptr<const TH1>& correHistLG);
+  void disableFADCCorrections() {m_haveFADCCorrections = false;}
+  
   bool LoadAndAnalyzeData(const std::vector<float>& ADCSamplesHG, const std::vector<float>& ADCSamplesLG);
 
   bool LoadAndAnalyzeData(const std::vector<float>& ADCSamplesHG, const std::vector<float>& ADCSamplesLG,

ForwardDetectors/ZDC/ZdcAnalysis/ZdcAnalysis/ZdcAnalysisTool.h

@@ -50,7 +50,8 @@ public:
 
   void setEnergyCalibrations(unsigned int runNumber);
   void setTimeCalibrations(unsigned int runNumber);
-
+  void setFADCCorrections(unsigned int runNumber = 0);
+  
   float getModuleSum(int side);
 
   float getCalibModuleSum(int side);
@@ -65,7 +66,6 @@ public:
 
   double getTriggerEfficiency(int side);
   double getTriggerEfficiencyUncertainty(int side);
-  bool m_doTimeCalib;
 
   const ZDCDataAnalyzer* getDataAnalyzer() {return m_zdcDataAnalyzer.get();}
 
@@ -151,9 +151,12 @@ private:
   const xAOD::ZdcModuleContainer* m_zdcSums {nullptr};
   bool m_flipEMDelay;
   unsigned int m_lowGainMode;
-  bool m_combineDelay;
-  bool m_doCalib;
-  bool m_doTrigEff;
+  bool m_combineDelay{false};
+  bool m_doCalib{false};
+  bool m_doTrigEff{false};
+  bool m_doTimeCalib{false};
+  bool m_doFADCCorr{false};
+  bool m_doFADCCorrPerSample{false};
   int m_forceCalibRun;
   int m_forceCalibLB;
   std::string m_calibVersion;