#include <EcalCoder.h>

Public Types
enum	{ NBITS = 12, MAXADC = 4095, ADCGAINSWITCH = 4079, NGAINS = 3 }

typedef CaloTSamples< float, 10 >	EcalSamples

typedef CorrelatedNoisifier< EcalCorrMatrix >	Noisifier

Public Member Functions
virtual void	analogToDigital (CLHEP::HepRandomEngine *, const EcalSamples &clf, EcalDataFrame &df) const
	from EcalSamples to EcalDataFrame More...

	EcalCoder (bool addNoise, bool PreMix1, Noisifier ebCorrNoise0, Noisifier eeCorrNoise0=nullptr, Noisifier ebCorrNoise1=nullptr, Noisifier eeCorrNoise1=nullptr, Noisifier ebCorrNoise2=nullptr, Noisifier eeCorrNoise2=nullptr)
	ctor More...

void	setFullScaleEnergy (double EBscale, double EEscale)

void	setGainRatios (const EcalGainRatios *gainRatios)

void	setIntercalibConstants (const EcalIntercalibConstantsMC *ical)

void	setPedestals (const EcalPedestals *pedestals)
	can be fetched every event from the EventSetup More...

virtual	~EcalCoder ()
	dtor More...

Private Member Functions
void	encode (const EcalSamples &ecalSamples, EcalDataFrame &df, CLHEP::HepRandomEngine *) const
	produce the pulse-shape More...

void	findGains (const DetId &detId, double theGains[]) const

void	findIntercalibConstant (const DetId &detId, double &icalconst) const

void	findPedestal (const DetId &detId, int gainId, double &pedestal, double &width) const
	not yet implemented More...

double	fullScaleEnergy (const DetId &did) const
	limit on the energy scale due to the electronics range More...

Private Attributes
bool	m_addNoise

const Noisifier *	m_ebCorrNoise [3]

const Noisifier *	m_eeCorrNoise [3]

const EcalGainRatios *	m_gainRatios

const EcalIntercalibConstantsMC *	m_intercals

double	m_maxEneEB

double	m_maxEneEE

const EcalPedestals *	m_peds

bool	m_PreMix1

Detailed Description

Definition at line 27 of file EcalCoder.h.

Member Typedef Documentation

◆ EcalSamples

typedef CaloTSamples<float, 10> EcalCoder::EcalSamples

Definition at line 29 of file EcalCoder.h.

◆ Noisifier

Definition at line 31 of file EcalCoder.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum

Enumerator
NBITS
MAXADC
ADCGAINSWITCH
NGAINS

Definition at line 33 of file EcalCoder.h.

        {
     NBITS = 12,            // number of available bits
     MAXADC = 4095,         // 2^12 -1,  adc max range
     ADCGAINSWITCH = 4079,  // adc gain switch
     NGAINS = 3             // number of electronic gains
   };

Constructor & Destructor Documentation

◆ EcalCoder()

EcalCoder::EcalCoder	(	bool	addNoise,
		bool	PreMix1,
		EcalCoder::Noisifier *	ebCorrNoise0,
		EcalCoder::Noisifier *	eeCorrNoise0 = `nullptr`,
		EcalCoder::Noisifier *	ebCorrNoise1 = `nullptr`,
		EcalCoder::Noisifier *	eeCorrNoise1 = `nullptr`,
		EcalCoder::Noisifier *	ebCorrNoise2 = `nullptr`,
		EcalCoder::Noisifier *	eeCorrNoise2 = `nullptr`
	)

ctor

Definition at line 12 of file EcalCoder.cc.

     : m_peds(nullptr),
       m_gainRatios(nullptr),
       m_intercals(nullptr),
       m_maxEneEB(1668.3),  // 4095(MAXADC)*12(gain 2)*0.035(GeVtoADC)*0.97
       m_maxEneEE(2859.9),  // 4095(MAXADC)*12(gain 2)*0.060(GeVtoADC)*0.97
       m_addNoise(addNoise),
       m_PreMix1(PreMix1) {
   m_ebCorrNoise[0] = ebCorrNoise0;
   assert(nullptr != m_ebCorrNoise[0]);
   m_eeCorrNoise[0] = eeCorrNoise0;
   m_ebCorrNoise[1] = ebCorrNoise1;
   m_eeCorrNoise[1] = eeCorrNoise1;
   m_ebCorrNoise[2] = ebCorrNoise2;
   m_eeCorrNoise[2] = eeCorrNoise2;
 }

References cms::cuda::assert(), m_ebCorrNoise, and m_eeCorrNoise.

◆ ~EcalCoder()

EcalCoder::~EcalCoder ( )

virtual

dtor

Definition at line 36 of file EcalCoder.cc.

36 {}

Member Function Documentation

◆ analogToDigital()

void EcalCoder::analogToDigital	(	CLHEP::HepRandomEngine *	engine,
		const EcalSamples &	clf,
		EcalDataFrame &	df
	)		const

virtual

from EcalSamples to EcalDataFrame

Definition at line 53 of file EcalCoder.cc.

                                                                                                              {
   df.setSize(clf.size());
   encode(clf, df, engine);
   /*   std::cout<<"   **Id=" ;
    if( CaloGenericDetId( clf.id() ).isEB() )
    {
       std::cout<<EBDetId( clf.id() ) ;
    }
    else
    {
       std::cout<<EEDetId( clf.id() ) ;
    }
    std::cout<<", size="<<df.size();
    for( int i ( 0 ) ; i != df.size() ; ++i )
    {
       std::cout<<", "<<df[i];
    }
    std::cout<<std::endl ;*/
 }

References hgcalPerformanceValidation::df, encode(), and CaloTSamplesBase< Ttype >::size().

◆ encode()

void EcalCoder::encode	(	const EcalSamples &	ecalSamples,
		EcalDataFrame &	df,
		CLHEP::HepRandomEngine *	engine
	)		const

private

produce the pulse-shape

Definition at line 73 of file EcalCoder.cc.

                                                                                                             {
   assert(nullptr != m_peds);
  
   const unsigned int csize(ecalSamples.size());
  
   DetId detId = ecalSamples.id();
   double Emax = fullScaleEnergy(detId);
  
   //....initialisation
   if (ecalSamples[5] > 0.)
     LogDebug("EcalCoder") << "Input caloSample"
                           << "\n"
                           << ecalSamples;
  
   double LSB[NGAINS + 1];
   double pedestals[NGAINS + 1];
   double widths[NGAINS + 1];
   double gains[NGAINS + 1];
   int maxADC[NGAINS + 1];
   double trueRMS[NGAINS + 1];
  
   double icalconst = 1.;
   findIntercalibConstant(detId, icalconst);
  
   for (unsigned int igain(0); igain <= NGAINS; ++igain) {
     // fill in the pedestal and width
     findPedestal(detId, igain, pedestals[igain], widths[igain]);
  
     if (0 < igain)
       trueRMS[igain] = std::sqrt(widths[igain] * widths[igain] - 1. / 12.);
  
     // set nominal value first
     findGains(detId, gains);
  
     LSB[igain] = 0.;
     if (igain > 0)
       LSB[igain] = Emax / (MAXADC * gains[igain]);
     maxADC[igain] = ADCGAINSWITCH;  // saturation at 4080 for middle and high gains x6 & x12
     if (igain == NGAINS)
       maxADC[igain] = MAXADC;  // saturation at 4095 for low gain x1
   }
  
   CaloSamples noiseframe[] = {CaloSamples(detId, csize), CaloSamples(detId, csize), CaloSamples(detId, csize)};
  
   const Noisifier* noisy[3] = {
       (nullptr == m_eeCorrNoise[0] || EcalBarrel == detId.subdetId() ? m_ebCorrNoise[0] : m_eeCorrNoise[0]),
       (EcalBarrel == detId.subdetId() ? m_ebCorrNoise[1] : m_eeCorrNoise[1]),
       (EcalBarrel == detId.subdetId() ? m_ebCorrNoise[2] : m_eeCorrNoise[2])};
  
   if (m_addNoise) {
     noisy[0]->noisify(noiseframe[0], engine);  // high gain
     if (nullptr == noisy[1])
       noisy[0]->noisify(noiseframe[1], engine,
                         &noisy[0]->vecgau());  // med
     if (nullptr == noisy[2])
       noisy[0]->noisify(noiseframe[2], engine,
                         &noisy[0]->vecgau());  // low
   }
  
   //   std::cout << " intercal, LSBs, gains " << icalconst << " " << LSB[0] << " " << LSB[1] << " " << gains[0] << " " << gains[1] << " " << Emax <<  std::endl;
  
   int wait = 0;
   int gainId = 0;
   bool isSaturated = false;
  
   for (unsigned int i(0); i != csize; ++i) {
     bool done(false);
     int adc = MAXADC;
     if (0 == wait)
       gainId = 1;
  
     // see which gain bin it fits in
     int igain(gainId - 1);
     do {
       ++igain;
  
       //     if( igain != gainId ) std::cout <<"$$$$ Gain switch from " << gainId
       //                     <<" to "<< igain << std::endl ;
  
       if (1 != igain &&                     // not high gain
           m_addNoise &&                     // want to add noise
           nullptr != noisy[igain - 1] &&    // exists
           noiseframe[igain - 1].isBlank())  // not already done
       {
         noisy[igain - 1]->noisify(noiseframe[igain - 1], engine, &noisy[0]->vecgau());
         //std::cout<<"....noisifying gain level = "<<igain<<std::endl ;
       }
  
       double signal;
  
       if (!m_PreMix1) {
         // noiseframe filled with zeros if !m_addNoise
         const double asignal(pedestals[igain] + ecalSamples[i] / (LSB[igain] * icalconst) +
                              trueRMS[igain] * noiseframe[igain - 1][i]);
         signal = asignal;
       } else {  // Any changes made here must be reverse-engineered in EcalSignalGenerator!
  
         if (igain == 1) {
           const double asignal(ecalSamples[i] * 1000.);  // save low level info
           signal = asignal;
         } else if (igain == 2) {
           const double asignal(ecalSamples[i] / (LSB[1] * icalconst));
           signal = asignal;
         } else if (igain == 3) {  // bet that no pileup hit has an energy over Emax/2
           const double asignal(ecalSamples[i] / (LSB[2] * icalconst));
           signal = asignal;
         } else {  //not sure we ever get here at gain=0, but hit wil be saturated anyway
           const double asignal(ecalSamples[i] / (LSB[3] * icalconst));  // just calculate something
           signal = asignal;
         }
         // old version
         //const double asignal ( // no pedestals for pre-mixing
         //           ecalSamples[i] /( LSB[igain]*icalconst ) );
         //signal = asignal;
       }
  
       //     std::cout << " " << ecalSamples[i] << " " << noiseframe[igain-1][i] << std::endl;
  
       const int isignal(signal);
       const int tmpadc(signal - (double)isignal < 0.5 ? isignal : isignal + 1);
       // LogDebug("EcalCoder") << "DetId " << detId.rawId() << " gain " << igain << " caloSample "
       //                       <<  ecalSamples[i] << " pededstal " << pedestals[igain]
       //                       << " noise " << widths[igain] << " conversion factor " << LSB[igain]
       //                       << " result (ped,tmpadc)= " << ped << " " << tmpadc;
  
       if (tmpadc <= maxADC[igain]) {
         adc = tmpadc;
         done = true;
       }
     } while (!done && igain < 3);
  
     if (igain == 1) {
       wait = 0;
       gainId = igain;
     } else {
       if (igain == gainId)
         --wait;
       else {
         wait = 5;
         gainId = igain;
       }
     }
  
     // change the gain for saturation
     int storeGainId = gainId;
     if (gainId == 3 && adc == MAXADC) {
       storeGainId = 0;
       isSaturated = true;
     }
     // LogDebug("EcalCoder") << " Writing out frame " << i << " ADC = " << adc << " gainId = " << gainId << " storeGainId = " << storeGainId ;
  
     df.setSample(i, EcalMGPASample(adc, storeGainId));
   }
   // handle saturation properly according to IN-2007/056
   // N.B. - FIXME
   // still missing the possibility for a signal to pass the saturation threshold
   // again within the 5 subsequent samples (higher order effect).
  
   if (isSaturated) {
     for (unsigned int i = 0; i < ecalSamples.size(); ++i) {
       if (df.sample(i).gainId() == 0) {
         unsigned int hyst = i + 1 + 2;
         for (unsigned int j = i + 1; j < hyst && j < ecalSamples.size(); ++j) {
           df.setSample(j, EcalMGPASample(MAXADC, 0));
         }
       }
     }
   }
 }

References ecalLiteDTU::adc(), ADCGAINSWITCH, cms::cuda::assert(), hgcalPerformanceValidation::df, fileCollector::done, EcalBarrel, findGains(), findIntercalibConstant(), findPedestal(), fullScaleEnergy(), ecalLiteDTU::gainId(), mps_fire::i, CaloTSamplesBase< Ttype >::id(), CastorSimpleRecAlgoImpl::isSaturated(), dqmiolumiharvest::j, LogDebug, m_addNoise, m_ebCorrNoise, m_eeCorrNoise, m_peds, m_PreMix1, MAXADC, NGAINS, CorrelatedNoisifier< M >::noisify(), CaloTSamplesBase< Ttype >::size(), mathSSE::sqrt(), and SequenceTypes::wait().

Referenced by analogToDigital().

◆ findGains()

void EcalCoder::findGains	(	const DetId &	detId,
		double	theGains[]
	)		const

private

Definition at line 296 of file EcalCoder.cc.

                                                                   {
   /*
     EcalGainRatios::EcalGainRatioMap::const_iterator grit=m_gainRatios->getMap().find(detId.rawId());
     EcalMGPAGainRatio mgpa;
     if( grit!=m_gainRatios->getMap().end() ){
     mgpa = grit->second;
     Gains[0] = 0.;
     Gains[3] = 1.;
     Gains[2] = mgpa.gain6Over1() ;
     Gains[1] = Gains[2]*(mgpa.gain12Over6()) ;
     LogDebug("EcalCoder") << "Gains for " << detId.rawId() << "\n" << " 1 = " << Gains[1] << "\n" << " 2 = " << Gains[2] << "\n" << " 3 = " << Gains[3];
     } else {
     edm::LogError("EcalCoder") << "Could not find gain ratios for " << detId.rawId() << " among the " << m_gainRatios->getMap().size();
     }
   */
  
   EcalGainRatioMap::const_iterator grit = m_gainRatios->getMap().find(detId);
   Gains[0] = 0.;
   Gains[3] = 1.;
   Gains[2] = (*grit).gain6Over1();
   Gains[1] = Gains[2] * ((*grit).gain12Over6());
  
   if ((detId.subdetId() != EcalBarrel) && (detId.subdetId() != EcalEndcap)) {
     edm::LogError("EcalCoder") << "Could not find gain ratios for " << detId.rawId() << " among the "
                                << m_gainRatios->getMap().size();
   }
 }

References EcalBarrel, EcalEndcap, EcalCondObjectContainer< T >::find(), EcalCondObjectContainer< T >::getMap(), m_gainRatios, DetId::rawId(), EcalCondObjectContainer< T >::size(), and DetId::subdetId().

Referenced by encode().

◆ findIntercalibConstant()

void EcalCoder::findIntercalibConstant	(	const DetId &	detId,
		double &	icalconst
	)		const

private

Definition at line 324 of file EcalCoder.cc.

                                                                                   {
   EcalIntercalibConstantMC thisconst = 1.;
   // find intercalib constant for this xtal
   const EcalIntercalibConstantMCMap& icalMap = m_intercals->getMap();
   EcalIntercalibConstantMCMap::const_iterator icalit = icalMap.find(detId);
   if (icalit != icalMap.end()) {
     thisconst = (*icalit);
     if (icalconst == 0.) {
       thisconst = 1.;
     }
   } else {
     edm::LogError("EcalCoder") << "No intercalib const found for xtal " << detId.rawId()
                                << "! something wrong with EcalIntercalibConstants in your DB? ";
   }
   icalconst = thisconst;
 }

References EcalCondObjectContainer< T >::end(), EcalCondObjectContainer< T >::find(), EcalCondObjectContainer< T >::getMap(), m_intercals, and DetId::rawId().

Referenced by encode().

◆ findPedestal()

void EcalCoder::findPedestal	(	const DetId &	detId,
		int	gainId,
		double &	pedestal,
		double &	width
	)		const

private

not yet implemented

Definition at line 243 of file EcalCoder.cc.

                                                                                              {
   /*
      EcalPedestalsMapIterator mapItr 
      = m_peds->m_pedestals.find(detId.rawId());
      // should I care if it doesn't get found?
      if(mapItr == m_peds->m_pedestals.end()) {
      edm::LogError("EcalCoder") << "Could not find pedestal for " << detId.rawId() << " among the " << m_peds->m_pedestals.size();
      } else {
      EcalPedestals::Item item = mapItr->second;
    */
  
   /*   
     EcalPedestals::Item const & item = (*m_peds)(detId);
     ped = item.mean(gainId);
     width = item.rms(gainId);
    */
  
   EcalPedestalsMap::const_iterator itped = m_peds->getMap().find(detId);
   ped = (*itped).mean(gainId);
   width = (*itped).rms(gainId);
  
   if ((detId.subdetId() != EcalBarrel) && (detId.subdetId() != EcalEndcap)) {
     edm::LogError("EcalCoder") << "Could not find pedestal for " << detId.rawId() << " among the "
                                << m_peds->getMap().size();
   }
  
   /*
     switch(gainId) {
     case 0:
       ped = 0.;
       width = 0.;
     case 1:
       ped = item.mean_x12;
       width = item.rms_x12;
       break;
     case 2:
       ped = item.mean_x6;
       width = item.rms_x6;
       break;
     case 3:
       ped = item.mean_x1;
       width = item.rms_x1;
       break;
     default:
       edm::LogError("EcalCoder") << "Bad Pedestal " << gainId;
       break;
     }
   */
  
   LogDebug("EcalCoder") << "Pedestals for " << detId.rawId() << " gain range " << gainId << " : \n"
                         << "Mean = " << ped << " rms = " << width;
 }

References EcalBarrel, EcalEndcap, EcalCondObjectContainer< T >::find(), ecalLiteDTU::gainId(), EcalCondObjectContainer< T >::getMap(), LogDebug, m_peds, DetId::rawId(), EcalCondObjectContainer< T >::size(), DetId::subdetId(), and ApeEstimator_cff::width.

Referenced by encode().

◆ fullScaleEnergy()

double EcalCoder::fullScaleEnergy ( const DetId & did ) const

private

limit on the energy scale due to the electronics range

Definition at line 49 of file EcalCoder.cc.

                                                           {
   return detId.subdetId() == EcalBarrel ? m_maxEneEB : m_maxEneEE;
 }

References EcalBarrel, m_maxEneEB, m_maxEneEE, and DetId::subdetId().

Referenced by encode().

◆ setFullScaleEnergy()

void EcalCoder::setFullScaleEnergy	(	double	EBscale,
		double	EEscale
	)

Definition at line 38 of file EcalCoder.cc.

                                                                  {
   m_maxEneEB = EBscale;
   m_maxEneEE = EEscale;
 }

References m_maxEneEB, and m_maxEneEE.

◆ setGainRatios()

void EcalCoder::setGainRatios ( const EcalGainRatios * gainRatios )

Definition at line 45 of file EcalCoder.cc.

45 { m_gainRatios = gainRatios; }

References m_gainRatios.

◆ setIntercalibConstants()

void EcalCoder::setIntercalibConstants ( const EcalIntercalibConstantsMC * ical )

Definition at line 47 of file EcalCoder.cc.

47 { m_intercals = ical; }

References m_intercals.

◆ setPedestals()

void EcalCoder::setPedestals ( const EcalPedestals * pedestals )

can be fetched every event from the EventSetup

Definition at line 43 of file EcalCoder.cc.

43 { m_peds = pedestals; }

References m_peds.

Member Data Documentation

◆ m_addNoise

bool EcalCoder::m_addNoise

private

Definition at line 93 of file EcalCoder.h.

Referenced by encode().

◆ m_ebCorrNoise

const Noisifier* EcalCoder::m_ebCorrNoise[3]

private

Definition at line 96 of file EcalCoder.h.

Referenced by EcalCoder(), and encode().

◆ m_eeCorrNoise

const Noisifier* EcalCoder::m_eeCorrNoise[3]

private

Definition at line 97 of file EcalCoder.h.

Referenced by EcalCoder(), and encode().

◆ m_gainRatios

const EcalGainRatios* EcalCoder::m_gainRatios

private

Definition at line 86 of file EcalCoder.h.

Referenced by findGains(), and setGainRatios().

◆ m_intercals

const EcalIntercalibConstantsMC* EcalCoder::m_intercals

private

Definition at line 88 of file EcalCoder.h.

Referenced by findIntercalibConstant(), and setIntercalibConstants().

◆ m_maxEneEB

double EcalCoder::m_maxEneEB

private

Definition at line 90 of file EcalCoder.h.

Referenced by fullScaleEnergy(), and setFullScaleEnergy().

◆ m_maxEneEE

double EcalCoder::m_maxEneEE

private

Definition at line 91 of file EcalCoder.h.

Referenced by fullScaleEnergy(), and setFullScaleEnergy().

◆ m_peds

const EcalPedestals* EcalCoder::m_peds

private

Definition at line 84 of file EcalCoder.h.

Referenced by encode(), findPedestal(), and setPedestals().

◆ m_PreMix1

bool EcalCoder::m_PreMix1

private

Definition at line 94 of file EcalCoder.h.

Referenced by encode().

Public Types

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Member Typedef Documentation

◆ EcalSamples

◆ Noisifier

Member Enumeration Documentation

◆ anonymous enum

Constructor & Destructor Documentation

◆ EcalCoder()

◆ ~EcalCoder()

Member Function Documentation

◆ analogToDigital()

◆ encode()

◆ findGains()

◆ findIntercalibConstant()

◆ findPedestal()

◆ fullScaleEnergy()

◆ setFullScaleEnergy()

◆ setGainRatios()

◆ setIntercalibConstants()

◆ setPedestals()

Member Data Documentation

◆ m_addNoise

◆ m_ebCorrNoise

◆ m_eeCorrNoise

◆ m_gainRatios

◆ m_intercals

◆ m_maxEneEB

◆ m_maxEneEE

◆ m_peds

◆ m_PreMix1