vsr_cga3D_op.h

Go to the documentation of this file.

#ifndef CGA3D_OPERATIONS_H_INCLUDED
#define CGA3D_OPERATIONS_H_INCLUDED  

#include <math.h> 
#include <vector>  

#include "vsr_cga3D_types.h"  
#include "vsr_cga3D_round.h"

#include "detail/vsr_generic_op.h"
#include "util/vsr_constants.h" 
#include "util/vsr_math.h"

namespace vsr{ namespace cga{               

    struct Op{
              
      static Rot AA( const Vec& s);   
      static Rot AA( const Dlp& s);   
      static Rot AA( const Cir& s);   
      static Rot AA( const Biv& s);   
      static Vec Pos( const Dlp& s);  
      static Pnt Pos( const Cir& s);  

      template<class A>  
      static auto dual( const A& a ) RETURNS (
        a.dual()
      ) 
    
      template<class A>  
      static auto undual( const A& a ) RETURNS (
        a.undual()
      )
    
      template<class A>  
      static auto duale( const A& a ) RETURNS (
        a.duale()
      )
       
      template<class A>  
      static auto unduale( const A& a ) RETURNS (
        a.unduale()
      )
      
      template<class T> static auto dl( const T& t ) RETURNS ( dual(t) )
      template<class T> static auto udl( const T& t) RETURNS ( udual(t) )
      template<class T> static auto dle( const T& t) RETURNS (duale(t))
      template<class T> static auto udle( const T& t) RETURNS (unduale(t))
     
      
      template<class A>
      static constexpr bool sign(const A& a, const A& b) {
        return (a / b)[0] > 0 ? 1 : 0;
      } 
    
      template<class A> 
      static constexpr bool sn(const A& a, const A& b){
        return sign(a,b);
      }
    
      template< class A, class B>
      static constexpr auto project(const A& a, const B& b) RETURNS ( 
        (a <= b ) / b
      ) 
    
      
      template< class A, class B>
      static constexpr auto reject(const A& a, const B& b) RETURNS ( 
        (a ^ b ) / b
      )
    
      template< class A, class B>
      static constexpr auto pj(const A& a, const B& b) RETURNS ( project(a,b) )
      template< class A, class B>
      static constexpr auto rj(const A& a, const B& b) RETURNS ( reject(a,b) )

      
    };
    

   struct Gen{  

      static Rot rot( const Biv& b);
      static Rot rotor( const Biv& b);
      static Bst bst(const Pair& p);
      static Bst boost(const Pair& p);
      static Tsd dil(const Pnt& p, VSR_PRECISION t);
      static Tsd dilator(const Pnt& p, VSR_PRECISION t);

      template<class A>
      static Trs trs(const A& a) { return nga::Gen::trs(a); }
      template<class A>
      static Trs translator(const A& a) { return nga::Gen::trs(a); }
      template<class A>
      static Trv trv(const A& a){ return nga::Gen::trv(a); }
      template<class A>
      static Trv transversor(const A& a){ return nga::Gen::trv(a); }


      static Rot ratio(const Vec& v, const Vec& v2);
      static Biv log(const Rot& r);


      static Rot rot(double theta, double phi);                

      static Rot rot( double yaw, double pitch, double roll);

   /*-----------------------------------------------------------------------------
    *  TWISTS (Motors, Plücker, etc)
    *-----------------------------------------------------------------------------*/
    static Mot mot( const Dll& dll);
    
    static Mot motor( const Dll& dll);

    static Dll log( const Mot& m);

    static Dll log(const Dll& a, const Dll& b, VSR_PRECISION t = 1.0);

    static Mot ratio( const vsr::cga::Dll& a, const vsr::cga::Dll& b, VSR_PRECISION t = 1.0);

    // Due to overloading, it is also possible to use Motors as Arguments

    static Mot ratio( const vsr::cga::Mot& a, const vsr::cga::Mot& b, VSR_PRECISION t);

    /*-----------------------------------------------------------------------------
     *  BOOSTS (Transversions, Conformal Rotors)
     *-----------------------------------------------------------------------------*/
 
      template <class A, class T>
      static Bst bst(const A& tnv, const Vec& vec, T t){
          Par s = Par( tnv.template copy<Tnv>() ).sp( nga::Gen::trs(vec) ); 
          return Gen::bst(s * t);
      }
 
      

   static Bst ratio( const DualSphere& a, const DualSphere& b, bool bFlip = true);
   
   static Pair atanh2( const Pair& p, VSR_PRECISION cs, bool bCW);
       
   static Pair log( const Bst& b, bool bCW = false);

   static Con ratio( const Circle& a, const Circle& b, bool bFlip = false, float theta = 0);
   
   static Con ratio( const Pair& a, const Pair& b,bool bFlip = false, float theta = 0);//{ return ratio( a.dual(), b.dual() ); }
   
   static vector<Pair> split(const Pair& par);   

   static vector<Pair> split(const Con& con);   

   
   static vector<Pair> log( const Con& rot);

   static vector<Pair> log( const Circle& ca, const Circle& cb, bool bFlip = false, VSR_PRECISION theta = 0);
   
   static vector<Pair> log( const Pair& ca, const Pair& cb, bool bFlip = false, VSR_PRECISION theta = 0);
   
   
   static Con con( const vector<Pair>& log, VSR_PRECISION amt);

   static Con con( const vector<Pair>& log, VSR_PRECISION amtA, VSR_PRECISION amtB);
 

  /* General Conformal Transformation from two circles */
   static Con con(  const Circle& ca, const Circle& cb, VSR_PRECISION amt);   

  /* General Conformal Transformation from two circles and two weights */
   static Con con(  const Circle& ca, const Circle& cb, VSR_PRECISION amtA, VSR_PRECISION amtB);   


     static Rotor xf(const Biv& b);
     static Motor xf(const DualLine& dll);
     static Dilator xf(const FlatPoint& flp);

     static Bst xf(const Pair& p);

};


    /*-----------------------------------------------------------------------------
     *  ROTORS
     *-----------------------------------------------------------------------------*/


    struct Construct {
   
      /*-----------------------------------------------------------------------------
       *  PAIRS
       *-----------------------------------------------------------------------------*/

      static Pair pair(const DualSphere& s, const Vec& v);     
      
      static Pair pair(VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z, Vec vec = Vec::y, VSR_PRECISION r=1.0);
     

      /*-----------------------------------------------------------------------------
       *  POINTS
       *-----------------------------------------------------------------------------*/
     static Point pointA(const Pair& pp);

     static Point pointB(const Pair& pp);
       
     static Point point(const Circle& c, VSR_PRECISION t);
     static Point point(const DualSphere& s, const Vec& v);
     static Point point(VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z);
     static Point point(const Vec& v);
     static Point point(const Line&, const Point& );
     static Point point(const DualLine&, const Point& );
     
      
      /*-----------------------------------------------------------------------------
       *  CIRCLES
       *-----------------------------------------------------------------------------*/
     
      static Circle circle(const Point& a, const Point& b, const Point& c);
 
      static Circle circle(const Point& p, VSR_PRECISION r, const Biv& biv = Biv::xy);            
      static Circle circle(const Biv& B);
     
      //circle Facing v
      static Circle circle(const Vec& v, VSR_PRECISION r=1.0);
      
      //Circle at x,y,z facing in biv
      static Circle circle(VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z, Biv biv = Biv::xy, VSR_PRECISION r=1.0);
 


      /*-----------------------------------------------------------------------------
       *  SPHERES
       *-----------------------------------------------------------------------------*/
      static Sphere sphere(const Pnt& a, const Pnt& b, const Pnt& c, const Pnt& d);
      static DualSphere sphere(VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z, VSR_PRECISION r=1.0);
      static DualSphere sphere(const Point& p, VSR_PRECISION r=1.0);
                             


      /*-----------------------------------------------------------------------------
       *  PLANES
       *-----------------------------------------------------------------------------*/

      static DualPlane plane( VSR_PRECISION a, VSR_PRECISION b, VSR_PRECISION c, VSR_PRECISION d=0.0);
      static DualPlane plane( const Vec& v, VSR_PRECISION d=0.0);
      static Plane plane(const Pnt& a, const Pnt& b, const Pnt& c);
 
      /*-----------------------------------------------------------------------------
       *  LINES
       *-----------------------------------------------------------------------------*/
 
      static DualLine axis(const Cir& c);
   
      static Line line( const Vec& a, const Vec& b);
     
      static Line line( VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z );
      
      static Line dualLine( VSR_PRECISION x, VSR_PRECISION y, VSR_PRECISION z );
      
      static Line line( const Point& a, const Point& b);
      static Line line( const Point& a, const Vec& b);
     
     
     
      static Circle hline( const Point& a, const Point& b);
      static Circle sline( const Point& a, const Point& b);


      
       static VSR_PRECISION distance( const Lin& lin, const Pnt& pnt);
    

      #pragma mark COINCIDENCE_FUNCTIONS

       static Circle meet( const Dls& s, const Dls& d);
       static Circle meet( const Dls& s, const Dlp& d);
       static Circle meet( const Dls& s, const Pln& d);
       static Circle meet( const Sphere& s, const DualPlane& d);
       static Circle meet( const Sphere& s, const Plane& d);


       static Point meet( const Line& lin, const DualPlane&);
       static Point meet( const DualLine&, const DualPlane&);    
       static Point meet( const Line& la, const Line& lb); 

       static Pair meet( const Cir& cir, const Dlp& dlp);  
       static Pair meet( const Cir& cir, const Dls& s);

       #pragma mark HIT_TESTS
     
       static bool hit(const Point&, const Pair&);
  
      static bool hit(const Point&, const Circle&);

      static double squaredDistance(const Point& a, const Point& b);


      #pragma mark HYPERBOLIC_FUNCTIONS
      /*-----------------------------------------------------------------------------
       *  hyperbolic functions (see alan cortzen's document on this)
       *-----------------------------------------------------------------------------*/

      static Point hnorm(const Pnt& p);

      
      static double hdist(const Pnt& pa, const Pnt& pb);


      static Pair hgen(const Pnt& pa, const Pnt& pb, double amt);

      static Point hspin(const Pnt& pa, const Pnt& pb, double amt);
     
    
      #pragma mark AFFINE combinations
      
      template<class A, class B>
      static Point affine( const A& a, const B& b, VSR_PRECISION t){
        return (a + (b-a)*t).null();
      }
};
     
     /*-----------------------------------------------------------------------------
      *  EVALUATION LAMBDAS
      *-----------------------------------------------------------------------------*/
 
    auto pointOnLine = []( const Line& lin, const Point& v){
      return Round::null( Flat::loc( lin, v, false) );
    };    

    auto pointOnCircle = [](const Circle& c, VSR_PRECISION t){
     return Construct::point(c,t);
    };
    auto pointsOnCircle = [](const Circle& c, int num){
     vector<Point> out;
     for (int i=0;i<=num;++i){
       out.push_back(pointOnCircle(c, TWOPI*(float)i/num) );
     }
     return out;
    };
    auto pairOnSphere = [](const DualSphere& s, VSR_PRECISION t, VSR_PRECISION p){
     return Construct::pair(s, Vec::x.sp( Gen::rot(t,p) ) );
    };
    auto pointOnSphere = [](const DualSphere& s, VSR_PRECISION t, VSR_PRECISION p){
      return Construct::pointA( pairOnSphere(s,t,p) ).null(); 
    };
    auto pointsOnSphere = [](const DualSphere& s, int u, int v){
      vector<Point> out;
      for (int i = 0; i < u; ++i){
       for (int j = 0; j < v; ++j){

        float tu= TWOPI * i/u;//-1 + 2.0 * i/num;
        float tv = -PIOVERTWO  + PI * j/v;

        out.push_back( pointOnSphere(s,tu,tv) );

      }
     }
     return out; 
    };


} //cga:: 


template<class Algebra, class B> template <class A>
Multivector<Algebra,B> Multivector<Algebra,B>::mot( const Multivector<Algebra,A>& t) const{
     return this -> sp ( cga::Gen::mot(t) );  
} 
template<class Algebra, class B> template <class A>
Multivector<Algebra,B> Multivector<Algebra,B>::motor( const Multivector<Algebra,A>& t) const{
     return this -> sp ( cga::Gen::mot(t) );  
} 
template<class Algebra, class B> template <class A>
Multivector<Algebra,B> Multivector<Algebra,B>::twist( const Multivector<Algebra,A>& t) const{
     return this -> sp ( cga::Gen::mot(t) );  
} 
    

#define E1 e1(1)
#define E2 e2(1)
#define E3 e3(1)

#define PT(x,y,z) vsr::cga::Round::null(vsr::cga::Vec(x,y,z))
#define DLS(r) vsr::cga::Round::dls(0,0,0,r)

#define PV(v) vsr::cga::Round::null(v)
#define PX(f) vsr::cga::Round::null(vsr::cga::Vec(f,0,0))
#define PY(f) vsr::cga::Round::null(vsr::cga::Vec(0,f,0))
#define PZ(f) vsr::cga::Round::null(vsr::cga::Vec(0,0,f))

#define PAIR(x,y,z) (PT(x,y,z)^PT(-x,-y,-z))
#define CXY(f) (PX(f)^PY(f)^PX(-f)).unit()
#define CXZ(f) (PX(f)^PZ(f)^PX(-f)).unit()
#define CYZ(f) (PY(f)^PY(-f)^PZ(f)).unit()
#define F2S(f) f*1000.0
#define S2F(f) f/1000.0

#define LN(x,y,z) ( vsr::cga::Point(0,0,0,1,.5)^PT(x,y,z)^vsr::cga::Inf(1) )
#define DLN(x,y,z) ( vsr::Op::dl(LN(x,y,z)) )
#define PAO vsr::cga::Point(0,0,0,1,0)   
#define EP vsr::cga::Dls(0,0,0,1,-.5)    
#define EM vsr::cga::Dls(0,0,0,1,.5)     
#define INFTY vsr::cga::Inf(1)           
#define HYPERBOLIC_INF EP
#define SPHERICAL_INF EM
#define EUCLIDEAN_INF INFTY
#define HLN(x,y,z) (vsr::cga::Ori(1)^PT(x,y,z)^EP) //hyperbolic line (circle)
#define HDLN(x,y,z) (vsr::Op::dl(HLN(x,y,z)))

}    //vsr::                                       

#endif 


//  template<bool _a, bool _b>
//  struct meet_impl{
//    template<class A, class B>
//    auto operator()(const A& a, const B& b) RETURNS
//    ( (a^b).dual() )
//  };


//     namespace topo{
//
//
//
//       /*!
//        *  \brief spin product in range [0,RANGE) of some type t and some bivector generator p
//        */
//      template<class T, class G>
//      vector<T> spin( const T& s, const G& p, int num, float range=PI){
//        vector<T> res;
//        for (int i=0;i<num;++i){
//          float t=range*(float)i/num;
//          res.push_back( s.spin( gen(p*t) ) ) ;
//        }
//        return res;
//      }
//
//      
//
//    } //topo::
  

//} //cga::