HcalTDC.cc

Go to the documentation of this file.

#include "SimCalorimetry/HcalSimAlgos/interface/HcalTDC.h"
#include "CalibFormats/HcalObjects/interface/HcalCalibrations.h"
#include "CalibFormats/HcalObjects/interface/HcalDbService.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"

#include "CLHEP/Random/RandGaussQ.h"

HcalTDC::HcalTDC(unsigned int thresholdDAC) : theTDCParameters(), 
                          theDbService(nullptr),
                          theDAC(thresholdDAC),
                          lsb(3.74) {}

HcalTDC::~HcalTDC() {
}

template <class Digi>
void HcalTDC::timing(const CaloSamples& lf, Digi& digi, CLHEP::HepRandomEngine* engine) const {

  float const TDC_Threshold(getThreshold(digi.id(), engine));
  float const TDC_Threshold_hyst(TDC_Threshold);
  bool risingReady(true), fallingReady(false);
  int tdcBins = theTDCParameters.nbins();
  // start with a loop over 10 samples
  bool hasTDCValues=true;
  if (lf.preciseSize()==0 ) hasTDCValues=false;

  // "AC" hysteresis from Tom:  There cannot be a second transition
  // within the QIE within 3ns
  int const rising_ac_hysteresis(5);
  int lastRisingEdge(0);

  // if (hasTDCValues)
  //   std::cout << digi.id() 
  //          << " threshold: " << TDC_Threshold
  //          << '\n';
  for (int ibin = 0; ibin < lf.size(); ++ibin) {
    /*
    If in a given 25ns bunch/time sample, the pulse is above
    TDC_Thresh already, then TDC_RisingEdge=0 if it was low in the
    last precision bin on the previous bunch crossing, otherwise,
    TDC_RisingEdge=63 if the pulse never crosses the threshold
    having started off, then the special code is 62 and then
    one can still have a TDC_FallingEdge that is valid.  If the pulse
    never falls below threshold having started above threshold (or
    goes above threshold in the bunch crossing and doesn't come down),
    then TDC_FallingEdge=.  If the pulse never went above
    threshold, then TDC_RisingEdge=63 and
    TDC_FallingEdge=62.
    */
    // special codes
    int TDC_RisingEdge = (risingReady) ? 
      theTDCParameters.noTransitionCode() :
      theTDCParameters.alreadyTransitionCode();
    int TDC_FallingEdge = (fallingReady) ? 
      theTDCParameters.alreadyTransitionCode() :
      theTDCParameters.noTransitionCode();
    int preciseBegin = ibin * tdcBins;
    int preciseEnd = preciseBegin + tdcBins;
    if ( hasTDCValues) {
      // std::cout << " alreadyOn: " << alreadyOn << '\n';
      for(int i = preciseBegin; i < preciseEnd; ++i) { 
      //    std::cout << " preciseBin: " << i
      //                << " preciseAt(i): " << lf.preciseAt(i);
    if( (fallingReady) && (i%3 == 0) && (lf.preciseAt(i) < TDC_Threshold)) {
      fallingReady = false;
      TDC_FallingEdge = i-preciseBegin;
      // std::cout << " falling ";
    }
    if ((risingReady) && (lf.preciseAt(i) > TDC_Threshold)) {
      risingReady = false;
      TDC_RisingEdge = i-preciseBegin;
      lastRisingEdge = 0;
      // std::cout << " rising ";
    }
    // std::cout << '\n';

    //This is the place for hysteresis code.  Need to to arm the
    //tdc by setting the risingReady or fallingReady to true.

    if ((!fallingReady) && (lf.preciseAt(i) > TDC_Threshold)) {
      fallingReady = true;
      if (TDC_FallingEdge == theTDCParameters.alreadyTransitionCode())
        TDC_FallingEdge = theTDCParameters.noTransitionCode();
    }
    if ((!risingReady) && (lastRisingEdge > rising_ac_hysteresis) &&
        (lf.preciseAt(i) < TDC_Threshold_hyst)) {
      risingReady = true;
      if (TDC_RisingEdge == theTDCParameters.alreadyTransitionCode())
        TDC_RisingEdge = theTDCParameters.noTransitionCode();
    }
    ++lastRisingEdge;
      }
    }
    // change packing to allow for special codes
    int packedTDC = TDC_RisingEdge + (tdcBins*2) * TDC_FallingEdge;
    digi.setSample(ibin, digi.adc(ibin), packedTDC, true);
    // if ( hasTDCValues) {
    //   std::cout << " sample: " << ibin 
    //      << " adc: " << digi.adc(ibin)
    //      << " fC: " << digi[ibin].nominal_fC()
    //      << " risingEdge: " << TDC_RisingEdge
    //      << " fallingEdge: " << TDC_FallingEdge
    //      << " packedTDC: " << packedTDC
    //      << std::endl;
    // }
  } // loop over bunch crossing bins
}

double HcalTDC::getThreshold(const HcalGenericDetId & detId, CLHEP::HepRandomEngine* engine) const {
  // subtract off pedestal and noise once
  double pedestal = theDbService->getHcalCalibrations(detId).pedestal(0);
  double pedestalWidth = theDbService->getHcalCalibrationWidths(detId).pedestal(0);
  // here the TDCthreshold_ is a multiple of the least significant bit 
  // for the TDC portion of the QIE.  The nominal reference is 40 uA which is
  // divided by an 8 bit DAC.  This means the least significant bit is 0.135 uA
  // in the TDC circuit or 3.74 uA at the input current.

  // the pedestal is assumed to be evenly distributed in time with some
  // random variation.

  return lsb*theDAC - CLHEP::RandGaussQ::shoot(engine, pedestal,  pedestalWidth)/theTDCParameters.deltaT()/theTDCParameters.nbins();
}

double HcalTDC::getHysteresisThreshold(double nominal) const {
  // the TDC will "re-arm" when it crosses the threshold minus 2.5 mV.
  // the lsb is 44kOhm * (3.74 uA / 24) = 6.86 mV.  2.5 mV is 0.36 lsb.
  // with the this means that the hysteresis it isn't far from the nominal 
  // threshold

  return nominal - 0.365*lsb;
}

void HcalTDC::setDbService(const HcalDbService * service) {
  theDbService = service;
}