EcalCoder Class Reference

#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

typedef CorrelatedNoisifier<EcalCorrMatrix> EcalCoder::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.

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

    : 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;
}

~EcalCoder()

EcalCoder::~EcalCoder ( )

virtual

dtor

Definition at line 36 of file EcalCoder.cc.

{}

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.

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

                                                                                                             {
  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 ;*/
}

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.

References gpuClustering::adc, ADCGAINSWITCH, cms::cuda::assert(), hcalRecHitTable_cff::detId, hgcalPerformanceValidation::df, fileCollector::done, EcalBarrel, findGains(), findIntercalibConstant(), findPedestal(), fullScaleEnergy(), ecalLiteDTU::gainId(), ecalph2::gains, 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().

                                                                                                            {
  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));
        }
      }
    }
  }
}

findGains()

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

private

Definition at line 296 of file EcalCoder.cc.

References hcalRecHitTable_cff::detId, EcalBarrel, EcalEndcap, EcalCondObjectContainer< T >::find(), EcalCondObjectContainer< T >::getMap(), m_gainRatios, and EcalCondObjectContainer< T >::size().

Referenced by encode().

                                                                  {
  /*
    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();
  }
}

findIntercalibConstant()

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

private

Definition at line 324 of file EcalCoder.cc.

References hcalRecHitTable_cff::detId, EcalCondObjectContainer< T >::end(), EcalCondObjectContainer< T >::find(), EcalCondObjectContainer< T >::getMap(), and m_intercals.

Referenced by encode().

                                                                                  {
  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;
}

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.

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

Referenced by encode().

                                                                                             {
  /*
     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;
}

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.

References hcalRecHitTable_cff::detId, EcalBarrel, m_maxEneEB, and m_maxEneEE.

Referenced by encode().

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

setFullScaleEnergy()

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

Definition at line 38 of file EcalCoder.cc.

References m_maxEneEB, and m_maxEneEE.

                                                                 {
  m_maxEneEB = EBscale;
  m_maxEneEE = EEscale;
}

setGainRatios()

void EcalCoder::setGainRatios

(

const EcalGainRatios *

gainRatios

)

Definition at line 45 of file EcalCoder.cc.

References m_gainRatios.

{ m_gainRatios = gainRatios; }

setIntercalibConstants()

void EcalCoder::setIntercalibConstants

(

const EcalIntercalibConstantsMC *

ical

)

Definition at line 47 of file EcalCoder.cc.

References m_intercals.

{ m_intercals = ical; }

setPedestals()

void EcalCoder::setPedestals

(

const EcalPedestals *

pedestals

)

can be fetched every event from the EventSetup

Definition at line 43 of file EcalCoder.cc.

References m_peds.

{ m_peds = pedestals; }

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().