CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/Fireworks/ParticleFlow/src/FWPFMaths.cc

Go to the documentation of this file.
00001 #include <math.h>
00002 #include "TMath.h"
00003 #include "TEveVector.h"
00004 
00005 namespace FWPFMaths
00006 {
00007 //______________________________________________________________________________
00008 float
00009 sgn( float val )
00010 {
00011    return ( val < 0 ) ? -1 : 1;
00012 }
00013 
00014 //______________________________________________________________________________
00015 TEveVector
00016 cross( const TEveVector &v1, const TEveVector &v2 )
00017 {
00018    TEveVector vec;
00019 
00020    vec.fX = ( v1.fY * v2.fZ ) - ( v1.fZ * v2.fY );
00021    vec.fY = ( v1.fZ * v2.fX ) - ( v1.fX * v2.fZ );
00022    vec.fZ = ( v1.fX * v2.fY ) - ( v1.fY * v2.fX );
00023 
00024    return vec;
00025 }
00026 
00027 //______________________________________________________________________________
00028 float
00029 dot( const TEveVector &v1, const TEveVector &v2 )
00030 {
00031    float result = ( v1.fX * v2.fX ) + ( v1.fY * v2.fY ) + ( v1.fZ * v1.fZ );
00032 
00033    return result;
00034 }
00035 
00036 //______________________________________________________________________________
00037 TEveVector
00038 lineCircleIntersect( const TEveVector &v1, const TEveVector &v2, float r )
00039 {
00040    // Definitions
00041    float x, y;
00042    float dx = v1.fX - v2.fX;
00043    float dy = v1.fY - v2.fY;
00044    float dr = sqrt( ( dx * dx ) + ( dy * dy ) );
00045    float D = ( ( v2.fX * v1.fY ) - ( v1.fX * v2.fY ) );
00046 
00047    float rtDescrim = sqrt( ( ( r * r ) * ( dr * dr ) ) - ( D * D ) );
00048 
00049    if( dy < 0 )   // Going down
00050    {
00051       x = ( D * dy ) - ( ( sgn(dy) * dx ) * rtDescrim );
00052       x /= ( dr * dr );
00053 
00054       y = ( -D * dx ) - ( fabs( dy ) * rtDescrim );
00055       y /= ( dr * dr );
00056    }
00057    else
00058    {
00059 
00060       x = ( D * dy ) + ( ( sgn(dy) * dx ) * rtDescrim );
00061       x /= ( dr * dr );
00062 
00063       y = ( -D * dx ) + ( fabs( dy ) * rtDescrim );
00064       y /= ( dr * dr );
00065    }
00066 
00067    TEveVector result = TEveVector( x, y, 0.001 );
00068    return result;
00069 }
00070 
00071 //______________________________________________________________________________
00072 TEveVector
00073 lineLineIntersect( const TEveVector &p1, const TEveVector &p2, const TEveVector &p3, const TEveVector &p4 )
00074 {
00075    TEveVector a = p2 - p1;
00076    TEveVector b = p4 - p3;
00077    TEveVector c = p3 - p1;
00078    TEveVector result;
00079    float s, val;
00080 
00081    s = dot( cross( c, b ), cross( a, b ) );  
00082    val = dot( cross( a, b ), cross( a, b ) );
00083    s /= val;
00084    
00085    result = p1 + ( a * s );
00086    return result; 
00087 }
00088 
00089 //______________________________________________________________________________
00090 float
00091 linearInterpolation( const TEveVector &p1, const TEveVector &p2, float z )
00092 {
00093    float y;
00094 
00095    y = ( ( z - fabs( p1.fZ ) ) * p2.fY ) + ( ( fabs( p2.fZ ) - z ) * p1.fY );
00096    y /= ( fabs( p2.fZ) - fabs( p1.fZ ) );
00097 
00098    return y;
00099 }
00100 
00101 //______________________________________________________________________________
00102 bool
00103 checkIntersect( const TEveVector &p, float r )
00104 {
00105    float h = sqrt( ( p.fX * p.fX ) + ( p.fY * p.fY ) );
00106    
00107    if( h >= r )
00108       return true;
00109 
00110    return false;
00111 }
00112 
00113 //______________________________________________________________________________
00114 float
00115 calculateEt( const TEveVector &centre, float e )
00116 {
00117    TEveVector vec = centre;
00118    float et;
00119 
00120    vec.Normalize();
00121    vec *= e;
00122    et = vec.Perp();
00123 
00124    return et;
00125 }
00126 }