HGCalRecAlgos/src/RecHitTools.cc

0001 #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
0002
0003 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0004 #include "DataFormats/ForwardDetId/interface/HGCScintillatorDetId.h"
0005 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0006 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0007 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0008 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0009 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
0010 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0011 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0012 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0013
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "FWCore/Framework/interface/EventSetup.h"
0018
0019 using namespace hgcal;
0020
0021 namespace {
0022   template <typename DDD>
0023   inline void check_ddd(const DDD* ddd) {
0024     if (nullptr == ddd) {
0025       throw cms::Exception("hgcal::RecHitTools") << "DDDConstants not accessibl to hgcal::RecHitTools!";
0026     }
0027   }
0028
0029   template <typename GEOM>
0030   inline void check_geom(const GEOM* geom) {
0031     if (nullptr == geom) {
0032       throw cms::Exception("hgcal::RecHitTools") << "Geometry not provided yet to hgcal::RecHitTools!";
0033     }
0034   }
0035
0036   inline const HGCalDDDConstants* get_ddd(const CaloSubdetectorGeometry* geom, const HGCalDetId& detid) {
0037     const HGCalGeometry* hg = static_cast<const HGCalGeometry*>(geom);
0038     const HGCalDDDConstants* ddd = &(hg->topology().dddConstants());
0039     check_ddd(ddd);
0040     return ddd;
0041   }
0042
0043   inline const HGCalDDDConstants* get_ddd(const CaloSubdetectorGeometry* geom, const HGCSiliconDetId& detid) {
0044     const HGCalGeometry* hg = static_cast<const HGCalGeometry*>(geom);
0045     const HGCalDDDConstants* ddd = &(hg->topology().dddConstants());
0046     check_ddd(ddd);
0047     return ddd;
0048   }
0049
0050   inline const HGCalDDDConstants* get_ddd(const CaloSubdetectorGeometry* geom, const HFNoseDetId& detid) {
0051     const HGCalGeometry* hg = static_cast<const HGCalGeometry*>(geom);
0052     const HGCalDDDConstants* ddd = &(hg->topology().dddConstants());
0053     check_ddd(ddd);
0054     return ddd;
0055   }
0056
0057   inline const HGCalDDDConstants* get_ddd(const CaloGeometry* geom, int type, int maxLayerEE, int layer) {
0058     DetId::Detector det = ((type == 0) ? DetId::Forward : ((layer > maxLayerEE) ? DetId::HGCalHSi : DetId::HGCalEE));
0059     int subdet = ((type != 0) ? ForwardSubdetector::ForwardEmpty
0060                               : ((layer > maxLayerEE) ? ForwardSubdetector::HGCEE : ForwardSubdetector::HGCHEF));
0061     const HGCalGeometry* hg = static_cast<const HGCalGeometry*>(geom->getSubdetectorGeometry(det, subdet));
0062     const HGCalDDDConstants* ddd = &(hg->topology().dddConstants());
0063     check_ddd(ddd);
0064     return ddd;
0065   }
0066
0067   enum CellType {
0068     CE_E_120 = 0,
0069     CE_E_200 = 1,
0070     CE_E_300 = 2,
0071     CE_H_120 = 3,
0072     CE_H_200 = 4,
0073     CE_H_300 = 5,
0074     CE_H_SCINT = 6,
0075     EnumSize = 7
0076   };
0077
0078 }  // namespace
0079
0080 void RecHitTools::setGeometry(const CaloGeometry& geom) {
0081   geom_ = &geom;
0082   unsigned int wmaxEE(0), wmaxFH(0);
0083   auto geomEE = static_cast<const HGCalGeometry*>(
0084       geom_->getSubdetectorGeometry(DetId::HGCalEE, ForwardSubdetector::ForwardEmpty));
0085   //check if it's the new geometry
0086   if (geomEE) {
0087     geometryType_ = 1;
0088     eeOffset_ = (geomEE->topology().dddConstants()).getLayerOffset();
0089     wmaxEE = (geomEE->topology().dddConstants()).waferCount(0);
0090     auto geomFH = static_cast<const HGCalGeometry*>(
0091         geom_->getSubdetectorGeometry(DetId::HGCalHSi, ForwardSubdetector::ForwardEmpty));
0092     fhOffset_ = (geomFH->topology().dddConstants()).getLayerOffset();
0093     wmaxFH = (geomFH->topology().dddConstants()).waferCount(0);
0094     fhLastLayer_ = fhOffset_ + (geomFH->topology().dddConstants()).lastLayer(true);
0095     auto geomBH = static_cast<const HGCalGeometry*>(
0096         geom_->getSubdetectorGeometry(DetId::HGCalHSc, ForwardSubdetector::ForwardEmpty));
0097     bhOffset_ = (geomBH->topology().dddConstants()).getLayerOffset();
0098     bhFirstLayer_ = bhOffset_ + (geomBH->topology().dddConstants()).firstLayer();
0099     bhLastLayer_ = bhOffset_ + (geomBH->topology().dddConstants()).lastLayer(true);
0100     bhMaxIphi_ = geomBH->topology().dddConstants().maxCells(true);
0101   } else {
0102     geometryType_ = 0;
0103     geomEE =
0104         static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCEE));
0105     eeOffset_ = (geomEE->topology().dddConstants()).getLayerOffset();
0106     wmaxEE = 1 + (geomEE->topology().dddConstants()).waferMax();
0107     auto geomFH =
0108         static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCHEF));
0109     fhOffset_ = (geomFH->topology().dddConstants()).getLayerOffset();
0110     fhLastLayer_ = fhOffset_ + (geomFH->topology().dddConstants()).layers(true);
0111     bhOffset_ = fhLastLayer_;
0112     bhFirstLayer_ = bhOffset_ + 1;
0113     wmaxFH = 1 + (geomFH->topology().dddConstants()).waferMax();
0114     auto geomBH =
0115         static_cast<const HcalGeometry*>(geom_->getSubdetectorGeometry(DetId::Hcal, HcalSubdetector::HcalEndcap));
0116     bhLastLayer_ = bhOffset_ + (geomBH->topology().dddConstants())->getMaxDepth(1);
0117   }
0118   maxNumberOfWafersPerLayer_ = std::max(wmaxEE, wmaxFH);
0119   // For nose geometry
0120   auto geomNose =
0121       static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HFNose));
0122   if (geomNose) {
0123     maxNumberOfWafersNose_ = (geomNose->topology().dddConstants()).waferCount(0);
0124     noseLastLayer_ = (geomNose->topology().dddConstants()).layers(true);
0125   } else {
0126     maxNumberOfWafersNose_ = 0;
0127     noseLastLayer_ = 0;
0128   }
0129 }
0130
0131 const CaloSubdetectorGeometry* RecHitTools::getSubdetectorGeometry(const DetId& id) const {
0132   DetId::Detector det = id.det();
0133   int subdet = (det == DetId::HGCalEE || det == DetId::HGCalHSi || det == DetId::HGCalHSc)
0134                    ? ForwardSubdetector::ForwardEmpty
0135                    : id.subdetId();
0136   auto geom = geom_->getSubdetectorGeometry(det, subdet);
0137   check_geom(geom);
0138   return geom;
0139 }
0140
0141 GlobalPoint RecHitTools::getPosition(const DetId& id) const {
0142   auto geom = getSubdetectorGeometry(id);
0143   GlobalPoint position;
0144   if (id.det() == DetId::Hcal || id.det() == DetId::Ecal) {
0145     position = geom->getGeometry(id)->getPosition();
0146   } else {
0147     auto hg = static_cast<const HGCalGeometry*>(geom);
0148     position = hg->getPosition(id);
0149   }
0150   return position;
0151 }
0152
0153 GlobalPoint RecHitTools::getPositionLayer(int layer, bool nose, bool barrel) const {
0154   unsigned int lay = std::abs(layer);
0155   double x(0), y(0), z(0);
0156   if (nose) {
0157     auto geomNose =
0158         static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HFNose));
0159     if (geomNose) {
0160       const HGCalDDDConstants* ddd = &(geomNose->topology().dddConstants());
0161       if (ddd)
0162         z = (layer > 0) ? ddd->waferZ(lay, true) : -ddd->waferZ(lay, true);
0163     }
0164   } else if (barrel) {
0165     if (layer == 0) {
0166       DetId firstValidDetId = (geom_->getValidDetIds(DetId::Ecal, EcalBarrel))[0];
0167       x = getPosition(firstValidDetId).x();
0168       y = getPosition(firstValidDetId).y();
0169     } else {
0170       const std::vector<DetId>& validDetIds = geom_->getValidDetIds(DetId::Hcal, HcalBarrel);
0171       auto firstValidDetId =
0172           std::find_if(validDetIds.begin(), validDetIds.end(), [&](auto id) { return HcalDetId(id).depth() == layer; });
0173       assert(firstValidDetId != validDetIds.end());
0174       x = getPosition(*firstValidDetId).x();
0175       y = getPosition(*firstValidDetId).y();
0176     }
0177   } else {
0178     const HGCalDDDConstants* ddd = get_ddd(geom_, geometryType_, fhOffset_, lay);
0179     if (geometryType_ == 1) {
0180       if (lay > fhOffset_)
0181         lay -= fhOffset_;
0182     }
0183     z = (layer > 0) ? ddd->waferZ(lay, true) : -ddd->waferZ(lay, true);
0184   }
0185   return GlobalPoint(x, y, z);
0186 }
0187
0188 int RecHitTools::zside(const DetId& id) const {
0189   int zside = 0;
0190   if (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi) {
0191     zside = HGCSiliconDetId(id).zside();
0192   } else if (id.det() == DetId::HGCalHSc) {
0193     zside = HGCScintillatorDetId(id).zside();
0194   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0195     zside = HFNoseDetId(id).zside();
0196   } else if (id.det() == DetId::Forward) {
0197     zside = HGCalDetId(id).zside();
0198   } else if (id.det() == DetId::Hcal && id.subdetId() == HcalEndcap) {
0199     zside = HcalDetId(id).zside();
0200   }
0201   return zside;
0202 }
0203
0204 std::float_t RecHitTools::getSiThickness(const DetId& id) const {
0205   auto geom = getSubdetectorGeometry(id);
0206   std::float_t thick(0.37);
0207   if (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi) {
0208     const HGCSiliconDetId hid(id);
0209     auto ddd = get_ddd(geom, hid);
0210     thick = ddd->cellThickness(hid.layer(), hid.waferU(), hid.waferV());
0211   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0212     const HFNoseDetId hid(id);
0213     auto ddd = get_ddd(geom, hid);
0214     thick = ddd->cellThickness(hid.layer(), hid.waferU(), hid.waferV());
0215   } else if (id.det() == DetId::Forward) {
0216     const HGCalDetId hid(id);
0217     auto ddd = get_ddd(geom, hid);
0218     thick = ddd->cellThickness(hid.layer(), hid.wafer(), 0);
0219   } else {
0220     LogDebug("getSiThickness::InvalidSiliconDetid")
0221         << "det id: " << std::hex << id.rawId() << std::dec << ":" << id.det() << " is not HGCal silicon!";
0222     // It does not make any sense to return 0.37 as thickness for a detector
0223     // that is not Silicon(does it?). Mark it as 0. to be able to distinguish
0224     // the case.
0225     thick = 0.f;
0226   }
0227   return thick;
0228 }
0229
0230 int RecHitTools::getSiThickIndex(const DetId& id) const {
0231   int thickIndex = -1;
0232   if (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi) {
0233     thickIndex = HGCSiliconDetId(id).type();
0234   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0235     thickIndex = HFNoseDetId(id).type();
0236   } else if (id.det() == DetId::Forward) {
0237     float thickness = getSiThickness(id);
0238     if (thickness > 99. && thickness < 121.)
0239       thickIndex = 0;
0240     else if (thickness > 199. && thickness < 201.)
0241       thickIndex = 1;
0242     else if (thickness > 299. && thickness < 301.)
0243       thickIndex = 2;
0244     else
0245       assert(thickIndex > 0 && "ERROR - silicon thickness has a nonsensical value");
0246   }
0247   return thickIndex;
0248 }
0249
0250 std::pair<float, float> RecHitTools::getScintDEtaDPhi(const DetId& id) const {
0251   if (!isScintillator(id)) {
0252     LogDebug("getScintDEtaDPhi::InvalidScintDetid")
0253         << "det id: " << std::hex << id.rawId() << std::dec << ":" << id.det() << " is not HGCal scintillator!";
0254     return {0.f, 0.f};
0255   }
0256   auto cellGeom = getSubdetectorGeometry(id)->getGeometry(id);
0257   return {cellGeom->etaSpan(), cellGeom->phiSpan()};
0258 }
0259
0260 std::float_t RecHitTools::getRadiusToSide(const DetId& id) const {
0261   auto geom = getSubdetectorGeometry(id);
0262   std::float_t size(std::numeric_limits<std::float_t>::max());
0263   if (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi) {
0264     const HGCSiliconDetId hid(id);
0265     auto ddd = get_ddd(geom, hid);
0266     size = ddd->cellSizeHex(hid.type());
0267   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0268     const HFNoseDetId hid(id);
0269     auto ddd = get_ddd(geom, hid);
0270     size = ddd->cellSizeHex(hid.type());
0271   } else if (id.det() == DetId::Forward) {
0272     const HGCalDetId hid(id);
0273     auto ddd = get_ddd(geom, hid);
0274     size = ddd->cellSizeHex(hid.waferType());
0275   } else {
0276     edm::LogError("getRadiusToSide::InvalidSiliconDetid")
0277         << "det id: " << std::hex << id.rawId() << std::dec << ":" << id.det() << " is not HGCal silicon!";
0278   }
0279   return size;
0280 }
0281
0282 unsigned int RecHitTools::getLayer(const ForwardSubdetector type) const {
0283   int layer(0);
0284   switch (type) {
0285     case (ForwardSubdetector::HGCEE): {
0286       auto geomEE =
0287           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCEE));
0288       layer = (geomEE->topology().dddConstants()).layers(true);
0289       break;
0290     }
0291     case (ForwardSubdetector::HGCHEF): {
0292       auto geomFH =
0293           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCHEF));
0294       layer = (geomFH->topology().dddConstants()).layers(true);
0295       break;
0296     }
0297     case (ForwardSubdetector::HGCHEB): {
0298       auto geomBH =
0299           static_cast<const HcalGeometry*>(geom_->getSubdetectorGeometry(DetId::Hcal, HcalSubdetector::HcalEndcap));
0300       layer = (geomBH->topology().dddConstants())->getMaxDepth(1);
0301       break;
0302     }
0303     case (ForwardSubdetector::ForwardEmpty): {
0304       auto geomEE =
0305           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCEE));
0306       layer = (geomEE->topology().dddConstants()).layers(true);
0307       auto geomFH =
0308           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCHEF));
0309       layer += (geomFH->topology().dddConstants()).layers(true);
0310       auto geomBH =
0311           static_cast<const HcalGeometry*>(geom_->getSubdetectorGeometry(DetId::Hcal, HcalSubdetector::HcalEndcap));
0312       if (geomBH)
0313         layer += (geomBH->topology().dddConstants())->getMaxDepth(1);
0314       auto geomBH2 =
0315           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HGCHEB));
0316       if (geomBH2)
0317         layer += (geomBH2->topology().dddConstants()).layers(true);
0318       break;
0319     }
0320     case (ForwardSubdetector::HFNose): {
0321       auto geomHFN =
0322           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HFNose));
0323       layer = (geomHFN->topology().dddConstants()).layers(true);
0324       break;
0325     }
0326     default:
0327       layer = 0;
0328   }
0329   return (unsigned int)(layer);
0330 }
0331
0332 unsigned int RecHitTools::getLayer(const DetId::Detector type, bool nose) const {
0333   int layer;
0334   switch (type) {
0335     case (DetId::HGCalEE): {
0336       auto geomEE =
0337           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(type, ForwardSubdetector::ForwardEmpty));
0338       layer = (geomEE->topology().dddConstants()).layers(true);
0339       break;
0340     }
0341     case (DetId::HGCalHSi): {
0342       auto geomFH =
0343           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(type, ForwardSubdetector::ForwardEmpty));
0344       layer = (geomFH->topology().dddConstants()).layers(true);
0345       break;
0346     }
0347     case (DetId::HGCalHSc): {
0348       auto geomBH =
0349           static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(type, ForwardSubdetector::ForwardEmpty));
0350       layer = (geomBH->topology().dddConstants()).layers(true);
0351       break;
0352     }
0353     case (DetId::Forward): {
0354       if (nose) {
0355         auto geomHFN = static_cast<const HGCalGeometry*>(
0356             geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HFNose));
0357         layer = (geomHFN->topology().dddConstants()).layers(true);
0358       } else {
0359         auto geomEE = static_cast<const HGCalGeometry*>(
0360             geom_->getSubdetectorGeometry(DetId::HGCalEE, ForwardSubdetector::ForwardEmpty));
0361         layer = (geomEE->topology().dddConstants()).layers(true);
0362         auto geomFH = static_cast<const HGCalGeometry*>(
0363             geom_->getSubdetectorGeometry(DetId::HGCalHSi, ForwardSubdetector::ForwardEmpty));
0364         layer += (geomFH->topology().dddConstants()).layers(true);
0365       }
0366       break;
0367     }
0368     default:
0369       layer = 0;
0370   }
0371   return (unsigned int)(layer);
0372 }
0373
0374 unsigned int RecHitTools::getLayer(const DetId& id) const {
0375   unsigned int layer = std::numeric_limits<unsigned int>::max();
0376   if (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi) {
0377     layer = HGCSiliconDetId(id).layer();
0378   } else if (id.det() == DetId::HGCalHSc) {
0379     layer = HGCScintillatorDetId(id).layer();
0380   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0381     layer = HFNoseDetId(id).layer();
0382   } else if (id.det() == DetId::Forward) {
0383     layer = HGCalDetId(id).layer();
0384   } else if (id.det() == DetId::Hcal && id.subdetId() != HcalEmpty) {
0385     if (id.subdetId() == HcalBarrel)
0386       layer = HcalDetId(id).depth();
0387     else if (id.subdetId() == HcalOuter)
0388       layer = HcalDetId(id).depth() + 1;
0389   } else if (id.det() == DetId::Ecal) {
0390     layer = 0;
0391   }
0392   return layer;
0393 }
0394
0395 unsigned int RecHitTools::getLayerWithOffset(const DetId& id) const {
0396   unsigned int layer = getLayer(id);
0397   if (id.det() == DetId::Forward && id.subdetId() == HGCHEF) {
0398     layer += fhOffset_;
0399   } else if (id.det() == DetId::HGCalHSi || id.det() == DetId::HGCalHSc) {
0400     // DetId::HGCalHSc hits include the offset w.r.t. EE already
0401     layer += fhOffset_;
0402   } else if (id.det() == DetId::Hcal && id.subdetId() == HcalEndcap) {
0403     layer += bhOffset_;
0404   }
0405   // no need to add offset for HFnose
0406   return layer;
0407 }
0408
0409 std::pair<int, int> RecHitTools::getWafer(const DetId& id) const {
0410   int waferU = std::numeric_limits<int>::max();
0411   int waferV = 0;
0412   if ((id.det() == DetId::HGCalEE) || (id.det() == DetId::HGCalHSi)) {
0413     waferU = HGCSiliconDetId(id).waferU();
0414     waferV = HGCSiliconDetId(id).waferV();
0415   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0416     waferU = HFNoseDetId(id).waferU();
0417     waferV = HFNoseDetId(id).waferV();
0418   } else if (id.det() == DetId::Forward) {
0419     waferU = HGCalDetId(id).wafer();
0420   } else {
0421     edm::LogError("getWafer::InvalidSiliconDetid")
0422         << "det id: " << std::hex << id.rawId() << std::dec << ":" << id.det() << " is not HGCal silicon!";
0423   }
0424   return std::pair<int, int>(waferU, waferV);
0425 }
0426
0427 std::pair<int, int> RecHitTools::getCell(const DetId& id) const {
0428   int cellU = std::numeric_limits<int>::max();
0429   int cellV = 0;
0430   if ((id.det() == DetId::HGCalEE) || (id.det() == DetId::HGCalHSi)) {
0431     cellU = HGCSiliconDetId(id).cellU();
0432     cellV = HGCSiliconDetId(id).cellV();
0433   } else if (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)) {
0434     cellU = HFNoseDetId(id).cellU();
0435     cellV = HFNoseDetId(id).cellV();
0436   } else if (id.det() == DetId::Forward) {
0437     cellU = HGCalDetId(id).cell();
0438   } else {
0439     edm::LogError("getCell::InvalidSiliconDetid")
0440         << "det id: " << std::hex << id.rawId() << std::dec << ":" << id.det() << " is not HGCal silicon!";
0441   }
0442   return std::pair<int, int>(cellU, cellV);
0443 }
0444
0445 bool RecHitTools::isHalfCell(const DetId& id) const {
0446   bool ishalf = false;
0447   if (id.det() == DetId::Forward) {
0448     HGCalDetId hid(id);
0449     auto geom = getSubdetectorGeometry(hid);
0450     auto ddd = get_ddd(geom, hid);
0451     const int waferType = ddd->waferTypeT(hid.waferType());
0452     return ddd->isHalfCell(waferType, hid.cell());
0453   }
0454   //new geometry is always false
0455   return ishalf;
0456 }
0457
0458 int RecHitTools::getCellType(const DetId& id) const {
0459   auto layer_number = getLayerWithOffset(id);
0460   auto thickness = getSiThickIndex(id);
0461   auto geomNose =
0462       static_cast<const HGCalGeometry*>(geom_->getSubdetectorGeometry(DetId::Forward, ForwardSubdetector::HFNose));
0463   auto isNose = geomNose ? true : false;
0464   auto isEELayer = (layer_number <= lastLayerEE(isNose));
0465   auto isScint = isScintillator(id);
0466   int layerType = -1;
0467
0468   if (isScint) {
0469     layerType = CE_H_SCINT;
0470   }
0471   if (isEELayer) {
0472     if (thickness == 0) {
0473       layerType = CE_E_120;
0474     } else if (thickness == 1) {
0475       layerType = CE_E_200;
0476     } else if (thickness == 2) {
0477       layerType = CE_E_300;
0478     }
0479   } else {
0480     if (thickness == 0) {
0481       layerType = CE_H_120;
0482     } else if (thickness == 1) {
0483       layerType = CE_H_200;
0484     } else if (thickness == 2) {
0485       layerType = CE_H_300;
0486     }
0487   }
0488   assert(layerType != -1);
0489   return layerType;
0490 }
0491
0492 bool RecHitTools::isSilicon(const DetId& id) const {
0493   return (id.det() == DetId::HGCalEE || id.det() == DetId::HGCalHSi ||
0494           (id.det() == DetId::Forward && id.subdetId() == static_cast<int>(HFNose)));
0495 }
0496
0497 bool RecHitTools::isScintillator(const DetId& id) const { return (id.det() == DetId::HGCalHSc); }
0498
0499 bool RecHitTools::isBarrel(const DetId& id) const { return (id.det() == DetId::Ecal || id.det() == DetId::Hcal); }
0500 bool RecHitTools::isOnlySilicon(const unsigned int layer) const {
0501   // HFnose TODO
0502   bool isonlysilicon = (layer % bhLastLayer_) < bhOffset_;
0503   return isonlysilicon;
0504 }
0505
0506 float RecHitTools::getEta(const GlobalPoint& position, const float& vertex_z) const {
0507   GlobalPoint corrected_position = GlobalPoint(position.x(), position.y(), position.z() - vertex_z);
0508   return corrected_position.eta();
0509 }
0510
0511 float RecHitTools::getEta(const DetId& id, const float& vertex_z) const {
0512   GlobalPoint position = getPosition(id);
0513   float eta = getEta(position, vertex_z);
0514   return eta;
0515 }
0516
0517 float RecHitTools::getPhi(const GlobalPoint& position) const {
0518   float phi = atan2(position.y(), position.x());
0519   return phi;
0520 }
0521
0522 float RecHitTools::getPhi(const DetId& id) const {
0523   GlobalPoint position = getPosition(id);
0524   float phi = atan2(position.y(), position.x());
0525   return phi;
0526 }
0527
0528 float RecHitTools::getPt(const GlobalPoint& position, const float& hitEnergy, const float& vertex_z) const {
0529   float eta = getEta(position, vertex_z);
0530   float pt = hitEnergy / cosh(eta);
0531   return pt;
0532 }
0533
0534 float RecHitTools::getPt(const DetId& id, const float& hitEnergy, const float& vertex_z) const {
0535   GlobalPoint position = getPosition(id);
0536   float eta = getEta(position, vertex_z);
0537   float pt = hitEnergy / cosh(eta);
0538   return pt;
0539 }
0540
0541 std::pair<uint32_t, uint32_t> RecHitTools::firstAndLastLayer(DetId::Detector det, int subdet) const {
0542   if ((det == DetId::HGCalEE) || ((det == DetId::Forward) && (subdet == HGCEE))) {
0543     return std::make_pair(eeOffset_ + 1, fhOffset_);
0544   } else if ((det == DetId::HGCalHSi) || ((det == DetId::Forward) && (subdet == HGCHEF))) {
0545     return std::make_pair(fhOffset_ + 1, fhLastLayer_);
0546   } else if ((det == DetId::Forward) && (subdet == HFNose)) {
0547     return std::make_pair(1, noseLastLayer_);
0548   } else {
0549     return std::make_pair(bhFirstLayer_, bhLastLayer_);
0550   }
0551 }
0552
0553 bool RecHitTools::maskCell(const DetId& id, int corners) const {
0554   if (id.det() == DetId::Hcal) {
0555     return false;
0556   } else {
0557     auto hg = static_cast<const HGCalGeometry*>(getSubdetectorGeometry(id));
0558     return hg->topology().dddConstants().maskCell(id, corners);
0559   }
0560 }