CastorCalib/src/CastorPulseShapes.cc

0001 #include "CalibCalorimetry/CastorCalib/interface/CastorPulseShapes.h"
0002 #include <cmath>
0003
0004 CastorPulseShapes::CastorPulseShapes() {
0005   //  computeHPDShape(hpdShape_);
0006   computeCastorShape(castorShape_);
0007 }
0008
0009 void CastorPulseShapes::computeCastorShape(CastorPulseShapes::Shape& sh) {
0010   //  cout << endl << " ===== computeShapeHF  !!! " << endl << endl;
0011
0012   const float ts = 3.0;  // time constant in   t * exp(-(t/ts)**2)
0013
0014   // first create pulse shape over a range of time 0 ns to 255 ns in 1 ns steps
0015   int nbin = 256;
0016   sh.setNBin(nbin);
0017   std::vector<float> ntmp(nbin, 0.0);  //
0018
0019   int j;
0020   float norm;
0021
0022   // CASTOR SHAPE
0023   norm = 0.0;
0024   for (j = 0; j < 3 * ts && j < nbin; j++) {
0025     ntmp[j] = ((float)j) * exp(-((float)(j * j)) / (ts * ts));
0026     norm += ntmp[j];
0027   }
0028   // normalize pulse area to 1.0
0029   for (j = 0; j < 3 * ts && j < nbin; j++) {
0030     ntmp[j] /= norm;
0031
0032     //    cout << " nt [" << j << "] = " <<  ntmp[j] << endl;
0033     sh.setShapeBin(j, ntmp[j]);
0034   }
0035 }
0036
0037 CastorPulseShapes::Shape::Shape() {
0038   nbin_ = 0;
0039   tpeak_ = 0;
0040 }
0041
0042 void CastorPulseShapes::Shape::setNBin(int n) {
0043   nbin_ = n;
0044   shape_ = std::vector<float>(n, 0.0f);
0045 }
0046
0047 void CastorPulseShapes::Shape::setShapeBin(int i, float f) {
0048   if (i >= 0 && i < nbin_)
0049     shape_[i] = f;
0050 }
0051
0052 float CastorPulseShapes::Shape::operator()(double t) const {
0053   // shape is in 1 ns steps
0054   int i = (int)(t + 0.5);
0055   float rv = 0;
0056   if (i >= 0 && i < nbin_)
0057     rv = shape_[i];
0058   return rv;
0059 }
0060
0061 float CastorPulseShapes::Shape::at(double t) const {
0062   // shape is in 1 ns steps
0063   int i = (int)(t + 0.5);
0064   float rv = 0;
0065   if (i >= 0 && i < nbin_)
0066     rv = shape_[i];
0067   return rv;
0068 }
0069
0070 float CastorPulseShapes::Shape::integrate(double t1, double t2) const {
0071   static const float int_delta_ns = 0.05f;
0072   double intval = 0.0;
0073
0074   for (double t = t1; t < t2; t += int_delta_ns) {
0075     float loedge = at(t);
0076     float hiedge = at(t + int_delta_ns);
0077     intval += (loedge + hiedge) * int_delta_ns / 2.0;
0078   }
0079
0080   return (float)intval;
0081 }