Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-04 00:15:38

0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "Geometry/HGCalGeometry/interface/FastTimeGeometryLoader.h"
0003 #include "Geometry/HGCalGeometry/interface/FastTimeGeometry.h"
0004 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0005 #include "Geometry/CaloGeometry/interface/FlatTrd.h"
0006 #include "Geometry/HGCalCommonData/interface/FastTimeDDDConstants.h"
0007 #include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
0008 
0009 #include <sstream>
0010 
0011 //#define EDM_ML_DEBUG
0012 
0013 typedef CaloCellGeometry::CCGFloat CCGFloat;
0014 typedef std::vector<float> ParmVec;
0015 
0016 FastTimeGeometryLoader::FastTimeGeometryLoader() {}
0017 FastTimeGeometryLoader::~FastTimeGeometryLoader() {}
0018 
0019 FastTimeGeometry* FastTimeGeometryLoader::build(const FastTimeTopology& topology) {
0020   // allocate geometry
0021   FastTimeGeometry* geom = new FastTimeGeometry(topology);
0022   unsigned int numberOfCells = topology.totalGeomModules();  // both sides
0023   int detType = topology.detectorType();
0024 #ifdef EDM_ML_DEBUG
0025   edm::LogVerbatim("FastTimeGeom") << "Number of Cells " << numberOfCells << " for type " << detType
0026                                    << " of sub-detector " << topology.subDetector() << " Shape parameters "
0027                                    << FastTimeGeometry::k_NumberOfShapes << ":"
0028                                    << FastTimeGeometry::k_NumberOfParametersPerShape;
0029 #endif
0030   geom->allocateCorners(numberOfCells);
0031   geom->allocatePar(FastTimeGeometry::k_NumberOfShapes, FastTimeGeometry::k_NumberOfParametersPerShape);
0032 
0033   // loop over modules
0034   ParmVec params(FastTimeGeometry::k_NumberOfParametersPerShape, 0);
0035   unsigned int counter(0);
0036 #ifdef EDM_ML_DEBUG
0037   edm::LogVerbatim("FastTimeGeom") << "FastTimeGeometryLoader with # of transformation matrices " << numberOfCells;
0038 #endif
0039   for (unsigned itr = 0; itr < numberOfCells; ++itr) {
0040     int zside = (itr == 0) ? 1 : -1;
0041 #ifdef EDM_ML_DEBUG
0042     edm::LogVerbatim("FastTimeGeom") << "FastTimeGeometryLoader:: Z:Layer:Type " << zside << ":" << detType;
0043 #endif
0044     double zv = zside * (topology.dddConstants().getZPos(detType));
0045     const CLHEP::HepRep3x3 rotation =
0046         (zside > 0) ? CLHEP::HepRep3x3(1, 0, 0, 0, 1, 0, 0, 0, 1) : CLHEP::HepRep3x3(-1, 0, 0, 0, 1, 0, 0, 0, -1);
0047     const CLHEP::HepRotation hr(rotation);
0048     const CLHEP::Hep3Vector h3v(0, 0, zv);
0049     const HepGeom::Transform3D ht3d(hr, h3v);
0050     DetId detId = (DetId)(FastTimeDetId(detType, 0, 0, zside));
0051 #ifdef EDM_ML_DEBUG
0052     edm::LogVerbatim("FastTimeGeom") << "FastTimeGeometryLoader:: transf " << ht3d.getTranslation() << " and "
0053                                      << ht3d.getRotation();
0054 #endif
0055     params[0] = topology.dddConstants().getZHalf(detType);
0056     params[1] = params[2] = 0;
0057     params[3] = params[7] = topology.dddConstants().getRin(detType);
0058     params[4] = params[8] = topology.dddConstants().getRout(detType);
0059     params[5] = params[9] = topology.dddConstants().getRout(detType);
0060     params[6] = params[10] = 0;
0061     params[11] = zside;
0062     buildGeom(params, ht3d, detId, topology, geom);
0063     counter++;
0064   }
0065 
0066   geom->sortDetIds();
0067 
0068   if (counter != numberOfCells) {
0069     edm::LogWarning("FastTimeGeom") << "inconsistent # of cells: expected " << numberOfCells << " , inited " << counter;
0070     assert(counter == numberOfCells);
0071   }
0072 
0073   return geom;
0074 }
0075 
0076 void FastTimeGeometryLoader::buildGeom(const ParmVec& params,
0077                                        const HepGeom::Transform3D& ht3d,
0078                                        const DetId& detId,
0079                                        const FastTimeTopology& topology,
0080                                        FastTimeGeometry* geom) {
0081 #ifdef EDM_ML_DEBUG
0082   std::ostringstream st1;
0083   st1 << "Volume Parameters";
0084   for (unsigned int i = 0; i < 12; ++i)
0085     st1 << " : " << params[i];
0086   edm::LogVerbatim("FastTimeGeom") << st1.str();
0087 #endif
0088   FastTimeDetId id = FastTimeDetId(detId);
0089   std::vector<GlobalPoint> corners = topology.dddConstants().getCorners(id.type(), 1, 1, id.zside());
0090 
0091   FlatTrd::createCorners(params, ht3d, corners);
0092 
0093   const CCGFloat* parmPtr(CaloCellGeometry::getParmPtr(params, geom->parMgr(), geom->parVecVec()));
0094 
0095   GlobalPoint front(0.25 * (corners[0].x() + corners[1].x() + corners[2].x() + corners[3].x()),
0096                     0.25 * (corners[0].y() + corners[1].y() + corners[2].y() + corners[3].y()),
0097                     0.25 * (corners[0].z() + corners[1].z() + corners[2].z() + corners[3].z()));
0098 
0099   GlobalPoint back(0.25 * (corners[4].x() + corners[5].x() + corners[6].x() + corners[7].x()),
0100                    0.25 * (corners[4].y() + corners[5].y() + corners[6].y() + corners[7].y()),
0101                    0.25 * (corners[4].z() + corners[5].z() + corners[6].z() + corners[7].z()));
0102 
0103   if (front.mag2() > back.mag2()) {  // front should always point to the center, so swap front and back
0104     std::swap(front, back);
0105     std::swap_ranges(corners.begin(), corners.begin() + 4, corners.begin() + 4);
0106   }
0107 
0108   geom->newCell(front, back, corners[0], parmPtr, detId);
0109 }