CaloGeometry/src/IdealObliquePrism.cc

0001 #include "Geometry/CaloGeometry/interface/IdealObliquePrism.h"
0002 #include <cmath>
0003
0004 typedef IdealObliquePrism::CCGFloat CCGFloat;
0005 typedef IdealObliquePrism::Pt3D Pt3D;
0006 typedef IdealObliquePrism::Pt3DVec Pt3DVec;
0007
0008 IdealObliquePrism::IdealObliquePrism() : CaloCellGeometry() {}
0009
0010 IdealObliquePrism::IdealObliquePrism(const IdealObliquePrism& idop) : CaloCellGeometry(idop) { *this = idop; }
0011
0012 IdealObliquePrism& IdealObliquePrism::operator=(const IdealObliquePrism& idop) {
0013   if (&idop != this)
0014     CaloCellGeometry::operator=(idop);
0015   return *this;
0016 }
0017
0018 IdealObliquePrism::IdealObliquePrism(const GlobalPoint& faceCenter, CornersMgr* mgr, const CCGFloat* parm)
0019     : CaloCellGeometry(faceCenter, mgr, parm) {
0020   initSpan();
0021 }
0022
0023 IdealObliquePrism::~IdealObliquePrism() {}
0024
0025 CCGFloat IdealObliquePrism::dEta() const { return param()[IdealObliquePrism::k_dEta]; }
0026
0027 CCGFloat IdealObliquePrism::dPhi() const { return param()[IdealObliquePrism::k_dPhi]; }
0028
0029 CCGFloat IdealObliquePrism::dz() const { return param()[IdealObliquePrism::k_dZ]; }
0030
0031 CCGFloat IdealObliquePrism::eta() const { return param()[IdealObliquePrism::k_Eta]; }
0032
0033 CCGFloat IdealObliquePrism::z() const { return param()[IdealObliquePrism::k_Z]; }
0034
0035 void IdealObliquePrism::vocalCorners(Pt3DVec& vec, const CCGFloat* pv, Pt3D& ref) const { localCorners(vec, pv, ref); }
0036
0037 GlobalPoint IdealObliquePrism::etaPhiPerp(float eta, float phi, float perp) {
0038   return GlobalPoint(perp * cosf(phi), perp * sinf(phi), perp * sinhf(eta));
0039 }
0040
0041 GlobalPoint IdealObliquePrism::etaPhiZ(float eta, float phi, float z) {
0042   return GlobalPoint(z * cosf(phi) / sinhf(eta), z * sinf(phi) / sinhf(eta), z);
0043 }
0044
0045 void IdealObliquePrism::localCorners(Pt3DVec& lc, const CCGFloat* pv, Pt3D& ref) {
0046   assert(8 == lc.size());
0047   assert(nullptr != pv);
0048
0049   const CCGFloat dEta(pv[IdealObliquePrism::k_dEta]);
0050   const CCGFloat dPhi(pv[IdealObliquePrism::k_dPhi]);
0051   const CCGFloat dz(pv[IdealObliquePrism::k_dZ]);
0052   const CCGFloat eta(pv[IdealObliquePrism::k_Eta]);
0053   const CCGFloat z(pv[IdealObliquePrism::k_Z]);
0054
0055   std::vector<GlobalPoint> gc(8, GlobalPoint(0, 0, 0));
0056
0057   const GlobalPoint p(etaPhiZ(eta, 0, z));
0058
0059   if (0 < dz) {
0060     const CCGFloat r_near(p.perp() / cos(dPhi));
0061     const CCGFloat r_far(r_near * ((p.mag() + 2 * dz) / p.mag()));
0062     gc[0] = etaPhiPerp(eta + dEta, +dPhi, r_near);  // (+,+,near)
0063     gc[1] = etaPhiPerp(eta + dEta, -dPhi, r_near);  // (+,-,near)
0064     gc[2] = etaPhiPerp(eta - dEta, -dPhi, r_near);  // (-,-,near)
0065     gc[3] = etaPhiPerp(eta - dEta, +dPhi, r_near);  // (-,+,near)
0066     gc[4] = etaPhiPerp(eta + dEta, +dPhi, r_far);   // (+,+,far)
0067     gc[5] = etaPhiPerp(eta + dEta, -dPhi, r_far);   // (+,-,far)
0068     gc[6] = etaPhiPerp(eta - dEta, -dPhi, r_far);   // (-,-,far)
0069     gc[7] = etaPhiPerp(eta - dEta, +dPhi, r_far);   // (-,+,far)
0070   } else {
0071     const CCGFloat z_near(z);
0072     const CCGFloat z_far(z * (1 - 2 * dz / p.mag()));
0073     gc[0] = etaPhiZ(eta + dEta, +dPhi, z_near);  // (+,+,near)
0074     gc[1] = etaPhiZ(eta + dEta, -dPhi, z_near);  // (+,-,near)
0075     gc[2] = etaPhiZ(eta - dEta, -dPhi, z_near);  // (-,-,near)
0076     gc[3] = etaPhiZ(eta - dEta, +dPhi, z_near);  // (-,+,near)
0077     gc[4] = etaPhiZ(eta + dEta, +dPhi, z_far);   // (+,+,far)
0078     gc[5] = etaPhiZ(eta + dEta, -dPhi, z_far);   // (+,-,far)
0079     gc[6] = etaPhiZ(eta - dEta, -dPhi, z_far);   // (-,-,far)
0080     gc[7] = etaPhiZ(eta - dEta, +dPhi, z_far);   // (-,+,far)
0081   }
0082   for (unsigned int i(0); i != 8; ++i) {
0083     lc[i] = Pt3D(gc[i].x(), gc[i].y(), gc[i].z());
0084   }
0085
0086   ref = 0.25 * (lc[0] + lc[1] + lc[2] + lc[3]);
0087 }
0088
0089 void IdealObliquePrism::initCorners(CaloCellGeometry::CornersVec& co) {
0090   if (co.uninitialized()) {
0091     CornersVec& corners(co);
0092     if (dz() > 0) {
0093       /* In this case, the faces are parallel to the zaxis.
0094         This implies that all corners will have the same
0095         cylindrical radius.
0096      */
0097       const GlobalPoint p(getPosition());
0098       const CCGFloat r_near(p.perp() / cos(dPhi()));
0099       const CCGFloat r_far(r_near * ((p.mag() + 2 * dz()) / p.mag()));
0100       const CCGFloat eta(p.eta());
0101       const CCGFloat phi(p.phi());
0102       corners[0] = etaPhiPerp(eta + dEta(), phi + dPhi(), r_near);  // (+,+,near)
0103       corners[1] = etaPhiPerp(eta + dEta(), phi - dPhi(), r_near);  // (+,-,near)
0104       corners[2] = etaPhiPerp(eta - dEta(), phi - dPhi(), r_near);  // (-,-,near)
0105       corners[3] = etaPhiPerp(eta - dEta(), phi + dPhi(), r_near);  // (-,+,near)
0106       corners[4] = etaPhiPerp(eta + dEta(), phi + dPhi(), r_far);   // (+,+,far)
0107       corners[5] = etaPhiPerp(eta + dEta(), phi - dPhi(), r_far);   // (+,-,far)
0108       corners[6] = etaPhiPerp(eta - dEta(), phi - dPhi(), r_far);   // (-,-,far)
0109       corners[7] = etaPhiPerp(eta - dEta(), phi + dPhi(), r_far);   // (-,+,far)
0110     } else {
0111       /* In this case, the faces are perpendicular to the zaxis.
0112         This implies that all corners will have the same
0113         z-dimension.
0114      */
0115       const GlobalPoint p(getPosition());
0116       const CCGFloat z_near(p.z());
0117       const CCGFloat mag(p.mag());
0118       const CCGFloat z_far(z_near * (1 - 2 * dz() / mag));  // negative to correct sign
0119       const CCGFloat eta(p.eta());
0120       const CCGFloat phi(p.phi());
0121
0122       corners[0] = etaPhiZ(eta + dEta(), phi + dPhi(), z_near);  // (+,+,near)
0123       corners[1] = etaPhiZ(eta + dEta(), phi - dPhi(), z_near);  // (+,-,near)
0124       corners[2] = etaPhiZ(eta - dEta(), phi - dPhi(), z_near);  // (-,-,near)
0125       corners[3] = etaPhiZ(eta - dEta(), phi + dPhi(), z_near);  // (-,+,near)
0126       corners[4] = etaPhiZ(eta + dEta(), phi + dPhi(), z_far);   // (+,+,far)
0127       corners[5] = etaPhiZ(eta + dEta(), phi - dPhi(), z_far);   // (+,-,far)
0128       corners[6] = etaPhiZ(eta - dEta(), phi - dPhi(), z_far);   // (-,-,far)
0129       corners[7] = etaPhiZ(eta - dEta(), phi + dPhi(), z_far);   // (-,+,far)
0130     }
0131   }
0132 }
0133
0134 std::ostream& operator<<(std::ostream& s, const IdealObliquePrism& cell) {
0135   s << "Center: " << cell.getPosition() << std::endl;
0136   s << "dEta = " << cell.dEta() << ", dPhi = " << cell.dPhi() << ", dz = " << cell.dz() << std::endl;
0137   return s;
0138 }