Allow RPD reconstruction geometry to be instantiation using a Modifier flag

Author: Sunanda

Configuration/Eras/python/Modifier_zdcAddRPD_cff.py

@@ -0,0 +1,3 @@
+import FWCore.ParameterSet.Config as cms
+
+zdcAddRPD =  cms.Modifier()

DataFormats/HcalDetId/interface/HcalElectronicsId.h

@@ -71,44 +71,44 @@ class HcalElectronicsId {
   }
 
   /// get subtype for this channel (valid for uTCA only)
-  constexpr int subtype() const { return (isUTCAid()) ? ((hcalElectronicsId_ >> 21) & 0x1F) : (-1); }
+  constexpr int32_t subtype() const { return (isUTCAid()) ? static_cast<int32_t>((hcalElectronicsId_ >> 21) & 0x1F) : (-1); }
   /// get the fiber channel id (which of channels on a fiber)
-  constexpr int fiberChanId() const { return (isVMEid()) ? (hcalElectronicsId_ & 0x3) : (hcalElectronicsId_ & 0xF); }
+  constexpr int32_t fiberChanId() const { return (isVMEid()) ? static_cast<int32_t>(hcalElectronicsId_ & 0x3) : (hcalElectronicsId_ & 0xF); }
   /// get the fiber index.  For VME 1-8 (which of eight fibers carried by a spigot), for uTCA fibers are zero-based
-  constexpr int fiberIndex() const {
+  constexpr int32_t fiberIndex() const {
     return (isVMEid()) ? (((hcalElectronicsId_ >> 2) & 0x7) + 1) : ((hcalElectronicsId_ >> 4) & 0x1F);
   }
   /// get the SLB channel index  (valid only for VME trigger-chain ids)
-  constexpr int slbChannelIndex() const { return hcalElectronicsId_ & 0x3; }
+  constexpr int32_t slbChannelIndex() const { return hcalElectronicsId_ & 0x3; }
   /// get the SLB site number (valid only for VME trigger-chain ids)
-  constexpr int slbSiteNumber() const { return ((hcalElectronicsId_ >> 2) & 0x7); }
+  constexpr int32_t slbSiteNumber() const { return ((hcalElectronicsId_ >> 2) & 0x7); }
 
   /// get the HTR channel id (1-24)
-  constexpr int htrChanId() const { return isVMEid() ? ((fiberChanId() + 1) + ((fiberIndex() - 1) * 3)) : (0); }
+  constexpr int32_t htrChanId() const { return isVMEid() ? ((fiberChanId() + 1) + ((fiberIndex() - 1) * 3)) : (0); }
 
   /// get the HTR-wide slb channel code (letter plus number)
   std::string slbChannelCode() const;
 
   /// get the spigot (input number on DCC, AMC card number for uTCA)
-  constexpr int spigot() const { return (isVMEid()) ? ((hcalElectronicsId_ >> 5) & 0xF) : slot(); }
+  constexpr int32_t spigot() const { return (isVMEid()) ? ((hcalElectronicsId_ >> 5) & 0xF) : slot(); }
   /// get the (Hcal local) DCC id for VME, crate number for uTCA
-  constexpr int dccid() const { return (isVMEid()) ? ((hcalElectronicsId_ >> 9) & 0x1F) : crateId(); }
+  constexpr int32_t dccid() const { return (isVMEid()) ? ((hcalElectronicsId_ >> 9) & 0x1F) : crateId(); }
   /// get the htr slot
-  constexpr int htrSlot() const { return slot(); }
+  constexpr int32_t htrSlot() const { return slot(); }
   /// get the htr or uHTR slot
-  constexpr int slot() const {
+  constexpr int32_t slot() const {
     return (isVMEid()) ? ((hcalElectronicsId_ >> 14) & 0x1F) : ((hcalElectronicsId_ >> 9) & 0xF);
   }
   /// get the htr top/bottom (1=top/0=bottom), valid for VME
-  constexpr int htrTopBottom() const { return (isVMEid()) ? ((hcalElectronicsId_ >> 19) & 0x1) : (-1); }
+  constexpr int32_t htrTopBottom() const { return (isVMEid()) ? static_cast<int32_t>((hcalElectronicsId_ >> 19) & 0x1) : (-1); }
   /// get the readout VME crate number
-  constexpr int readoutVMECrateId() const { return crateId(); }
+  constexpr int32_t readoutVMECrateId() const { return crateId(); }
   /// get the readout VME crate number
-  constexpr int crateId() const {
+  constexpr int32_t crateId() const {
     return (isVMEid()) ? ((hcalElectronicsId_ >> 20) & 0x1F) : ((hcalElectronicsId_ >> 13) & 0x3F);
   }
   /// get a fast, compact, unique index for linear lookups
-  constexpr int linearIndex() const {
+  constexpr int32_t linearIndex() const {
     return (isVMEid()) ? ((hcalElectronicsId_)&0x3FFF) : ((hcalElectronicsId_)&0x7FFFF);
   }
 

DataFormats/HcalDetId/interface/HcalZDCDetId.h

@@ -34,6 +34,14 @@ class HcalZDCDetId : public DetId {
 
   static constexpr int32_t SubdetectorId = 2;
 
+  static constexpr int32_t kDepEM = 5;
+  static constexpr int32_t kDepHAD = 4;
+  static constexpr int32_t kDepLUM = 2;
+  static constexpr int32_t kDepRPD = 16;
+  static constexpr int32_t kDepRun1 = (kDepEM + kDepHAD + kDepLUM);
+  static constexpr int32_t kDepTot = (kDepRun1 + kDepRPD);
+  static constexpr int32_t kDepRun3 = kDepTot;
+
   /** Create a null cellid*/
   constexpr HcalZDCDetId() : DetId() {}
   /** Create cellid from raw id (0=invalid tower id) */
@@ -162,16 +170,6 @@ class HcalZDCDetId : public DetId {
   }
 
 private:
-  enum {
-    kDepEM = 5,
-    kDepHAD = 4,
-    kDepLUM = 2,
-    kDepRPD = 16,
-    kDepRun1 = kDepEM + kDepHAD + kDepLUM,
-    kDepTot = kDepRun1 + kDepRPD,
-    kDepRun3 = kDepTot
-  };
-
   constexpr static uint32_t packHcalZDCDetId(const Section& se, const bool& zside, const int32_t& channel) {
     uint32_t id = DetId(DetId::Calo, SubdetectorId);
     id |= kZDCnewFormat;

Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h

@@ -1,6 +1,7 @@
 #ifndef GEOMETRY_CALOGEOMETRY_CALOSUBDETECTORGEOMETRY_H
 #define GEOMETRY_CALOGEOMETRY_CALOSUBDETECTORGEOMETRY_H 1
 
+#include <algorithm>
 #include <memory>
 #include <vector>
 #include <set>
@@ -99,6 +100,8 @@ class CaloSubdetectorGeometry {
 
   virtual void initializeParms() { return; }
 
+  virtual bool valid(const DetId& id) const { return (std::find(m_validIds.begin(), m_validIds.end(), id) != m_validIds.end()); }
+
 protected:
   virtual unsigned int indexFor(const DetId& id) const;
   virtual unsigned int sizeForDenseIndex(const DetId& id) const;

Geometry/ForwardGeometry/interface/ZdcGeometry.h

@@ -57,6 +57,11 @@ class ZdcGeometry : public CaloSubdetectorGeometry {
                const CCGFloat* parm,
                const DetId& detId) override;
 
+  void getSummary(CaloSubdetectorGeometry::TrVec& tVec,
+		  CaloSubdetectorGeometry::IVec& iVec,
+		  CaloSubdetectorGeometry::DimVec& dVec,
+		  CaloSubdetectorGeometry::IVec& dins) const override;
+
 protected:
   unsigned int indexFor(const DetId& id) const override { return HcalZDCDetId(id).denseIndex(); }
 

Geometry/ForwardGeometry/interface/ZdcHardcodeGeometryLoader.h

@@ -23,6 +23,7 @@ class ZdcHardcodeGeometryLoader {
 
   virtual ReturnType load(DetId::Detector det, int subdet);
   ReturnType load();
+  void setAddRPD(bool flag) { m_zdcAddRPD = flag; }
 
 private:
   void init();
@@ -31,7 +32,7 @@ class ZdcHardcodeGeometryLoader {
 
   ZdcTopology* theTopology;
   const ZdcTopology* extTopology;
-
+  bool m_zdcAddRPD;
   float theEMSectiondX;
   float theEMSectiondY;
   float theEMSectiondZ;

Geometry/ForwardGeometry/plugins/ZdcHardcodeGeometryEP.cc

@@ -19,8 +19,7 @@
 #include "Geometry/ForwardGeometry/interface/ZdcGeometry.h"
 #include "FWCore/MessageLogger/interface/MessageLogger.h"
 
-ZdcHardcodeGeometryEP::ZdcHardcodeGeometryEP(const edm::ParameterSet& ps)
-    : m_loader(nullptr), m_topology(), m_applyAlignment(ps.getUntrackedParameter<bool>("applyAlignment", false)) {
+ZdcHardcodeGeometryEP::ZdcHardcodeGeometryEP(const edm::ParameterSet& ps) : m_loader(nullptr), m_topology(), m_applyAlignment(ps.getParameter<bool>("applyAlignment")), m_zdcAddRPD(ps.getParameter<bool>("zdcAddRPD")) {
   //the following line is needed to tell the framework what
   // data is being produced
   setWhatProduced(this, ZdcGeometry::producerTag());
@@ -31,7 +30,7 @@ ZdcHardcodeGeometryEP::ZdcHardcodeGeometryEP(const edm::ParameterSet& ps)
   //		    edm::es::Label( "ZDC" ) );
 }
 
-ZdcHardcodeGeometryEP::~ZdcHardcodeGeometryEP() { delete m_loader; }
+ZdcHardcodeGeometryEP::~ZdcHardcodeGeometryEP() { }
 
 //
 // member functions
@@ -41,7 +40,14 @@ ZdcHardcodeGeometryEP::~ZdcHardcodeGeometryEP() { delete m_loader; }
 
 ZdcHardcodeGeometryEP::ReturnType ZdcHardcodeGeometryEP::produce(const ZDCGeometryRecord& iRecord) {
   //   ZdcHardcodeGeometryLoader loader ( m_topology ) ;
-  m_loader = new ZdcHardcodeGeometryLoader(m_topology);
-
+  m_loader = std::make_unique<ZdcHardcodeGeometryLoader>(m_topology);
+  m_loader->setAddRPD(m_zdcAddRPD);
   return ReturnType(m_loader->load());
 }
+
+void ZdcHardcodeGeometryEP::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
+  edm::ParameterSetDescription desc;
+  desc.add<bool>("applyAlignment", false);
+  desc.add<bool>("zdcAddRPD", false);
+  descriptions.addWithDefaultLabel(desc);
+}

Geometry/ForwardGeometry/plugins/ZdcHardcodeGeometryEP.h

@@ -21,17 +21,20 @@ class ZdcHardcodeGeometryEP : public edm::ESProducer {
   ~ZdcHardcodeGeometryEP() override;
 
   using ReturnType = std::unique_ptr<CaloSubdetectorGeometry>;
-
+  
   ReturnType produce(const ZDCGeometryRecord&);
 
+  static void fillDescriptions(edm::ConfigurationDescriptions&);
+
 private:
   // ----------member data ---------------------------
 
-  ZdcHardcodeGeometryLoader* m_loader;
+  std::unique_ptr<ZdcHardcodeGeometryLoader> m_loader;
 
   ZdcTopology m_topology;
 
   bool m_applyAlignment;
+  bool m_zdcAddRPD;
 };
 
 #endif

Geometry/ForwardGeometry/python/ZdcGeometry_cfi.py

@@ -1,3 +1,4 @@
-import FWCore.ParameterSet.Config as cms
+from Geometry.ForwardGeometry.zdcHardcodeGeometryEP_cfi import zdcHardcodeGeometryEP
+from Configuration.Eras.Modifier_zdcAddRPD_cff import zdcAddRPD
 
-ZdcHardcodeGeometryEP = cms.ESProducer( "ZdcHardcodeGeometryEP" )
+zdcAddRPD.toModify(zdcHardcodeGeometryEP, zdcAddRPD = True)

Geometry/ForwardGeometry/src/ZdcGeometry.cc

@@ -6,9 +6,18 @@
 #include "Geometry/ForwardGeometry/interface/ZdcHardcodeGeometryData.h"
 #include <algorithm>
 
-typedef CaloCellGeometry::CCGFloat CCGFloat;
+#include <Math/Transform3D.h>
+#include <Math/EulerAngles.h>
+
+#include <algorithm>
+
 typedef CaloCellGeometry::Pt3D Pt3D;
 typedef CaloCellGeometry::Pt3DVec Pt3DVec;
+typedef CaloCellGeometry::Tr3D Tr3D;
+
+typedef CaloSubdetectorGeometry::CCGFloat CCGFloat;
+
+//#define EDM_ML_DEBUG
 
 ZdcGeometry::ZdcGeometry()
     : theTopology(new ZdcTopology),
@@ -66,6 +75,9 @@ void ZdcGeometry::newCell(const GlobalPoint& f1,
                           const CCGFloat* parm,
                           const DetId& detId) {
   const CaloGenericDetId cgid(detId);
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("ZDCGeom") << "ZDCGeometry " << HcalZDCDetId(detId) << " Generic ID " << std::hex << cgid.rawId() << std::dec << " ZDC? " << cgid.isZDC();
+#endif
   assert(cgid.isZDC());
 
   const unsigned int di(cgid.denseIndex());
@@ -81,3 +93,70 @@ const CaloCellGeometry* ZdcGeometry::getGeometryRawPtr(uint32_t index) const {
   const CaloCellGeometry* cell(&m_cellVec[index]);
   return (((cell == nullptr) || (nullptr == cell->param())) ? nullptr : cell);
 }
+
+void ZdcGeometry::getSummary(CaloSubdetectorGeometry::TrVec& tVec,
+			     CaloSubdetectorGeometry::IVec& iVec,
+			     CaloSubdetectorGeometry::DimVec& dVec,
+			     CaloSubdetectorGeometry::IVec& /*dins*/) const {
+  tVec.reserve(m_validIds.size() * numberOfTransformParms());
+  iVec.reserve(numberOfShapes() == 1 ? 1 : m_validIds.size());
+  dVec.reserve(numberOfShapes() * numberOfParametersPerShape());
+
+  for (const auto& pv : parVecVec()) {
+    for (float iv : pv) {
+      dVec.emplace_back(iv);
+    }
+  }
+
+  for (uint32_t i(0); i != m_validIds.size(); ++i) {
+    Tr3D tr;
+    std::shared_ptr<const CaloCellGeometry> ptr(cellGeomPtr(i));
+#ifdef EDM_ML_DEBUG
+    edm::LogVerbatim("ZDCGeom") << "ZDCGeometry:Summary " << i << ":" << HcalZDCDetId::kSizeForDenseIndexingRun1 << " Pointer " << ptr << ":" << (nullptr != ptr);
+#endif
+    if (i < static_cast<int32_t>(HcalZDCDetId::kSizeForDenseIndexingRun1))
+      assert(nullptr != ptr);
+    if (ptr != nullptr) {
+      ptr->getTransform(tr, (Pt3DVec*)nullptr);
+
+      if (Tr3D() == tr) {  // for preshower there is no rotation
+	const GlobalPoint& gp(ptr->getPosition());
+	tr = HepGeom::Translate3D(gp.x(), gp.y(), gp.z());
+      }
+
+      const CLHEP::Hep3Vector tt(tr.getTranslation());
+      tVec.emplace_back(tt.x());
+      tVec.emplace_back(tt.y());
+      tVec.emplace_back(tt.z());
+      if (6 == numberOfTransformParms()) {
+	const CLHEP::HepRotation rr(tr.getRotation());
+	const ROOT::Math::Transform3D rtr(rr.xx(), rr.xy(), rr.xz(), tt.x(), rr.yx(), rr.yy(), rr.yz(), tt.y(), rr.zx(), rr.zy(), rr.zz(), tt.z());
+	ROOT::Math::EulerAngles ea;
+	rtr.GetRotation(ea);
+	tVec.emplace_back(ea.Phi());
+	tVec.emplace_back(ea.Theta());
+	tVec.emplace_back(ea.Psi());
+      }
+
+      const CCGFloat* par(ptr->param());
+
+      unsigned int ishape(9999);
+      for (unsigned int ivv(0); ivv != parVecVec().size(); ++ivv) {
+	bool ok(true);
+	const CCGFloat* pv(&(*parVecVec()[ivv].begin()));
+	for (unsigned int k(0); k != numberOfParametersPerShape(); ++k) {
+	  ok = ok && (fabs(par[k] - pv[k]) < 1.e-6);
+	}
+	if (ok) {
+	  ishape = ivv;
+	  break;
+	}
+      }
+      assert(9999 != ishape);
+
+      const unsigned int nn((numberOfShapes() == 1) ? (unsigned int)1 : m_validIds.size());
+      if (iVec.size() < nn)
+	iVec.emplace_back(ishape);
+    }
+  }
+}