EcalLaserAnalyzer/src/TAPDPulse.cc

0001 /*
0002  *  \class TAPDPulse
0003  *
0004  *  \author: Julie Malcles - CEA/Saclay
0005  */
0006
0007 #include <CalibCalorimetry/EcalLaserAnalyzer/interface/TAPDPulse.h>
0008 #include <TMath.h>
0009 #include <iostream>
0010 #include <cassert>
0011 using namespace std;
0012
0013 //ClassImp(TAPDPulse)
0014
0015 // Default Constructor...
0016 TAPDPulse::TAPDPulse() { init(10, 3, 1, 2, 2, 9, 3, 8, 0.4, 0.95, 0.8); }
0017
0018 // Constructor...
0019 TAPDPulse::TAPDPulse(int nsamples,
0020                      int presample,
0021                      int firstsample,
0022                      int lastsample,
0023                      int timingcutlow,
0024                      int timingcuthigh,
0025                      int timingquallow,
0026                      int timingqualhigh,
0027                      double ratiomincutlow,
0028                      double ratiomincuthigh,
0029                      double ratiomaxcutlow) {
0030   init(nsamples,
0031        presample,
0032        firstsample,
0033        lastsample,
0034        timingcutlow,
0035        timingcuthigh,
0036        timingquallow,
0037        timingqualhigh,
0038        ratiomincutlow,
0039        ratiomincuthigh,
0040        ratiomaxcutlow);
0041 }
0042
0043 // Destructor
0044 TAPDPulse::~TAPDPulse() {}
0045
0046 void TAPDPulse::init(int nsamples,
0047                      int presample,
0048                      int firstsample,
0049                      int lastsample,
0050                      int timingcutlow,
0051                      int timingcuthigh,
0052                      int timingquallow,
0053                      int timingqualhigh,
0054                      double ratiomincutlow,
0055                      double ratiomincuthigh,
0056                      double ratiomaxcutlow) {
0057   _nsamples = 10;
0058   assert(nsamples == _nsamples);
0059   assert(presample != 0);
0060   adc_ = new double[10];
0061
0062   _presample = presample;
0063   _firstsample = firstsample;
0064   _lastsample = lastsample;
0065
0066   _timingcutlow = timingcutlow;
0067   _timingcuthigh = timingcuthigh;
0068   _timingquallow = timingquallow;
0069   _timingqualhigh = timingqualhigh;
0070   _ratiomincutlow = ratiomincutlow;
0071   _ratiomincuthigh = ratiomincuthigh;
0072   _ratiomaxcutlow = ratiomaxcutlow;
0073
0074   for (int i = 0; i < _nsamples; i++) {
0075     adc_[i] = 0.0;
0076   }
0077
0078   adcMax_ = 0;
0079   iadcMax_ = 0;
0080   pedestal_ = 0;
0081
0082   isMaxFound_ = false;
0083   isPedCalc_ = false;
0084 }
0085
0086 bool TAPDPulse::setPulse(double *adc) {
0087   bool done = false;
0088   adc_ = adc;
0089   done = true;
0090   isMaxFound_ = false;
0091   isPedCalc_ = false;
0092   return done;
0093 }
0094 double TAPDPulse::getMax() {
0095   if (isMaxFound_)
0096     return adcMax_;
0097
0098   int iadcmax = 0;
0099   double adcmax = 0.0;
0100   for (int i = 0; i < _nsamples; i++) {
0101     if (adc_[i] > adcmax) {
0102       iadcmax = i;
0103       adcmax = adc_[i];
0104     }
0105   }
0106   iadcMax_ = iadcmax;
0107   adcMax_ = adcmax;
0108   return adcMax_;
0109 }
0110
0111 int TAPDPulse::getMaxSample() {
0112   if (!isMaxFound_)
0113     getMax();
0114   return iadcMax_;
0115 }
0116 double TAPDPulse::getDelta(int n1, int n2) {
0117   assert(n1 < _nsamples && n1 >= 0);
0118   assert(n2 < _nsamples && n2 >= 0);
0119
0120   double delta = adc_[n1] - adc_[n2];
0121   return delta;
0122 }
0123 double TAPDPulse::getRatio(int n1, int n2) {
0124   assert(n1 < _nsamples && n1 >= 0);
0125   assert(n2 < _nsamples && n2 >= 0);
0126
0127   double ped = 0;
0128   if (isPedCalc_)
0129     ped = pedestal_;
0130   else
0131     ped = adc_[0];
0132
0133   double ratio = (adc_[n1] - ped) / (adc_[n2] - ped);
0134   return ratio;
0135 }
0136
0137 bool TAPDPulse::isTimingOK() {
0138   bool ok = true;
0139   if (!isMaxFound_)
0140     getMax();
0141   if (iadcMax_ <= _timingcutlow || iadcMax_ >= _timingcuthigh)
0142     ok = false;
0143   return ok;
0144 }
0145 bool TAPDPulse::isTimingQualOK() {
0146   bool ok = true;
0147   if (!isMaxFound_)
0148     getMax();
0149   if (iadcMax_ <= _timingquallow || iadcMax_ >= _timingqualhigh)
0150     ok = false;
0151   return ok;
0152 }
0153
0154 bool TAPDPulse::areFitSamplesOK() {
0155   bool ok = true;
0156   if (!isMaxFound_)
0157     getMax();
0158   if ((iadcMax_ - _firstsample) < _presample || (iadcMax_ + _lastsample) > _nsamples - 1)
0159     ok = false;
0160   return ok;
0161 }
0162 bool TAPDPulse::isPulseOK() {
0163   bool okSamples = areFitSamplesOK();
0164   bool okTiming = isTimingOK();
0165   bool okPulse = arePulseRatioOK();
0166
0167   bool ok = (okSamples && okTiming && okPulse);
0168
0169   return ok;
0170 }
0171 bool TAPDPulse::arePulseRatioOK() {
0172   bool ok = true;
0173
0174   if (!isMaxFound_)
0175     getMax();
0176   if (iadcMax_ < 1 || iadcMax_ >= _nsamples - 1)
0177     return false;
0178
0179   double ratioNm1 = getRatio(iadcMax_ - 1, iadcMax_);
0180   double ratioNp1 = getRatio(iadcMax_ + 1, iadcMax_);
0181   double ratioMax = TMath::Max(ratioNm1, ratioNp1);
0182   double ratioMin = TMath::Min(ratioNm1, ratioNp1);
0183
0184   if (ratioMax < _ratiomaxcutlow)
0185     ok = false;
0186   if (ratioMin < _ratiomincutlow || ratioMin > _ratiomincuthigh)
0187     ok = false;
0188
0189   return ok;
0190 }
0191 bool TAPDPulse::isPulseRatioMaxOK() {
0192   bool ok = true;
0193
0194   if (!isMaxFound_)
0195     getMax();
0196   if (iadcMax_ < 1 || iadcMax_ >= _nsamples - 1)
0197     return false;
0198
0199   double ratioNm1 = getRatio(iadcMax_ - 1, iadcMax_);
0200   double ratioNp1 = getRatio(iadcMax_ + 1, iadcMax_);
0201   double ratioMax = TMath::Max(ratioNm1, ratioNp1);
0202
0203   if (ratioMax < _ratiomaxcutlow)
0204     ok = false;
0205   return ok;
0206 }
0207 bool TAPDPulse::isPulseRatioMinOK() {
0208   bool ok = true;
0209
0210   if (!isMaxFound_)
0211     getMax();
0212   if (iadcMax_ < 1 || iadcMax_ >= _nsamples - 1)
0213     return false;
0214
0215   double ratioNm1 = getRatio(iadcMax_ - 1, iadcMax_);
0216   double ratioNp1 = getRatio(iadcMax_ + 1, iadcMax_);
0217   double ratioMin = TMath::Min(ratioNm1, ratioNp1);
0218
0219   if (ratioMin < _ratiomincutlow || ratioMin > _ratiomincuthigh)
0220     ok = false;
0221   return ok;
0222 }
0223
0224 double TAPDPulse::getPedestal() {
0225   if (isPedCalc_)
0226     return pedestal_;
0227   double ped = 0;
0228   for (int i = 0; i < _presample; i++) {
0229     ped += adc_[i];
0230   }
0231   ped /= double(_presample);
0232   pedestal_ = ped;
0233   isPedCalc_ = true;
0234   return pedestal_;
0235 }
0236
0237 double *TAPDPulse::getAdcWithoutPedestal() {
0238   double ped;
0239   if (!isPedCalc_)
0240     ped = getPedestal();
0241   else
0242     ped = pedestal_;
0243
0244   double *adcNoPed = new double[10];
0245   for (int i = 0; i < _nsamples; i++) {
0246     adcNoPed[i] = adc_[i] - ped;
0247   }
0248   return adcNoPed;
0249 }
0250
0251 void TAPDPulse::setPresamples(int presample) {
0252   isPedCalc_ = false;
0253   _presample = presample;
0254 }