Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Run4-hgx306 Backport the bugfix to HGCal geometry from the PR's 37330 and 37460 #37466

Merged
merged 2 commits into from Apr 8, 2022
Merged

Run4-hgx306 Backport the bugfix to HGCal geometry from the PR's 37330 and 37460 #37466

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
250803c
Backport the bugfix to HGCal geometry from the PR's 37330 and 37460
Apr 5, 2022
7aa146a
Add test code as suggested
Apr 7, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
11 changes: 9 additions & 2 deletions Geometry/HGCalCommonData/interface/HGCalDDDConstants.h
View file
Open in desktop
Expand Up @@ -84,8 +84,15 @@ class HGCalDDDConstants {
std::pair<float, float> localToGlobal8(
int lay, int waferU, int waferV, double localX, double localY, bool reco, bool debug) const;
std::pair<float, float> locateCell(int cell, int lay, int type, bool reco) const;
std::pair<float, float> locateCell(
int lay, int waferU, int waferV, int cellU, int cellV, bool reco, bool all, bool debug = false) const;
std::pair<float, float> locateCell(int lay,
int waferU,
int waferV,
int cellU,
int cellV,
bool reco,
bool all,
bool norot = false,
bool debug = false) const;
std::pair<float, float> locateCell(const HGCSiliconDetId&, bool debug = false) const;
std::pair<float, float> locateCell(const HGCScintillatorDetId&, bool debug = false) const;
std::pair<float, float> locateCellHex(int cell, int wafer, bool reco) const;
Expand Down
23 changes: 8 additions & 15 deletions Geometry/HGCalCommonData/src/HGCalDDDConstants.cc
View file
Open in desktop
Expand Up @@ -213,7 +213,7 @@ bool HGCalDDDConstants::cellInLayer(int waferU, int waferV, int cellU, int cellV
const auto& indx = getIndex(lay, true);
if (indx.first >= 0) {
if (waferHexagon8() || waferHexagon6()) {
const auto& xy = ((waferHexagon8()) ? locateCell(lay, waferU, waferV, cellU, cellV, reco, true, true)
const auto& xy = ((waferHexagon8()) ? locateCell(lay, waferU, waferV, cellU, cellV, reco, true, false, false)
: locateCell(cellU, lay, waferU, reco));
double rpos = sqrt(xy.first * xy.first + xy.second * xy.second);
return ((rpos >= hgpar_->rMinLayHex_[indx.first]) && (rpos <= hgpar_->rMaxLayHex_[indx.first]));
Expand Down Expand Up @@ -621,7 +621,7 @@ std::pair<float, float> HGCalDDDConstants::localToGlobal8(
bool rotx =
((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[lay - hgpar_->firstLayer_] == HGCalTypes::WaferCenterR));
if (debug)
edm::LogVerbatim("HGCalGeom") << "LocalToGlobal " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << rotx
edm::LogVerbatim("HGCalGeom") << "LocalToGlobal8 " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << rotx
<< " Local (" << x << ":" << y << ") Reco " << reco;
if (!reco) {
x *= HGCalParameters::k_ScaleToDDD;
Expand All @@ -631,7 +631,7 @@ std::pair<float, float> HGCalDDDConstants::localToGlobal8(
x += xy.first;
y += xy.second;
if (debug)
edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by addong " << xy.first << ":" << xy.second;
edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
}

Expand Down Expand Up @@ -667,44 +667,39 @@ std::pair<float, float> HGCalDDDConstants::locateCell(int cell, int lay, int typ
}

std::pair<float, float> HGCalDDDConstants::locateCell(
int lay, int waferU, int waferV, int cellU, int cellV, bool reco, bool all, bool debug) const {
int lay, int waferU, int waferV, int cellU, int cellV, bool reco, bool all, bool norot, bool debug) const {
double x(0), y(0);
int indx = HGCalWaferIndex::waferIndex(lay, waferU, waferV);
auto itr = hgpar_->typesInLayers_.find(indx);
int type = ((itr == hgpar_->typesInLayers_.end()) ? 2 : hgpar_->waferTypeL_[itr->second]);
bool rotx =
((!hgpar_->layerType_.empty()) && (hgpar_->layerType_[lay - hgpar_->firstLayer_] == HGCalTypes::WaferCenterR));
#ifdef EDM_ML_DEBUG
bool rotx = (norot) ? false
: ((!hgpar_->layerType_.empty()) &&
(hgpar_->layerType_[lay - hgpar_->firstLayer_] == HGCalTypes::WaferCenterR));
if (debug) {
edm::LogVerbatim("HGCalGeom") << "LocateCell " << lay << ":" << (lay - hgpar_->firstLayer_) << ":" << rotx << ":"
<< waferU << ":" << waferV << ":" << indx << ":"
<< (itr == hgpar_->typesInLayers_.end()) << ":" << type << " Flags " << reco << ":"
<< all;
}
#endif
int kndx = cellV * 100 + cellU;
if (type == 0) {
auto ktr = hgpar_->cellFineIndex_.find(kndx);
if (ktr != hgpar_->cellFineIndex_.end()) {
x = hgpar_->cellFineX_[ktr->second];
y = hgpar_->cellFineY_[ktr->second];
}
#ifdef EDM_ML_DEBUG
if (debug)
edm::LogVerbatim("HGCalGeom") << "Fine " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
<< (ktr != hgpar_->cellFineIndex_.end());
#endif
} else {
auto ktr = hgpar_->cellCoarseIndex_.find(kndx);
if (ktr != hgpar_->cellCoarseIndex_.end()) {
x = hgpar_->cellCoarseX_[ktr->second];
y = hgpar_->cellCoarseY_[ktr->second];
}
#ifdef EDM_ML_DEBUG
if (debug)
edm::LogVerbatim("HGCalGeom") << "Coarse " << cellU << ":" << cellV << ":" << kndx << ":" << x << ":" << y << ":"
<< (ktr != hgpar_->cellCoarseIndex_.end());
#endif
}
if (!reco) {
x *= HGCalParameters::k_ScaleToDDD;
Expand All @@ -714,10 +709,8 @@ std::pair<float, float> HGCalDDDConstants::locateCell(
const auto& xy = waferPositionNoRot(lay, waferU, waferV, reco, debug);
x += xy.first;
y += xy.second;
#ifdef EDM_ML_DEBUG
if (debug)
edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by addong " << xy.first << ":" << xy.second;
#endif
edm::LogVerbatim("HGCalGeom") << "With wafer " << x << ":" << y << " by adding " << xy.first << ":" << xy.second;
}
return (rotx ? getXY(lay, x, y, false) : std::make_pair(x, y));
}
Expand Down
6 changes: 3 additions & 3 deletions Geometry/HGCalGeometry/interface/HGCalGeometry.h
View file
Open in desktop
Expand Up @@ -50,7 +50,7 @@ class HGCalGeometry final : public CaloSubdetectorGeometry {

HGCalGeometry(const HGCalTopology& topology);

~HGCalGeometry() override;
~HGCalGeometry() override = default;

void localCorners(Pt3DVec& lc, const CCGFloat* pv, unsigned int i, Pt3D& ref);

Expand All @@ -70,7 +70,7 @@ class HGCalGeometry final : public CaloSubdetectorGeometry {
CaloSubdetectorGeometry::DimVec& dimVector,
CaloSubdetectorGeometry::IVec& dinsVector) const override;

GlobalPoint getPosition(const DetId& id) const;
GlobalPoint getPosition(const DetId& id, bool debug = false) const;
GlobalPoint getWaferPosition(const DetId& id) const;

/// Returns area of a cell
Expand All @@ -79,7 +79,7 @@ class HGCalGeometry final : public CaloSubdetectorGeometry {
/// Returns the corner points of this cell's volume.
CornersVec getCorners(const DetId& id) const;
CornersVec get8Corners(const DetId& id) const;
CornersVec getNewCorners(const DetId& id) const;
CornersVec getNewCorners(const DetId& id, bool debug = false) const;

// Get neighbor in z along a direction
DetId neighborZ(const DetId& idin, const GlobalVector& p) const;
Expand Down
83 changes: 49 additions & 34 deletions Geometry/HGCalGeometry/src/HGCalGeometry.cc
View file
Open in desktop
Expand Up @@ -23,6 +23,8 @@ typedef std::vector<float> ParmVec;

//#define EDM_ML_DEBUG

const bool debugLocate = false;

HGCalGeometry::HGCalGeometry(const HGCalTopology& topology_)
: m_topology(topology_),
m_validGeomIds(topology_.totalGeomModules()),
Expand All @@ -40,8 +42,6 @@ HGCalGeometry::HGCalGeometry(const HGCalTopology& topology_)
#endif
}

HGCalGeometry::~HGCalGeometry() {}

void HGCalGeometry::fillNamedParams(DDFilteredView fv) {}

void HGCalGeometry::initializeParms() {}
Expand Down Expand Up @@ -183,7 +183,7 @@ std::shared_ptr<const CaloCellGeometry> HGCalGeometry::getGeometry(const DetId&
return nullptr; // nothing to get
DetId geomId = getGeometryDetId(detId);
const uint32_t cellIndex(m_topology.detId2denseGeomId(geomId));
const GlobalPoint pos = (detId != geomId) ? getPosition(detId) : GlobalPoint();
const GlobalPoint pos = (detId != geomId) ? getPosition(detId, false) : GlobalPoint();
return cellGeomPtr(cellIndex, pos);
}

Expand All @@ -195,7 +195,7 @@ bool HGCalGeometry::present(const DetId& detId) const {
return (nullptr != getGeometryRawPtr(index));
}

GlobalPoint HGCalGeometry::getPosition(const DetId& detid) const {
GlobalPoint HGCalGeometry::getPosition(const DetId& detid, bool debug) const {
unsigned int cellIndex = indexFor(detid);
GlobalPoint glob;
unsigned int maxSize = (m_topology.tileTrapezoid() ? m_cellVec2.size() : m_cellVec.size());
Expand All @@ -206,32 +206,30 @@ GlobalPoint HGCalGeometry::getPosition(const DetId& detid) const {
xy = m_topology.dddConstants().locateCellHex(id.iCell1, id.iSec1, true);
const HepGeom::Point3D<float> lcoord(xy.first, xy.second, 0);
glob = m_cellVec[cellIndex].getPosition(lcoord);
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("HGCalGeom") << "getPosition:: index " << cellIndex << " Local " << lcoord.x() << ":"
<< lcoord.y() << " ID " << id.iCell1 << ":" << id.iSec1 << " Global " << glob;
#endif
if (debug)
edm::LogVerbatim("HGCalGeom") << "getPosition:: index " << cellIndex << " Local " << lcoord.x() << ":"
<< lcoord.y() << " ID " << id.iCell1 << ":" << id.iSec1 << " Global " << glob;
} else if (m_topology.tileTrapezoid()) {
const HepGeom::Point3D<float> lcoord(0, 0, 0);
glob = m_cellVec2[cellIndex].getPosition(lcoord);
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("HGCalGeom") << "getPositionTrap:: index " << cellIndex << " Local " << lcoord.x() << ":"
<< lcoord.y() << " ID " << id.iLay << ":" << id.iSec1 << ":" << id.iCell1
<< " Global " << glob;
#endif
if (debug)
edm::LogVerbatim("HGCalGeom") << "getPositionTrap:: index " << cellIndex << " Local " << lcoord.x() << ":"
<< lcoord.y() << " ID " << id.iLay << ":" << id.iSec1 << ":" << id.iCell1
<< " Global " << glob;
} else {
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("HGCalGeom") << "getPosition for " << HGCSiliconDetId(detid) << " Layer " << id.iLay << " Wafer "
<< id.iSec1 << ":" << id.iSec2 << " Cell " << id.iCell1 << ":" << id.iCell2;
#endif
xy = m_topology.dddConstants().locateCell(id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, true, true);
if (debug)
edm::LogVerbatim("HGCalGeom") << "getPosition for " << HGCSiliconDetId(detid) << " Layer " << id.iLay
<< " Wafer " << id.iSec1 << ":" << id.iSec2 << " Cell " << id.iCell1 << ":"
<< id.iCell2;
xy = m_topology.dddConstants().locateCell(
id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, true, false, debug);
double xx = id.zSide * xy.first;
double zz = id.zSide * m_topology.dddConstants().waferZ(id.iLay, true);
glob = GlobalPoint(xx, xy.second, zz);
#ifdef EDM_ML_DEBUG
edm::LogVerbatim("HGCalGeom") << "getPositionWafer:: index " << cellIndex << " Local " << xy.first << ":"
<< xy.second << " ID " << id.iLay << ":" << id.iSec1 << ":" << id.iSec2 << ":"
<< id.iCell1 << ":" << id.iCell2 << " Global " << glob;
#endif
if (debug)
edm::LogVerbatim("HGCalGeom") << "getPositionWafer:: index " << cellIndex << " Local " << xy.first << ":"
<< xy.second << " ID " << id.iLay << ":" << id.iSec1 << ":" << id.iSec2 << ":"
<< id.iCell1 << ":" << id.iCell2 << " Global " << glob;
}
}
return glob;
Expand Down Expand Up @@ -276,7 +274,7 @@ HGCalGeometry::CornersVec HGCalGeometry::getCorners(const DetId& detid) const {
unsigned int cellIndex = indexFor(detid);
HGCalTopology::DecodedDetId id = m_topology.decode(detid);
if (cellIndex < m_cellVec2.size() && m_det == DetId::HGCalHSc) {
GlobalPoint v = getPosition(detid);
GlobalPoint v = getPosition(detid, false);
int type = std::min(id.iType, 1);
std::pair<double, double> rr = m_topology.dddConstants().cellSizeTrap(type, id.iSec1);
float dr = k_half * (rr.second - rr.first);
Expand Down Expand Up @@ -307,7 +305,8 @@ HGCalGeometry::CornersVec HGCalGeometry::getCorners(const DetId& detid) const {
co[i] = m_cellVec[cellIndex].getPosition(lcoord);
}
} else {
xy = m_topology.dddConstants().locateCell(id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false);
xy = m_topology.dddConstants().locateCell(
id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false, true, debugLocate);
float zz = m_topology.dddConstants().waferZ(id.iLay, true);
float dx = k_fac2 * m_cellVec[cellIndex].param()[FlatHexagon::k_r];
float dy = k_fac1 * m_cellVec[cellIndex].param()[FlatHexagon::k_R];
Expand All @@ -332,7 +331,7 @@ HGCalGeometry::CornersVec HGCalGeometry::get8Corners(const DetId& detid) const {
unsigned int cellIndex = indexFor(detid);
HGCalTopology::DecodedDetId id = m_topology.decode(detid);
if (cellIndex < m_cellVec2.size() && m_det == DetId::HGCalHSc) {
GlobalPoint v = getPosition(detid);
GlobalPoint v = getPosition(detid, false);
int type = std::min(id.iType, 1);
std::pair<double, double> rr = m_topology.dddConstants().cellSizeTrap(type, id.iSec1);
float dr = k_half * (rr.second - rr.first);
Expand Down Expand Up @@ -363,7 +362,8 @@ HGCalGeometry::CornersVec HGCalGeometry::get8Corners(const DetId& detid) const {
co[i] = m_cellVec[cellIndex].getPosition(lcoord);
}
} else {
xy = m_topology.dddConstants().locateCell(id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false);
xy = m_topology.dddConstants().locateCell(
id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false, true, debugLocate);
dx = k_fac2 * m_cellVec[cellIndex].param()[FlatHexagon::k_r];
float dy = k_fac1 * m_cellVec[cellIndex].param()[FlatHexagon::k_R];
float dz = -id.zSide * m_cellVec[cellIndex].param()[FlatHexagon::k_dZ];
Expand All @@ -379,13 +379,16 @@ HGCalGeometry::CornersVec HGCalGeometry::get8Corners(const DetId& detid) const {
return co;
}

HGCalGeometry::CornersVec HGCalGeometry::getNewCorners(const DetId& detid) const {
HGCalGeometry::CornersVec HGCalGeometry::getNewCorners(const DetId& detid, bool debug) const {
unsigned int ncorner = (m_det == DetId::HGCalHSc) ? 5 : 7;
HGCalGeometry::CornersVec co(ncorner, GlobalPoint(0, 0, 0));
unsigned int cellIndex = indexFor(detid);
HGCalTopology::DecodedDetId id = m_topology.decode(detid);
if (debug)
edm::LogVerbatim("HGCalGeom") << "NewCorners for Layer " << id.iLay << " Wafer " << id.iSec1 << ":" << id.iSec2
<< " Cell " << id.iCell1 << ":" << id.iCell2;
if (cellIndex < m_cellVec2.size() && m_det == DetId::HGCalHSc) {
GlobalPoint v = getPosition(detid);
GlobalPoint v = getPosition(detid, false);
int type = std::min(id.iType, 1);
std::pair<double, double> rr = m_topology.dddConstants().cellSizeTrap(type, id.iSec1);
float dr = k_half * (rr.second - rr.first);
Expand All @@ -407,18 +410,30 @@ HGCalGeometry::CornersVec HGCalGeometry::getNewCorners(const DetId& detid) const
float dz = -id.zSide * m_cellVec[cellIndex].param()[FlatHexagon::k_dZ];
static const int signx[] = {1, -1, -2, -1, 1, 2};
static const int signy[] = {1, 1, 0, -1, -1, 0};
#ifdef EDM_ML_DEBUG
if (debug)
edm::LogVerbatim("HGCalGeom") << "kfac " << k_fac1 << ":" << k_fac2 << " dx:dy:dz " << dx << ":" << dy << ":"
<< dz;
#endif
if (m_topology.waferHexagon6()) {
xy = m_topology.dddConstants().locateCellHex(id.iCell1, id.iSec1, true);
for (unsigned int i = 0; i < ncorner - 1; ++i) {
const HepGeom::Point3D<float> lcoord(xy.first + signx[i] * dx, xy.second + signy[i] * dy, dz);
co[i] = m_cellVec[cellIndex].getPosition(lcoord);
}
} else {
xy = m_topology.dddConstants().locateCell(id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false);
xy = m_topology.dddConstants().locateCell(
id.iLay, id.iSec1, id.iSec2, id.iCell1, id.iCell2, true, false, true, debug);
float zz = m_topology.dddConstants().waferZ(id.iLay, true);
for (unsigned int i = 0; i < ncorner; ++i) {
auto xyglob = m_topology.dddConstants().localToGlobal8(
id.iLay, id.iSec1, id.iSec2, (xy.first + signx[i] * dx), (xy.second + signy[i] * dy), true, false);
double xloc = xy.first + signx[i] * dx;
double yloc = xy.second + signy[i] * dy;
#ifdef EDM_ML_DEBUG
if (debug)
edm::LogVerbatim("HGCalGeom") << "Corner " << i << " x " << xy.first << ":" << xloc << " y " << xy.second
<< ":" << yloc << " z " << zz << ":" << id.zSide * (zz + dz);
#endif
auto xyglob = m_topology.dddConstants().localToGlobal8(id.iLay, id.iSec1, id.iSec2, xloc, yloc, true, debug);
double xx = id.zSide * xyglob.first;
co[i] = GlobalPoint(xx, xyglob.second, id.zSide * (zz + dz));
}
Expand All @@ -440,7 +455,7 @@ DetId HGCalGeometry::neighborZ(const DetId& idin, const GlobalVector& momentum)
#endif
if ((lay >= m_topology.dddConstants().firstLayer()) && (lay <= m_topology.dddConstants().lastLayer(true)) &&
(momentum.z() != 0.0)) {
GlobalPoint v = getPosition(idin);
GlobalPoint v = getPosition(idin, false);
double z = id.zSide * m_topology.dddConstants().waferZ(lay, true);
double grad = (z - v.z()) / momentum.z();
GlobalPoint p(v.x() + grad * momentum.x(), v.y() + grad * momentum.y(), z);
Expand Down Expand Up @@ -471,7 +486,7 @@ DetId HGCalGeometry::neighborZ(const DetId& idin,
#endif
if ((lay >= m_topology.dddConstants().firstLayer()) && (lay <= m_topology.dddConstants().lastLayer(true)) &&
(momentum.z() != 0.0)) {
GlobalPoint v = getPosition(idin);
GlobalPoint v = getPosition(idin, false);
double z = id.zSide * m_topology.dddConstants().waferZ(lay, true);
FreeTrajectoryState fts(v, momentum, charge, bField);
Plane::PlanePointer nPlane = Plane::build(Plane::PositionType(0, 0, z), Plane::RotationType());
Expand Down