MTDGeometryBuilder/src/MTDGeomUtil.cc

0001 #include "Geometry/MTDGeometryBuilder/interface/MTDGeomUtil.h"
0002
0003 #include "Geometry/MTDCommonData/interface/MTDTopologyMode.h"
0004
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007
0008 using namespace mtd;
0009
0010 void MTDGeomUtil::setGeometry(const MTDGeometry* geom) { geom_ = geom; }
0011
0012 void MTDGeomUtil::setTopology(const MTDTopology* topo) { topology_ = topo; }
0013
0014 bool MTDGeomUtil::isETL(const DetId& id) const {
0015   MTDDetId hid{id};
0016   const auto& subDet = hid.mtdSubDetector();
0017   if (subDet == 0)
0018     throw cms::Exception("mtd::MTDGeomUtil") << "DetId " << hid.rawId() << " not in MTD!" << std::endl;
0019   if (subDet == MTDDetId::MTDType::ETL)
0020     return true;
0021   return false;
0022 }
0023
0024 bool MTDGeomUtil::isBTL(const DetId& id) const { return !(isETL(id)); }
0025
0026 // row and column set to 0 by default since they are not needed for BTL
0027 std::pair<LocalPoint, GlobalPoint> MTDGeomUtil::position(const DetId& id, int row, int column) const {
0028   LocalPoint local_point(0., 0., 0.);
0029   GlobalPoint global_point(0., 0., 0.);
0030   if (isBTL(id)) {
0031     BTLDetId detId{id};
0032     DetId geoId = detId.geographicalId(MTDTopologyMode::crysLayoutFromTopoMode(topology_->getMTDTopologyMode()));
0033     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0034     if (thedet == nullptr)
0035       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0036                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0037     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0038     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0039
0040     local_point = topo.pixelToModuleLocalPoint(local_point, detId.row(topo.nrows()), detId.column(topo.nrows()));
0041     global_point = thedet->toGlobal(local_point);
0042   } else {
0043     ETLDetId detId{id};
0044     DetId geoId = detId.geographicalId();
0045     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0046     if (thedet == nullptr)
0047       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0048                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0049     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0050     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0051
0052     local_point = LocalPoint(topo.localX(row), topo.localY(column), 0.);
0053     global_point = thedet->toGlobal(local_point);
0054   }
0055   return {local_point, global_point};
0056 }
0057
0058 GlobalPoint MTDGeomUtil::globalPosition(const DetId& id, const LocalPoint& local_point) const {
0059   auto global_point = GlobalPoint(0., 0., 0.);
0060   if (isBTL(id)) {
0061     BTLDetId detId{id};
0062     DetId geoId = detId.geographicalId(MTDTopologyMode::crysLayoutFromTopoMode(topology_->getMTDTopologyMode()));
0063     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0064     if (thedet == nullptr)
0065       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0066                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0067     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0068     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0069     auto local_point_sim =
0070         topo.pixelToModuleLocalPoint(local_point, detId.row(topo.nrows()), detId.column(topo.nrows()));
0071     global_point = thedet->toGlobal(local_point_sim);
0072   } else {
0073     ETLDetId detId{id};
0074     DetId geoId = detId.geographicalId();
0075     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0076     if (thedet == nullptr)
0077       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0078                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0079     global_point = thedet->toGlobal(local_point);
0080   }
0081   return global_point;
0082 }
0083
0084 int MTDGeomUtil::zside(const DetId& id) const {
0085   const MTDDetId hid(id);
0086   return hid.zside();
0087 }
0088
0089 unsigned int MTDGeomUtil::layer(const DetId& id) const {
0090   unsigned int layer(0);
0091   if (isETL(id)) {
0092     ETLDetId hid(id);
0093     layer = hid.nDisc();
0094   }
0095   return layer;
0096 }
0097
0098 int MTDGeomUtil::module(const DetId& id) const {
0099   int module = -1;
0100   if (isETL(id)) {
0101     ETLDetId hid(id);
0102     module = hid.module();
0103   } else {
0104     BTLDetId hid(id);
0105     module = hid.module();
0106   }
0107   return module;
0108 }
0109
0110 // returns the local position as a pair (x, y) - for ETL
0111 std::pair<float, float> MTDGeomUtil::pixelInModule(const DetId& id, const int row, const int column) const {
0112   if (isBTL(id))
0113     throw cms::Exception("mtd::MTDGeomUtil")
0114         << "ID: " << std::hex << id.rawId() << " from BTL. This method is for ETL only." << std::endl;
0115   ETLDetId detId(id);
0116   DetId geoId = detId.geographicalId();
0117   const MTDGeomDet* thedet = geom_->idToDet(geoId);
0118   if (thedet == nullptr)
0119     throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " (" << detId.rawId()
0120                                              << ") is invalid!" << std::dec << std::endl;
0121   const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0122   const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0123   const Local3DPoint local_point(topo.localX(row), topo.localY(column), 0.);
0124   return topo.pixel(local_point);
0125 }
0126
0127 // returns row and column as a pair (row, col)
0128 std::pair<uint8_t, uint8_t> MTDGeomUtil::pixelInModule(const DetId& id, const LocalPoint& local_point) const {
0129   if (isETL(id)) {
0130     ETLDetId detId(id);
0131     DetId geoId = detId.geographicalId();
0132     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0133     if (thedet == nullptr)
0134       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0135                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0136     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0137     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0138     const auto& thepixel = topo.pixelIndex(local_point);
0139     uint8_t row = static_cast<uint8_t>(thepixel.first);
0140     uint8_t col = static_cast<uint8_t>(thepixel.second);
0141     return std::pair<uint8_t, uint8_t>(row, col);
0142   } else {
0143     BTLDetId detId(id);
0144     DetId geoId = detId.geographicalId(MTDTopologyMode::crysLayoutFromTopoMode(topology_->getMTDTopologyMode()));
0145     const MTDGeomDet* thedet = geom_->idToDet(geoId);
0146     if (thedet == nullptr)
0147       throw cms::Exception("mtd::MTDGeomUtil") << "GeographicalID: " << std::hex << geoId.rawId() << " ("
0148                                                << detId.rawId() << ") is invalid!" << std::dec << std::endl;
0149     const ProxyMTDTopology& topoproxy = static_cast<const ProxyMTDTopology&>(thedet->topology());
0150     const RectangularMTDTopology& topo = static_cast<const RectangularMTDTopology&>(topoproxy.specificTopology());
0151     auto topo_point = topo.pixelToModuleLocalPoint(local_point, detId.row(topo.nrows()), detId.column(topo.nrows()));
0152     const auto& thepixel = topo.pixelIndex(topo_point);
0153     uint8_t row = static_cast<uint8_t>(thepixel.first);
0154     uint8_t col = static_cast<uint8_t>(thepixel.second);
0155     return std::pair<uint8_t, uint8_t>(row, col);
0156   }
0157 }
0158
0159 int MTDGeomUtil::crystalInModule(const DetId& id) const {
0160   BTLDetId hid(id);
0161   return hid.crystal();
0162 }
0163
0164 // returns the sensor module id
0165 uint32_t MTDGeomUtil::sensorModuleId(const DetId& id) const {
0166   if (isBTL(id)) {
0167     BTLDetId detId(id);
0168     DetId geoId = detId.geographicalId(MTDTopologyMode::crysLayoutFromTopoMode(topology_->getMTDTopologyMode()));
0169     return (geoId.rawId());
0170   } else {
0171     ETLDetId detId(id);
0172     DetId geoId = detId.geographicalId();
0173     return (geoId.rawId());
0174   }
0175 }
0176
0177 float MTDGeomUtil::eta(const GlobalPoint& position, const float& vertex_z) const {
0178   GlobalPoint corrected_position = GlobalPoint(position.x(), position.y(), position.z() - vertex_z);
0179   return corrected_position.eta();
0180 }
0181
0182 float MTDGeomUtil::eta(const DetId& id, const LocalPoint& local_point, const float& vertex_z) const {
0183   GlobalPoint position = globalPosition(id, local_point);
0184   float Eta = eta(position, vertex_z);
0185   return Eta;
0186 }
0187
0188 float MTDGeomUtil::phi(const GlobalPoint& position) const {
0189   float phi = (position.x() == 0 && position.y() == 0) ? 0 : atan2(position.y(), position.x());
0190   return phi;
0191 }
0192
0193 float MTDGeomUtil::phi(const DetId& id, const LocalPoint& local_point) const {
0194   GlobalPoint position = globalPosition(id, local_point);
0195   float phi = (position.x() == 0 && position.y() == 0) ? 0 : atan2(position.y(), position.x());
0196   return phi;
0197 }
0198
0199 float MTDGeomUtil::pt(const GlobalPoint& position, const float& hitEnergy, const float& vertex_z) const {
0200   float Eta = eta(position, vertex_z);
0201   float pt = hitEnergy / cosh(Eta);
0202   return pt;
0203 }
0204
0205 float MTDGeomUtil::pt(const DetId& id,
0206                       const LocalPoint& local_point,
0207                       const float& hitEnergy,
0208                       const float& vertex_z) const {
0209   GlobalPoint position = globalPosition(id, local_point);
0210   float Eta = eta(position, vertex_z);
0211   float pt = hitEnergy / cosh(Eta);
0212   return pt;
0213 }