File indexing completed on 2024-04-06 12:25:45
0001 #include "ESRecHitWorker.h"
0002
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/EcalDigi/interface/EcalDigiCollections.h"
0005 #include "DataFormats/EcalRecHit/interface/EcalUncalibratedRecHit.h"
0006 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009
0010 #include <cmath>
0011 #include <iomanip>
0012 #include <iostream>
0013
0014 ESRecHitWorker::ESRecHitWorker(const edm::ParameterSet &ps, edm::ConsumesCollector cc) : ESRecHitWorkerBaseClass(ps) {
0015 recoAlgo_ = ps.getParameter<int>("ESRecoAlgo");
0016 esgainToken_ = cc.esConsumes<ESGain, ESGainRcd>();
0017 esMIPToGeVToken_ = cc.esConsumes<ESMIPToGeVConstant, ESMIPToGeVConstantRcd>();
0018 esWeightsToken_ = cc.esConsumes<ESTimeSampleWeights, ESTimeSampleWeightsRcd>();
0019 esPedestalsToken_ = cc.esConsumes<ESPedestals, ESPedestalsRcd>();
0020 esMIPsToken_ = cc.esConsumes<ESIntercalibConstants, ESIntercalibConstantsRcd>();
0021 esChannelStatusToken_ = cc.esConsumes<ESChannelStatus, ESChannelStatusRcd>();
0022 esRatioCutsToken_ = cc.esConsumes<ESRecHitRatioCuts, ESRecHitRatioCutsRcd>();
0023 esAngleCorrFactorsToken_ = cc.esConsumes<ESAngleCorrectionFactors, ESAngleCorrectionFactorsRcd>();
0024
0025 if (recoAlgo_ == 0)
0026 algoW_ = new ESRecHitSimAlgo();
0027 else if (recoAlgo_ == 1)
0028 algoF_ = new ESRecHitFitAlgo();
0029 else
0030 algoA_ = new ESRecHitAnalyticAlgo();
0031 }
0032
0033 ESRecHitWorker::~ESRecHitWorker() {
0034 if (recoAlgo_ == 0)
0035 delete algoW_;
0036 else if (recoAlgo_ == 1)
0037 delete algoF_;
0038 else
0039 delete algoA_;
0040 }
0041
0042 void ESRecHitWorker::set(const edm::EventSetup &es) {
0043 esgain_ = es.getHandle(esgainToken_);
0044 const ESGain *gain = esgain_.product();
0045
0046 esMIPToGeV_ = es.getHandle(esMIPToGeVToken_);
0047 const ESMIPToGeVConstant *mipToGeV = esMIPToGeV_.product();
0048
0049 double ESGain = gain->getESGain();
0050 double ESMIPToGeV = (ESGain == 1) ? mipToGeV->getESValueLow() : mipToGeV->getESValueHigh();
0051
0052 esWeights_ = es.getHandle(esWeightsToken_);
0053 const ESTimeSampleWeights *wgts = esWeights_.product();
0054
0055 float w0 = wgts->getWeightForTS0();
0056 float w1 = wgts->getWeightForTS1();
0057 float w2 = wgts->getWeightForTS2();
0058
0059 esPedestals_ = es.getHandle(esPedestalsToken_);
0060 const ESPedestals *peds = esPedestals_.product();
0061
0062 esMIPs_ = es.getHandle(esMIPsToken_);
0063 const ESIntercalibConstants *mips = esMIPs_.product();
0064
0065 esAngleCorrFactors_ = es.getHandle(esAngleCorrFactorsToken_);
0066 const ESAngleCorrectionFactors *ang = esAngleCorrFactors_.product();
0067
0068 esChannelStatus_ = es.getHandle(esChannelStatusToken_);
0069 const ESChannelStatus *channelStatus = esChannelStatus_.product();
0070
0071 esRatioCuts_ = es.getHandle(esRatioCutsToken_);
0072 const ESRecHitRatioCuts *ratioCuts = esRatioCuts_.product();
0073
0074 if (recoAlgo_ == 0) {
0075 algoW_->setESGain(ESGain);
0076 algoW_->setMIPGeV(ESMIPToGeV);
0077 algoW_->setW0(w0);
0078 algoW_->setW1(w1);
0079 algoW_->setW2(w2);
0080 algoW_->setPedestals(peds);
0081 algoW_->setIntercalibConstants(mips);
0082 algoW_->setChannelStatus(channelStatus);
0083 algoW_->setRatioCuts(ratioCuts);
0084 algoW_->setAngleCorrectionFactors(ang);
0085 } else if (recoAlgo_ == 1) {
0086 algoF_->setESGain(ESGain);
0087 algoF_->setMIPGeV(ESMIPToGeV);
0088 algoF_->setPedestals(peds);
0089 algoF_->setIntercalibConstants(mips);
0090 algoF_->setChannelStatus(channelStatus);
0091 algoF_->setRatioCuts(ratioCuts);
0092 algoF_->setAngleCorrectionFactors(ang);
0093 } else {
0094 algoA_->setESGain(ESGain);
0095 algoA_->setMIPGeV(ESMIPToGeV);
0096 algoA_->setPedestals(peds);
0097 algoA_->setIntercalibConstants(mips);
0098 algoA_->setChannelStatus(channelStatus);
0099 algoA_->setRatioCuts(ratioCuts);
0100 algoA_->setAngleCorrectionFactors(ang);
0101 }
0102 }
0103
0104 bool ESRecHitWorker::run(const ESDigiCollection::const_iterator &itdg, ESRecHitCollection &result) {
0105 if (recoAlgo_ == 0)
0106 result.push_back(algoW_->reconstruct(*itdg));
0107 else if (recoAlgo_ == 1)
0108 result.push_back(algoF_->reconstruct(*itdg));
0109 else
0110 result.push_back(algoA_->reconstruct(*itdg));
0111 return true;
0112 }
0113
0114 #include "FWCore/Framework/interface/MakerMacros.h"
0115 #include "RecoLocalCalo/EcalRecProducers/interface/ESRecHitWorkerFactory.h"
0116 DEFINE_EDM_PLUGIN(ESRecHitWorkerFactory, ESRecHitWorker, "ESRecHitWorker");