EgammaCoreTools/plugins/EcalClusterEnergyUncertainty.cc

0001 /** \class EcalClusterEnergyUncertainty
0002   *  Function that provides uncertainty on supercluster energy measurement
0003   *  Available numbers: total effective uncertainty (in GeV)
0004   *                     assymetric uncertainties (positive and negative)
0005   *
0006   *  $Id: EcalClusterEnergyUncertainty.h
0007   *  $Date:
0008   *  $Revision:
0009   *  \author Yurii Maravin, KSU, March 2009
0010   */
0011
0012 #include "CondFormats/DataRecord/interface/EcalClusterEnergyUncertaintyParametersRcd.h"
0013 #include "CondFormats/EcalObjects/interface/EcalClusterEnergyUncertaintyParameters.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/Framework/interface/ConsumesCollector.h"
0016 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionBaseClass.h"
0017
0018 class EcalClusterEnergyUncertainty : public EcalClusterFunctionBaseClass {
0019 public:
0020   EcalClusterEnergyUncertainty(const edm::ParameterSet &, edm::ConsumesCollector iC) : paramsToken_(iC.esConsumes()) {}
0021
0022   // get/set explicit methods for parameters
0023   const EcalClusterEnergyUncertaintyParameters *getParameters() const { return params_; }
0024   // check initialization
0025   void checkInit() const;
0026
0027   // compute the correction
0028   float getValue(const reco::SuperCluster &, const int mode) const override;
0029   float getValue(const reco::BasicCluster &, const EcalRecHitCollection &) const override { return 0.; };
0030
0031   // set parameters
0032   void init(const edm::EventSetup &es) override;
0033
0034 private:
0035   const edm::ESGetToken<EcalClusterEnergyUncertaintyParameters, EcalClusterEnergyUncertaintyParametersRcd> paramsToken_;
0036   const EcalClusterEnergyUncertaintyParameters *params_ = nullptr;
0037 };
0038
0039 void EcalClusterEnergyUncertainty::init(const edm::EventSetup &es) { params_ = &es.getData(paramsToken_); }
0040
0041 void EcalClusterEnergyUncertainty::checkInit() const {
0042   if (!params_) {
0043     // non-initialized function parameters: throw exception
0044     throw cms::Exception("EcalClusterEnergyUncertainty::checkInit()")
0045         << "Trying to access an uninitialized crack correction function.\n"
0046            "Please call `init( edm::EventSetup &)' before any use of the function.\n";
0047   }
0048 }
0049
0050 float EcalClusterEnergyUncertainty::getValue(const reco::SuperCluster &superCluster, const int mode) const {
0051   checkInit();
0052   // mode  = -1 returns negative energy uncertainty
0053   //       = +1 returns positive energy uncertainty
0054   //       =  0 (default) returns overall  energy uncertainty
0055   float en = superCluster.energy();
0056   float eta = fabs(superCluster.eta());
0057   float et = en / cosh(eta);
0058   //fixing divide by zero issue for brem varible, this is the case for single crystal superclusters
0059   //as these "superclusters" are likely noise or spikes so setting value to 0 as the uncertainties
0060   //will be incorrect regardless so doesnt matter what it is
0061   float brem = superCluster.etaWidth() != 0 ? superCluster.phiWidth() / superCluster.etaWidth() : 0;
0062
0063   int offset = 0;
0064
0065   //if ( superCluster.algoID() == reco::CaloCluster::hybrid ) offset = 0;
0066   //else if ( superCluster.algoID() == reco::CaloCluster::multi5x5 ) offset = 36;
0067   // TEMPORARY FIX!!
0068   if (eta < 1.5)
0069     offset = 0;
0070   else if (eta >= 1.5)
0071     offset = 36;
0072   else {
0073     // not supported now
0074     //std::cout << "Not supported value " << superCluster.algoID() << std::endl;
0075     //std::cout << "eta = " << superCluster.eta() << std::endl;
0076     //std::cout << "phi = " << superCluster.phi() << std::endl;
0077     //std::cout << "En  = " << superCluster.energy() << std::endl;
0078     return -1;
0079   }
0080   if (mode == 0)
0081     offset += 0;  // total effective uncertainty
0082   else if (mode == -1)
0083     offset += 12;  // negative energy uncertainty
0084   else if (mode == 1)
0085     offset += 24;  // positive energy uncertainty
0086   else {
0087     // wrong input
0088     return 0;
0089   }
0090
0091   float br0_0 = (params_->params())[offset + 0];
0092   float br0_1 = (params_->params())[offset + 1];
0093   float br0_2 = (params_->params())[offset + 2];
0094   float br0_3 = (params_->params())[offset + 3];
0095
0096   float br1_0 = (params_->params())[offset + 4];
0097   float br1_1 = (params_->params())[offset + 5];
0098   float br1_2 = (params_->params())[offset + 6];
0099   float br1_3 = (params_->params())[offset + 7];
0100
0101   float br2_0 = (params_->params())[offset + 8];
0102   float br2_1 = (params_->params())[offset + 9];
0103   float br2_2 = (params_->params())[offset + 10];
0104   float br2_3 = (params_->params())[offset + 11];
0105
0106   float p0 = (br0_0 + br0_1 * brem) + (br0_2 + br0_3 * brem) / et;
0107   float p1 = (br1_0 + br1_1 * brem) + (br1_2 + br1_3 * brem) / et;
0108   float p2 = (br2_0 + br2_1 * brem) + (br2_2 + br2_3 * brem) / et;
0109   //std::cout << "====================================================" << std::endl;
0110   //std::cout << "et = " << et << "\t eta = " << eta << std::endl;
0111   //std::cout << "p0 = " << p0 << "\t p1 = " << p1 << "\t p2 = " << p2 << std::endl;
0112   float uncertainty = en * (p0 + p1 * fabs(eta) + p2 * eta * eta);
0113   //std::cout << uncertainty << std::endl;
0114   //std::cout << std::endl;
0115   return uncertainty;
0116 }
0117
0118 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionFactory.h"
0119 DEFINE_EDM_PLUGIN(EcalClusterFunctionFactory, EcalClusterEnergyUncertainty, "EcalClusterEnergyUncertainty");