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