EcalCorrections/src/EcalShowerContainmentCorrections.cc

0001 // Implementation of class EcalShowerContainmentCorrections
0002 // Author: Stefano Argiro'
0003 // $Id: EcalShowerContainmentCorrections.cc,v 1.1 2007/05/15 20:37:22 argiro Exp $
0004
0005 #include "CondFormats/EcalCorrections/interface/EcalShowerContainmentCorrections.h"
0006 #include <DataFormats/EcalDetId/interface/EBDetId.h>
0007 #include <FWCore/MessageLogger/interface/MessageLogger.h>
0008 //#include <iostream>
0009
0010 #include "CondCore/CondDB/interface/Serialization.h"
0011 #include "CondFormats/External/interface/EcalDetID.h"
0012
0013 const EcalShowerContainmentCorrections::Coefficients EcalShowerContainmentCorrections::correctionCoefficients(
0014     const EBDetId& centerxtal) const {
0015   GroupMap::const_iterator iter = groupmap_.find(centerxtal.rawId());
0016
0017   if (iter != groupmap_.end()) {
0018     int group = iter->second;
0019     return coefficients_[group - 1];
0020   }
0021
0022   edm::LogError("ShowerContaiment Correction not found");
0023   return Coefficients();
0024 }
0025
0026 void EcalShowerContainmentCorrections::fillCorrectionCoefficients(const EBDetId& xtal,
0027                                                                   int group,
0028                                                                   const Coefficients& coefficients) {
0029   // do not replace if we already have the xtal
0030   if (groupmap_.find(xtal) != groupmap_.end())
0031     return;
0032   groupmap_[xtal] = group;
0033
0034   if (coefficients_.size() < (unsigned int)(group)) {
0035     coefficients_.resize(group);
0036     coefficients_[group - 1] = coefficients;
0037   }
0038
0039   // we don't need to fill coefficients if the group has already been inserted
0040 }
0041
0042 void EcalShowerContainmentCorrections::fillCorrectionCoefficients(const int supermodule,
0043                                                                   const int module,
0044                                                                   const Coefficients& coefficients) {
0045   if (module > EBDetId::kModulesPerSM) {
0046     edm::LogError("Invalid Module Number");
0047     return;
0048   }
0049
0050   // what is EBDetID::kModuleBoundaries ? we better redefine them here ...
0051   const int kModuleLow[] = {1, 501, 901, 1301};
0052   const int kModuleHigh[] = {500, 900, 1300, 1700};
0053
0054   for (int xtal = kModuleLow[module - 1]; xtal <= kModuleHigh[module - 1]; ++xtal) {
0055     EBDetId detid(supermodule, xtal, EBDetId::SMCRYSTALMODE);
0056     fillCorrectionCoefficients(detid, module, coefficients);
0057   }
0058 }
0059
0060 /** Calculate the correction for the given direction and type  */
0061 const double EcalShowerContainmentCorrections::correctionXY(const EBDetId& xtal,
0062                                                             double position,
0063                                                             EcalShowerContainmentCorrections::Direction dir,
0064                                                             EcalShowerContainmentCorrections::Type type) const {
0065   GroupMap::const_iterator iter = groupmap_.find(xtal);
0066   if (iter == groupmap_.end())
0067     return -1;
0068
0069   int group = iter->second;
0070   EcalShowerContainmentCorrections::Coefficients coeff = coefficients_[group - 1];
0071
0072   int offset = 0;
0073
0074   if (dir == eY)
0075     offset += 2 * Coefficients::kPolynomialDegree;
0076   if (position < 0)
0077     offset += Coefficients::kPolynomialDegree;
0078   if (type == e5x5)
0079     offset += 4 * Coefficients::kPolynomialDegree;
0080
0081   double corr = 0;
0082
0083   for (int i = offset; i < offset + Coefficients::kPolynomialDegree; ++i) {
0084     corr += coeff.data[i] * pow(position, i - offset);
0085   }
0086
0087   return corr;
0088 }
0089
0090 const double EcalShowerContainmentCorrections::correction3x3(const EBDetId& xtal, const math::XYZPoint& pos) const {
0091   double x = pos.X() * 10;  // correction functions use mm
0092   double y = pos.Y() * 10;
0093
0094   double corrx = correctionXY(xtal, x, eX, e3x3);
0095   double corry = correctionXY(xtal, y, eY, e3x3);
0096
0097   return corrx * corry;
0098 }
0099
0100 const double EcalShowerContainmentCorrections::correction5x5(const EBDetId& xtal, const math::XYZPoint& pos) const {
0101   double x = pos.X() * 10;  // correction functions use mm
0102   double y = pos.Y() * 10;
0103
0104   double corrx = correctionXY(xtal, x, eX, e5x5);
0105   double corry = correctionXY(xtal, y, eY, e5x5);
0106
0107   return corrx * corry;
0108 }