几何中点（三维向量）、线、面的类

AssPoint.h // 版权所有(C) 梁意, 2004 // 最后修改: 2004.6.杭州 #ifndef AssPoint_HEADER #define AssPoint_HEADER #include<cmath> #include<limits> #include<vector> #include<deque> #include<iostream> //#define Container vector #define Container deque using namespace std; #define PI 3.1415926535897932384626 //M_PI #define ASS_DELTA 1e-6 //or FLT_EPSILON (1.192092896e-07F) or DBL_EPSILON in float.h template class AssPoint { public: typedef AssPoint _Myt; AssPoint()//构造函数，原点 { X=0; Y=0; Z=0; } AssPoint(T _X,T _Y,T _Z=0)//构造函数，X,Y,Z { X=_X; Y=_Y; Z=_Z; } AssPoint(const _Myt &P)//构造函数，X,Y,Z { X=P.X; Y=P.Y; Z=P.Z; } bool operator == ( _Myt & P) //比较两点是否相等,坐标误差DELTA { return(fabs(X-P.X) } _Myt operator -( _Myt & P)//两点相减 { return _Myt(X-P.X,Y-P.Y,Z-P.Z); } _Myt operator +( _Myt & P)//两点相加 { return _Myt(X+P.X,Y+P.Y,Z+P.Z); } _Myt operator -( T d)//三方向同时减去d { return _Myt(X-d,Y-d,Z-d); } _Myt operator +( T d)//三方向同时加上d { return _Myt(X+d,Y+d,Z+d); } _Myt operator *( T d)//三方向同时乘上d { return _Myt(X*d,Y*d,Z*d); } _Myt operator /( T d)//三方向同时除以d { return _Myt(X/d,Y/d,Z/d); } T operator *( _Myt & P)//点积 { return X*P.X+Y*P.Y+Z*P.Z; } T Dis( _Myt &P)//到另一点的距离 { return sqrt(pow(X-P.X,2)+pow(Y-P.Y,2)+pow(Z-P.Z,2)); } _Myt _X( _Myt &P)//叉积 { return _Myt(Y*P.Z-Z*P.Y,Z*P.X-X*P.Z,X*P.Y-Y*P.X); } T Modul()//求模 { return sqrt( (*this) * (* this) ); } bool HasTheSameDirection( _Myt & _P)//向量是否方向相同 { return (_P/_P.Modul())==((*this)/Modul()); } bool IsNearPt( _Myt & _P,T _dis)//是否在_P附近 { return Dis( _P )<=fabs(_dis); } _Myt Scale(T Sx,T Sy,T Sz, _Myt &BasePoint=_Myt() )//以BasePoint为基点,三方向分别缩放Sx,Sy,Sz { _Myt tmp=(*this)-BasePoint; return _Myt(tmp.X*Sx,tmp.Y*Sy,tmp.Z*Sz)+BasePoint; } _Myt Transfrom(T tx,T ty,T tz )//移动tx,ty,tz { return _Myt(X+tx,Y+ty,Z+tz); } T AngleYZ()//与YZ面夹角 { return asin(X/Modul()); } T AngleZX()//与ZX面夹角 { return asin(Y/Modul()); } T AngleXY()//与XY面夹角 { return asin(Z/Modul()); } T AngleX_XY()//XY面的投影与X轴正向的夹角 { if(Y>=0) return acos(X/sqrt(pow(X,2)+pow(Y,2))); return PI*2.0-acos(X/sqrt(pow(X,2)+pow(Y,2))); } T AngleY_XY()//XY面的投影与Y轴正向的夹角 { if(X<=0) return acos(Y/sqrt(pow(X,2)+pow(Y,2))); return PI*2.0-acos(Y/sqrt(pow(X,2)+pow(Y,2))); } T AngleY_YZ()//YZ面的投影与Y轴正向的夹角 { if(Z>=0) return acos(Y/sqrt(pow(Y,2)+pow(Z,2))); return PI*2.0-acos(Y/sqrt(pow(Y,2)+pow(Z,2))); } T AngleZ_YZ()//YZ面的投影与Z轴正向的夹角 { if(Y<=0) return acos(Z/sqrt(pow(Y,2)+pow(Z,2))); return PI*2.0-acos(Z/sqrt(pow(Y,2)+pow(Z,2))); } T AngleZ_ZX()//ZX面的投影与Z轴正向的夹角 { if(X>=0) return acos(Z/sqrt(pow(Z,2)+pow(X,2))); return PI*2.0-acos(Z/sqrt(pow(Z,2)+pow(X,2))); } T AngleX_ZX()//ZX面的投影与X轴正向的夹角 { if(Z<=0) return acos(X/sqrt(pow(Z,2)+pow(X,2))); return PI*2.0-acos(X/sqrt(pow(Z,2)+pow(X,2))); } _Myt RotateX(T ceta)//绕X轴旋转ceta { return( _Myt( X,Y*cos(ceta)-Z*sin(ceta),Z*cos(ceta)+Y*sin(ceta) ) ); } _Myt RotateY(T ceta)//绕Y轴旋转ceta { return( _Myt( X*cos(ceta)+Z*sin(ceta) ,Y ,Z*cos(ceta)-X*sin(ceta) ) ); } _Myt RotateZ(T ceta)//绕Z轴旋转ceta { return( _Myt( X*cos(ceta)-Y*sin(ceta) ,Y*cos(ceta)+X*sin(ceta),Z) ); } T Angle( _Myt &other)//与另一向量的夹角 { return acos( (*this)*other/this->Modul()/other.Modul() ); } _Myt Rotate(T ceta, _Myt &Dir)//绕一原点为起点方向为Dir的轴旋转ceta { T a=-Dir.AngleZ_YZ(); T b=-Dir.AngleYZ(); return RotateX(a).RotateY(b).RotateZ(ceta).RotateY(-b).RotateX(-a); } T X; T Y; T Z; }; #endif AssLine.h // 版权所有(C) 梁意, 2004 // 最后修改: 2004.6.杭州 #ifndef AssLine_HEADER #define AssLine_HEADER #include"AssPoint.h" template class AssLine { public: enum Input_type { Point_Slope, TwoPoints }; typedef AssPoint _MyPoint; typedef AssLine _Myt; AssLine()//两点一线 { } AssLine(const _Myt &other)//两点一线 { P1=other.P1; P2=other.P2; } AssLine(_MyPoint & _P1,_MyPoint & _P2)//两点一线 { P1=_P1; P2=_P2; } AssLine(_MyPoint &P,T A,T B,T C,Input_type type=_Myt::Point_Slope) { if(type==_Myt::Point_Slope) { P2=_MyPoint(P.X+A,P.Y+B,P.Z+C); } if(type==_Myt::TwoPoints) { P2=_MyPoint(A,B,C); } P1=P; } AssLine(T X,T Y,T Z,T X2,T Y2,T Z2)//两点一线 { P1=_MyPoint(X,Y,Z); P2=_MyPoint(X2,Y2,Z2); } _Myt & operator =( const _Myt &other) { P1=other.P1; P2=other.P2; return (*this); } T GetLength()//返回长度 { return sqrt(pow(P1.X-P2.X,2)+pow(P1.Y-P2.Y,2)+pow(P1.Z-P2.Z,2)); } T AngleX()//对X轴夹角 { return acos( (P2-P1).X/(P2-P1).Modul() ); } T AngleY()//对Y轴夹角 { return acos( (P2-P1).Y/(P2-P1).Modul() ); } T AngleZ()//对Z轴夹角 { return acos( (P2-P1).Z/(P2-P1).Modul() ); } T AngleXY()//对XY面夹角 { return (P2-P1).AngleXY(); } T AngleZX()//对ZX面夹角 { return (P2-P1).AngleZX(); } T AngleYZ()//对YZ面夹角 { return (P2-P1).AngleYZ(); } bool IsLine()//长度太小就不是线,DELTA控制 { return GetLength()>=ASS_DELTA; } bool operator == (_Myt & _L)//两条线是否重合 { return (*this || _L) && IsOnLine(_L.P1); } bool IsOnLine(_MyPoint &P)//点是否在线上 { return (P-P1).HasTheSameDirection(P2-P1); } bool IsInLineSegment(_MyPoint &P)//点是否在线段内 { return IsOnLine(P) && ( (P-P1).HasTheSameDirection(P2-P) ); } T Angle(AssLine & _L)//与另一线的夹角 { return acos( (P2-P1)*(_L.P2-_L.P1)/( (P2-P1).Modul() * (_L.P2-_L.P1).Modul() ) ); } bool IsPerpendicular(_Myt & _L)//是否与另一线垂直 { return (P2-P1)*(_L.P2-_L.P1)< ASS_DELTA; } bool operator || (_Myt & _L)//是否与另一线平行 { _MyPoint tmp1=P2-P1; _MyPoint tmp2=_L.P2-_L.P1; return ((tmp1/tmp1.Modul())._X(tmp2/tmp2.Modul())).Modul() } _Myt ParallelLine(_MyPoint & P)//过点的平行线 { return _Myt(P,P2-P1+P); } _Myt PerpendicularLine(_MyPoint &P)//过点的垂线 { return _Myt(P,(P2-P1) * ( (P-P1)*(P2-P1)/pow((P2-P1).Modul(),2) ) + P1 ); } T Dis(_MyPoint & P)//到点的距离 { return (P-P1)._X(P2-P1).Modul()/(P2-P1).Modul(); } T Dis(_Myt & _L)//到线的距离 { _MyPoint P=(P2-P1)._X(_L.P2-_L.P1); return fabs((P1-_L.P1)*P)/P.Modul(); } _MyPoint Intersection(_Myt & _L)//与线的交点 { if( Dis(_L)>ASS_DELTA ) return _MyPoint( DBL_MAX,0); _MyPoint _tmp1(P2-P1); _MyPoint _tmp2(_L.P2-_L.P1); _MyPoint _tmp4=_tmp1._X(_tmp2)._X( _tmp2 ); _MyPoint _tmp5=_L.P1-P1; T t= (_tmp4*_tmp5)/(_tmp4*_tmp1); return ( P1+( _tmp1*t) ); } _Myt PerpendicularLine( AssLine & _L)//同时与另一线垂直的线 { if( ( (*this) || _L ) ) return _Myt(0,0,0,DBL_MAX,DBL_MAX,DBL_MAX); _MyPoint _tmp1(P2-P1); _MyPoint _tmp2(_L.P2-_L.P1); _MyPoint _tmp3=_tmp1._X(_tmp2); _MyPoint _tmp4=_tmp3._X(_tmp2); _MyPoint _tmp5=_L.P1-P1; T t= (_tmp4*_tmp5)/(_tmp4*_tmp1); return _Myt(P1+( _tmp1*t),P1+( _tmp1*t)+_tmp3); } bool IsNearPt(_MyPoint & _P,T _dis)//是否在离线_dis距离内 { return Dis( _P )<=fabs(_dis) ; } void GetNDivide( Container< _MyPoint > & Points,size_t n)//n等分点 { Points.clear(); // if(n<=0) return; if(n==1) { Points.push_back(P1); Points.push_back(P2); return; } T tmp_stepX=(P2.X-P1.X)/((T)n); T tmp_stepY=(P2.Y-P1.Y)/((T)n); T tmp_stepZ=(P2.Z-P1.Z)/((T)n); for(int i=0;i<=n;i++) { Points.push_back(P1+ ( _MyPoint(tmp_stepX,tmp_stepY,tmp_stepZ)*((T)i)) ); } } _MyPoint Rotate(_MyPoint P,T ceta)//旋转 { return (P-P1).Rotate(ceta,P2-P1)+P1; } _MyPoint P1,P2; }; #endif AssPlane.h // 版权所有(C) 梁意, 2004 // 最后修改: 2004.6.杭州 #ifndef AssPlane_HEADER #define AssPlane_HEADER #include"AssLine.h" template class AssPlane { public: typedef AssPoint _MyPoint; typedef AssLine _MyLine; typedef AssPlane _Myt; AssPlane()//默认构造函数,方向(1,1,1),过原点的平面 { A=1; B=1; C=1; } AssPlane(const _Myt &other) { A=other.A; B=other.B; C=other.C; P=other.P; } AssPlane(T _A,T _B,T _C,_MyPoint & _P=_MyPoint() )//方向(_A,_B,_C),过点_P的平面 { A=_A; B=_B; C=_C; P=_P; } AssPlane(_MyPoint & vector,_MyPoint _P=_MyPoint() )//方向vector,过点_P的平面 { A=vector.X ; B=vector.Y; C=vector.Z; P=_P; } AssPlane(_MyPoint &P1,_MyPoint &P2,_MyPoint &P3)//过3点的平面 { _MyPoint tmp1=P2-P1; _MyPoint tmp2=P3-P1; _MyPoint tmp3=tmp1._X(tmp2); A=tmp3.X; B=tmp3.Y; C=tmp3.Z; P=P1; } AssPlane(_MyLine &_L,_MyPoint &_P)//过点和线的平面 { _MyPoint tmp1=_L.P1-_P; _MyPoint tmp2=_L.P2-_P; _MyPoint tmp3=tmp1._X(tmp2); A=tmp3.X; B=tmp3.Y; C=tmp3.Z; P=_P; } _MyLine PerpendicularLine(_MyPoint & _P)//过点_P与面垂直的线 { return _MyLine(_P,A,B,C,_MyLine::Point_Slope); } _Myt ParallelPlane(_MyPoint & _P)//过点与面平行的平面 { return _Myt(A,B,C,_P); } T Dis(_MyPoint & _P)//面到点距离 { _MyPoint _tmp1(A,B,C); return fabs( _tmp1*( _P-P ) )/_tmp1.Modul(); } T Dis(_Myt & _PL)//面到面的距离 { _MyPoint _tmp1(A,B,C); _MyPoint _tmp2(_PL.A,_PL.B,_PL.C); _MyPoint _tmp3=_tmp1._X(_tmp2); if(_tmp3.Modul()>ASS_DELTA) return 0; return Dis( _PL.P ); } _Myt PerpendicularPlane(_MyLine & _L)//过线与面垂直的面 { _MyPoint _tmp1=_L.P2-_L.P1; _MyPoint _tmp2=_tmp1._X(_MyPoint(A,B,C)); return _Myt(_tmp2.X,_tmp2.Y,_tmp2.Z,_L.P1); } _Myt ParallelPlane(_MyLine & _L)//过线与面平行的面 { if( ! ( (*this) || _L ) ) return AssPlane(0,0,0,Point(0,0,0)); return _Myt(A,B,C,_L.P1); } _MyPoint Intersection(_MyLine &_L)//面与线交点 { if( (*this) || _L ) return _MyPoint(DBL_MAX,DBL_MAX,DBL_MAX); _MyPoint _tmp1=_L.P2-_L.P1; _MyPoint _tmp2=P- _L.P1; _MyPoint _tmp3(A,B,C); return ( _tmp1* ( ( _tmp3*_tmp2 )/( _tmp3*_tmp1) ) ) + _L.P1; } bool operator || (_MyLine & _L)//与线平行判断 { return fabs( (_L.P2-_L.P1)*_MyPoint(A,B,C) ) } bool operator || (_Myt & _PL)//与面平行判断 { _MyPoint _tmp1(A,B,C); _MyPoint _tmp2(_PL.A,_PL.B,_PL.C); return ((_tmp1/_tmp1.Modul())._X(_tmp2/_tmp2.Modul())).Modul() } T Angle(_MyLine & _L)//面与线的夹角 { _MyPoint _tmp1=_L.P2-_L.P1; _MyPoint _tmp2(A,B,C); return acos( _tmp1*_tmp2 /( _tmp1.Modul()* _tmp2.Modul() ) ); } T Angle(_Myt & _PL)//与面的夹角 { _MyPoint _tmp1(A,B,C); _MyPoint _tmp2(_PL.A,_PL.B,_PL.C); return acos( (_tmp1 * _tmp2 )/( _tmp1.Modul() *_tmp2.Modul() ) ); } _MyLine Intersection(AssPlane & _PL)//面与面的交 { if( (*this) || _PL ) return _MyLine (0,0,0,0,0,0); _MyPoint _tmp1(A,B,C); _MyPoint _tmp2(_PL.A,_PL.B,_PL.C); _MyPoint _tmp3=_tmp1._X(_tmp2); _MyLine _tmpL1(_PL.P,_PL.P+_tmp3._X(_tmp2)); _MyLine _tmpL2(P,P+_tmp3._X(_tmp1)); if(_tmpL1.Dis(_tmpL2) return _MyLine(_tmpL1.Intersection(_tmpL2),_tmp3.X,_tmp3.Y,_tmp3.Z,_MyLine::Point_Slope); return _tmpL1.PerpendicularLine(_tmpL2); } bool IsInPlane(_MyPoint &_P)//点在平面内判断 { return fabs((_P-P)*_MyPoint(A,B,C)) } bool IsInPlane(_MyLine &_L)//线在平面内判断 { return IsInPlane(_L.P1) && IsInPlane(_L.P2); } _MyPoint GetNormal()//单位方向 { _MyPoint _tmp(A,B,C); return _tmp/_tmp.Modul(); } bool IsNearPt(_MyPoint & _P,T _dis)//是否在离面_dis距离内 { return Dis( _P )<=fabs(_dis) ; } _Myt Offset(T _h)//偏移_h得到的平面 { _MyPoint tmpP(A,B,C); T t=_h/tmpP.Modul(); return _Myt(A,B,C,P+(_MyPoint(A,B,C)*t)); } T A,B,C; _MyPoint P; }; #endif