Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-02-04 00:11:15

0001 /* Implementation of the  GEMGeometryParsFromDD Class
0002  *  Build the GEMGeometry from the DDD and DD4hep description
0003  *  
0004  *  DD4hep part added to the original old file (DD version) made by M. Maggi (INFN Bari)
0005  *  Author:  Sergio Lo Meo (sergio.lo.meo@cern.ch) 
0006  *  Created:  Mon, 15 Feb 2021 
0007  *
0008  */
0009 #include "Geometry/GEMGeometryBuilder/interface/GEMGeometryParsFromDD.h"
0010 #include "DataFormats/MuonDetId/interface/GEMDetId.h"
0011 
0012 #include "CondFormats/GeometryObjects/interface/RecoIdealGeometry.h"
0013 #include "DetectorDescription/Core/interface/DDFilter.h"
0014 #include "DetectorDescription/Core/interface/DDFilteredView.h"
0015 #include "DetectorDescription/Core/interface/DDSolid.h"
0016 
0017 #include "Geometry/MuonNumbering/interface/MuonGeometryNumbering.h"
0018 #include "Geometry/MuonNumbering/interface/MuonBaseNumber.h"
0019 #include "Geometry/MuonNumbering/interface/GEMNumberingScheme.h"
0020 #include "Geometry/MuonNumbering/interface/MuonGeometryConstants.h"
0021 
0022 #include "DataFormats/GeometryVector/interface/Basic3DVector.h"
0023 
0024 #include "DetectorDescription/DDCMS/interface/DDFilteredView.h"
0025 #include "DetectorDescription/DDCMS/interface/DDCompactView.h"
0026 #include "DetectorDescription/DDCMS/interface/DDSpecParRegistry.h"
0027 
0028 #include <iostream>
0029 #include <algorithm>
0030 
0031 GEMGeometryParsFromDD::GEMGeometryParsFromDD() {}
0032 
0033 GEMGeometryParsFromDD::~GEMGeometryParsFromDD() {}
0034 
0035 // DDD
0036 
0037 void GEMGeometryParsFromDD::build(const DDCompactView* cview,
0038                                   const MuonGeometryConstants& muonConstants,
0039                                   RecoIdealGeometry& rgeo) {
0040   std::string attribute = "MuStructure";
0041   std::string value = "MuonEndCapGEM";
0042 
0043   // Asking only for the MuonGEM's
0044   DDSpecificsMatchesValueFilter filter{DDValue(attribute, value, 0.0)};
0045   DDFilteredView fv(*cview, filter);
0046   DDFilteredView fv2(*cview, filter);
0047 
0048   this->buildGeometry(fv, fv2, muonConstants, rgeo);
0049 }
0050 
0051 void GEMGeometryParsFromDD::buildGeometry(DDFilteredView& fv,
0052                                           DDFilteredView& fvGE2,
0053                                           const MuonGeometryConstants& muonConstants,
0054                                           RecoIdealGeometry& rgeo) {
0055   LogDebug("GEMGeometryParsFromDD") << "Building the geometry service";
0056   LogDebug("GEMGeometryParsFromDD") << "About to run through the GEM structure\n"
0057                                     << " First logical part " << fv.logicalPart().name().name();
0058 
0059   edm::LogVerbatim("GEMGeometryParsFromDD") << "(0) GEMGeometryParsFromDD - DDD ";
0060   MuonGeometryNumbering muonDDDNumbering(muonConstants);
0061   GEMNumberingScheme gemNumbering(muonConstants);
0062 
0063   // Check for the demonstrator geometry (only 1 chamber of GE2/1)
0064   int nGE21 = 0;
0065   bool doSuper = fvGE2.firstChild();
0066   while (doSuper) {
0067     // getting chamber id from eta partitions
0068     fvGE2.firstChild();
0069     doSuper = fvGE2.firstChild();
0070     if (doSuper) {
0071       int rawidCh = gemNumbering.baseNumberToUnitNumber(muonDDDNumbering.geoHistoryToBaseNumber(fvGE2.geoHistory()));
0072       GEMDetId detIdCh = GEMDetId(rawidCh);
0073       if (detIdCh.station() == 2)
0074         nGE21++;
0075 
0076       // back to chambers
0077       fvGE2.parent();
0078       fvGE2.parent();
0079       doSuper = (nGE21 < 2 && fvGE2.nextSibling());
0080     } else {
0081       edm::LogError("GEMGeometryParsFromDD") << "Failed to find next child volume. Cannot determine presence of GE 2/1";
0082     }
0083   }
0084   bool demonstratorGeometry = nGE21 == 1;
0085 
0086 #ifdef EDM_ML_DEBUG
0087   edm::LogVerbatim("Geometry") << "Found " << nGE21 << " GE2/1 chambers. Demonstrator geometry on? "
0088                                << demonstratorGeometry;
0089 #endif
0090 
0091   doSuper = fv.firstChild();
0092 
0093   LogDebug("GEMGeometryParsFromDD") << "doSuperChamber = " << doSuper;
0094   // loop over superchambers
0095   while (doSuper) {
0096     // getting chamber id from eta partitions
0097     fv.firstChild();
0098     doSuper = fv.firstChild();
0099     if (doSuper) {
0100       GEMDetId detIdCh =
0101           GEMDetId(gemNumbering.baseNumberToUnitNumber(muonDDDNumbering.geoHistoryToBaseNumber(fv.geoHistory())));
0102       // back to chambers
0103       fv.parent();
0104       fv.parent();
0105 
0106       // currently there is no superchamber in the geometry
0107       // only 2 chambers are present separated by a gap.
0108       // making superchamber out of the first chamber layer including the gap between chambers
0109 
0110       // In Run 3 we also have a single GE2/1 chamber at layer 2. We
0111       // make sure the superchamber gets built but also we build on the
0112       // first layer for the other stations so the superchamber is in
0113       // the right position there.
0114       if ((detIdCh.layer() == 1) || (detIdCh.layer() == 2 and detIdCh.station() == 2 and demonstratorGeometry)) {
0115         buildSuperChamber(fv, detIdCh, rgeo);
0116       }
0117       buildChamber(fv, detIdCh, rgeo);
0118 
0119       // loop over chambers
0120       // only 1 chamber
0121       bool doChambers = fv.firstChild();
0122       while (doChambers) {
0123         // loop over GEMEtaPartitions
0124         bool doEtaPart = fv.firstChild();
0125         while (doEtaPart) {
0126           GEMDetId detId =
0127               GEMDetId(gemNumbering.baseNumberToUnitNumber(muonDDDNumbering.geoHistoryToBaseNumber(fv.geoHistory())));
0128           buildEtaPartition(fv, detId, rgeo);
0129 
0130           doEtaPart = fv.nextSibling();
0131         }
0132         fv.parent();
0133         doChambers = fv.nextSibling();
0134       }
0135       fv.parent();
0136       doSuper = fv.nextSibling();
0137     } else {
0138       edm::LogError("GEMGeometryParsFromDD") << "Failed to find next child volume. Cannot build GEM chambers.";
0139     }
0140   }
0141 }
0142 
0143 void GEMGeometryParsFromDD::buildSuperChamber(DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0144   LogDebug("GEMGeometryParsFromDD") << "buildSuperChamber " << fv.logicalPart().name().name() << " " << detId
0145                                     << std::endl;
0146 
0147   DDBooleanSolid solid = (DDBooleanSolid)(fv.logicalPart().solid());
0148   std::vector<double> dpar = solid.solidA().parameters();
0149 
0150   GEMDetId gemid = detId.superChamberId();
0151 
0152   double dy = dpar[0];   //length is along local Y
0153   double dz = dpar[3];   // thickness is long local Z
0154   double dx1 = dpar[4];  // bottom width is along local X
0155   double dx2 = dpar[8];  // top width is along local X
0156   dpar = solid.solidB().parameters();
0157 
0158   dz += dpar[3];  // chamber thickness
0159   dz *= 2;        // 2 chambers in superchamber
0160   dz += 2.105;    // gap between chambers
0161 
0162   std::vector<double> pars{dx1, dx2, dy, dz};
0163   std::vector<double> vtra = getTranslation(fv);
0164   std::vector<double> vrot = getRotation(fv);
0165 
0166   LogDebug("GEMGeometryParsFromDD") << "dimension dx1 " << dx1 << ", dx2 " << dx2 << ", dy " << dy << ", dz " << dz;
0167   edm::LogVerbatim("GEMGeometryParsFromDD")
0168       << "(3) DDD, SuperChamber DetID " << gemid.rawId() << " Name " << fv.logicalPart().name().name() << " dx1 " << dx1
0169       << " dx2 " << dx2 << " dy " << dy << " dz " << dz;
0170   rgeo.insert(gemid.rawId(), vtra, vrot, pars, {fv.logicalPart().name().name()});
0171 }
0172 
0173 void GEMGeometryParsFromDD::buildChamber(DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0174   LogDebug("GEMGeometryParsFromDD") << "buildChamber " << fv.logicalPart().name().name() << " " << detId << std::endl;
0175 
0176   DDBooleanSolid solid = (DDBooleanSolid)(fv.logicalPart().solid());
0177   std::vector<double> dpar = solid.solidA().parameters();
0178 
0179   double dy = dpar[0];   //length is along local Y
0180   double dz = dpar[3];   // thickness is long local Z
0181   double dx1 = dpar[4];  // bottom width is along local X
0182   double dx2 = dpar[8];  // top width is along local X
0183   dpar = solid.solidB().parameters();
0184   dz += dpar[3];  // chamber thickness
0185 
0186   GEMDetId gemid = detId.chamberId();
0187 
0188   std::vector<double> pars{dx1, dx2, dy, dz};
0189   std::vector<double> vtra = getTranslation(fv);
0190   std::vector<double> vrot = getRotation(fv);
0191 
0192   LogDebug("GEMGeometryParsFromDD") << "dimension dx1 " << dx1 << ", dx2 " << dx2 << ", dy " << dy << ", dz " << dz;
0193   edm::LogVerbatim("GEMGeometryParsFromDD")
0194       << "(4) DDD, Chamber DetID " << gemid.rawId() << " Name " << fv.logicalPart().name().name() << " dx1 " << dx1
0195       << " dx2 " << dx2 << " dy " << dy << " dz " << dz;
0196   rgeo.insert(gemid.rawId(), vtra, vrot, pars, {fv.logicalPart().name().name()});
0197 }
0198 
0199 void GEMGeometryParsFromDD::buildEtaPartition(DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0200   LogDebug("GEMGeometryParsFromDD") << "buildEtaPartition " << fv.logicalPart().name().name() << " " << detId
0201                                     << std::endl;
0202 
0203   // EtaPartition specific parameter (nstrips and npads)
0204   DDValue numbOfStrips("nStrips");
0205   DDValue numbOfPads("nPads");
0206   DDValue delPhi("dPhi");
0207   const std::vector<const DDsvalues_type*>& specs = fv.specifics();
0208   double nStrips = 0., nPads = 0., dPhi = 0.;
0209   for (auto const& is : specs) {
0210     if (DDfetch(is, numbOfStrips))
0211       nStrips = numbOfStrips.doubles()[0];
0212     if (DDfetch(is, numbOfPads))
0213       nPads = numbOfPads.doubles()[0];
0214     if (DDfetch(is, delPhi))
0215       dPhi = delPhi.doubles()[0];
0216   }
0217   LogDebug("GEMGeometryParsFromDD") << ((nStrips == 0.) ? ("No nStrips found!!")
0218                                                         : ("Number of strips: " + std::to_string(nStrips)));
0219   LogDebug("GEMGeometryParsFromDD") << ((nPads == 0.) ? ("No nPads found!!")
0220                                                       : ("Number of pads: " + std::to_string(nPads)));
0221 
0222   // EtaPartition specific parameter (size)
0223   std::vector<double> dpar = fv.logicalPart().solid().parameters();
0224 
0225   double dy = dpar[0];   //length is along local Y
0226   double dz = dpar[3];   //0.4;// thickness is long local Z
0227   double dx1 = dpar[4];  // bottom width is along local X
0228   double dx2 = dpar[8];  // top width is along local X
0229 
0230   std::vector<double> pars{dx1, dx2, dy, dz, nStrips, nPads, dPhi};
0231   std::vector<double> vtra = getTranslation(fv);
0232   std::vector<double> vrot = getRotation(fv);
0233 
0234   LogDebug("GEMGeometryParsFromDD") << " dx1 " << dx1 << " dx2 " << dx2 << " dy " << dy << " dz " << dz << " nStrips "
0235                                     << nStrips << " nPads " << nPads << " dPhi " << dPhi;
0236 
0237   edm::LogVerbatim("GEMGeometryParsFromDD")
0238       << "(5) DDD, Eta Partion DetID " << detId.rawId() << " Name " << fv.logicalPart().name().name() << " dx1 " << dx1
0239       << " dx2 " << dx2 << " dy " << dy << " dz " << dz << " nStrips " << nStrips << " nPads " << nPads << " dPhi "
0240       << dPhi;
0241   rgeo.insert(detId.rawId(), vtra, vrot, pars, {fv.logicalPart().name().name()});
0242 }
0243 
0244 std::vector<double> GEMGeometryParsFromDD::getTranslation(DDFilteredView& fv) {
0245   const DDTranslation& tran = fv.translation();
0246   edm::LogVerbatim("GEMGeometryParsFromDD")
0247       << "(1) DDD, tran vector " << tran.x() << "  " << tran.y() << "  " << tran.z();
0248   return {tran.x(), tran.y(), tran.z()};
0249 }
0250 
0251 std::vector<double> GEMGeometryParsFromDD::getRotation(DDFilteredView& fv) {
0252   const DDRotationMatrix& rota = fv.rotation();  //.Inverse();
0253   DD3Vector x, y, z;
0254   rota.GetComponents(x, y, z);
0255   edm::LogVerbatim("GEMGeometryParsFromDD")
0256       << "(2) DDD, rot matrix " << x.X() << "  " << x.Y() << "  " << x.Z() << " " << y.X() << "  " << y.Y() << "  "
0257       << y.Z() << " " << z.X() << "  " << z.Y() << "  " << z.Z();
0258   return {x.X(), x.Y(), x.Z(), y.X(), y.Y(), y.Z(), z.X(), z.Y(), z.Z()};
0259 }
0260 
0261 // DD4hep
0262 
0263 void GEMGeometryParsFromDD::build(const cms::DDCompactView* cview,
0264                                   const MuonGeometryConstants& muonConstants,
0265                                   RecoIdealGeometry& rgeo) {
0266   std::string attribute = "MuStructure";
0267   std::string value = "MuonEndCapGEM";
0268 
0269   const cms::DDFilter filter(attribute, value);
0270   cms::DDFilteredView fv(*cview, filter);
0271 
0272   this->buildGeometry(fv, muonConstants, rgeo);
0273 }
0274 
0275 void GEMGeometryParsFromDD::buildGeometry(cms::DDFilteredView& fv,
0276                                           const MuonGeometryConstants& muonConstants,
0277                                           RecoIdealGeometry& rgeo) {
0278   edm::LogVerbatim("GEMGeometryParsFromDD") << "(0) GEMGeometryParsFromDD - DD4hep ";
0279 
0280   MuonGeometryNumbering mdddnum(muonConstants);
0281   GEMNumberingScheme gemNum(muonConstants);
0282   static constexpr uint32_t levelChamb = 7;
0283   int chamb(0), region(0);
0284   int theLevelPart = muonConstants.getValue("level");
0285   int theRingLevel = muonConstants.getValue("mg_ring") / theLevelPart;
0286   int theSectorLevel = muonConstants.getValue("mg_sector") / theLevelPart;
0287 
0288   // Check for the demonstrator geometry (only 1 chamber of GE2/1)
0289   auto start = fv.copyNos();
0290   int nGE21 = 0;
0291   while (nGE21 < 2 && fv.firstChild()) {
0292     const auto& history = fv.history();
0293     MuonBaseNumber num(mdddnum.geoHistoryToBaseNumber(history));
0294     GEMDetId detId(gemNum.baseNumberToUnitNumber(num));
0295     if (fv.level() == levelChamb && detId.station() == 2) {
0296       nGE21++;
0297     }
0298   }
0299   bool demonstratorGeometry = nGE21 == 1;
0300 #ifdef EDM_ML_DEBUG
0301   edm::LogVerbatim("Geometry") << "Found " << nGE21 << " GE2/1 chambers. Demonstrator geometry on? "
0302                                << demonstratorGeometry;
0303 #endif
0304 
0305   fv.goTo(start);
0306   while (fv.firstChild()) {
0307     const auto& history = fv.history();
0308     MuonBaseNumber num(mdddnum.geoHistoryToBaseNumber(history));
0309     GEMDetId detId(gemNum.baseNumberToUnitNumber(num));
0310 
0311     if (detId.station() == GEMDetId::minStationId0) {
0312       if (num.getLevels() == theRingLevel) {
0313         if (detId.region() != region) {
0314           region = detId.region();
0315           chamb = 0;
0316         }
0317         ++chamb;
0318         detId = GEMDetId(detId.region(), detId.ring(), detId.station(), detId.layer(), chamb, 0);
0319         buildSuperChamber(fv, detId, rgeo);
0320       } else if (num.getLevels() == theSectorLevel) {
0321         buildChamber(fv, detId, rgeo);
0322       } else {
0323         buildEtaPartition(fv, detId, rgeo);
0324       }
0325     } else {
0326       if (fv.level() == levelChamb) {
0327         if ((detId.layer() == 1) || (detId.layer() == 2 and detId.station() == 2 and demonstratorGeometry)) {
0328           buildSuperChamber(fv, detId, rgeo);
0329         }
0330         buildChamber(fv, detId, rgeo);
0331       } else if (num.getLevels() > theSectorLevel) {
0332         buildEtaPartition(fv, detId, rgeo);
0333       }
0334     }
0335   }
0336 }
0337 
0338 void GEMGeometryParsFromDD::buildSuperChamber(cms::DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0339   cms::DDSolid solid(fv.solid());
0340   auto solidA = solid.solidA();
0341   std::vector<double> dpar = solidA.dimensions();
0342 
0343   double dy = dpar[3] / dd4hep::mm;   //length is along local Y
0344   double dz = dpar[2] / dd4hep::mm;   // thickness is long local Z
0345   double dx1 = dpar[0] / dd4hep::mm;  // bottom width is along local X
0346   double dx2 = dpar[1] / dd4hep::mm;  // top width is along loc
0347 
0348   auto solidB = solid.solidB();
0349   dpar = solidB.dimensions();
0350   const int nch = 2;
0351   const double chgap = 2.105;
0352 
0353   GEMDetId gemid = detId.superChamberId();
0354   std::string_view name = fv.name();
0355 
0356   dz += (dpar[2] / dd4hep::mm);  // chamber thickness
0357   dz *= nch;                     // 2 chambers in superchamber
0358   dz += chgap;                   // gap between chambers
0359 
0360   std::vector<double> pars{dx1, dx2, dy, dz};
0361   std::vector<double> vtra = getTranslation(fv);
0362   std::vector<double> vrot = getRotation(fv);
0363 
0364   edm::LogVerbatim("GEMGeometryParsFromDD")
0365       << "(3) DD4hep, SuperChamber DetID " << gemid.rawId() << " Name " << std::string(name) << " dx1 " << dx1
0366       << " dx2 " << dx2 << " dy " << dy << " dz " << dz;
0367   rgeo.insert(gemid.rawId(), vtra, vrot, pars, {std::string(name)});
0368 }
0369 
0370 void GEMGeometryParsFromDD::buildChamber(cms::DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0371   cms::DDSolid solid(fv.solid());
0372   auto solidA = solid.solidA();
0373   std::vector<double> dpar = solidA.dimensions();
0374 
0375   double dy = dpar[3] / dd4hep::mm;   //length is along local Y
0376   double dz = dpar[2] / dd4hep::mm;   // thickness is long local Z
0377   double dx1 = dpar[0] / dd4hep::mm;  // bottom width is along local X
0378   double dx2 = dpar[1] / dd4hep::mm;  // top width is along local X
0379 
0380   auto solidB = solid.solidB();
0381   dpar = solidB.dimensions();
0382 
0383   dz += (dpar[2] / dd4hep::mm);  // chamber thickness
0384 
0385   GEMDetId gemid = detId.chamberId();
0386   std::string_view name = fv.name();
0387 
0388   std::vector<double> pars{dx1, dx2, dy, dz};
0389   std::vector<double> vtra = getTranslation(fv);
0390   std::vector<double> vrot = getRotation(fv);
0391 
0392   edm::LogVerbatim("GEMGeometryParsFromDD")
0393       << "(4) DD4hep, Chamber DetID " << gemid.rawId() << " Name " << std::string(name) << " dx1 " << dx1 << " dx2 "
0394       << dx2 << " dy " << dy << " dz " << dz;
0395   rgeo.insert(gemid.rawId(), vtra, vrot, pars, {std::string(name)});
0396 }
0397 
0398 void GEMGeometryParsFromDD::buildEtaPartition(cms::DDFilteredView& fv, GEMDetId detId, RecoIdealGeometry& rgeo) {
0399   auto nStrips = fv.get<double>("nStrips");
0400   auto nPads = fv.get<double>("nPads");
0401   auto dPhi = fv.get<double>("dPhi");
0402 
0403   std::vector<double> dpar = fv.parameters();
0404   std::string_view name = fv.name();
0405 
0406   double dx1 = dpar[0] / dd4hep::mm;
0407   double dx2 = dpar[1] / dd4hep::mm;
0408   double dy = dpar[3] / dd4hep::mm;
0409   double dz = dpar[2] / dd4hep::mm;
0410 
0411   std::vector<double> pars{dx1, dx2, dy, dz, nStrips, nPads, dPhi};
0412   std::vector<double> vtra = getTranslation(fv);
0413   std::vector<double> vrot = getRotation(fv);
0414 
0415   edm::LogVerbatim("GEMGeometryParsFromDD")
0416       << "(5) DD4hep, Eta Partion DetID " << detId.rawId() << " Name " << std::string(name) << " dx1 " << dx1 << " dx2 "
0417       << dx2 << " dy " << dy << " dz " << dz << " nStrips " << nStrips << " nPads " << nPads << " dPhi " << dPhi;
0418   rgeo.insert(detId.rawId(), vtra, vrot, pars, {std::string(name)});
0419 }
0420 
0421 std::vector<double> GEMGeometryParsFromDD::getTranslation(cms::DDFilteredView& fv) {
0422   std::vector<double> tran(3);
0423   tran[0] = static_cast<double>(fv.translation().X()) / dd4hep::mm;
0424   tran[1] = static_cast<double>(fv.translation().Y()) / dd4hep::mm;
0425   tran[2] = static_cast<double>(fv.translation().Z()) / dd4hep::mm;
0426 
0427   edm::LogVerbatim("GEMGeometryParsFromDD")
0428       << "(1) DD4hep, tran vector " << tran[0] << "  " << tran[1] << "  " << tran[2];
0429   return {tran[0], tran[1], tran[2]};
0430 }
0431 
0432 std::vector<double> GEMGeometryParsFromDD::getRotation(cms::DDFilteredView& fv) {
0433   DDRotationMatrix rota;
0434   fv.rot(rota);
0435   DD3Vector x, y, z;
0436   rota.GetComponents(x, y, z);
0437   const std::vector<double> rot = {x.X(), x.Y(), x.Z(), y.X(), y.Y(), y.Z(), z.X(), z.Y(), z.Z()};
0438   edm::LogVerbatim("GEMGeometryParsFromDD")
0439       << "(2) DD4hep, rot matrix " << rot[0] << "  " << rot[1] << "  " << rot[2] << " " << rot[3] << "  " << rot[4]
0440       << "  " << rot[5] << " " << rot[6] << "  " << rot[7] << "  " << rot[8];
0441   return {rot[0], rot[1], rot[2], rot[3], rot[4], rot[5], rot[6], rot[7], rot[8]};
0442 }