Go to the documentation of this file.00001 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
00002 #include <CLHEP/Geometry/Plane3D.h>
00003
00004 const float CaloCellGeometry::k_ScaleFromDDDtoGeant ( 0.1 ) ;
00005
00006 const CaloCellGeometry::CornersVec&
00007 CaloCellGeometry::getCorners() const
00008 {
00009 return m_corners ;
00010 }
00011
00012 std::ostream& operator<<( std::ostream& s, const CaloCellGeometry& cell )
00013 {
00014 s << ", Center: " << cell.getPosition() << std::endl;
00015
00016 if( cell.emptyCorners() )
00017 {
00018 s << "Corners vector is empty." << std::endl ;
00019 }
00020 else
00021 {
00022 const CaloCellGeometry::CornersVec& corners ( cell.getCorners() ) ;
00023 for ( unsigned int i ( 0 ) ; i != corners.size() ; ++i )
00024 {
00025 s << "Corner: " << corners[ i ] << std::endl;
00026 }
00027 }
00028 return s ;
00029 }
00030
00031 HepGeom::Transform3D
00032 CaloCellGeometry::getTransform( std::vector<HepGeom::Point3D<double> >* lptr ) const
00033 {
00034 const GlobalPoint& p ( CaloCellGeometry::getPosition() ) ;
00035 const HepGeom::Point3D<double> gFront ( p.x(), p.y(), p.z() ) ;
00036
00037 HepGeom::Point3D<double> lFront ;
00038 assert( 0 != param() ) ;
00039 std::vector<HepGeom::Point3D<double> > lc ( vocalCorners( param(), lFront ) ) ;
00040
00041 HepGeom::Point3D<double> lBack ( 0.25*(lc[4]+lc[5]+lc[6]+lc[7]) ) ;
00042
00043 const HepGeom::Point3D<double> lOne ( lc[0] ) ;
00044
00045 const CornersVec& cor ( getCorners() ) ;
00046 std::vector<HepGeom::Point3D<double> > kor ( 8, HepGeom::Point3D<double> (0,0,0) ) ;
00047 for( unsigned int i ( 0 ) ; i != 8 ; ++i )
00048 {
00049 kor[i] = HepGeom::Point3D<double> ( cor[i].x(), cor[i].y(), cor[i].z() ) ;
00050 }
00051
00052 HepGeom::Point3D<double> gBack ( 0.25*( kor[4]+kor[5]+kor[6]+kor[7] ) ) ;
00053
00054 const HepGeom::Vector3D<double> gAxis ( (gBack-gFront).unit() ) ;
00055
00056 gBack = ( gFront + (lBack-lFront).mag()*gAxis ) ;
00057 const HepGeom::Point3D<double> gOneT ( gFront + ( lOne - lFront ).mag()*( kor[0] - gFront ).unit() ) ;
00058
00059 const double langle ( ( lBack - lFront).angle( lOne - lFront ) ) ;
00060 const double gangle ( ( gBack - gFront).angle( gOneT- gFront ) ) ;
00061 const double dangle ( langle - gangle ) ;
00062
00063 const HepGeom::Plane3D<double> gPl ( gFront, gOneT, gBack ) ;
00064 const HepGeom::Point3D<double> p2 ( gFront + gPl.normal().unit() ) ;
00065
00066 const HepGeom::Point3D<double> gOne ( gFront + HepGeom::Rotate3D( -dangle, gFront, p2 )*
00067 HepGeom::Vector3D<double> ( gOneT - gFront ) ) ;
00068
00069 const HepGeom::Transform3D tr ( lFront , lBack , lOne ,
00070 gFront , gBack , gOne ) ;
00071
00072 if( 0 != lptr ) (*lptr) = lc ;
00073
00074 return tr ;
00075 }
00076
00077 const double*
00078 CaloCellGeometry::checkParmPtr(
00079 const std::vector<double>& vv ,
00080 CaloCellGeometry::ParVecVec& pvv )
00081 {
00082 const double* pP ( 0 ) ;
00083
00084 for( unsigned int ii ( 0 ) ; ii != pvv.size() ; ++ii )
00085 {
00086 const ParVec& v ( pvv[ii] ) ;
00087 assert( v.size() == vv.size() ) ;
00088
00089 bool same ( true ) ;
00090 for( unsigned int j ( 0 ) ; j != vv.size() ; ++j )
00091 {
00092 same = same && ( fabs( vv[j] - v[j] )<1.e-6 ) ;
00093 if( !same ) break ;
00094 }
00095 if( same )
00096 {
00097 pP = &(*v.begin()) ;
00098 break ;
00099 }
00100 }
00101 return pP ;
00102 }
00103
00104 const double*
00105 CaloCellGeometry::getParmPtr(
00106 const std::vector<double>& vv ,
00107 CaloCellGeometry::ParMgr* mgr ,
00108 CaloCellGeometry::ParVecVec& pvv )
00109 {
00110 const double* pP ( checkParmPtr( vv, pvv ) ) ;
00111
00112 if( 0 == pP )
00113 {
00114 pvv.push_back( ParVec( mgr ) ) ;
00115 ParVec& back ( pvv.back() ) ;
00116 for( unsigned int i ( 0 ) ; i != vv.size() ; ++i )
00117 {
00118 back[i] = vv[i] ;
00119 }
00120 pP = &(*back.begin()) ;
00121 }
00122 return pP ;
00123 }
00124
00125 bool
00126 CaloCellGeometry::inside( const GlobalPoint& point ) const
00127 {
00128 bool ans ( false ) ;
00129 const HepGeom::Point3D<double> p ( point.x(), point.y(), point.z() ) ;
00130 const CornersVec& cog ( getCorners() ) ;
00131 HepGeom::Point3D<double> co[8] ;
00132 for( unsigned int i ( 0 ) ; i != 8 ; ++i )
00133 {
00134 co[i] = HepGeom::Point3D<double> ( cog[i].x(), cog[i].y(), cog[i].z() ) ;
00135 }
00136
00137 const HepGeom::Plane3D<double> AA ( co[0], co[1], co[2] ) ;
00138 const HepGeom::Plane3D<double> BB ( co[6], co[5], co[4] ) ;
00139
00140 if( AA.distance(p)*BB.distance(p) >= 0 )
00141 {
00142 const HepGeom::Plane3D<double> CC ( co[0], co[4], co[5] ) ;
00143 const HepGeom::Plane3D<double> DD ( co[2], co[6], co[7] ) ;
00144 if( CC.distance(p)*DD.distance(p) >= 0 )
00145 {
00146 const HepGeom::Plane3D<double> EE ( co[3], co[7], co[4] ) ;
00147 const HepGeom::Plane3D<double> FF ( co[1], co[5], co[6] ) ;
00148 if( EE.distance(p)*FF.distance(p) >= 0 )
00149 {
00150 ans = true ;
00151 }
00152 }
00153 }
00154 return ans ;
00155 }