00001 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
00002 #include "Geometry/CaloGeometry/interface/CaloGenericDetId.h"
00003
00004 #include <Math/Transform3D.h>
00005 #include <Math/EulerAngles.h>
00006
00007 typedef CaloCellGeometry::Pt3D Pt3D ;
00008 typedef CaloCellGeometry::Pt3DVec Pt3DVec ;
00009 typedef CaloCellGeometry::Tr3D Tr3D ;
00010
00011 typedef CaloSubdetectorGeometry::CCGFloat CCGFloat ;
00012
00013 CaloSubdetectorGeometry::CaloSubdetectorGeometry() :
00014 m_parMgr ( 0 ) ,
00015 m_cmgr ( 0 ) ,
00016 m_sortedIds (false) ,
00017 m_deltaPhi ( 0 ) ,
00018 m_deltaEta ( 0 )
00019 {
00020 }
00021
00022
00023 CaloSubdetectorGeometry::~CaloSubdetectorGeometry()
00024 {
00025 delete m_cmgr ;
00026 delete m_parMgr ;
00027 delete m_deltaPhi ;
00028 delete m_deltaEta ;
00029 }
00030
00031 const std::vector<DetId>&
00032 CaloSubdetectorGeometry::getValidDetIds( DetId::Detector det ,
00033 int subdet ) const
00034 {
00035 if( !m_sortedIds )
00036 {
00037 m_sortedIds = true ;
00038 std::sort( m_validIds.begin(), m_validIds.end() ) ;
00039 }
00040 return m_validIds ;
00041 }
00042
00043 const CaloCellGeometry*
00044 CaloSubdetectorGeometry::getGeometry( const DetId& id ) const
00045 {
00046 return cellGeomPtr( CaloGenericDetId( id ).denseIndex() ) ;
00047 }
00048
00049 bool
00050 CaloSubdetectorGeometry::present( const DetId& id ) const
00051 {
00052 return ( 0 != getGeometry( id ) ) ;
00053 }
00054
00055 DetId
00056 CaloSubdetectorGeometry::getClosestCell( const GlobalPoint& r ) const
00057 {
00058 const CCGFloat eta ( r.eta() ) ;
00059 const CCGFloat phi ( r.phi() ) ;
00060 uint32_t index ( ~0 ) ;
00061 CCGFloat closest ( 1e9 ) ;
00062
00063 for( uint32_t i ( 0 ); i != m_validIds.size() ; ++i )
00064 {
00065 const CaloCellGeometry* cell ( getGeometry( m_validIds[ i ] ) ) ;
00066 if( 0 != cell )
00067 {
00068 const GlobalPoint& p ( cell->getPosition() ) ;
00069 const CCGFloat eta0 ( p.eta() ) ;
00070 const CCGFloat phi0 ( p.phi() ) ;
00071 const CCGFloat dR2 ( reco::deltaR2( eta0, phi0, eta, phi ) ) ;
00072 if( dR2 < closest )
00073 {
00074 closest = dR2 ;
00075 index = i ;
00076 }
00077 }
00078 }
00079 return ( closest > 0.9e9 ||
00080 (uint32_t)(~0) == index ? DetId(0) :
00081 m_validIds[index] ) ;
00082 }
00083
00084 CaloSubdetectorGeometry::DetIdSet
00085 CaloSubdetectorGeometry::getCells( const GlobalPoint& r,
00086 double dR ) const
00087 {
00088 const double dR2 ( dR*dR ) ;
00089 const double eta ( r.eta() ) ;
00090 const double phi ( r.phi() ) ;
00091
00092 DetIdSet dss;
00093
00094 if( 0.000001 < dR )
00095 {
00096 for( uint32_t i ( 0 ); i != m_validIds.size() ; ++i )
00097 {
00098 const CaloCellGeometry* cell ( getGeometry( m_validIds[i] ) ) ;
00099 if( 0 != cell )
00100 {
00101 const GlobalPoint& p ( cell->getPosition() ) ;
00102 const CCGFloat eta0 ( p.eta() ) ;
00103 if( fabs( eta - eta0 ) < dR )
00104 {
00105 const CCGFloat phi0 ( p.phi() ) ;
00106 CCGFloat delp ( fabs( phi - phi0 ) ) ;
00107 if( delp > M_PI ) delp = 2*M_PI - delp ;
00108 if( delp < dR )
00109 {
00110 const CCGFloat dist2 ( reco::deltaR2( eta0, phi0, eta, phi ) ) ;
00111 if( dist2 < dR2 ) dss.insert( m_validIds[i] ) ;
00112 }
00113 }
00114 }
00115 }
00116 }
00117 return dss;
00118 }
00119
00120 void
00121 CaloSubdetectorGeometry::allocateCorners( CaloCellGeometry::CornersVec::size_type n )
00122 {
00123 assert( 0 == m_cmgr ) ;
00124 m_cmgr = new CaloCellGeometry::CornersMgr( n*( CaloCellGeometry::k_cornerSize ),
00125 CaloCellGeometry::k_cornerSize ) ;
00126
00127 m_validIds.reserve( n ) ;
00128 }
00129
00130 void
00131 CaloSubdetectorGeometry::allocatePar( ParVec::size_type n,
00132 unsigned int m )
00133 {
00134 assert( 0 == m_parMgr ) ;
00135 m_parMgr = new ParMgr( n*m, m ) ;
00136 }
00137
00138 void
00139 CaloSubdetectorGeometry::getSummary( CaloSubdetectorGeometry::TrVec& tVec ,
00140 CaloSubdetectorGeometry::IVec& iVec ,
00141 CaloSubdetectorGeometry::DimVec& dVec ) const
00142 {
00143 tVec.reserve( m_validIds.size()*numberOfTransformParms() ) ;
00144 iVec.reserve( numberOfShapes()==1 ? 1 : m_validIds.size() ) ;
00145 dVec.reserve( numberOfShapes()*numberOfParametersPerShape() ) ;
00146
00147 for( ParVecVec::const_iterator ivv ( parVecVec().begin() ) ; ivv != parVecVec().end() ; ++ivv )
00148 {
00149 const ParVec& pv ( *ivv ) ;
00150 for( ParVec::const_iterator iv ( pv.begin() ) ; iv != pv.end() ; ++iv )
00151 {
00152 dVec.push_back( *iv ) ;
00153 }
00154 }
00155
00156 for( uint32_t i ( 0 ) ; i != m_validIds.size() ; ++i )
00157 {
00158 Tr3D tr ;
00159 const CaloCellGeometry* ptr ( cellGeomPtr( i ) ) ;
00160 assert( 0 != ptr ) ;
00161 ptr->getTransform( tr, ( Pt3DVec* ) 0 ) ;
00162
00163 if( Tr3D() == tr )
00164 {
00165 const GlobalPoint& gp ( ptr->getPosition() ) ;
00166 tr = HepGeom::Translate3D( gp.x(), gp.y(), gp.z() ) ;
00167 }
00168
00169 const CLHEP::Hep3Vector tt ( tr.getTranslation() ) ;
00170 tVec.push_back( tt.x() ) ;
00171 tVec.push_back( tt.y() ) ;
00172 tVec.push_back( tt.z() ) ;
00173 if( 6 == numberOfTransformParms() )
00174 {
00175 const CLHEP::HepRotation rr ( tr.getRotation() ) ;
00176 const ROOT::Math::Transform3D rtr ( rr.xx(), rr.xy(), rr.xz(), tt.x(),
00177 rr.yx(), rr.yy(), rr.yz(), tt.y(),
00178 rr.zx(), rr.zy(), rr.zz(), tt.z() ) ;
00179 ROOT::Math::EulerAngles ea ;
00180 rtr.GetRotation( ea ) ;
00181 tVec.push_back( ea.Phi() ) ;
00182 tVec.push_back( ea.Theta() ) ;
00183 tVec.push_back( ea.Psi() ) ;
00184 }
00185
00186 const CCGFloat* par ( ptr->param() ) ;
00187
00188 unsigned int ishape ( 9999 ) ;
00189 for( unsigned int ivv ( 0 ) ; ivv != parVecVec().size() ; ++ivv )
00190 {
00191 bool ok ( true ) ;
00192 const CCGFloat* pv ( &(*parVecVec()[ivv].begin() ) ) ;
00193 for( unsigned int k ( 0 ) ; k != numberOfParametersPerShape() ; ++k )
00194 {
00195 ok = ok && ( fabs( par[k] - pv[k] ) < 1.e-6 ) ;
00196 }
00197 if( ok )
00198 {
00199 ishape = ivv ;
00200 break ;
00201 }
00202 }
00203 assert( 9999 != ishape ) ;
00204
00205 const unsigned int nn (( numberOfShapes()==1) ? (unsigned int)1 : m_validIds.size() ) ;
00206 if( iVec.size() < nn ) iVec.push_back( ishape ) ;
00207 }
00208 }
00209
00210 CCGFloat
00211 CaloSubdetectorGeometry::deltaPhi( const DetId& detId ) const
00212 {
00213 const CaloGenericDetId cgId ( detId ) ;
00214
00215 if( 0 == m_deltaPhi )
00216 {
00217 const uint32_t kSize ( cgId.sizeForDenseIndexing() ) ;
00218 m_deltaPhi = new std::vector<CCGFloat> ( kSize ) ;
00219 for( uint32_t i ( 0 ) ; i != kSize ; ++i )
00220 {
00221 const CaloCellGeometry* cellPtr ( cellGeomPtr( i ) ) ;
00222 if( 0 != cellPtr )
00223 {
00224 const CaloCellGeometry& cell ( *cellPtr ) ;
00225 CCGFloat dPhi1 ( fabs(
00226 GlobalPoint( ( cell.getCorners()[0].x() +
00227 cell.getCorners()[1].x() )/2. ,
00228 ( cell.getCorners()[0].y() +
00229 cell.getCorners()[1].y() )/2. ,
00230 ( cell.getCorners()[0].z() +
00231 cell.getCorners()[1].z() )/2. ).phi() -
00232 GlobalPoint( ( cell.getCorners()[2].x() +
00233 cell.getCorners()[3].x() )/2. ,
00234 ( cell.getCorners()[2].y() +
00235 cell.getCorners()[3].y() )/2. ,
00236 ( cell.getCorners()[2].z() +
00237 cell.getCorners()[3].z() )/2. ).phi() ) ) ;
00238 CCGFloat dPhi2 ( fabs(
00239 GlobalPoint( ( cell.getCorners()[0].x() +
00240 cell.getCorners()[3].x() )/2. ,
00241 ( cell.getCorners()[0].y() +
00242 cell.getCorners()[3].y() )/2. ,
00243 ( cell.getCorners()[0].z() +
00244 cell.getCorners()[3].z() )/2. ).phi() -
00245 GlobalPoint( ( cell.getCorners()[2].x() +
00246 cell.getCorners()[1].x() )/2. ,
00247 ( cell.getCorners()[2].y() +
00248 cell.getCorners()[1].y() )/2. ,
00249 ( cell.getCorners()[2].z() +
00250 cell.getCorners()[1].z() )/2. ).phi() ) ) ;
00251 if( M_PI < dPhi1 ) dPhi1 = fabs( dPhi1 - 2.*M_PI ) ;
00252 if( M_PI < dPhi2 ) dPhi2 = fabs( dPhi2 - 2.*M_PI ) ;
00253 (*m_deltaPhi)[i] = dPhi1>dPhi2 ? dPhi1 : dPhi2 ;
00254 }
00255 }
00256 }
00257 return (*m_deltaPhi)[ cgId.denseIndex() ] ;
00258 }
00259
00260 CCGFloat
00261 CaloSubdetectorGeometry::deltaEta( const DetId& detId ) const
00262 {
00263 const CaloGenericDetId cgId ( detId ) ;
00264
00265 if( 0 == m_deltaEta )
00266 {
00267 const uint32_t kSize ( cgId.sizeForDenseIndexing() ) ;
00268 m_deltaEta = new std::vector<CCGFloat> ( kSize ) ;
00269 for( uint32_t i ( 0 ) ; i != kSize ; ++i )
00270 {
00271 const CaloCellGeometry* cellPtr ( cellGeomPtr( i ) ) ;
00272 if( 0 != cellPtr )
00273 {
00274 const CaloCellGeometry& cell ( *cellPtr ) ;
00275 const CCGFloat dEta1 ( fabs(
00276 GlobalPoint( ( cell.getCorners()[0].x() +
00277 cell.getCorners()[1].x() )/2. ,
00278 ( cell.getCorners()[0].y() +
00279 cell.getCorners()[1].y() )/2. ,
00280 ( cell.getCorners()[0].z() +
00281 cell.getCorners()[1].z() )/2. ).eta() -
00282 GlobalPoint( ( cell.getCorners()[2].x() +
00283 cell.getCorners()[3].x() )/2. ,
00284 ( cell.getCorners()[2].y() +
00285 cell.getCorners()[3].y() )/2. ,
00286 ( cell.getCorners()[2].z() +
00287 cell.getCorners()[3].z() )/2. ).eta() ) ) ;
00288 const CCGFloat dEta2 ( fabs(
00289 GlobalPoint( ( cell.getCorners()[0].x() +
00290 cell.getCorners()[3].x() )/2. ,
00291 ( cell.getCorners()[0].y() +
00292 cell.getCorners()[3].y() )/2. ,
00293 ( cell.getCorners()[0].z() +
00294 cell.getCorners()[3].z() )/2. ).eta() -
00295 GlobalPoint( ( cell.getCorners()[2].x() +
00296 cell.getCorners()[1].x() )/2. ,
00297 ( cell.getCorners()[2].y() +
00298 cell.getCorners()[1].y() )/2. ,
00299 ( cell.getCorners()[2].z() +
00300 cell.getCorners()[1].z() )/2. ).eta() ) ) ;
00301 (*m_deltaEta)[i] = dEta1>dEta2 ? dEta1 : dEta2 ;
00302 }
00303 }
00304 }
00305 return (*m_deltaEta)[ cgId.denseIndex() ] ;
00306 }
1.7.3