In This Package:

SolidMath.h

Go to the documentation of this file.

// $Id: SolidMath.h,v 1.8 2007/09/20 15:13:21 wouter Exp $
// ============================================================================
// CVS tag $Name: v18r12 $
// ============================================================================
// $Log: SolidMath.h,v $
// Revision 1.8  2007/09/20 15:13:21  wouter
// changed LineIntersThePhi such that it only intersects with half the
// phi plane. important for optimization in Tubs and Cons
//
// Revision 1.7  2005/12/08 19:20:01  jpalac
//
// ! 2005-12-08 - Juan Palacios
//  - Provide inInside and intersectionTicks methods to solids that accept
//    points and vectors in coordinates other than cartesian. New methods pass
//    call to templated implementations.
//
// Revision 1.6  2005/12/05 16:18:43  jpalac
//
// ! 2005-12-05 - Juan Palacios
//  - Add class Gaudi::Plane3D as stop-gap while MathCore equivalent is implemented.
//    Supports only double precision cartesian representation.
//  - Fix remaining MathCore-related compilation problems.
//
// Revision 1.5  2005/12/02 18:36:55  jpalac
//
// ! 2005-02-12 - Juan Palacios
//  - Add Kernel/LHCbDefinitions to requirements
//  - Change CLHEP geometry classes to LHCbDefinitions typedefs:
//             * These typedefs point to MathCore classes with the
//               exception of Gaudi::Plane3D, which pointe to HepPoint3D.
//               Awaiting implementation of MathCore Plane3D class.
//  - Make changes to all code to adapt to MathCore syntax
//  - TO-DO: Not compiling due to Plane3D operaitons with MathCore classes
//
// Revision 1.4  2003/11/24 14:43:44  cattanem
// fix doxygen comments
//
// Revision 1.3  2002/10/30 15:02:44  ibelyaev
//  fix a bug in DetDesc/SolidMath.h
//
// Revision 1.2  2002/05/13 18:29:53  ibelyaev
//  see $DETDESCROOT/doc/release.notes 13 May 2002
//
// Revision 1.1  2002/05/11 18:25:46  ibelyaev
//  see $DETDESCROOT/doc/release.notes 11 May 2002
// 
// ============================================================================
#ifndef DETDESC_SOLIDMATH_H 
#define DETDESC_SOLIDMATH_H 1
// Include files
#include "DetDesc/ISolid.h" 
#include "DetDesc/SolidTicks.h"

namespace SolidTicks
{

  template < class OUTPUTTYPE > 
  inline unsigned int 
  SolveQuadraticEquation
  ( const double a   , 
    const double b   , 
    const double c   ,
    OUTPUTTYPE   out ) 
  {
    if( 0 == a ) // it is indeed  a linear equation:  b*x + c = 0 
      {
        // no solution!
        if( b == 0 ) { return 0 ; }   // RETURN !!! 
        // 1 solution!
        *out++ = -1.0 * c / b ; 
        *out++ = -1.0 * c / b ;       // double the solutions 
        return 1;                     // RETURN !!!   
      }
    double d = b * b - 4.0 * a * c ; 
    // no solutions 
    if(  d < 0  )   { return 0; }     // RETURN !!!
    // 1 or 2 solution
    d = sqrt( d )                  ;   
    *out++ = 0.5 * ( -b - d ) / a  ; 
    *out++ = 0.5 * ( -b + d ) / a  ; 
    // return number of solutions;
    return 0 == d ? 1 : 2 ;           // RETURN !!! 
  };

  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheSphere
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        radius ,
    OUTPUTTYPE          out    )
  {
    // sphere with non-positive radius is not able to intersect the line! 
    if( radius <= 0 ) { return 0 ; } 
    // line with null direction vector does not  intersect the sphere! 
    double v2 = vect.mag2(); 
    if( v2 <= 0     ) { return 0 ; }
    double p2 = point.mag2()    ; 
    double pv = point.Dot(vect) ; 
    const double a = v2                   ; 
    const double b = 2.0*pv               ; 
    const double c = p2 - radius * radius ;
    // solve the equation! 
    return SolidTicks::SolveQuadraticEquation( a , b , c , out );
  };
  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheSphere2
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        r2     ,
    OUTPUTTYPE          out    )
  {
    // sphere with non-positive radius is not able to intersect the line! 
    if( r2 <= 0     ) { return 0 ; } 
    // line with null direction vector does not  intersect the sphere! 
    double v2 = vect.mag2(); 
    if( v2 <= 0     ) { return 0 ; }
    double p2 = point.mag2()    ; 
    double pv = point.Dot(vect) ; 
    const double a = v2      ; 
    const double b = 2.0*pv  ; 
    const double c = p2 - r2 ;
    // solve the equation! 
    return SolidTicks::SolveQuadraticEquation( a , b , c , out );
  };
  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheCylinder
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        radius ,
    OUTPUTTYPE          out    )
    
  {
    // Cylinder with non-positive radius is not able to intersect the line! 
    if( radius <= 0 ) { return 0 ; }     
    // line with null direction vector is not able to intersect the cylinder! 
    const double v2 = vect.x()*vect.x() + vect.y()*vect.y() ; 
    if( v2 <= 0     ) { return 0 ; }
    //
    const double p2 = point.x() * point.x() + point.y() * point.y() ; 
    const double pv = point.x() * vect.x()  + point.y() * vect.y()  ; 
    const double a = v2                   ; 
    const double b = 2.0*pv               ; 
    const double c = p2 - radius * radius ;
    // solve the equation! 
    return SolidTicks::SolveQuadraticEquation( a , b , c , out );  
  };
  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheX
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        X      ,
    OUTPUTTYPE          out   )
    
  {
    if( 0 == vect.x() ) { return 0; } 
    *out++ = ( X - point.x() ) / vect.x() ; 
    return 1;      
  };

  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheY
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        Y      ,
    OUTPUTTYPE          out    )
    
  {    
    if( 0 == vect.y() ) { return 0; } 
    *out++ = ( Y - point.y() ) / vect.y() ; 
    return 1;      
  };

  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheZ
  ( const aPoint & point  , 
    const aVector& vect   , 
    const double       Z      ,
    OUTPUTTYPE         out    )
    
  {    
    if( 0 == vect.z() ) { return 0; } 
    *out++ = ( Z - point.z() ) / vect.z() ; 
    return 1;      
  };
  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsThePhi
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        Phi    ,
    OUTPUTTYPE          out    )
    
  {
    const double sinphi = sin( Phi ) ; 
    const double cosphi = cos( Phi ) ; 
    const double d      = vect.x() * sinphi - vect.y() * cosphi ; 
    if( 0 == d ) { return 0; } 
    // only accept half the phi plane !
    const double e      = vect.y() * point.x() - vect.x() * point.y() ;
    if( e * d > 0 ) { return 0 ; }
    *out++ = ( point.y() * cosphi - point.x() * sinphi ) / d ; 
    return 1; 
  };
  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheTheta
  ( const aPoint & point  , 
    const aVector & vect   , 
    const double       Theta  ,
    OUTPUTTYPE         out    )
    
  { 
    const double sinthe = sin( Theta )    ; 
    const double costhe = cos( Theta )    ;
    
    const double c2     = costhe * costhe ; 
    const double s2     = sinthe * sinthe ; 
    
    const double a = 
      c2 *  vect.x() *  vect.x() + 
      c2 *  vect.y() *  vect.y() - 
      s2 *  vect.z() *  vect.z() ; 
    double       b = 
      c2 *  vect.x() * point.x() + 
      c2 *  vect.y() * point.y() - 
      s2 *  vect.z() * point.z() ;
    const double c = 
      c2 * point.x() * point.x() + 
      c2 * point.y() * point.y() - 
      s2 * point.z() * point.z() ;
    
    b *= 2.0;

    // return the result
    return SolidTicks::SolveQuadraticEquation( a , b, c, out );
  };

  
  template < class OUTPUTTYPE, class aPoint, class aVector > 
  inline unsigned int  
  LineIntersectsTheCone
  ( const aPoint  & point  , 
    const aVector & vect   , 
    const double        r1     ,
    const double        r2     ,
    const double        z1     ,
    const double        z2     ,
    OUTPUTTYPE          out    )
    
  {
    const double drdz = (r2-r1)/(z2-z1)               ;
    const double p1   = r1 - z1*drdz + drdz*point.z() ;
    const double p2   = drdz * vect.z()               ;
    
    double a = vect.x () * vect.x () + vect.y () * vect.y () ;
    a -= p2*p2 ; 
    double b = vect.x () * point.x() + vect.y () * point.y() ;
    b -= p2*p1 ;
    b *= 2.0   ; 
    double c = point.x() * point.x() + point.y() * point.y() ; 
    c -= p1*p1 ;
    
    // return result
    return SolidTicks::SolveQuadraticEquation( a , b, c, out );
  };

}; // end of namespace 

// ============================================================================
// The END
// ============================================================================
#endif // DETDESC_SOLIDMATH_H
// ============================================================================

---