Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-06-22 02:23:47 

0001 #include "Geometry/HGCalTBCommonData/interface/HGCalTBDDDConstants.h"
0002 
0003 #include "DataFormats/Math/interface/GeantUnits.h"
0004 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/Utilities/interface/Exception.h"
0007 #include "Geometry/HGCalCommonData/interface/HGCalGeometryMode.h"
0008 #include "Geometry/HGCalCommonData/interface/HGCalTypes.h"
0009 #include "Geometry/HGCalTBCommonData/interface/HGCalTBGeomParameters.h"
0010 
0011 #include <algorithm>
0012 #include <bitset>
0013 #include <iterator>
0014 #include <functional>
0015 #include <numeric>
0016 
0017 //#define EDM_ML_DEBUG
0018 using namespace geant_units::operators;
0019 
0020 HGCalTBDDDConstants::HGCalTBDDDConstants(const HGCalTBParameters* hp, const std::string& name)
0021     : hgpar_(hp), sqrt3_(std::sqrt(3.0)), mode_(hgpar_->mode_) {
0022 #ifdef EDM_ML_DEBUG
0023   edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants::Mode " << mode_;
0024 #endif
0025   if (waferHexagon6()) {
0026     rmax_ = (HGCalTBParameters::k_ScaleFromDDD * (hgpar_->waferR_) * std::cos(30._deg));
0027     rmaxT_ = rmax_ + 0.5 * hgpar_->sensorSeparation_;
0028     hexside_ = 2.0 * rmax_ * tan30deg_;
0029     hexsideT_ = 2.0 * rmaxT_ * tan30deg_;
0030 #ifdef EDM_ML_DEBUG
0031     edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants::rmax_ " << rmax_ << ":" << rmaxT_ << ":" << hexside_ << ":"
0032                                   << hexsideT_ << " CellSize "
0033                                   << 0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[0] << ":"
0034                                   << 0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[1];
0035 #endif
0036   }
0037   // init maps and constants
0038   modHalf_ = 0;
0039   maxWafersPerLayer_ = 0;
0040   for (int simreco = 0; simreco < 2; ++simreco) {
0041     tot_layers_[simreco] = layersInit((bool)simreco);
0042     max_modules_layer_[simreco].resize(tot_layers_[simreco] + 1);
0043     for (unsigned int layer = 1; layer <= tot_layers_[simreco]; ++layer) {
0044       max_modules_layer_[simreco][layer] = modulesInit(layer, (bool)simreco);
0045       if (simreco == 1) {
0046         modHalf_ += max_modules_layer_[simreco][layer];
0047         maxWafersPerLayer_ = std::max(maxWafersPerLayer_, max_modules_layer_[simreco][layer]);
0048 #ifdef EDM_ML_DEBUG
0049         edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants::Layer " << layer << " with "
0050                                       << max_modules_layer_[simreco][layer] << ":" << modHalf_ << " modules in RECO";
0051       } else {
0052         edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants::Layer " << layer << " with "
0053                                       << max_modules_layer_[simreco][layer] << " modules in SIM";
0054 #endif
0055       }
0056     }
0057 #ifdef EDM_ML_DEBUG
0058     edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants::SimReco " << simreco << " with " << tot_layers_[simreco]
0059                                   << " Layers";
0060 #endif
0061   }
0062   tot_wafers_ = wafers();
0063 
0064 #ifdef EDM_ML_DEBUG
0065   edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants initialized for " << name << " with " << layers(false) << ":"
0066                                 << layers(true) << " layers, " << wafers() << ":" << 2 * modHalf_
0067                                 << " wafers with maximum " << maxWafersPerLayer_ << " per layer and "
0068                                 << "maximum of " << maxCells(false) << ":" << maxCells(true) << " cells";
0069 #endif
0070   if (waferHexagon6()) {
0071     int wminT(9999999), wmaxT(-9999999), kount1(0), kount2(0);
0072     for (unsigned int i = 0; i < getTrFormN(); ++i) {
0073       int lay0 = getTrForm(i).lay;
0074       int wmin(9999999), wmax(-9999999), kount(0);
0075       for (int wafer = 0; wafer < sectors(); ++wafer) {
0076         bool waferIn = waferInLayer(wafer, lay0, true);
0077         if (waferIn) {
0078           if (wafer < wmin)
0079             wmin = wafer;
0080           if (wafer > wmax)
0081             wmax = wafer;
0082           ++kount;
0083         }
0084       }
0085       if (wminT > wmin)
0086         wminT = wmin;
0087       if (wmaxT < wmax)
0088         wmaxT = wmax;
0089       if (kount1 < kount)
0090         kount1 = kount;
0091       kount2 += kount;
0092 #ifdef EDM_ML_DEBUG
0093       int lay1 = getIndex(lay0, true).first;
0094       edm::LogVerbatim("HGCalGeom") << "Index " << i << " Layer " << lay0 << ":" << lay1 << " Wafer " << wmin << ":"
0095                                     << wmax << ":" << kount;
0096 #endif
0097       HGCWaferParam a1{{wmin, wmax, kount}};
0098       waferLayer_[lay0] = a1;
0099     }
0100     waferMax_ = std::array<int, 4>{{wminT, wmaxT, kount1, kount2}};
0101 #ifdef EDM_ML_DEBUG
0102     edm::LogVerbatim("HGCalGeom") << "Overall wafer statistics: " << wminT << ":" << wmaxT << ":" << kount1 << ":"
0103                                   << kount2;
0104 #endif
0105   }
0106 }
0107 
0108 std::pair<int, int> HGCalTBDDDConstants::assignCell(float x, float y, int lay, int subSec, bool reco) const {
0109   const auto& index = getIndex(lay, reco);
0110   if (index.first < 0)
0111     return std::make_pair(-1, -1);
0112   if (waferHexagon6()) {
0113     float xx = (reco) ? x : HGCalTBParameters::k_ScaleFromDDD * x;
0114     float yy = (reco) ? y : HGCalTBParameters::k_ScaleFromDDD * y;
0115 
0116     // First the wafer
0117     int wafer = cellHex(xx, yy, rmax_, hgpar_->waferPosX_, hgpar_->waferPosY_);
0118     if (wafer < 0 || wafer >= static_cast<int>(hgpar_->waferTypeT_.size())) {
0119       edm::LogWarning("HGCalGeom") << "Wafer no. out of bound for " << wafer << ":" << (hgpar_->waferTypeT_).size()
0120                                    << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0121                                    << " ***** ERROR *****";
0122       return std::make_pair(-1, -1);
0123     } else {
0124       // Now the cell
0125       xx -= hgpar_->waferPosX_[wafer];
0126       yy -= hgpar_->waferPosY_[wafer];
0127       if (hgpar_->waferTypeT_[wafer] == 1)
0128         return std::make_pair(wafer,
0129                               cellHex(xx,
0130                                       yy,
0131                                       0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[0],
0132                                       hgpar_->cellFineX_,
0133                                       hgpar_->cellFineY_));
0134       else
0135         return std::make_pair(wafer,
0136                               cellHex(xx,
0137                                       yy,
0138                                       0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
0139                                       hgpar_->cellCoarseX_,
0140                                       hgpar_->cellCoarseY_));
0141     }
0142   } else {
0143     return std::make_pair(-1, -1);
0144   }
0145 }
0146 
0147 double HGCalTBDDDConstants::cellSizeHex(int type) const {
0148   int indx = (type == 1) ? 1 : 0;
0149   double cell = 0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[indx];
0150   return cell;
0151 }
0152 
0153 double HGCalTBDDDConstants::cellThickness(int layer, int wafer) const {
0154   double thick(-1);
0155   int type = waferType(layer, wafer);
0156   if (type >= 0)
0157     thick = 100.0 * (type + 1);  // type = 1,2,3 for 100,200,300 micron
0158   return thick;
0159 }
0160 
0161 double HGCalTBDDDConstants::distFromEdgeHex(double x, double y, double z) const {
0162   // Assming the point is within a hexagonal plane of the wafer, calculate
0163   // the shortest distance from the edge
0164   if (z < 0)
0165     x = -x;
0166   double dist(0);
0167   // Input x, y in Geant4 unit and transformed to CMSSW standard
0168   double xx = HGCalTBParameters::k_ScaleFromDDD * x;
0169   double yy = HGCalTBParameters::k_ScaleFromDDD * y;
0170   int sizew = static_cast<int>(hgpar_->waferPosX_.size());
0171   int wafer = sizew;
0172   // Transform to the local coordinate frame of the wafer first
0173   for (int k = 0; k < sizew; ++k) {
0174     double dx = std::abs(xx - hgpar_->waferPosX_[k]);
0175     double dy = std::abs(yy - hgpar_->waferPosY_[k]);
0176     if ((dx <= rmax_) && (dy <= hexside_) && ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy)))) {
0177       wafer = k;
0178       xx -= hgpar_->waferPosX_[k];
0179       yy -= hgpar_->waferPosY_[k];
0180       break;
0181     }
0182   }
0183   // Look at only one quarter (both x,y are positive)
0184   if (wafer < sizew) {
0185     if (std::abs(yy) < 0.5 * hexside_) {
0186       dist = rmax_ - std::abs(xx);
0187     } else {
0188       dist = 0.5 * ((rmax_ - std::abs(xx)) - sqrt3_ * (std::abs(yy) - 0.5 * hexside_));
0189     }
0190   } else {
0191     dist = 0;
0192   }
0193   dist *= HGCalTBParameters::k_ScaleToDDD;
0194 #ifdef EDM_ML_DEBUG
0195   edm::LogVerbatim("HGCalGeom") << "DistFromEdgeHex: Local " << xx << ":" << yy << " wafer " << wafer << " flag "
0196                                 << (wafer < sizew) << " Distance " << rmax_ << ":" << (rmax_ - std::abs(xx)) << ":"
0197                                 << (std::abs(yy) - 0.5 * hexside_) << ":" << 0.5 * hexside_ << ":" << dist;
0198 #endif
0199   return dist;
0200 }
0201 
0202 int HGCalTBDDDConstants::getLayer(double z, bool reco) const {
0203   // Get the layer # from the gloabl z coordinate
0204   unsigned int k = 0;
0205   double zz = (reco ? std::abs(z) : HGCalTBParameters::k_ScaleFromDDD * std::abs(z));
0206   const auto& zLayerHex = hgpar_->zLayerHex_;
0207   auto itr = std::find_if(zLayerHex.begin() + 1, zLayerHex.end(), [&k, &zz, &zLayerHex](double zLayer) {
0208     ++k;
0209     return zz < 0.5 * (zLayerHex[k - 1] + zLayerHex[k]);
0210   });
0211   int lay = (itr == zLayerHex.end()) ? static_cast<int>(zLayerHex.size()) : k;
0212   if (waferHexagon6() && reco) {
0213     int indx = layerIndex(lay, false);
0214     if (indx >= 0)
0215       lay = hgpar_->layerGroupO_[indx];
0216   } else {
0217     lay += (hgpar_->firstLayer_ - 1);
0218   }
0219   return lay;
0220 }
0221 
0222 HGCalTBParameters::hgtrap HGCalTBDDDConstants::getModule(unsigned int indx, bool hexType, bool reco) const {
0223   HGCalTBParameters::hgtrap mytr;
0224   if (hexType) {
0225     if (indx >= hgpar_->waferTypeL_.size())
0226       edm::LogWarning("HGCalGeom") << "Wafer no. out bound for index " << indx << ":" << (hgpar_->waferTypeL_).size()
0227                                    << ":" << (hgpar_->waferPosX_).size() << ":" << (hgpar_->waferPosY_).size()
0228                                    << " ***** ERROR *****";
0229     unsigned int type =
0230         ((indx < hgpar_->waferTypeL_.size()) ? hgpar_->waferTypeL_[indx] - 1 : HGCalTBParameters::HGCalCoarseThick);
0231     mytr = hgpar_->getModule(type, reco);
0232   } else {
0233     mytr = hgpar_->getModule(indx, reco);
0234   }
0235   return mytr;
0236 }
0237 
0238 std::vector<HGCalTBParameters::hgtrap> HGCalTBDDDConstants::getModules() const {
0239   std::vector<HGCalTBParameters::hgtrap> mytrs;
0240   for (unsigned int k = 0; k < hgpar_->moduleLayR_.size(); ++k)
0241     mytrs.emplace_back(hgpar_->getModule(k, true));
0242   return mytrs;
0243 }
0244 
0245 std::vector<HGCalTBParameters::hgtrform> HGCalTBDDDConstants::getTrForms() const {
0246   std::vector<HGCalTBParameters::hgtrform> mytrs;
0247   for (unsigned int k = 0; k < hgpar_->trformIndex_.size(); ++k)
0248     mytrs.emplace_back(hgpar_->getTrForm(k));
0249   return mytrs;
0250 }
0251 
0252 bool HGCalTBDDDConstants::isHalfCell(int waferType, int cell) const {
0253   if (waferType < 1 || cell < 0)
0254     return false;
0255   return waferType == 2 ? hgpar_->cellCoarseHalf_[cell] : hgpar_->cellFineHalf_[cell];
0256 }
0257 
0258 bool HGCalTBDDDConstants::isValidHex(int lay, int mod, int cell, bool reco) const {
0259   // Check validity for a layer|wafer|cell of pre-TDR version
0260   bool result(false), resultMod(false);
0261   int cellmax(0);
0262   if (waferHexagon6()) {
0263     int32_t copyNumber = hgpar_->waferCopy_[mod];
0264     result = ((lay > 0 && lay <= static_cast<int>(layers(reco))));
0265     if (result) {
0266       const int32_t lay_idx = reco ? (lay - 1) * 3 + 1 : lay;
0267       const auto& the_modules = hgpar_->copiesInLayers_[lay_idx];
0268       auto moditr = the_modules.find(copyNumber);
0269       result = resultMod = (moditr != the_modules.end());
0270 #ifdef EDM_ML_DEBUG
0271       if (!result)
0272         edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants: Layer " << lay << ":" << lay_idx << " Copy "
0273                                       << copyNumber << ":" << mod << " Flag " << result;
0274 #endif
0275       if (result) {
0276         if (moditr->second >= 0) {
0277           if (mod >= static_cast<int>(hgpar_->waferTypeT_.size()))
0278             edm::LogWarning("HGCalGeom") << "Module no. out of bound for " << mod << " to be compared with "
0279                                          << (hgpar_->waferTypeT_).size() << " ***** ERROR *****";
0280           cellmax = ((hgpar_->waferTypeT_[mod] - 1 == HGCalTBParameters::HGCalFine)
0281                          ? static_cast<int>(hgpar_->cellFineX_.size())
0282                          : static_cast<int>(hgpar_->cellCoarseX_.size()));
0283           result = (cell >= 0 && cell <= cellmax);
0284         } else {
0285           result = isValidCell(lay_idx, mod, cell);
0286         }
0287       }
0288     }
0289   }
0290 
0291 #ifdef EDM_ML_DEBUG
0292   if (!result)
0293     edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants: Layer " << lay << ":"
0294                                   << (lay > 0 && (lay <= static_cast<int>(layers(reco)))) << " Module " << mod << ":"
0295                                   << resultMod << " Cell " << cell << ":" << cellmax << ":"
0296                                   << (cell >= 0 && cell <= cellmax) << ":" << maxCells(reco);
0297 #endif
0298   return result;
0299 }
0300 
0301 int HGCalTBDDDConstants::lastLayer(bool reco) const {
0302   return (hgpar_->firstLayer_ + tot_layers_[static_cast<int>(reco)] - 1);
0303 }
0304 
0305 int HGCalTBDDDConstants::layerIndex(int lay, bool reco) const {
0306   int ll = lay - hgpar_->firstLayer_;
0307   if (ll < 0 || ll >= static_cast<int>(hgpar_->layerIndex_.size()))
0308     return -1;
0309   if (waferHexagon6()) {
0310     if (reco && ll >= static_cast<int>(hgpar_->depthIndex_.size()))
0311       return -1;
0312     return (reco ? hgpar_->depthLayerF_[ll] : hgpar_->layerIndex_[ll]);
0313   }
0314   return -1;
0315 }
0316 
0317 unsigned int HGCalTBDDDConstants::layers(bool reco) const { return tot_layers_[static_cast<int>(reco)]; }
0318 
0319 unsigned int HGCalTBDDDConstants::layersInit(bool reco) const {
0320   return (reco ? hgpar_->depthIndex_.size() : hgpar_->layerIndex_.size());
0321 }
0322 
0323 std::pair<float, float> HGCalTBDDDConstants::locateCell(int cell, int lay, int type, bool reco) const {
0324   // type refers to wafer # for hexagon cell
0325   float x(999999.), y(999999.);
0326   const auto& index = getIndex(lay, reco);
0327   int i = index.first;
0328   if (i < 0)
0329     return std::make_pair(x, y);
0330   if (waferHexagon6()) {
0331     x = hgpar_->waferPosX_[type];
0332     y = hgpar_->waferPosY_[type];
0333 #ifdef EDM_ML_DEBUG
0334     float x0(x), y0(y);
0335 #endif
0336     if (hgpar_->waferTypeT_[type] - 1 == HGCalTBParameters::HGCalFine) {
0337       x += hgpar_->cellFineX_[cell];
0338       y += hgpar_->cellFineY_[cell];
0339     } else {
0340       x += hgpar_->cellCoarseX_[cell];
0341       y += hgpar_->cellCoarseY_[cell];
0342     }
0343 #ifdef EDM_ML_DEBUG
0344     edm::LogVerbatim("HGCalGeom") << "LocateCell (Wafer) " << x0 << ":" << y0 << " Final " << x << ":" << y;
0345 #endif
0346     if (!reco) {
0347       x *= HGCalTBParameters::k_ScaleToDDD;
0348       y *= HGCalTBParameters::k_ScaleToDDD;
0349     }
0350   }
0351   return std::make_pair(x, y);
0352 }
0353 
0354 std::pair<float, float> HGCalTBDDDConstants::locateCellHex(int cell, int wafer, bool reco) const {
0355   float x(0), y(0);
0356   if (hgpar_->waferTypeT_[wafer] - 1 == HGCalTBParameters::HGCalFine) {
0357     x = hgpar_->cellFineX_[cell];
0358     y = hgpar_->cellFineY_[cell];
0359   } else {
0360     x = hgpar_->cellCoarseX_[cell];
0361     y = hgpar_->cellCoarseY_[cell];
0362   }
0363   if (!reco) {
0364     x *= HGCalTBParameters::k_ScaleToDDD;
0365     y *= HGCalTBParameters::k_ScaleToDDD;
0366   }
0367   return std::make_pair(x, y);
0368 }
0369 
0370 int HGCalTBDDDConstants::maxCells(bool reco) const {
0371   int cells(0);
0372   for (unsigned int i = 0; i < layers(reco); ++i) {
0373     int lay = reco ? hgpar_->depth_[i] : hgpar_->layer_[i];
0374     if (cells < maxCells(lay, reco))
0375       cells = maxCells(lay, reco);
0376   }
0377   return cells;
0378 }
0379 
0380 int HGCalTBDDDConstants::maxCells(int lay, bool reco) const {
0381   const auto& index = getIndex(lay, reco);
0382   if ((index.first < 0) || (!waferHexagon6()))
0383     return 0;
0384   unsigned int cells(0);
0385   for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
0386     if (waferInLayerTest(k, index.first)) {
0387       unsigned int cell = (hgpar_->waferTypeT_[k] - 1 == HGCalTBParameters::HGCalFine) ? (hgpar_->cellFineX_.size())
0388                                                                                        : (hgpar_->cellCoarseX_.size());
0389       if (cell > cells)
0390         cells = cell;
0391     }
0392   }
0393   return static_cast<int>(cells);
0394 }
0395 
0396 int HGCalTBDDDConstants::maxRows(int lay, bool reco) const {
0397   int kymax(0);
0398   const auto& index = getIndex(lay, reco);
0399   int i = index.first;
0400   if ((i >= 0) && waferHexagon6()) {
0401     for (unsigned int k = 0; k < hgpar_->waferCopy_.size(); ++k) {
0402       if (waferInLayerTest(k, i)) {
0403         int ky = ((hgpar_->waferCopy_[k]) / 100) % 100;
0404         if (ky > kymax)
0405           kymax = ky;
0406       }
0407     }
0408   }
0409   return kymax;
0410 }
0411 
0412 int HGCalTBDDDConstants::modifyUV(int uv, int type1, int type2) const {
0413   // Modify u/v for transition of type1 to type2
0414   return (((type1 == type2) || (type1 * type2 != 0)) ? uv : ((type1 == 0) ? (2 * uv + 1) / 3 : (3 * uv) / 2));
0415 }
0416 
0417 int HGCalTBDDDConstants::modules(int lay, bool reco) const {
0418   if (getIndex(lay, reco).first < 0)
0419     return 0;
0420   else
0421     return max_modules_layer_[static_cast<int>(reco)][lay];
0422 }
0423 
0424 int HGCalTBDDDConstants::modulesInit(int lay, bool reco) const {
0425   int nmod(0);
0426   const auto& index = getIndex(lay, reco);
0427   if (index.first < 0)
0428     return nmod;
0429   for (unsigned int k = 0; k < hgpar_->waferPosX_.size(); ++k) {
0430     if (waferInLayerTest(k, index.first))
0431       ++nmod;
0432   }
0433   return nmod;
0434 }
0435 
0436 double HGCalTBDDDConstants::mouseBite(bool reco) const {
0437   return (reco ? hgpar_->mouseBite_ : HGCalTBParameters::k_ScaleToDDD * hgpar_->mouseBite_);
0438 }
0439 
0440 int HGCalTBDDDConstants::numberCells(bool reco) const {
0441   int cells(0);
0442   unsigned int nlayer = (reco) ? hgpar_->depth_.size() : hgpar_->layer_.size();
0443   for (unsigned k = 0; k < nlayer; ++k) {
0444     std::vector<int> ncells = numberCells(((reco) ? hgpar_->depth_[k] : hgpar_->layer_[k]), reco);
0445     cells = std::accumulate(ncells.begin(), ncells.end(), cells);
0446   }
0447   return cells;
0448 }
0449 
0450 std::vector<int> HGCalTBDDDConstants::numberCells(int lay, bool reco) const {
0451   const auto& index = getIndex(lay, reco);
0452   int i = index.first;
0453   std::vector<int> ncell;
0454   if ((i >= 0) && (waferHexagon6())) {
0455     for (unsigned int k = 0; k < hgpar_->waferTypeT_.size(); ++k) {
0456       if (waferInLayerTest(k, i)) {
0457         unsigned int cell = (hgpar_->waferTypeT_[k] - 1 == HGCalTBParameters::HGCalFine)
0458                                 ? (hgpar_->cellFineX_.size())
0459                                 : (hgpar_->cellCoarseX_.size());
0460         ncell.emplace_back(static_cast<int>(cell));
0461       }
0462     }
0463   }
0464   return ncell;
0465 }
0466 
0467 int HGCalTBDDDConstants::numberCellsHexagon(int wafer) const {
0468   if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferTypeT_.size())) {
0469     if (hgpar_->waferTypeT_[wafer] - 1 == 0)
0470       return static_cast<int>(hgpar_->cellFineX_.size());
0471     else
0472       return static_cast<int>(hgpar_->cellCoarseX_.size());
0473   } else {
0474     return 0;
0475   }
0476 }
0477 
0478 std::pair<double, double> HGCalTBDDDConstants::rangeR(double z, bool reco) const {
0479   double rmin(0), rmax(0);
0480   if (!reco) {
0481     rmin *= HGCalTBParameters::k_ScaleToDDD;
0482     rmax *= HGCalTBParameters::k_ScaleToDDD;
0483   }
0484 #ifdef EDM_ML_DEBUG
0485   edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants:rangeR: " << z << ":0"
0486                                 << " R " << rmin << ":" << rmax;
0487 #endif
0488   return std::make_pair(rmin, rmax);
0489 }
0490 
0491 std::pair<double, double> HGCalTBDDDConstants::rangeRLayer(int lay, bool reco) const {
0492   double rmin(0), rmax(0);
0493   const auto& index = getIndex(lay, reco);
0494   if (index.first >= 0 && index.first < static_cast<int>(hgpar_->rMinLayHex_.size())) {
0495     rmin = hgpar_->rMinLayHex_[index.first];
0496     rmax = hgpar_->rMaxLayHex_[index.first];
0497   }
0498   if (!reco) {
0499     rmin *= HGCalTBParameters::k_ScaleToDDD;
0500     rmax *= HGCalTBParameters::k_ScaleToDDD;
0501   }
0502 #ifdef EDM_ML_DEBUG
0503   edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants:rangeR: " << lay << ":" << index.first << " R " << rmin << ":"
0504                                 << rmax;
0505 #endif
0506   return std::make_pair(rmin, rmax);
0507 }
0508 
0509 std::pair<double, double> HGCalTBDDDConstants::rangeZ(bool reco) const {
0510   double zmin = (hgpar_->zLayerHex_[0] - hgpar_->waferThick_);
0511   double zmax = (hgpar_->zLayerHex_[hgpar_->zLayerHex_.size() - 1] + hgpar_->waferThick_);
0512 #ifdef EDM_ML_DEBUG
0513   edm::LogVerbatim("HGCalGeom") << "HGCalTBDDDConstants:rangeZ: " << zmin << ":" << zmax << ":" << hgpar_->waferThick_;
0514 #endif
0515   if (!reco) {
0516     zmin *= HGCalTBParameters::k_ScaleToDDD;
0517     zmax *= HGCalTBParameters::k_ScaleToDDD;
0518   }
0519   return std::make_pair(zmin, zmax);
0520 }
0521 
0522 std::pair<int, int> HGCalTBDDDConstants::rowColumnWafer(int wafer) const {
0523   int row(0), col(0);
0524   if (wafer < static_cast<int>(hgpar_->waferCopy_.size())) {
0525     int copy = hgpar_->waferCopy_[wafer];
0526     col = HGCalTypes::getUnpackedU(copy);
0527     row = HGCalTypes::getUnpackedV(copy);
0528   }
0529   return std::make_pair(row, col);
0530 }
0531 
0532 std::pair<int, int> HGCalTBDDDConstants::simToReco(int cell, int lay, int mod, bool half) const {
0533   return std::make_pair(cell, lay);
0534 }
0535 
0536 int HGCalTBDDDConstants::waferFromCopy(int copy) const {
0537   const int ncopies = hgpar_->waferCopy_.size();
0538   int wafer(ncopies);
0539   bool result(false);
0540   for (int k = 0; k < ncopies; ++k) {
0541     if (copy == hgpar_->waferCopy_[k]) {
0542       wafer = k;
0543       result = true;
0544       break;
0545     }
0546   }
0547   if (!result) {
0548     wafer = -1;
0549 #ifdef EDM_ML_DEBUG
0550     edm::LogVerbatim("HGCalGeom") << "Cannot find " << copy << " in a list of " << ncopies << " members";
0551     for (int k = 0; k < ncopies; ++k)
0552       edm::LogVerbatim("HGCalGeom") << "[" << k << "] " << hgpar_->waferCopy_[k];
0553 #endif
0554   }
0555 #ifdef EDM_ML_DEBUG
0556   edm::LogVerbatim("HGCalGeom") << "WaferFromCopy " << copy << ":" << wafer << ":" << result;
0557 #endif
0558   return wafer;
0559 }
0560 
0561 void HGCalTBDDDConstants::waferFromPosition(const double x, const double y, int& wafer, int& icell, int& celltyp) const {
0562   // Input x, y in Geant4 unit and transformed to CMSSW standard
0563   double xx = HGCalTBParameters::k_ScaleFromDDD * x;
0564   double yy = HGCalTBParameters::k_ScaleFromDDD * y;
0565   int size_ = static_cast<int>(hgpar_->waferCopy_.size());
0566   wafer = size_;
0567   for (int k = 0; k < size_; ++k) {
0568     double dx = std::abs(xx - hgpar_->waferPosX_[k]);
0569     double dy = std::abs(yy - hgpar_->waferPosY_[k]);
0570     if (dx <= rmax_ && dy <= hexside_) {
0571       if ((dy <= 0.5 * hexside_) || (dx * tan30deg_ <= (hexside_ - dy))) {
0572         wafer = k;
0573         celltyp = hgpar_->waferTypeT_[k];
0574         xx -= hgpar_->waferPosX_[k];
0575         yy -= hgpar_->waferPosY_[k];
0576         break;
0577       }
0578     }
0579   }
0580   if (wafer < size_) {
0581     if (celltyp - 1 == HGCalTBParameters::HGCalFine)
0582       icell = cellHex(xx,
0583                       yy,
0584                       0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[0],
0585                       hgpar_->cellFineX_,
0586                       hgpar_->cellFineY_);
0587     else
0588       icell = cellHex(xx,
0589                       yy,
0590                       0.5 * HGCalTBParameters::k_ScaleFromDDD * hgpar_->cellSize_[1],
0591                       hgpar_->cellCoarseX_,
0592                       hgpar_->cellCoarseY_);
0593   } else {
0594     wafer = -1;
0595 #ifdef EDM_ML_DEBUG
0596     edm::LogWarning("HGCalGeom") << "Cannot get wafer type corresponding to " << x << ":" << y << "    " << xx << ":"
0597                                  << yy;
0598 #endif
0599   }
0600 #ifdef EDM_ML_DEBUG
0601   edm::LogVerbatim("HGCalGeom") << "Position " << x << ":" << y << " Wafer " << wafer << ":" << size_ << " XX " << xx
0602                                 << ":" << yy << " Cell " << icell << " Type " << celltyp;
0603 #endif
0604 }
0605 
0606 bool HGCalTBDDDConstants::waferInLayer(int wafer, int lay, bool reco) const {
0607   const auto& indx = getIndex(lay, reco);
0608   if (indx.first < 0)
0609     return false;
0610   return waferInLayerTest(wafer, indx.first);
0611 }
0612 
0613 bool HGCalTBDDDConstants::waferFullInLayer(int wafer, int lay, bool reco) const {
0614   const auto& indx = getIndex(lay, reco);
0615   if (indx.first < 0)
0616     return false;
0617   return waferInLayerTest(wafer, indx.first);
0618 }
0619 
0620 std::pair<double, double> HGCalTBDDDConstants::waferParameters(bool reco) const {
0621   if (reco)
0622     return std::make_pair(rmax_, hexside_);
0623   else
0624     return std::make_pair(HGCalTBParameters::k_ScaleToDDD * rmax_, HGCalTBParameters::k_ScaleToDDD * hexside_);
0625 }
0626 
0627 std::pair<double, double> HGCalTBDDDConstants::waferPosition(int wafer, bool reco) const {
0628   double xx(0), yy(0);
0629   if (wafer >= 0 && wafer < static_cast<int>(hgpar_->waferPosX_.size())) {
0630     xx = hgpar_->waferPosX_[wafer];
0631     yy = hgpar_->waferPosY_[wafer];
0632   }
0633   if (!reco) {
0634     xx *= HGCalTBParameters::k_ScaleToDDD;
0635     yy *= HGCalTBParameters::k_ScaleToDDD;
0636   }
0637   return std::make_pair(xx, yy);
0638 }
0639 
0640 int HGCalTBDDDConstants::wafers() const { return static_cast<int>(hgpar_->moduleLayR_.size()); }
0641 
0642 int HGCalTBDDDConstants::wafers(int layer, int type) const {
0643   int wafer(0);
0644   auto itr = waferLayer_.find(layer);
0645   if (itr != waferLayer_.end()) {
0646     unsigned ity = (type > 0 && type <= 2) ? type : 0;
0647     wafer = (itr->second)[ity];
0648   }
0649   return wafer;
0650 }
0651 
0652 int HGCalTBDDDConstants::waferType(DetId const& id) const {
0653   return waferType(HGCalDetId(id).layer(), HGCalDetId(id).wafer());
0654 }
0655 
0656 int HGCalTBDDDConstants::waferType(int layer, int wafer) const {
0657   int type(HGCalTBParameters::HGCalCoarseThick);
0658   if ((wafer >= 0) && (wafer < static_cast<int>(hgpar_->waferTypeL_.size())))
0659     type = (hgpar_->waferTypeL_[wafer] - 1);
0660   return type;
0661 }
0662 
0663 bool HGCalTBDDDConstants::waferVirtual(int layer, int wafer) const {
0664   int wl = HGCalWaferIndex::waferIndex(layer, wafer, 0, true);
0665   return (hgpar_->waferTypes_.find(wl) != hgpar_->waferTypes_.end());
0666 }
0667 
0668 double HGCalTBDDDConstants::waferZ(int lay, bool reco) const {
0669   const auto& index = getIndex(lay, reco);
0670   if (index.first < 0)
0671     return 0;
0672   else
0673     return (reco ? hgpar_->zLayerHex_[index.first] : HGCalTBParameters::k_ScaleToDDD * hgpar_->zLayerHex_[index.first]);
0674 }
0675 
0676 int HGCalTBDDDConstants::cellHex(
0677     double xx, double yy, const double& cellR, const std::vector<double>& posX, const std::vector<double>& posY) const {
0678   int num(0);
0679   const double tol(0.00001);
0680   double cellY = 2.0 * cellR * tan30deg_;
0681   for (unsigned int k = 0; k < posX.size(); ++k) {
0682     double dx = std::abs(xx - posX[k]);
0683     double dy = std::abs(yy - posY[k]);
0684     if (dx <= (cellR + tol) && dy <= (cellY + tol)) {
0685       double xmax = (dy <= 0.5 * cellY) ? cellR : (cellR - (dy - 0.5 * cellY) / tan30deg_);
0686       if (dx <= (xmax + tol)) {
0687         num = k;
0688         break;
0689       }
0690     }
0691   }
0692   return num;
0693 }
0694 
0695 std::pair<int, float> HGCalTBDDDConstants::getIndex(int lay, bool reco) const {
0696   int indx = layerIndex(lay, reco);
0697   if (indx < 0)
0698     return std::make_pair(-1, 0);
0699   float cell(0);
0700   if (waferHexagon6()) {
0701     cell = (reco ? hgpar_->moduleCellR_[0] : hgpar_->moduleCellS_[0]);
0702   }
0703   return std::make_pair(indx, cell);
0704 }
0705 
0706 int HGCalTBDDDConstants::layerFromIndex(int index, bool reco) const {
0707   int ll(-1);
0708   if (waferHexagon6() && reco) {
0709     ll = static_cast<int>(std::find(hgpar_->depthLayerF_.begin(), hgpar_->depthLayerF_.end(), index) -
0710                           hgpar_->depthLayerF_.begin());
0711     if (ll == static_cast<int>(hgpar_->depthLayerF_.size()))
0712       ll = -1;
0713   }
0714 #ifdef EDM_ML_DEBUG
0715   edm::LogVerbatim("HGCalGeom") << "LayerFromIndex for " << index << ":" << reco << ":" << waferHexagon6() << " is"
0716                                 << ll << ":" << (ll + hgpar_->firstLayer_);
0717 #endif
0718   return ((ll < 0) ? ll : (ll + hgpar_->firstLayer_));
0719 }
0720 
0721 bool HGCalTBDDDConstants::isValidCell(int lay, int wafer, int cell) const {
0722   // Calculate the position of the cell
0723   // Works for options HGCalHexagon/HGCalHexagonFull
0724   double x = hgpar_->waferPosX_[wafer];
0725   double y = hgpar_->waferPosY_[wafer];
0726   if (hgpar_->waferTypeT_[wafer] - 1 == HGCalTBParameters::HGCalFine) {
0727     x += hgpar_->cellFineX_[cell];
0728     y += hgpar_->cellFineY_[cell];
0729   } else {
0730     x += hgpar_->cellCoarseX_[cell];
0731     y += hgpar_->cellCoarseY_[cell];
0732   }
0733   double rr = sqrt(x * x + y * y);
0734   bool result = ((rr >= hgpar_->rMinLayHex_[lay - 1]) && (rr <= hgpar_->rMaxLayHex_[lay - 1]) &&
0735                  (wafer < static_cast<int>(hgpar_->waferPosX_.size())));
0736 #ifdef EDM_ML_DEBUG
0737   if (!result)
0738     edm::LogVerbatim("HGCalGeom") << "Input " << lay << ":" << wafer << ":" << cell << " Position " << x << ":" << y
0739                                   << ":" << rr << " Compare Limits " << hgpar_->rMinLayHex_[lay - 1] << ":"
0740                                   << hgpar_->rMaxLayHex_[lay - 1] << " Flag " << result;
0741 #endif
0742   return result;
0743 }
0744 
0745 int32_t HGCalTBDDDConstants::waferIndex(int wafer, int index) const {
0746   int layer = layerFromIndex(index, true);
0747   int indx = HGCalWaferIndex::waferIndex(layer, wafer, 0);
0748 #ifdef EDM_ML_DEBUG
0749   edm::LogVerbatim("HGCalGeom") << "WaferIndex for " << wafer << ":" << index << " (" << layer << ") " << indx;
0750 #endif
0751   return indx;
0752 }
0753 
0754 #include "FWCore/Utilities/interface/typelookup.h"
0755 
0756 TYPELOOKUP_DATA_REG(HGCalTBDDDConstants);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!