CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/SimCalorimetry/HcalSimAlgos/src/HcalSiPMHitResponse.cc

Go to the documentation of this file.
00001 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSiPMHitResponse.h"
00002 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSiPM.h"
00003 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSiPMRecovery.h"
00004 #include "SimCalorimetry/HcalSimAlgos/interface/HcalSimParameters.h"
00005 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 #include "SimCalorimetry/CaloSimAlgos/interface/CaloVHitFilter.h"
00008 #include "SimCalorimetry/CaloSimAlgos/interface/CaloVSimParameterMap.h"
00009 #include "SimCalorimetry/CaloSimAlgos/interface/CaloSimParameters.h"
00010 #include "SimCalorimetry/CaloSimAlgos/interface/CaloVHitCorrection.h"
00011 #include "SimCalorimetry/CaloSimAlgos/interface/CaloVShape.h"
00012 
00013 #include <set>
00014 #include <map>
00015 
00016 class PCaloHitCompareTimes {
00017 public:
00018   bool operator()(const PCaloHit * a, 
00019                   const PCaloHit * b) const {
00020     return a->time()<b->time();
00021   }
00022 };
00023 
00024 HcalSiPMHitResponse::HcalSiPMHitResponse(const CaloVSimParameterMap * parameterMap,
00025                                          const CaloVShape * shape) :
00026   CaloHitResponse(parameterMap, shape), theSiPM(0), theRecoveryTime(250.) {
00027   theSiPM = new HcalSiPM(14000);
00028 }
00029 
00030 HcalSiPMHitResponse::~HcalSiPMHitResponse() {
00031   if (theSiPM)
00032     delete theSiPM;
00033 }
00034 
00035 void HcalSiPMHitResponse::run(MixCollection<PCaloHit> & hits) {
00036 
00037   typedef std::multiset <const PCaloHit *, PCaloHitCompareTimes> SortedHitSet;
00038 
00039   std::map< DetId, SortedHitSet > sortedhits;
00040   for (MixCollection<PCaloHit>::MixItr hitItr = hits.begin();
00041        hitItr != hits.end(); ++hitItr) {
00042     if (!((hitItr.bunch() < theMinBunch) || (hitItr.bunch() > theMaxBunch)) &&
00043         !(isnan(hitItr->time())) &&
00044         ((theHitFilter == 0) || (theHitFilter->accepts(*hitItr)))) {
00045       DetId id(hitItr->id());
00046       if (sortedhits.find(id)==sortedhits.end())
00047         sortedhits.insert(std::pair<DetId, SortedHitSet>(id, SortedHitSet()));
00048       sortedhits[id].insert(&(*hitItr));
00049     }
00050   }
00051   int pixelIntegral, oldIntegral;
00052   HcalSiPMRecovery pixelHistory(theRecoveryTime);
00053   const PCaloHit * hit;
00054   for (std::map<DetId, SortedHitSet>::iterator i = sortedhits.begin(); 
00055        i!=sortedhits.end(); ++i) {
00056     pixelHistory.clearHistory();
00057     for (SortedHitSet::iterator itr = i->second.begin(); 
00058          itr != i->second.end(); ++itr) {
00059       hit = *itr;
00060       pixelIntegral = pixelHistory.getIntegral(hit->time());
00061       oldIntegral = pixelIntegral;
00062       CaloSamples signal(makeSiPMSignal(*hit, pixelIntegral));
00063       pixelHistory.addToHistory(hit->time(), pixelIntegral-oldIntegral);
00064       add(signal);
00065     }
00066   }
00067 }
00068 
00069 
00070 void HcalSiPMHitResponse::setRandomEngine(CLHEP::HepRandomEngine & engine)
00071 {
00072   theSiPM->initRandomEngine(engine);
00073   CaloHitResponse::setRandomEngine(engine);
00074 }
00075 
00076 
00077 CaloSamples HcalSiPMHitResponse::makeSiPMSignal(const PCaloHit & inHit, 
00078                                                 int & integral ) const {
00079   PCaloHit hit = inHit;
00080   if (theHitCorrection != 0)
00081     theHitCorrection->correct(hit);
00082 
00083   DetId id(hit.id());
00084   const HcalSimParameters& pars = dynamic_cast<const HcalSimParameters&>(theParameterMap->simParameters(id));
00085   theSiPM->setNCells(pars.pixels());
00086 
00087   double signal = analogSignalAmplitude(hit, pars);
00088   int photons = static_cast<int>(signal + 0.5);
00089   int pixels = theSiPM->hitCells(photons, integral);
00090   integral += pixels;
00091   signal = double(pixels);
00092 
00093   CaloSamples result(makeBlankSignal(id));
00094 
00095   if(pixels > 0)
00096   {
00097     double jitter = hit.time() - timeOfFlight(id);
00098 
00099     const double tzero = pars.timePhase() - jitter -
00100       BUNCHSPACE*(pars.binOfMaximum() - thePhaseShift_);
00101     double binTime = tzero;
00102 
00103     for (int bin = 0; bin < result.size(); bin++) {
00104       result[bin] += (*theShape)(binTime)*signal;
00105       binTime += BUNCHSPACE;
00106     }
00107   }
00108 
00109   return result;
00110 }
00111