Phase2-hgx335J Try to add toplogy class for HGCal TB application

Geometry/CaloEventSetup/plugins/HGCalTBTopologyBuilder.cc

@@ -0,0 +1,82 @@
+// -*- C++ -*-
+//
+// Package:    CaloEventSetup
+// Class:      HGCalTBTopologyBuilder
+//
+/**\class HGCalTBTopologyBuilder HGCalTBTopologyBuilder.h 
+
+ Description: <one line class summary>
+
+ Implementation:
+     <Notes on implementation>
+*/
+//
+// Original Author:  Sunanda Banerjee
+//
+//
+
+// system include files
+#include <memory>
+
+// user include files
+#include <FWCore/Framework/interface/ModuleFactory.h>
+#include "FWCore/Framework/interface/ESProducer.h"
+
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+#include "FWCore/Utilities/interface/ESGetToken.h"
+
+#include "Geometry/CaloTopology/interface/HGCalTBTopology.h"
+#include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
+#include "Geometry/HGCalTBCommonData/interface/HGCalTBDDDConstants.h"
+#include "Geometry/Records/interface/IdealGeometryRecord.h"
+#include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
+#include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
+
+//#define EDM_ML_DEBUG
+
+//
+// class decleration
+//
+
+class HGCalTBTopologyBuilder : public edm::ESProducer {
+public:
+  HGCalTBTopologyBuilder(const edm::ParameterSet& iP);
+  ~HGCalTBTopologyBuilder() override;
+
+  using ReturnType = std::unique_ptr<HGCalTBTopology>;
+
+  ReturnType produce(const IdealGeometryRecord&);
+
+private:
+  // ----------member data ---------------------------
+  edm::ESGetToken<HGCalTBDDDConstants, IdealGeometryRecord> hgcToken_;
+  int det_;
+};
+
+HGCalTBTopologyBuilder::HGCalTBTopologyBuilder(const edm::ParameterSet& iConfig) {
+  auto name = iConfig.getParameter<std::string>("Name");
+  det_ = iConfig.getParameter<int>("Type");
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("HGCalGeom") << "constructing HGCalTBTopology for " << name << " and det " << det_;
+#endif
+  hgcToken_ = setWhatProduced(this, name).consumes<HGCalTBDDDConstants>(edm::ESInputTag{"", name});
+}
+
+HGCalTBTopologyBuilder::~HGCalTBTopologyBuilder() {}
+
+//
+// member functions
+//
+
+// ------------ method called to produce the data  ------------
+HGCalTBTopologyBuilder::ReturnType HGCalTBTopologyBuilder::produce(const IdealGeometryRecord& iRecord) {
+  const HGCalTBDDDConstants& hgdc = iRecord.get(hgcToken_);
+
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("HGCalGeom") << "Create HGCalTBTopology(hgdc,det)";
+#endif
+  return std::make_unique<HGCalTBTopology>(&hgdc, det_);
+}
+
+DEFINE_FWK_EVENTSETUP_MODULE(HGCalTBTopologyBuilder);

Geometry/CaloEventSetup/python/HGCalTBTopology_cfi.py

@@ -0,0 +1,14 @@
+import FWCore.ParameterSet.Config as cms
+
+#
+# This cfi should be included to build the HGCal TB Topologies
+#
+HGCalTBEETopologyBuilder = cms.ESProducer("HGCalTBTopologyBuilder",
+                                          Name = cms.string("HGCalEESensitive"),
+                                          Type = cms.int32(3) )
+
+
+HGCalTBHESilTopologyBuilder = cms.ESProducer("HGCalTBTopologyBuilder",
+                                             Name = cms.string("HGCalHESiliconSensitive"),
+                                             Type = cms.int32(4) )
+

Geometry/CaloTopology/BuildFile.xml

@@ -6,6 +6,7 @@
 <use name="Geometry/CaloGeometry"/>
 <use name="Geometry/HcalCommonData"/>
 <use name="Geometry/HGCalCommonData"/>
+<use name="Geometry/HGCalTBCommonData"/>
 <use name="boost"/>
 <export>
   <lib name="1"/>

Geometry/CaloTopology/interface/HGCalTBTopology.h

@@ -0,0 +1,133 @@
+#ifndef Geometry_CaloTopology_HGCalTBTopology_h
+#define Geometry_CaloTopology_HGCalTBTopology_h 1
+
+#include "DataFormats/ForwardDetId/interface/ForwardSubdetector.h"
+#include "DataFormats/DetId/interface/DetId.h"
+#include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
+#include "Geometry/HGCalCommonData/interface/HGCalGeometryMode.h"
+#include "Geometry/HGCalTBCommonData/interface/HGCalTBDDDConstants.h"
+#include <vector>
+#include <iostream>
+
+class HGCalTBTopology : public CaloSubdetectorTopology {
+public:
+  /// create a new Topology
+  HGCalTBTopology(const HGCalTBDDDConstants* hdcons, int subdet);
+
+  /// default destructor
+
+  /// move the Topology north (increment iy)
+  DetId goNorth(const DetId& id) const override { return changeXY(id, 0, +1); }
+  std::vector<DetId> north(const DetId& id) const override {
+    DetId nextId = goNorth(id);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  /// move the Topology south (decrement iy)
+  DetId goSouth(const DetId& id) const override { return changeXY(id, 0, -1); }
+  std::vector<DetId> south(const DetId& id) const override {
+    DetId nextId = goSouth(id);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  /// move the Topology east (positive ix)
+  DetId goEast(const DetId& id) const override { return changeXY(id, +1, 0); }
+  std::vector<DetId> east(const DetId& id) const override {
+    DetId nextId = goEast(id);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  /// move the Topology west (negative ix)
+  DetId goWest(const DetId& id) const override { return changeXY(id, -1, 0); }
+  std::vector<DetId> west(const DetId& id) const override {
+    DetId nextId = goWest(id);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  std::vector<DetId> up(const DetId& id) const override {
+    DetId nextId = changeZ(id, +1);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  std::vector<DetId> down(const DetId& id) const override {
+    DetId nextId = changeZ(id, -1);
+    std::vector<DetId> vNeighborsDetId;
+    if (!(nextId == DetId(0)))
+      vNeighborsDetId.emplace_back(nextId);
+    return vNeighborsDetId;
+  }
+
+  ///Geometry mode
+  HGCalGeometryMode::GeometryMode geomMode() const { return mode_; }
+
+  ///Dense indexing
+  uint32_t detId2denseId(const DetId& id) const override;
+  DetId denseId2detId(uint32_t denseId) const override;
+  virtual uint32_t detId2denseGeomId(const DetId& id) const;
+
+  ///Is this a valid cell id
+  bool valid(const DetId& id) const override;
+  bool validHashIndex(uint32_t ix) const { return (ix < kSizeForDenseIndexing); }
+
+  unsigned int totalModules() const { return kSizeForDenseIndexing; }
+  unsigned int totalGeomModules() const { return (unsigned int)(2 * kHGeomHalf_); }
+  unsigned int allGeomModules() const;
+
+  const HGCalTBDDDConstants& dddConstants() const { return *hdcons_; }
+
+  /** returns a new DetId offset by nrStepsX and nrStepsY (can be negative),
+   * returns DetId(0) if invalid */
+  DetId offsetBy(const DetId startId, int nrStepsX, int nrStepsY) const;
+  DetId switchZSide(const DetId startId) const;
+
+  /// Use subSector in square mode as wafer type in hexagon mode
+  static const int subSectors_ = 2;
+
+  struct DecodedDetId {
+    DecodedDetId() : iCell1(0), iCell2(0), iLay(0), iSec1(0), iSec2(0), iType(0), zSide(0), det(0) {}
+    int iCell1, iCell2, iLay, iSec1, iSec2, iType, zSide, det;
+  };
+
+  DecodedDetId geomDenseId2decId(const uint32_t& hi) const;
+  DecodedDetId decode(const DetId& id) const;
+  DetId encode(const DecodedDetId& id_) const;
+
+  DetId::Detector detector() const { return det_; }
+  ForwardSubdetector subDetector() const { return subdet_; }
+  bool detectorType() const { return false; }
+
+private:
+  HGCalTBTopology() : hdcons_(nullptr) {}
+
+  /// move the nagivator along x, y
+  DetId changeXY(const DetId& id, int nrStepsX, int nrStepsY) const;
+
+  /// move the nagivator along z
+  DetId changeZ(const DetId& id, int nrStepsZ) const;
+
+  const HGCalTBDDDConstants* hdcons_;
+  HGCalGeometryMode::GeometryMode mode_;
+
+  DetId::Detector det_;
+  ForwardSubdetector subdet_;
+  int sectors_, layers_, cells_, types_, firstLay_;
+  int kHGhalf_, kHGeomHalf_, kHGhalfType_;
+  unsigned int kSizeForDenseIndexing;
+};
+
+#endif

Geometry/CaloTopology/src/HGCalTBTopology.cc

@@ -0,0 +1,137 @@
+#include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
+#include "Geometry/CaloTopology/interface/HGCalTBTopology.h"
+#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
+
+#define EDM_ML_DEBUG
+
+HGCalTBTopology::HGCalTBTopology(const HGCalTBDDDConstants* hdcons, int det) : hdcons_(hdcons) {
+  sectors_ = hdcons_->sectors();
+  layers_ = hdcons_->layers(true);
+  cells_ = hdcons_->maxCells(true);
+  mode_ = hdcons_->geomMode();
+  kHGhalf_ = sectors_ * layers_ * cells_;
+  firstLay_ = hdcons_->firstLayer();
+  det_ = DetId::Forward;
+  subdet_ = (ForwardSubdetector)(det);
+  kHGeomHalf_ = sectors_ * layers_;
+  types_ = 2;
+  kHGhalfType_ = sectors_ * layers_ * cells_ * types_;
+  kSizeForDenseIndexing = static_cast<unsigned int>(2 * kHGhalf_);
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("HGCalGeom") << "HGCalTBTopology initialized for detector " << det << ":" << det_ << ":" << subdet_
+                                << " having " << sectors_ << " Sectors, " << layers_ << " Layers from " << firstLay_
+                                << ", " << cells_ << " cells and total channels " << kSizeForDenseIndexing << ":"
+                                << (2 * kHGeomHalf_);
+#endif
+}
+
+unsigned int HGCalTBTopology::allGeomModules() const { return (static_cast<unsigned int>(2 * hdcons_->wafers())); }
+
+DetId HGCalTBTopology::denseId2detId(uint32_t hi) const {
+  HGCalTBTopology::DecodedDetId id;
+  if (validHashIndex(hi)) {
+    id.zSide = ((int)(hi) < kHGhalfType_ ? -1 : 1);
+    int di = ((int)(hi) % kHGhalfType_);
+    int type = (di % types_);
+    id.iType = (type == 0 ? -1 : 1);
+    id.iSec1 = (((di - type) / types_) % sectors_);
+    id.iLay = (((((di - type) / types_) - id.iSec1 + 1) / sectors_) % layers_ + 1);
+    id.iCell1 = (((((di - type) / types_) - id.iSec1 + 1) / sectors_ - id.iLay + 1) / layers_ + 1);
+#ifdef EDM_ML_DEBUG
+    edm::LogVerbatim("HGCalGeom") << "Input Hex " << hi << " o/p " << id.zSide << ":" << id.iLay << ":" << id.iType
+                                  << ":" << id.iSec1 << ":" << id.iCell1;
+#endif
+  }
+  return encode(id);
+}
+
+uint32_t HGCalTBTopology::detId2denseGeomId(const DetId& idin) const {
+  HGCalTBTopology::DecodedDetId id = decode(idin);
+  uint32_t idx;
+  idx = (uint32_t)(((id.zSide > 0) ? kHGeomHalf_ : 0) + (id.iLay - 1) * sectors_ + id.iSec1);
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("HGCalGeom") << "Geom Hex I/P " << id.zSide << ":" << id.iLay << ":" << id.iSec1 << ":" << id.iType
+                                << " Constants " << kHGeomHalf_ << ":" << layers_ << ":" << sectors_ << " o/p " << idx;
+#endif
+  return idx;
+}
+
+bool HGCalTBTopology::valid(const DetId& idin) const {
+  HGCalTBTopology::DecodedDetId id = decode(idin);
+  bool flag;
+  flag = (idin.det() == det_ && idin.subdetId() == (int)(subdet_) && id.iCell1 >= 0 && id.iCell1 < cells_ &&
+          id.iLay > 0 && id.iLay <= layers_ && id.iSec1 >= 0 && id.iSec1 <= sectors_);
+  if (flag)
+    flag = hdcons_->isValidHex(id.iLay, id.iSec1, id.iCell1, true);
+  return flag;
+}
+
+DetId HGCalTBTopology::offsetBy(const DetId startId, int nrStepsX, int nrStepsY) const {
+  if (startId.det() == DetId::Forward && startId.subdetId() == static_cast<int>(subdet_)) {
+    DetId id = changeXY(startId, nrStepsX, nrStepsY);
+    if (valid(id))
+      return id;
+  }
+  return DetId(0);
+}
+
+DetId HGCalTBTopology::switchZSide(const DetId startId) const {
+  HGCalTBTopology::DecodedDetId id_ = decode(startId);
+  id_.zSide = -id_.zSide;
+  DetId id = encode(id_);
+  if (valid(id))
+    return id;
+  else
+    return DetId(0);
+}
+
+HGCalTBTopology::DecodedDetId HGCalTBTopology::geomDenseId2decId(const uint32_t& hi) const {
+  HGCalTBTopology::DecodedDetId id;
+  if (hi < totalGeomModules()) {
+    id.zSide = ((int)(hi) < kHGeomHalf_ ? -1 : 1);
+    int di = ((int)(hi) % kHGeomHalf_);
+    id.iSec1 = (di % sectors_);
+    di = (di - id.iSec1) / sectors_;
+    id.iLay = (di % layers_) + 1;
+    id.iType = ((di - id.iLay + 1) / layers_ == 0) ? -1 : 1;
+#ifdef EDM_ML_DEBUG
+    edm::LogVerbatim("HGCalGeom") << "Geom Hex I/P " << hi << " O/P " << id.zSide << ":" << id.iType << ":" << id.iLay
+                                  << ":" << id.iSec1;
+#endif
+  }
+  return id;
+}
+
+HGCalTBTopology::DecodedDetId HGCalTBTopology::decode(const DetId& startId) const {
+  HGCalTBTopology::DecodedDetId idx;
+  HGCalDetId id(startId);
+  idx.iCell1 = id.cell();
+  idx.iCell2 = 0;
+  idx.iLay = id.layer();
+  idx.iSec1 = id.wafer();
+  idx.iSec2 = 0;
+  idx.iType = id.waferType();
+  idx.zSide = id.zside();
+  idx.det = id.subdetId();
+  return idx;
+}
+
+DetId HGCalTBTopology::encode(const HGCalTBTopology::DecodedDetId& idx) const {
+  DetId id;
+#ifdef EDM_ML_DEBUG
+  edm::LogVerbatim("HGCalGeomX") << "Encode " << idx.det << ":" << idx.zSide << ":" << idx.iType << ":" << idx.iLay
+                                 << ":" << idx.iSec1 << ":" << idx.iSec2 << ":" << idx.iCell1 << ":" << idx.iCell2;
+#endif
+  id = HGCalDetId((ForwardSubdetector)(idx.det), idx.zSide, idx.iLay, ((idx.iType > 0) ? 1 : 0), idx.iSec1, idx.iCell1)
+           .rawId();
+  return id;
+}
+
+DetId HGCalTBTopology::changeXY(const DetId& id, int nrStepsX, int nrStepsY) const { return DetId(); }
+
+DetId HGCalTBTopology::changeZ(const DetId& id, int nrStepsZ) const { return DetId(); }
+
+#include "FWCore/Utilities/interface/typelookup.h"
+
+TYPELOOKUP_DATA_REG(HGCalTBTopology);