#include <ZdcRecChannelTool.h>

Inheritance diagram for ZdcRecChannelTool:

Collaboration diagram for ZdcRecChannelTool:

Public Member Functions
	ZdcRecChannelTool (const std::string &type, const std::string &name, const IInterface *parent)
virtual	~ZdcRecChannelTool ()
virtual StatusCode	initialize ()
virtual StatusCode	finalize ()
int	makeRawFromDigits (const ZdcDigitsCollection &data_collection, ZdcRawChannelCollection &raw_collection)
int	getCalibration (const ZdcDigitsCollection &data_collection, ZdcRawChannelCollection &raw_collection)
ServiceHandle< StoreGateSvc > &	evtStore ()
	The standard `StoreGateSvc` (event store) Returns (kind of) a pointer to the `StoreGateSvc`.
const ServiceHandle< StoreGateSvc > &	detStore () const
	The standard `StoreGateSvc/DetectorStore` Returns (kind of) a pointer to the `StoreGateSvc`.
virtual StatusCode	sysInitialize () override
	Perform system initialization for an algorithm.
virtual StatusCode	sysStart () override
	Handle START transition.
virtual std::vector< Gaudi::DataHandle * >	inputHandles () const override
	Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * >	outputHandles () const override
	Return this algorithm's output handles.
Gaudi::Details::PropertyBase &	declareProperty (Gaudi::Property< T, V, H > &t)
void	updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream &	msg () const
bool	msgLvl (const MSG::Level lvl) const

Static Public Member Functions
static const InterfaceID &	interfaceID ()
	AlgTool InterfaceID.

Protected Member Functions
void	renounceArray (SG::VarHandleKeyArray &handlesArray)
	remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void >	renounce (T &h)
void	extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
	Add StoreName to extra input/output deps as needed.

Private Types
typedef ServiceHandle< StoreGateSvc >	StoreGateSvc_t

Private Member Functions
int	getTimingCFD (const Identifier &id, const std::vector< std::vector< int > > &wfm)
int	getTimingSinc (const Identifier &id, const std::vector< std::vector< int > > &wfm)
int	getTimingSinc2 (const Identifier &id, const std::vector< std::vector< int > > &wfm)
Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
	specialization for handling Gaudi::Property<SG::VarHandleKey>

Static Private Member Functions
static double	fx (double x0, void *params)
static double	fxCallback (double d, void *v)

Private Attributes
gsl_interp_accel *	m_interp_acc {}
gsl_spline *	m_spline {}
std::vector< float >	m_cfd_result
int	m_bwl_time_resolution
std::vector< float >	m_wfm_bwl
std::vector< float >	m_bwl_vpeak
std::vector< float >	m_bwl_tpeak
std::vector< float >	m_bwl_vpeak2
std::vector< float >	m_bwl_tpeak2
unsigned int	m_nsamples
float	m_sample_time
int	m_delta
int	m_zeroSupress
float	m_cfd_fraction
float	m_cfd_delay
const ZdcID *	m_zdcId {}
StoreGateSvc_t	m_evtStore
	Pointer to StoreGate (event store by default)
StoreGateSvc_t	m_detStore
	Pointer to StoreGate (detector store by default)
std::vector< SG::VarHandleKeyArray * >	m_vhka
bool	m_varHandleArraysDeclared

Static Private Attributes
static const int	s_FADC_SATURATION = 1022
static const int	s_CFD_FRACTION = -3
static const int	s_SAMPLING_TIME = 25000
static const int	s_HALF_SAMPLING_TIME = 12500
static const int	s_GAIN_RATIO = 1

Detailed Description

Definition at line 71 of file ZdcRecChannelTool.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t

privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ ZdcRecChannelTool()

ZdcRecChannelTool::ZdcRecChannelTool	(	const std::string &	type,
		const std::string &	name,
		const IInterface *	parent )

Definition at line 52 of file ZdcRecChannelTool.cxx.

                                               :
                           AthAlgTool(type, name, parent),
                           m_bwl_time_resolution(100),
                           m_nsamples(7),
                           m_sample_time (25.),
                           m_cfd_fraction (0.35),
                           m_cfd_delay (20.)
 
{
    //Declare properties here...
 
    declareInterface<ZdcRecChannelTool>(this);
 
    declareProperty("ZeroSuppress", m_zeroSupress = 0,
                            "Supress channels with only 0");
 
    declareProperty("DeltaPeak", m_delta = 5,
                            "Minimum difference between min and max to be considered a signal");
 
 
}

◆ ~ZdcRecChannelTool()

ZdcRecChannelTool::~ZdcRecChannelTool ( )

virtual

Definition at line 78 of file ZdcRecChannelTool.cxx.

79{

80}

Member Function Documentation

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty	(	Gaudi::Property< T, V, H > &	hndl,
		const SG::VarHandleKeyType &	)

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

  {
    return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                       hndl.value(), 
                                                       hndl.documentation());
 
  }

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t )

inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

                                                                     {
    typedef typename SG::HandleClassifier<T>::type htype;
    return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
  }

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const

inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

95{ return m_detStore; }

AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore

StoreGateSvc_t m_detStore

Definition AthCommonDataStore.h:393

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )

inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

85{ return m_evtStore; }

AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore

StoreGateSvc_t m_evtStore

Definition AthCommonDataStore.h:390

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps )

protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ finalize()

StatusCode ZdcRecChannelTool::finalize ( )

virtual

Definition at line 138 of file ZdcRecChannelTool.cxx.

{
    msg(MSG::INFO) << "Finalizing " << name() << endmsg;
 
    gsl_spline_free (m_spline);
    gsl_interp_accel_free (m_interp_acc);
 
    return StatusCode::SUCCESS;
}

◆ fx()

double ZdcRecChannelTool::fx	(	double	x0,
		void *	params )

staticprivate

Definition at line 617 of file ZdcRecChannelTool.cxx.

    {
        ZdcRecChannelTool* p = static_cast<ZdcRecChannelTool*>(params);
        //std::cout << " MMMMMMMMMMMMAAAAAAAAA " << p->m_cfd_delay << std::endl;
        return gsl_spline_eval (p->m_spline, x0, p->m_interp_acc);
    }

◆ fxCallback()

double ZdcRecChannelTool::fxCallback	(	double	d,
		void *	v )

inlinestaticprivate

Definition at line 105 of file ZdcRecChannelTool.h.

      {
        CCallbackHolder* h = static_cast<CCallbackHolder*>(v);
        return h->cls->fx(d, h->cls);
      }

◆ getCalibration()

int ZdcRecChannelTool::getCalibration	(	const ZdcDigitsCollection &	data_collection,
		ZdcRawChannelCollection &	raw_collection )

Definition at line 457 of file ZdcRecChannelTool.cxx.

{
    int ncha = 0;
    for (const ZdcDigits* digits_p : mydata) {
        digits_p->identify();
        ncha++;
    }
    int ncal = ChannelCollection.size();
    msg(MSG::DEBUG) << " Zdc ---> Calibration for " << ncal << " channels " << endmsg;
    return  ncha;
}

◆ getTimingCFD()

int ZdcRecChannelTool::getTimingCFD	(	const Identifier &	id,
		const std::vector< std::vector< int > > &	wfm )

private

Definition at line 473 of file ZdcRecChannelTool.cxx.

{
 
    unsigned int i = 0;
    //unsigned int j = 0;
    unsigned int nsamples_local = 14;
    double v[14];
    double t[14];
    double td[14];
    double vd[14];
    double vf[14];
 
    std::vector<float> vaux;
    unsigned int tmin = 0;
    unsigned int tmax = 0;
 
    double x_lo = 0.0;
    double x_hi = 200.0;
    float r = 0;
    int status = 0;
 
    int iterations = 0;
    int max_iterations = 100;
 
    // Thread-safe according to https://www.gnu.org/software/gsl/doc/html/roots.html
    const gsl_root_fsolver_type *T ATLAS_THREAD_SAFE = gsl_root_fsolver_brent;
    gsl_root_fsolver *s = gsl_root_fsolver_alloc (T);
    gsl_function F;
 
    CCallbackHolder cc{};
 
 
    std::vector<int>  y;
 
    std::vector<std::vector<int> >::const_iterator vit;
    std::vector<int>::iterator it;
 
    std:: vector<float> result;
 
    m_cfd_result.clear();
 
    msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool::getTimingCFD: Id " << id.getString() ;
 
    for (i=0;i<14;i++) {
        t[i]  = 0.;
        td[i] = 0.;
        v[i]  = 0.;
        vd[i] = 0.;
        vf[i] = 0.;
    }
 
    for (vit = wfm.begin(); vit<wfm.end(); ++vit) {
        y = *vit;
        i = 4;
        for (it=y.begin();it != y.end();++it) {
            v[i] = *it;
            i++;
        }
        for (i=0;i<nsamples_local;i++) {
            t[i] = m_sample_time*i ; //TODO: Do it only once somewhere else
            td[i] = t[i]-m_cfd_delay;
        }
 
        gsl_spline_init (m_spline, t, v, nsamples_local);
 
 
 
        for (i=0;i<nsamples_local;i++) {
            vd[i] = gsl_spline_eval (m_spline, td[i], m_interp_acc);
            vf[i] = vd[i] - v[i]* m_cfd_fraction;
        }
 
        //gsl_spline_free (m_spline); // Not sure if this is necessary
 
        gsl_spline_init (m_spline, t, vf, nsamples_local);
 
        cc.cls = this;
        F.function = &ZdcRecChannelTool::fxCallback;
        F.params = &cc;
 
 
 
 
        //Get the range for root finding
        vaux.assign(vf,vf+14);
        tmin = std::min_element(vaux.begin(), vaux.end()) - vaux.begin();
        tmax = std::max_element(vaux.begin(), vaux.end()) - vaux.begin();
        if (tmin < tmax)
            {
              x_lo = t[tmin];
              x_hi = t[tmax];
            }
          else
            {
              x_lo = t[tmax];
              x_hi = t[tmin];
            }
        //std::cout << "::::::::::tmin,tmax " << x_lo << " " << x_hi << std::endl;
 
        /*
        std::cout << "::::::::::V  " << id.getString() << " "<< j << " " ;
        for (i=0;i<nsamples_local;i++) {
            std::cout << v[i] << " " ;
        }
        std::cout << std::endl;
 
        std::cout << "::::::::::VD  " << id.getString() << " "<< j << " " ;
        for (i=0;i<nsamples_local;i++) {
            std::cout << vd[i] << " " ;
        }
        std::cout << std::endl;
 
        std::cout << "::::::::::VF  " << id.getString() << " "<< j << " " ;
        for (i=0;i<nsamples_local;i++) {
                    std::cout << vf[i] << " " ;
        }
        std::cout << std::endl;
       */
        //j++;
 
        gsl_root_fsolver_set (s, &F, x_lo, x_hi);
 
        iterations = 0;
        do
            {
                iterations++;
                status = gsl_root_fsolver_iterate(s);
                r =      gsl_root_fsolver_root(s);
                x_lo =   gsl_root_fsolver_x_lower(s);
                x_hi =   gsl_root_fsolver_x_upper(s);
                status = gsl_root_test_interval (x_lo, x_hi, 0, 0.001);
                //if (status == GSL_SUCCESS) printf ("Converged:\n");
            }
 
        while (status == GSL_CONTINUE && iterations < max_iterations);
        //std::cout << " COnvergence ----> " << iterations << " " << r << " " << x_lo << " " << x_hi << std::endl;
        m_cfd_result.push_back(r-100); //subtraction accounts for zero padding
    }
    gsl_root_fsolver_free (s);
    return 0;
}

◆ getTimingSinc()

int ZdcRecChannelTool::getTimingSinc	(	const Identifier &	id,
		const std::vector< std::vector< int > > &	wfm )

private

Definition at line 633 of file ZdcRecChannelTool.cxx.

{
    unsigned int i = 0;
    unsigned int k = 0;
    int t = 0;
    int tpeak = 0;
    float tfpeak = 0.;
    float vpeak = 0.;
    float x = 0.;
    float z = 0.;
    float t_cor = 0. ;
    int mSide = 0;
    int mModule = 0;
    int mType = 0;
    int mChannel = 0;
    int do_tcor = false;
 
 
    std::vector<float> vt(7);
    std::vector<int> y;
 
    std::vector<float>::iterator vf_it;
    std::vector<std::vector<int> >::const_iterator vvi_it;
 
    //initialize the vectors
    m_bwl_vpeak.resize(4,0);
    m_bwl_tpeak.resize(4,0);
 
    //zId = id.get_identifier32().get_compact();
    msg(MSG::DEBUG) << "--> ZDC : START OF MODIFICATION 2 " << endmsg ;
        mSide    = m_zdcId->side(id);
        mModule  = m_zdcId->module(id);
        mType    = m_zdcId->type(id);
        mChannel = m_zdcId->channel(id);
        msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool::getTimingSinc: "
                        << " id;Side;Module;Type;Channel:  "
                        << id.getString() << ";"
                        << mSide    << ";"
                        << mModule  << ";"
                        << mType    << ";"
                        << mChannel << endmsg;
 
    for (vvi_it = wfm.begin(); vvi_it<wfm.end(); ++vvi_it) {
        //FIXME: Change to the method of ID identification
        //if ((zId == 0xec000000) || (zId == 0xed000000) ||
        //  (zId == 0xec200000) || (zId == 0xed200000) ||
        //  (zId == 0xec400000) || (zId == 0xed400000) ||
        //  (zId == 0xec600000) || (zId == 0xed600000) )
        //  {
            if (mType == 0) {
 
            y = *vvi_it;
 
            t = 0;
            //Fill the interpolated vector
            for (vf_it = m_wfm_bwl.begin(); vf_it != m_wfm_bwl.end(); ++vf_it ) {
                z = 0.;
                for (i=0;i<m_nsamples;i++) {
                    x =  (TMath::Pi() *(t*0.1 - i*25.))/25. ;
                    if (fabs(x) < 1e-8) {
                        z = z + y[i]*1.0;
                    }
                    else {
                        z =  z + y[i]*(sin(x)/(x));
                    }
                }
                *vf_it = z;
                //if (id.get_identifier32().get_compact() == 0xec400000)  {
                //std::cout << "========== " << t << "  " << *vf_it << std::endl;
                //}
                t++;
            }
 
            //Get The maximum and the position of the peak
            tpeak = std::max_element(m_wfm_bwl.begin(),m_wfm_bwl.end()) - m_wfm_bwl.begin();
            vpeak = m_wfm_bwl[tpeak];
 
            //Account for 100 ps step
            tfpeak = tpeak/10.;
 
 
            //The correction for non-linear behavior
            //Right now only two channels implemented:
            //First hadronic modules for A and C sides
            //zId = id.get_identifier32().get_compact();
 
            //Side C, Mod 0, Low Gain
            //if (id.get_compact() == 0xec400000) {
            if ( (mSide == -1) && (mType == 0) && (mModule == 0) && (mChannel == 0) ) {
                t_cor = -777.277 +
                    33.685      * tfpeak +
                    -0.483265   * pow(tfpeak,2) +
                    0.00234842 * pow(tfpeak,3);
            }
 
            //Side A, Mod 0, Low Gain
            //if (id.get_compact() == 0xed400000) {
            if ( (mSide == 1) && (mType == 0) && (mModule == 0) && (mChannel == 0) ) {
            //piecewise fit
                if (tfpeak <= 68) {
                    t_cor = -1749.27 +
                        80.5426 * tfpeak +
                        -1.229    * pow(tfpeak,2) +
                        0.006276 * pow(tfpeak,3);
                }
                else {
                    t_cor = -1965.52 +
                        84.0379     * tfpeak +
                        -1.19169    * pow(tfpeak,2) +
                        0.005659   * pow(tfpeak,3);
                }
            }
 
            if (do_tcor) m_bwl_tpeak[k] = t_cor;
            else m_bwl_tpeak[k] = tpeak;
            m_bwl_vpeak[k] = vpeak;
        }
        else {
            m_bwl_tpeak[k] = -999;
            m_bwl_vpeak[k] = -999;
        }
        k++;
    }
    return k;
}

◆ getTimingSinc2()

int ZdcRecChannelTool::getTimingSinc2	(	const Identifier &	id,
		const std::vector< std::vector< int > > &	wfm )

private

Definition at line 764 of file ZdcRecChannelTool.cxx.

{
      int i = 0;
      int wfmIndex = 0;
      //float energy = 0.;
      //float w = 0.;
      //float timing = 0.;
      int error = 0;
      int warning = 0;
 
      double Slices[10];
      int NSlice = 5;
      double gain=1.;
      double pedestal=0;
      double CFD_frac = 1.0;
      bool corr=0;
 
      std::vector<std::vector<int> >::const_iterator vit;
      std::vector<int>  y;
 
      int mType    = 0;
 
 
 
      //Identifier::value_type zId;
 
      /*
       Id:       Side       Module Number  Module Type
       0xec000000   0  (C)           0               0
       0xec400000   0  (C)           1               0
       0xec800000   0  (C)           2               0
       0xecc00000   0  (C)           3               0
       0xed000000   1  (A)           0               0
       0xed400000   1  (A)           1               0
       0xed800000   1  (A)           2               0
       0xedc00000   1  (A)           3               0
       */
      //Initialize the vectors
      m_bwl_vpeak2.resize(4,0);
      m_bwl_tpeak2.resize(4,0);
 
      mType    = m_zdcId->type(id);
 
 
 
      //This controls if we get timing from all channels or
      //from only a subset
      wfmIndex = 2;
      i = 0;
      for (vit = wfm.begin(); vit<wfm.end(); ++vit) {
          if ( (i < wfmIndex) && (mType == 0) )
            {
              y = *vit;
 
              // This is for the sinx/x interpolation
 
              for(int I=0;I<NSlice;I++)  {Slices[I]=y[I];}
              ZdcSignalSinc zdcSignalSinc(NSlice);
 
              zdcSignalSinc.process(Slices,gain,pedestal,CFD_frac,corr);
 
 
              m_bwl_vpeak2[i] = zdcSignalSinc.getAmp() ;
              m_bwl_tpeak2[i] = zdcSignalSinc.getTime() ;
 
              error   = zdcSignalSinc.getError();
              warning = zdcSignalSinc.getWarning();
 
              msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool::getTimingSinc2: "
                      << error << ";" << warning << endmsg;
 
            }
          else {
              m_bwl_vpeak2[i] = -999 ;
              m_bwl_tpeak2[i] = -999 ;
          }
          i++;
      }
      return 0;
 
}

◆ initialize()

StatusCode ZdcRecChannelTool::initialize ( )

virtual

Definition at line 85 of file ZdcRecChannelTool.cxx.

{
    msg(MSG::INFO) << "Initializing " << name() << endmsg;
 
    //Get the pedestal information for the channels.
    //For now, this is a file; later on it will be moved to a database
 
 
    // Bandwidth interpolation method initializations
    // making sure we use only 50 ps steps for the resolution
 
    if (m_bwl_time_resolution < 50)  m_bwl_time_resolution = 50;
    m_bwl_time_resolution = floor(m_bwl_time_resolution/50.) * 50;
 
    int s = m_nsamples*(s_SAMPLING_TIME/m_bwl_time_resolution);
 
    m_wfm_bwl.resize(s);
 
    msg(MSG::INFO) << "Using a time step of  "
                   <<  m_bwl_time_resolution
                   << "ps for the Bandwidth Limited Sin(x)/x Interpolation: "
                   << " Vector Size = "
                   << s
                   << endmsg;
 
 
 
 
    // GSL Interpolation functions initializations
    m_interp_acc = gsl_interp_accel_alloc ();
    // Thread-safe according to https://www.gnu.org/software/gsl/doc/html/roots.html
    const gsl_interp_type *interp_type ATLAS_THREAD_SAFE = gsl_interp_akima;
    m_spline = gsl_spline_alloc (interp_type, 14); //FIXME: TWICE THE SAMPLES
 
        //Load Mapping
    msg(MSG::DEBUG) << "--> ZDC : START OF MODIFICATION 0" << endmsg ;
 
    const ZdcID* zdcId = nullptr;
    if (detStore()->retrieve( zdcId ).isFailure() ) {
        msg(MSG::ERROR) << "execute: Could not retrieve ZdcID object from the detector store" << endmsg;
        return StatusCode::FAILURE;
    }
    else {
        msg(MSG::DEBUG) << "execute: retrieved ZdcID" << endmsg;
    }
    m_zdcId = zdcId;
 
    msg(MSG::DEBUG) << "--> ZDC : END OF MODIFICATION 0" << endmsg ;
    return StatusCode::SUCCESS;
}

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const

overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ interfaceID()

const InterfaceID & ZdcRecChannelTool::interfaceID ( )

static

AlgTool InterfaceID.

Definition at line 44 of file ZdcRecChannelTool.cxx.

{
    return IID_IZdcRecChannelTool;
 
}

◆ makeRawFromDigits()

int ZdcRecChannelTool::makeRawFromDigits	(	const ZdcDigitsCollection &	data_collection,
		ZdcRawChannelCollection &	raw_collection )

Traditional code bellow (Energy by peak or sum, time by ratio of samples

Definition at line 151 of file ZdcRecChannelTool.cxx.

{
    //ZdcRawChannel *z;
    Identifier id;
    float soma = 0.;
    float pico = 0.;
    float tzc  = 0.;
    float tsr  = 0.;
    float t0   = 0.;
    float t1   = 0.;
    float m    = 0.;
 
    //bool pp = false;
    int i = 0;
    //int j = 0;
    int k = 0;
    //int s = 0;
    int ped = 0;
    int imax = 0;
 
    int mSide    = 0;
    int mModule  = 0;
    int mType    = 0;
    int mChannel = 0;
 
    int ncha = 0;
 
 
    //TODO use BOOST containers here, more flexible and inteligent, and it will
    //      nevertless stay inside this method, so no fussing with complex collections
 
    std::vector<std::vector<int> > wfm;
    std::vector<std::vector<int> >::iterator vit;
    std::vector<int>::iterator it;
 
    std::vector<float> energy_sum;
    std::vector<float> energy_peak;
    std::vector<float>  time_cfd;
    std::vector<float>  time_sratio;
    std::vector<float>  chi;
 
    std::vector<int>  v;
    std::vector<int>  vd;
    std::vector<int>  vi;
 
    std::vector<float> vcfd1;
    std::vector<float> vcfd2;
    std::vector<float> vcfd3;
 
    //Identifier::value_type zId;
 
 
    // Main procedure: Reduce the data to energy and time, and store it into
    // ZdcRawDataCollection
 
    //digits_p = *mydata.begin();
    //nsamples = digits_p->m_nSamples;
    //vi.resize(2*nsamples);
    vi.resize(14);
 
    for (const ZdcDigits* digits_p : mydata) {
 
            msg(MSG::DEBUG) << "--> ZDC : START OF MODIFICATION " << endmsg ;
        id = digits_p->identify();
        mSide    = m_zdcId->side(id);
        mModule  = m_zdcId->module(id);
        mType    = m_zdcId->type(id);
        mChannel = m_zdcId->channel(id);
 
        msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool::makeRawFromDigits: "
                            << " id;Side;Module;Type;Channel:  "
                            << id.getString() << ";"
                            << mSide    << ";"
                            << mModule  << ";"
                            << mType    << ";"
                            << mChannel << endmsg;
 
 
        //1) Check the high gain. If it has saturation (max = FADC_SATURATION)
        //   drop and use the low gain
        //
        wfm.clear();
        wfm.resize(2);
        energy_sum.clear();
        energy_peak.clear();
        time_cfd.clear();
        time_sratio.clear();
        chi.clear();
 
        wfm[0] = digits_p->get_digits_gain0_delay0();
        if (wfm[0].empty()) {
            msg(MSG::DEBUG) << "ZERO SIZE g0d0 at ID  " << id.getString() << endmsg;
            wfm[0].resize(7);
        }
 
        wfm[1] = digits_p->get_digits_gain1_delay0();
        if (wfm[1].empty()) {
            msg(MSG::DEBUG) << "ZERO SIZE g1d0 at ID  " << id.getString() << endmsg;
            wfm[1].resize(7);
        }
 
 
        msg(MSG::DEBUG) << " ------- 1 " << endmsg;
 
        //2) Check to see if there are delayed information
        //        module type 0 -> (full Energy) delayed info
        //        module type 1 -> (segmented) do not have delayed info
        //        - use the identifiers here -
        //
        //if (((id.get_identifier32().get_compact() >> 21) & 1) == 0) {
 
        //zId = id.get_identifier32().get_compact();
 
        //Only type 0 (Total Module Energy) has delay information
        //if ( ( (zId == 0xec000000) || (zId == 0xed000000) ||
        //      (zId == 0xec200000) || (zId == 0xed200000) ||
        //      (zId == 0xec400000) || (zId == 0xed400000) ||
        //      (zId == 0xec600000) || (zId == 0xed600000) ) ) {
        if (mType == 0) {
            //std::cout << "*** Resize 4 at Id " << zId << std::endl ;
            wfm.resize(4);
            wfm[2] = digits_p->get_digits_gain0_delay1();
            wfm[3] = digits_p->get_digits_gain1_delay1();
        }
 
 
 
 
        msg(MSG::DEBUG) << " ------- 2 " << endmsg;
        //3) Subtratct baseline pedestal
        //We need to be carefull. Sometimes the vector size here is zero (PPM flaw) and
        //the code crashs if we do not treat this.
        i = 0;
        for (vit = wfm.begin(); vit<wfm.end(); ++vit) {
            if (vit->empty()) vit->resize(7);
            ped = *vit->begin();
            for (it=vit->begin(); it<vit->end();++it) {
                (*it) -= ped;
                //if ((i==1) || (i==3) ) (*it) = (*it) * s_GAIN_RATIO;
            }
            i++;
        }
        /*
        if (((id.get_identifier32().get_compact() >> 21) & 1) == 0) {
            if ( ( *(std::max_element(wfm[0].begin(),wfm[0].end()) ) > 300 ) ||
                 ( *(std::max_element(wfm[1].begin(),wfm[1].end()) ) > 300 ) ||
                 ( *(std::max_element(wfm[2].begin(),wfm[2].end()) ) > 300 ) ||
                 ( *(std::max_element(wfm[3].begin(),wfm[3].end()) ) > 300 ) ) {
                //std::cout << "&&&&& OVFLW &&&&&" << std::endl;
                msg(MSG::DEBUG) << "%%%%%%%% ID "   << id.getString() << endmsg;
                msg(MSG::DEBUG) << "%%%%%%%% g0d0 " << wfm[0] << endmsg;
                msg(MSG::DEBUG) << "%%%%%%%% g1d0 " << wfm[1] << endmsg;
                msg(MSG::DEBUG) << "%%%%%%%% g0d1 " << wfm[2] << endmsg;
                msg(MSG::DEBUG) << "%%%%%%%% g1d1 " << wfm[3] << endmsg;
 
            }
        }
        else {
            if ( ( *(std::max_element(wfm[0].begin(),wfm[0].end()) ) > 20 )  ||
                 (*(std::max_element(wfm[1].begin(),wfm[1].end()) ) > 20 ) )    {
                //std::cout << "**** g1d0 OVFLW ****" << std::endl;
                msg(MSG::DEBUG) << "******** ID "   << id.getString() << endmsg;
                msg(MSG::DEBUG) << "******** g0d0 " << wfm[0] << endmsg;
                msg(MSG::DEBUG) << "******** g1d0 " << wfm[1] << endmsg;
 
            }
        }
        */
        msg(MSG::DEBUG) << " ------- 3 " << endmsg;
        //4) Calculate the quantities of interest. For now it will be Energy (by SUM and Peak)
        //   and timing (using undelayed, delayed and combined)
        //   Include a quality factor chi which will follow design definitions for
        //   chi <  0 -> something very bad (for example saturation)
        //   0 >= chi >= 1 -> depends on the algorithm to be used; 1 is very good
 
        k = 0;
        // GSL for CFD Timing
 
        getTimingCFD(id, wfm);
        // Sinx/x Timing (GSL)
        // Not yet needed
        //getTimingSinc(id, wfm);
 
        // Sinx/x Timing (Andrei's coding)
 
        getTimingSinc2(id, wfm);
        for (vit = wfm.begin(); vit<wfm.end(); ++vit) {
            v = *vit;
            soma =   std::accumulate(v.begin(), v.end(), 0) ;
            pico =  *(std::max_element(v.begin(),v.end()) );
            energy_sum.push_back ((float) soma);
            energy_peak.push_back((float) pico);
            if (pico >= (s_FADC_SATURATION -42)) {
                chi.push_back(-1);
            }
            else {
                chi.push_back(0); //nothing to say here
            }
 
 
            // Now, lets also calculate the timing by the ratio of the
            // second to the third sample. Because the correction is module
            // dependent, we postpone the calibration to the next stage
            // Time will be
            // t = (r-R0)/(R1-R0) ns
            // where r =  v[imax -1]/v[imax], (imax = peak)
            // R0, R1 are constants that depends on the module
 
            // Get the position of the maximum
            // Take care of div by 0
            imax = std::max_element(v.begin(),v.end()) - v.begin();
            if ((v[imax] != 0) && (imax>0)) {
                tsr = ((float) v[imax-1])/v[imax];
            }
            else {
                tsr = v[imax];
            }
 
            time_sratio.push_back(tsr);
 
 
            msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool " << id.getString() <<
                           "  Type=" << k  <<
                           "  Energy Sum="  << soma  <<
                           "  Energy Peak=" << pico  <<
                           "  Chi=" << *(chi.end() -1) <<
                           "  t0="  << t0 << "  t1=" << t1 <<
                           "  m="   << m <<
                           "  Time by CFD=" << tzc <<
                           "  Peak at imax = " << imax <<
                                               "  v[imax-1] = " << (imax > 0 ? v[imax-1] : 0) <<
                           "  V[imax] = " << v[imax] <<
                           "  Sample Ratio A1/A2=" << tsr << endmsg;
            k++;
        }
 
        msg(MSG::DEBUG) << " ------- 5 " << endmsg;
        ZdcRawChannel *z = new ZdcRawChannel(id);
 
        msg(MSG::DEBUG) << "CFD EXACT***** ID "<< id.getString() << m_cfd_result << endmsg;
        msg(MSG::DEBUG) << "CFD APPRO***** ID "<< id.getString() << m_bwl_tpeak << endmsg;
 
        z->setSize(3*k); //FIXME very important, but not smart move ...
        /*
         *   What can go here:
         *
         *   energy_peak -> max sample peak
         *   energy_sum  -> sum of samples
         *   time_sratio -> from ratio of samples
         *
         *   m_bwl_vpeak -> Energy  peak from sinc interpolation
         *   m_bwl_tpeak -> Time at peak from sinc interpolation
         *
         *   m_bwl_vpeak2 -> Energy  peak from sinc interpolation (AP code)
         *   m_bwl_tpeak2 -> Time at peak from sinc interpolation (AP code)
         *
         *   m_cfd_result -> timing by software CFD
         */
        for (i=0;i<k;i++) {
            z->setEnergy (i, energy_peak[i]);
            z->setTime   (i, time_sratio[i]);
            z->setChi    (i, chi[i]);
 
            z->setEnergy (i+k, m_bwl_vpeak2[i]);
            z->setTime   (i+k, m_bwl_tpeak2[i]);
            z->setChi    (i+k, chi[i]);
 
            z->setEnergy (i+2*k, energy_sum[i]);
            z->setTime   (i+2*k, m_cfd_result[i]);
            z->setChi    (i+2*k, chi[i]);
 
            //z->setEnergy(i+k, energy_peak[i]);
            //z->setTime(i,time_cfd[i]);
 
 
        }
        ncha++;
        ChannelCollection.push_back(z);
        msg(MSG::DEBUG) << " ------- 6 " << endmsg;
    }
    msg(MSG::DEBUG) << "--> ZDC : ZdcRecChannelTool ChannelCollection size " << ChannelCollection.size() << endmsg ;
    return ncha;
}

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const

inlineinherited

Definition at line 24 of file AthCommonMsg.h.

                                {
    return this->msgStream();
  }

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl ) const

inlineinherited

Definition at line 30 of file AthCommonMsg.h.

                                               {
    return this->msgLevel(lvl);
  }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const

overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h )

inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

  {
    h.renounce();
    PBASE::renounce (h);
  }

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray )

inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

                                                          {
    handlesArray.renounce();
  }

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )

overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in asg::AsgMetadataTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and DerivationFramework::CfAthAlgTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )

overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )

inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

                                               {
    // debug() << "updateVHKA for property " << p.name() << " " << p.toString() 
    //         << "  size: " << m_vhka.size() << endmsg;
    for (auto &a : m_vhka) {
      std::vector<SG::VarHandleKey*> keys = a->keys();
      for (auto k : keys) {
        k->setOwner(this);
      }
    }
  }

Member Data Documentation

◆ m_bwl_time_resolution

int ZdcRecChannelTool::m_bwl_time_resolution

private

Definition at line 115 of file ZdcRecChannelTool.h.

◆ m_bwl_tpeak

std::vector<float> ZdcRecChannelTool::m_bwl_tpeak

private

Definition at line 118 of file ZdcRecChannelTool.h.

◆ m_bwl_tpeak2

std::vector<float> ZdcRecChannelTool::m_bwl_tpeak2

private

Definition at line 123 of file ZdcRecChannelTool.h.

◆ m_bwl_vpeak

std::vector<float> ZdcRecChannelTool::m_bwl_vpeak

private

Definition at line 117 of file ZdcRecChannelTool.h.

◆ m_bwl_vpeak2

std::vector<float> ZdcRecChannelTool::m_bwl_vpeak2

private

Definition at line 122 of file ZdcRecChannelTool.h.

◆ m_cfd_delay

float ZdcRecChannelTool::m_cfd_delay

private

Definition at line 142 of file ZdcRecChannelTool.h.

◆ m_cfd_fraction

float ZdcRecChannelTool::m_cfd_fraction

private

Definition at line 141 of file ZdcRecChannelTool.h.

◆ m_cfd_result

std::vector<float> ZdcRecChannelTool::m_cfd_result

private

Definition at line 111 of file ZdcRecChannelTool.h.

◆ m_delta

int ZdcRecChannelTool::m_delta

private

Definition at line 133 of file ZdcRecChannelTool.h.

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore

privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore

privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_interp_acc

gsl_interp_accel* ZdcRecChannelTool::m_interp_acc {}

private

Definition at line 100 of file ZdcRecChannelTool.h.

100{};

◆ m_nsamples

unsigned int ZdcRecChannelTool::m_nsamples

private

Definition at line 130 of file ZdcRecChannelTool.h.

◆ m_sample_time

float ZdcRecChannelTool::m_sample_time

private

Definition at line 131 of file ZdcRecChannelTool.h.

◆ m_spline

gsl_spline* ZdcRecChannelTool::m_spline {}

private

Definition at line 101 of file ZdcRecChannelTool.h.

101{};

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared

privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka

privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ m_wfm_bwl

std::vector<float> ZdcRecChannelTool::m_wfm_bwl

private

Definition at line 116 of file ZdcRecChannelTool.h.

◆ m_zdcId

const ZdcID* ZdcRecChannelTool::m_zdcId {}

private

Definition at line 156 of file ZdcRecChannelTool.h.

156{};

◆ m_zeroSupress

int ZdcRecChannelTool::m_zeroSupress

private

Definition at line 134 of file ZdcRecChannelTool.h.

◆ s_CFD_FRACTION

const int ZdcRecChannelTool::s_CFD_FRACTION = -3

staticprivate

Definition at line 144 of file ZdcRecChannelTool.h.

◆ s_FADC_SATURATION

const int ZdcRecChannelTool::s_FADC_SATURATION = 1022

staticprivate

Definition at line 137 of file ZdcRecChannelTool.h.

◆ s_GAIN_RATIO

const int ZdcRecChannelTool::s_GAIN_RATIO = 1

staticprivate

Definition at line 152 of file ZdcRecChannelTool.h.

◆ s_HALF_SAMPLING_TIME

const int ZdcRecChannelTool::s_HALF_SAMPLING_TIME = 12500

staticprivate

Definition at line 146 of file ZdcRecChannelTool.h.

◆ s_SAMPLING_TIME

const int ZdcRecChannelTool::s_SAMPLING_TIME = 25000

staticprivate

Definition at line 145 of file ZdcRecChannelTool.h.

The documentation for this class was generated from the following files:

Public Member Functions

Static Public Member Functions

Protected Member Functions

Private Types

Private Member Functions

Static Private Member Functions

Private Attributes

Static Private Attributes

Detailed Description

Member Typedef Documentation

◆ StoreGateSvc_t

Constructor & Destructor Documentation

◆ ZdcRecChannelTool()

◆ ~ZdcRecChannelTool()

Member Function Documentation

◆ declareGaudiProperty()

◆ declareProperty()

◆ detStore()

◆ evtStore()

◆ extraDeps_update_handler()

◆ finalize()

◆ fx()

◆ fxCallback()

◆ getCalibration()

◆ getTimingCFD()

◆ getTimingSinc()

◆ getTimingSinc2()

◆ initialize()

◆ inputHandles()

◆ interfaceID()

◆ makeRawFromDigits()

◆ msg()

◆ msgLvl()

◆ outputHandles()

◆ renounce()

◆ renounceArray()

◆ sysInitialize()

◆ sysStart()

◆ updateVHKA()

Member Data Documentation

◆ m_bwl_time_resolution

◆ m_bwl_tpeak

◆ m_bwl_tpeak2

◆ m_bwl_vpeak

◆ m_bwl_vpeak2

◆ m_cfd_delay

◆ m_cfd_fraction

◆ m_cfd_result

◆ m_delta

◆ m_detStore

◆ m_evtStore

◆ m_interp_acc

◆ m_nsamples

◆ m_sample_time

◆ m_spline

◆ m_varHandleArraysDeclared

◆ m_vhka

◆ m_wfm_bwl

◆ m_zdcId

◆ m_zeroSupress

◆ s_CFD_FRACTION

◆ s_FADC_SATURATION

◆ s_GAIN_RATIO

◆ s_HALF_SAMPLING_TIME

◆ s_SAMPLING_TIME