Go to the documentation of this file.00001 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
00002 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
00003 #include "Geometry/HcalTowerAlgo/interface/HcalDDDGeometry.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 #include "Geometry/CaloGeometry/interface/CaloGenericDetId.h"
00006
00007 #include <algorithm>
00008
00009
00010
00011 HcalDDDGeometry::HcalDDDGeometry() :
00012 lastReqDet_(DetId::Detector(0)),
00013 lastReqSubdet_(0),
00014 etaMax_(0),
00015 firstHFQuadRing_(40) ,
00016 m_hbCellVec ( HcalDetId::kHBSize ) ,
00017 m_heCellVec ( HcalDetId::kHESize ) ,
00018 m_hoCellVec ( HcalDetId::kHOSize ) ,
00019 m_hfCellVec ( HcalDetId::kHFSize )
00020 {
00021 twopi = M_PI + M_PI;
00022 deg = M_PI/180.;
00023 }
00024
00025
00026 HcalDDDGeometry::~HcalDDDGeometry()
00027 {
00028 }
00029
00030
00031 std::vector<DetId> const & HcalDDDGeometry::getValidDetIds(DetId::Detector det,
00032 int subdet) const {
00033
00034 const std::vector<DetId>& baseIds(CaloSubdetectorGeometry::getValidDetIds());
00035 if (det == DetId::Detector( 0 ) && subdet == 0) {
00036 return baseIds ;
00037 }
00038
00039 if (lastReqDet_ != det || lastReqSubdet_ != subdet ) {
00040 lastReqDet_ = det ;
00041 lastReqSubdet_ = subdet ;
00042 m_validIds.clear();
00043 m_validIds.reserve( baseIds.size() ) ;
00044 }
00045
00046 if (m_validIds.empty() ) {
00047 for (unsigned int i = 0 ; i != baseIds.size() ; ++i ) {
00048 const DetId id ( baseIds[i] );
00049 if (id.det() == det && id.subdetId() == subdet ) {
00050 m_validIds.push_back( id ) ;
00051 }
00052 }
00053 std::sort(m_validIds.begin(),m_validIds.end());
00054 }
00055
00056 #ifdef DebugLog
00057 LogDebug("HCalGeom") << "HcalDDDGeometry::getValidDetIds: "
00058 << m_validIds.size() << " valid IDs found for detector "
00059 << det << " Sub-detector " << subdet;
00060 #endif
00061 return m_validIds;
00062 }
00063
00064
00065 DetId HcalDDDGeometry::getClosestCell(const GlobalPoint& r) const {
00066
00067
00068
00069 double abseta = fabs(r.eta());
00070 double phi = r.phi();
00071 if (phi < 0) phi += twopi;
00072 double radius = r.mag();
00073 double z = fabs(r.z());
00074 #ifdef DebugLog
00075 LogDebug("HCalGeom") << "HcalDDDGeometry::getClosestCell for eta "
00076 << r.eta() << " phi " << phi/deg << " z " << r.z()
00077 << " radius " << radius;
00078 #endif
00079
00080 HcalDetId bestId;
00081 if (abseta <= etaMax_) {
00082 for (unsigned int i=0; i<hcalCells_.size(); i++) {
00083 if (abseta >=hcalCells_[i].etaMin() && abseta <=hcalCells_[i].etaMax()) {
00084 HcalSubdetector bc = hcalCells_[i].detType();
00085 int etaring = hcalCells_[i].etaBin();
00086 int phibin = 0;
00087 if (hcalCells_[i].unitPhi() == 4) {
00088
00089
00090
00091
00092 phibin = static_cast<int>(((phi/deg)+hcalCells_[i].phiOffset()+
00093 0.5*hcalCells_[i].phiBinWidth())/
00094 hcalCells_[i].phiBinWidth());
00095 if (phibin == 0) phibin = hcalCells_[i].nPhiBins();
00096 phibin = phibin*4 - 1;
00097 } else {
00098 phibin = static_cast<int>(((phi/deg)+hcalCells_[i].phiOffset())/
00099 hcalCells_[i].phiBinWidth()) + 1;
00100
00101 phibin = (phibin-1)*(hcalCells_[i].unitPhi()) + 1;
00102 }
00103
00104 int dbin = 1;
00105 int etabin = (r.z() > 0) ? etaring : -etaring;
00106 if (bc == HcalForward) {
00107 bestId = HcalDetId(bc, etabin, phibin, dbin);
00108 break;
00109 } else {
00110 double rz = z;
00111 if (hcalCells_[i].depthType()) rz = radius;
00112 if (rz < hcalCells_[i].depthMax()) {
00113 dbin = hcalCells_[i].depthSegment();
00114 bestId = HcalDetId(bc, etabin, phibin, dbin);
00115 break;
00116 }
00117 }
00118 }
00119 }
00120 }
00121 #ifdef DebugLog
00122 LogDebug("HCalGeom") << "HcalDDDGeometry::getClosestCell " << bestId;
00123 #endif
00124 return bestId;
00125 }
00126
00127
00128 int HcalDDDGeometry::insertCell(std::vector<HcalCellType> const & cells){
00129
00130 hcalCells_.insert(hcalCells_.end(), cells.begin(), cells.end());
00131 int num = static_cast<int>(hcalCells_.size());
00132 for (unsigned int i=0; i<cells.size(); i++) {
00133 if (cells[i].etaMax() > etaMax_ ) etaMax_ = cells[i].etaMax();
00134 }
00135 #ifdef DebugLog
00136 LogDebug("HCalGeom") << "HcalDDDGeometry::insertCell " << cells.size()
00137 << " cells inserted == Total " << num
00138 << " EtaMax = " << etaMax_;
00139 #endif
00140 return num;
00141 }
00142
00143 void
00144 HcalDDDGeometry::newCell( const GlobalPoint& f1 ,
00145 const GlobalPoint& f2 ,
00146 const GlobalPoint& f3 ,
00147 const CCGFloat* parm ,
00148 const DetId& detId )
00149 {
00150 const CaloGenericDetId cgid ( detId ) ;
00151
00152 const unsigned int din ( cgid.denseIndex() ) ;
00153
00154 assert( cgid.isHcal() ) ;
00155
00156 if( cgid.isHB() )
00157 {
00158 m_hbCellVec[ din ] = IdealObliquePrism( f1, cornersMgr(), parm ) ;
00159 }
00160 else
00161 {
00162 if( cgid.isHE() )
00163 {
00164 const unsigned int index ( din - m_hbCellVec.size() ) ;
00165 m_heCellVec[ index ] = IdealObliquePrism( f1, cornersMgr(), parm ) ;
00166 }
00167 else
00168 {
00169 if( cgid.isHO() )
00170 {
00171 const unsigned int index ( din
00172 - m_hbCellVec.size()
00173 - m_heCellVec.size() ) ;
00174 m_hoCellVec[ index ] = IdealObliquePrism( f1, cornersMgr(), parm ) ;
00175 }
00176 else
00177 {
00178 const unsigned int index ( din
00179 - m_hbCellVec.size()
00180 - m_heCellVec.size()
00181 - m_hoCellVec.size() ) ;
00182 m_hfCellVec[ index ] = IdealZPrism( f1, cornersMgr(), parm ) ;
00183 }
00184 }
00185 }
00186 m_validIds.push_back( detId ) ;
00187 }
00188
00189 const CaloCellGeometry*
00190 HcalDDDGeometry::cellGeomPtr( uint32_t din ) const
00191 {
00192 const CaloCellGeometry* cell ( 0 ) ;
00193 if( m_hbCellVec.size() > din )
00194 {
00195 cell = &m_hbCellVec[ din ] ;
00196 }
00197 else
00198 {
00199 if( m_hbCellVec.size() +
00200 m_heCellVec.size() > din )
00201 {
00202 const unsigned int index ( din - m_hbCellVec.size() ) ;
00203 cell = &m_heCellVec[ index ] ;
00204 }
00205 else
00206 {
00207 if( m_hbCellVec.size() +
00208 m_heCellVec.size() +
00209 m_hoCellVec.size() > din )
00210 {
00211 const unsigned int index ( din
00212 - m_hbCellVec.size()
00213 - m_heCellVec.size() ) ;
00214 cell = &m_hoCellVec[ index ] ;
00215 }
00216 else
00217 {
00218 if( m_hbCellVec.size() +
00219 m_heCellVec.size() +
00220 m_hoCellVec.size() +
00221 m_hfCellVec.size() > din )
00222 {
00223 const unsigned int index ( din
00224 - m_hbCellVec.size()
00225 - m_heCellVec.size()
00226 - m_hoCellVec.size() ) ;
00227 cell = &m_hfCellVec[ index ] ;
00228 }
00229 }
00230 }
00231 }
00232 return ( 0 == cell || 0 == cell->param() ? 0 : cell ) ;
00233 }