CSCFindPeakTime Class Reference

#include <CSCFindPeakTime.h>

Public Member Functions
float	averageTime (int tmax, const float *adc)
	Weighted average of time bins. More...
	CSCFindPeakTime (const edm::ParameterSet &ps)
void	fivePoleFitCharge (int tmax, const float *adc, const float &t_zero, const float &t_peak, std::vector< float > &adcsFit)
float	fivePoleFitTime (int tmax, const float *adc, float t_peak)
float	parabolaFitTime (int tmax, const float *adc)
	Parabolic fit to three time bins centered on maximum. More...
float	peakTime (int tmax, const float *adc, float t_peak)
	Basic result of this class. More...
	~CSCFindPeakTime ()

Private Attributes
bool	useAverageTime
bool	useFivePoleFit
bool	useParabolaFit

Detailed Description

This is CSCFindPeakTime

Used to provide improved estimate of SCA peak time.

Definition at line 12 of file CSCFindPeakTime.h.

Constructor & Destructor Documentation

CSCFindPeakTime()

CSCFindPeakTime::CSCFindPeakTime ( const edm::ParameterSet &  ps )

explicit

Definition at line 8 of file CSCFindPeakTime.cc.

References edm::ParameterSet::getParameter(), LogTrace, useAverageTime, useFivePoleFit, and useParabolaFit.

    : useAverageTime(false), useParabolaFit(false), useFivePoleFit(false) {
  useAverageTime = ps.getParameter<bool>("UseAverageTime");
  useParabolaFit = ps.getParameter<bool>("UseParabolaFit");
  useFivePoleFit = ps.getParameter<bool>("UseFivePoleFit");
  LogTrace("CSCRecHit") << "[CSCFindPeakTime] useAverageTime=" << useAverageTime
                        << ", useParabolaFit=" << useParabolaFit << ", useFivePoleFit=" << useFivePoleFit;
}

~CSCFindPeakTime()

CSCFindPeakTime::~CSCFindPeakTime ( )

inline

Definition at line 16 of file CSCFindPeakTime.h.

{}

Member Function Documentation

averageTime()

float CSCFindPeakTime::averageTime	(	int	tmax,
		const float *	adc
	)

Weighted average of time bins.

Definition at line 30 of file CSCFindPeakTime.cc.

References gpuClustering::adc, ALCARECOEcalPhiSym_cff::float, mps_fire::i, and tmax.

Referenced by peakTime().

                                                             {
  float sum = 0.;
  float sumt = 0.;
  for (size_t i = 0; i < 4; ++i) {
    sum += adc[i];
    sumt += adc[i] * float(tmax - 1 + i);
  }
  return sumt / sum * 50.;  //@@ in ns. May be some bin width offset things to handle here?
}

fivePoleFitCharge()

void CSCFindPeakTime::fivePoleFitCharge	(	int	tmax,
		const float *	adc,
		const float &	t_zero,
		const float &	t_peak,
		std::vector< float > &	adcsFit
	)

Integrated charge after fivePoleFitTime

Definition at line 147 of file CSCFindPeakTime.cc.

References gpuClustering::adc, JetChargeProducer_cfi::exp, ALCARECOEcalPhiSym_cff::float, mps_fire::i, dqmiolumiharvest::j, N, submitPVValidationJobs::t, tmax, x, and y.

                                                                                                     {
  //@@ This code can certainly be replaced by fivePoleFitTime above, but I haven't time to do that now (Tim).

  float p0 = 4.f / t_peak;
  float tt0 = t_zero;
  int n_fit = 4;
  if (tmax == 6)
    n_fit = 3;

  float tb[4], y[4];
  for (int t = 0; t < 4; ++t) {
    tb[t] = float(tmax + t - 1) * 50.f;
    y[t] = adc[t];
  }

  // Find the normalization factor for the function
  float x[4];
  float sx2 = 0.f;
  float sxy = 0.f;
  for (int j = 0; j < n_fit; ++j) {
    float t = tb[j];
    x[j] = ((t - tt0) * (t - tt0)) * ((t - tt0) * (t - tt0)) * std::exp(-p0 * (t - tt0));
    sx2 = sx2 + x[j] * x[j];
    sxy = sxy + x[j] * y[j];
  }
  float N = sxy / sx2;

  // Now compute charge for a given t  --> only need charges at: t_peak-50, t_peak and t_peak+50
  for (int i = 0; i < 3; ++i) {
    float t = t_peak + float(i - 1) * 50.f;
    float q_fitted = float(N) * ((t - tt0) * (t - tt0)) * ((t - tt0) * (t - tt0)) * std::exp(-p0 * (t - tt0));
    adcsFit.push_back(q_fitted);
  }
  return;
}

fivePoleFitTime()

float CSCFindPeakTime::fivePoleFitTime	(	int	tmax,
		const float *	adc,
		float	t_peak
	)

Based on RecoLocalMuon/CSCStandAlone/interface/PulseTime.h by S. Durkin, and ported by D. Fortin. Comments updated by Tim Cox Apr 2009.

The SCA pulse shape should be representable by a function
N*(p0^2/256/exp(-4)) * (t-t0)^4 * exp( -p0*(t-t0) )

Rather than do a full fit with varying peak time too, assume the peak time is fixed to 133 nsec w.r.t. start time, t0, and fit for t0. The fit uses a binary search in t0, and at each step calculates the overall normalization factor between the function and the SCA pulse height as a least-squares fit over the 4 time bins tmax -1, tmax, tmax+1, tmax+2

Note: t0peak =4/p0 = 133 nsec, and adc[0] is arbitrarily defined a time of 0.0 nsec.

Definition at line 64 of file CSCFindPeakTime.cc.

References gpuClustering::adc, isoTrack_cff::chi2, ALPAKA_ACCELERATOR_NAMESPACE::brokenline::constexpr(), JetChargeProducer_cfi::exp, f, hcaldqm::quantity::fN, dqmiolumiharvest::j, FrontierCondition_GT_autoExpress_cfi::t0, hcalRecHitTable_cff::time, tmax, and x.

Referenced by peakTime().

                                                                               {
  // Input is
  // tmax   = bin# 0-7 containing max SCA pulse height
  // adc    = 4-dim array containing SCA pulse heights in bins tmax-1 to tmax+2
  // t_peak = input estimate for SCA peak time

  // Returned value is improved (we hope) estimate for SCA peak time

  // Algorithm is to fit five-pole Semi-Gaussian function for start time of SCA pulse, t0
  // (The SCA peak is assumed to be 133 ns from t0.)
  // Note that t^4 in time domain corresponds to 1/t^5 in frequency domain (that's the 5 poles).

  // Initialize parameters to sensible (?) values

  float t0 = 0.f;
  constexpr float t0peak = 133.f;  // this is offset of peak from start time t0
  constexpr float p0 = 4.f / t0peak;

  // Require that tmax is in range 2-6 of bins the eight SCA time bins 0-7
  // (Bins 0, 1 used for dynamic ped)

  if (tmax < 2 || tmax > 6)
    return t_peak;  //@@ Just return the input value

  // Set up time bins to match adc[4] input

  float tb[4];
  for (int time = 0; time < 4; ++time) {
    tb[time] = (tmax + time - 1) * 50.f;
  }

  // How many time bins are we fitting?

  int n_fit = 4;
  if (tmax == 6)
    n_fit = 3;

  float chi_min = 1.e10f;
  float chi_last = 1.e10f;
  float tt0 = 0.f;
  float chi2 = 0.f;
  float del_t = 100.f;

  float x[4];
  float sx2 = 0.f;
  float sxy = 0.f;
  float fN = 0.f;

  while (del_t > 1.f) {
    sx2 = 0.f;
    sxy = 0.f;

    for (int j = 0; j < n_fit; ++j) {
      float tdif = tb[j] - tt0;
      x[j] = (tdif * tdif) * (tdif * tdif) * std::exp(-p0 * tdif);
      sx2 += x[j] * x[j];
      sxy += x[j] * adc[j];
    }
    fN = sxy / sx2;  // least squares fit over time bins i to adc[i] = fN * fivePoleFunction[i]

    // Compute chi^2
    chi2 = 0.0;
    for (int j = 0; j < n_fit; ++j)
      chi2 += (adc[j] - fN * x[j]) * (adc[j] - fN * x[j]);

    // Test on chi^2 to decide what to do
    if (chi_last > chi2) {
      if (chi2 < chi_min) {
        t0 = tt0;
      }
      chi_last = chi2;
      tt0 = tt0 + del_t;
    } else {
      tt0 = tt0 - 2.f * del_t;
      del_t = 0.5 * del_t;
      tt0 = tt0 + del_t;
      chi_last = 1.0e10f;
    }
  }

  return t0 + t0peak;
}

parabolaFitTime()

float CSCFindPeakTime::parabolaFitTime	(	int	tmax,
		const float *	adc
	)

Parabolic fit to three time bins centered on maximum.

Definition at line 40 of file CSCFindPeakTime.cc.

References gpuClustering::adc, f, LogTrace, tmax, testProducerWithPsetDescEmpty_cfi::y1, and testProducerWithPsetDescEmpty_cfi::y2.

Referenced by peakTime().

                                                                 {
  // 3-point parabolic fit, from Andy Kubik

  // We calculate offset to tmax by finding the peak of a parabola through three points
  float tpeak = tmax;
  float tcorr = 0;

  // By construction, input array adc is for bins tmax-1 to tmax+2
  float y1 = adc[0];
  float y2 = adc[1];
  float y3 = adc[2];

  // Checked and simplified... Tim Cox 08-Apr-2009
  // Denominator is not zero unless we fed in nonsense values with y2 not the peak!
  if ((y1 + y3) < 2.f * y2)
    tcorr = 0.5f * (y1 - y3) / (y1 - 2. * y2 + y3);
  tpeak += tcorr;

  LogTrace("CSCFindPeakTime") << "[CSCFindPeakTime] tmax=" << tmax << ", parabolic peak time is tmax+" << tcorr
                              << " bins, or " << tpeak * 50.f << " ns";

  return tpeak * 50.f;  // convert to ns.
}

peakTime()

float CSCFindPeakTime::peakTime	(	int	tmax,
		const float *	adc,
		float	t_peak
	)

Basic result of this class.

Definition at line 17 of file CSCFindPeakTime.cc.

References gpuClustering::adc, averageTime(), fivePoleFitTime(), parabolaFitTime(), tmax, useAverageTime, useFivePoleFit, and useParabolaFit.

                                                                        {
  if (useAverageTime) {
    return averageTime(tmax, adc);
  } else if (useParabolaFit) {
    return parabolaFitTime(tmax, adc);
  } else if (useFivePoleFit) {
    return fivePoleFitTime(tmax, adc, t_peak);
  } else {
    // return something, anyway.. may as well be average
    return averageTime(tmax, adc);
  }
}

Member Data Documentation

useAverageTime

bool CSCFindPeakTime::useAverageTime

private

Definition at line 53 of file CSCFindPeakTime.h.

Referenced by CSCFindPeakTime(), and peakTime().

useFivePoleFit

bool CSCFindPeakTime::useFivePoleFit

private

Definition at line 55 of file CSCFindPeakTime.h.

Referenced by CSCFindPeakTime(), and peakTime().

useParabolaFit

bool CSCFindPeakTime::useParabolaFit

private

Definition at line 54 of file CSCFindPeakTime.h.

Referenced by CSCFindPeakTime(), and peakTime().