CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch1/src/Geometry/HcalTowerAlgo/src/HcalDDDGeometry.cc

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 //#define DebugLog
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   // Now find the closest eta_bin, eta value of a bin i is average
00068   // of eta[i] and eta[i-1]
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           // rings 40 and 41 are offset wrt the other phi numbering
00089           //  1        1         1         2
00090           //  ------------------------------
00091           //  72       36        36        1
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           // convert to the convention of numbering 1,3,5, in 36 phi bins
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 }