CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/CalibCalorimetry/CastorCalib/src/CastorPulseContainmentCorrection.cc

Go to the documentation of this file.
00001 #include "CalibCalorimetry/CastorCalib/interface/CastorPulseContainmentCorrection.h"
00002 #include "CalibCalorimetry/CastorCalib/interface/CastorTimeSlew.h"
00003 #include "CalibCalorimetry/CastorCalib/interface/CastorPulseShapes.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 #include <algorithm> // for "max","min"
00006 #include <math.h>
00007 #include <iostream>
00008 
00009 #include "CalibCalorimetry/HcalAlgos/interface/genlkupmap.h"
00010 
00011 
00012 template<>
00013 RecoFCcorFactorAlgo<CastorPulseShapes::Shape>::RecoFCcorFactorAlgo(int num_samples, double fixedphase_ns)
00014 {
00015   fixedphasens_ = fixedphase_ns;
00016   CastorPulseShapes shapes;
00017   shape_=shapes.castorShape();
00018   const int binsize_ns = 25;
00019 
00020   // First set up controlling parameters for calculating the correction factor:
00021   // Integration window size...
00022   //
00023   integrationwindowns_ = (double)(binsize_ns*num_samples);
00024 
00025   // First find the point at which time bin "1" exceeds time bin "0",
00026   // and call that point "time 0".
00027   //
00028   for (int shift_ns=0; shift_ns<binsize_ns; shift_ns++) {
00029 
00030     // Digitize by integrating to find all time sample
00031     // bin values for this shift.
00032     //
00033     double tmin    = -(double)shift_ns;
00034     double bin0val = (double)shape_.integrate(tmin, tmin+binsize_ns);
00035     double bin1val = (double)shape_.integrate(tmin+binsize_ns, tmin+2*binsize_ns);
00036 
00037 #if 0
00038     char s[80];
00039     sprintf (s, "%7.3f %8.5f %8.5f\n", tmin, bin0val, bin1val);
00040     cout << s;
00041 #endif
00042 
00043     if (bin1val > bin0val) {
00044       time0shiftns_ = shift_ns;
00045       break;
00046     }
00047   }
00048 
00049 #if 0
00050   cout << "time0shiftns_ = " << time0shiftns_ << endl;
00051 #endif
00052 }
00053 
00054 template <>
00055 std::pair<double,double> RecoFCcorFactorAlgo<CastorPulseShapes::Shape>::calcpair(double truefc)
00056 {
00057   double timeslew_ns = CastorTimeSlew::delay(std::max(0.0,(double)truefc),
00058                                            CastorTimeSlew::Medium);
00059   double shift_ns  = fixedphasens_ - time0shiftns_ + timeslew_ns;
00060 
00061   double tmin      = -shift_ns;
00062   double tmax      = tmin+integrationwindowns_;
00063 
00064   double integral  = shape_.integrate( tmin, tmax );
00065   double corfactor = 1.0/integral;
00066   double recofc    = (double)truefc * integral;
00067 
00068 #if 0
00069   char s[80];
00070   sprintf (s, "%8.2f %8.4f %8.4f %8.5f %8.5f %8.5f ",
00071            truefc, tmin, tmax, integral, corfactor, recofc);
00072   cout << s;
00073 #endif
00074 
00075   std::pair<double,double> thepair(recofc,corfactor);
00076   return thepair;
00077 }
00078 
00080 //
00081 CastorPulseContainmentCorrection::CastorPulseContainmentCorrection(int num_samples,
00082                                                                float fixedphase_ns,
00083                                                                float max_fracerror)
00084 {
00085   RecoFCcorFactorAlgo<CastorPulseShapes::Shape> corFalgo(num_samples, (double)fixedphase_ns);
00086 
00087   // Generate lookup map for the correction function, never exceeding
00088   // a maximum fractional error for lookups.
00089   //
00090   genlkupmap< RecoFCcorFactorAlgo<CastorPulseShapes::Shape> > (1.0, 5000.0f,  // generation domain
00091                                      max_fracerror,     // maximum fractional error
00092                                      1.0,   // min_xstep = minimum true fC increment
00093                                      corFalgo,
00094                                      mCorFactors_);     // return lookup map
00095 }
00096 
00097 
00098 double CastorPulseContainmentCorrection::getCorrection(double fc_ampl) const
00099 {
00100   double correction;
00101 
00102   std::map<double,double>::const_iterator fcupper,fclower;
00103 
00104   fcupper = mCorFactors_.upper_bound(fc_ampl);
00105   fclower = fcupper;
00106   fclower--;
00107 
00108   if (fcupper == mCorFactors_.end()) {
00109     correction = fclower->second;
00110   }
00111   else if (fcupper == mCorFactors_.begin()) {
00112     correction = fcupper->second;
00113   }
00114   else {
00115     if (fabs(fclower->first - fc_ampl) <
00116         fabs(fcupper->first - fc_ampl) )
00117       correction = fclower->second;
00118     else
00119       correction = fcupper->second;
00120   }
00121 
00122 #if 0
00123   char s[80];
00124   sprintf (s, "%7.1f (%8.5f %8.5f) (%8.5f %8.5f) %8.5f",
00125            fc_ampl,
00126            fclower->first, fclower->second,
00127            fcupper->first, fcupper->second,
00128            correction);
00129   cout << s << endl;
00130 #endif
00131 
00132   return correction;
00133 }