1.向量Vector3d
using System; namespace RGeos.Geometry
{
/// <summary>
/// 3D向量类
/// </summary>
public class Vector3d
{
public double[] vector;
private const double E = 0.0000001f;
/// <summary>
///
/// </summary>
/// <param name="x"></param>
/// <param name="y"></param>
/// <param name="z"></param>
///
public Vector3d()
{
vector = new double[];
}
public Vector3d(double x, double y, double z)
{
vector = new double[] { x, y, z };
}
public Vector3d(Vector3d vct)
{
vector = new double[];
vector[] = vct.X;
vector[] = vct.Y;
vector[] = vct.Z;
}
#region 属性
/// <summary>
/// X向量
/// </summary>
public double X
{
get { return vector[]; }
set { vector[] = value; }
}
/// <summary>
/// Y向量
/// </summary>
public double Y
{
get { return vector[]; }
set { vector[] = value; }
}
/// <summary>
/// Z向量
/// </summary>
public double Z
{
get { return vector[]; }
set { vector[] = value; }
}
#endregion #region 向量操作
/// <summary>
/// /// <summary>
/// 向量加法+
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static Vector3d operator +(Vector3d lhs, Vector3d rhs)//向量加
{
Vector3d result = new Vector3d(lhs);
result.X += rhs.X;
result.Y += rhs.Y;
result.Z += rhs.Z;
return result;
}
/// <summary>
/// 向量减-
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static Vector3d operator -(Vector3d lhs, Vector3d rhs)//向量减法
{
Vector3d result = new Vector3d(lhs);
result.X -= rhs.X;
result.Y -= rhs.Y;
result.Z -= rhs.Z;
return result;
}
/// <summary>
/// 向量除
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static Vector3d operator /(Vector3d lhs, double rhs)//向量除以数量
{
if (rhs != )
return new Vector3d(lhs.X / rhs, lhs.Y / rhs, lhs.Z / rhs);
else
return new Vector3d(, , );
}
/// <summary>
/// 向量数乘*
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static Vector3d operator *(double lhs, Vector3d rhs)//左乘数量
{
return new Vector3d(lhs * rhs.X, lhs * rhs.Y, lhs * rhs.Z);
}
/// <summary>
/// 向量数乘
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static Vector3d operator *(Vector3d lhs, double rhs)//右乘数量
{
return new Vector3d(lhs.X * rhs, lhs.Y * rhs, lhs.Z * rhs);
} /// <summary>
/// 判断量向量是否相等
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns>True 或False</returns>
public static bool operator ==(Vector3d lhs, Vector3d rhs)
{
if (Math.Abs(lhs.X - rhs.X) < E && Math.Abs(lhs.Y - rhs.Y) < E && Math.Abs(lhs.Z - rhs.Z) < E)
return true;
else
return false;
}
public static bool operator !=(Vector3d lhs, Vector3d rhs)
{
return !(lhs == rhs);
}
public override bool Equals(object obj)
{
return base.Equals(obj);
}
public override int GetHashCode()
{
return base.GetHashCode();
}
public override string ToString()
{
return "(" + X + "," + Y + "," + Z + ")";
}
/// <summary>
/// 向量叉积，求与两向量垂直的向量
/// </summary>
public static Vector3d Cross(Vector3d v1, Vector3d v2)
{
Vector3d r = new Vector3d(, , );
r.X = (v1.Y * v2.Z) - (v1.Z * v2.Y);
r.Y = (v1.Z * v2.X) - (v1.X * v2.Z);
r.Z = (v1.X * v2.Y) - (v1.Y * v2.X);
return r;
}
/// <summary>
/// 向量数量积
/// </summary>
/// <param name="lhs"></param>
/// <param name="rhs"></param>
/// <returns></returns>
public static double operator *(Vector3d lhs, Vector3d rhs)//
{
return lhs.X * rhs.X + lhs.Y * rhs.Y + lhs.Z * rhs.Z;
}
/// <summary>
/// 内积
/// </summary>
/// <param name="v1"></param>
/// <param name="v2"></param>
/// <returns></returns>
public static double InnerMultiply(Vector3d v1, Vector3d v2)
{
double inner = 0.0;
inner = v1.X * v2.X + v1.Y * v2.Y + v1.Z * v2.Z;
return inner;
}
/// <summary>
/// 求向量长度，向量的模
/// </summary>
public static double Magnitude(Vector3d v1)
{
return (double)Math.Sqrt((v1.X * v1.X) + (v1.Y * v1.Y) + (v1.Z * v1.Z));
}
/// <summary>
/// 单位化向量
/// </summary>
public static Vector3d Normalize(Vector3d v1)
{
double magnitude = Magnitude(v1);
v1 = v1 / magnitude;
return v1;
}
#endregion
}
}
2. 计算基础
using System;
using RGeos.Geometry; namespace RGeos.Basic
{
public class RMath
{
public static double SMALL_NUM = 0.0000000001; // anything that avoids division overflow
/// <summary>
/// dot product (3D) which allows vector operations in arguments
/// </summary>
/// <param name="u"></param>
/// <param name="v"></param>
/// <returns></returns>
public static double dot(Vector3d u, Vector3d v)
{
return ((u).X * (v).X + (u).Y * (v).Y + (u).Z * (v).Z);
}
/// <summary>
/// 2D数量积，点乘
/// </summary>
/// <param name="u"></param>
/// <param name="v"></param>
/// <returns></returns>
public static double dot2(Vector3d u, Vector3d v)
{
return ((u).X * (v).X + (u).Y * (v).Y);
}
/// <summary>
/// 2D矢量叉积，定义为(0,0),P1,P2和P1P2包围四边形的带符号面积
/// </summary>
/// <param name="u"></param>
/// <param name="v"></param>
/// <returns></returns>
public static double perp(Vector3d u, Vector3d v)
{
return ((u).X * (v).Y - (u).Y * (v).X);
}
/// <summary>
/// 向量的模
/// </summary>
/// <param name="v"></param>
/// <returns></returns>
public static double norm(Vector3d v)
{
return Math.Sqrt(dot(v, v)); // norm = length of vector
}
/// <summary>
/// 两点间距离
/// </summary>
/// <param name="u"></param>
/// <param name="v"></param>
/// <returns></returns>
public static double d(Vector3d u, Vector3d v)
{
return norm(u - v); // distance = norm of difference
} public static double d(RPoint P1, RPoint P2)
{
return GetDistance(P1, P2); // distance = norm of difference
} // 判断点P2在直线P0P1的左边还是在右边，还是在直线上
//isLeft(): tests if a point is Left|On|Right of an infinite line.
// Input: three points P0, P1, and P2
// Return: >0 for P2 left of the line through P0 and P1
// =0 for P2 on the line
// <0 for P2 right of the line
public static int isLeft(RPoint P0, RPoint P1, RPoint P2)
{
double l = ((P1.X - P0.X) * (P2.Y - P0.Y) - (P2.X - P0.X) * (P1.Y - P0.Y));
return (int)l;
}
/// <summary>
/// 获取由两个点所形成的向量的象限角度
/// </summary>
/// <param name="preCoord">第一个点的坐标</param>
/// <param name="nextCoord">第二个点的坐标</param>
/// <returns></returns>
public static double GetQuadrantAngle(RPoint preCoord, RPoint nextCoord)
{
return GetQuadrantAngle(nextCoord.X - preCoord.X, nextCoord.Y - preCoord.Y);
}
/// <summary>
/// 由增量X和增量Y所形成的向量的象限角度
/// 区别方位角：方位角以正北方向顺时针
/// | |
/// | / |b /
/// | / a |^/
/// |/_)____象限角 |/______方位角
/// </summary>
/// <param name="x">增量X</param>
/// <param name="y">增量Y</param>
/// <returns>象限角</returns>
public static double GetQuadrantAngle(double x, double y)
{
double theta = Math.Atan(y / x);
if (x > && y == ) return ;
if (x == && y > ) return Math.PI / ;
if (x < && y == ) return Math.PI;
if (x == && y < ) return * Math.PI / ; if (x > && y > ) return theta;
if (x > && y < ) return Math.PI * + theta;
if (x < && y > ) return theta + Math.PI;
if (x < && y < ) return theta + Math.PI;
return theta;
}
/// <summary>
/// 获取由相邻的三个点A-B-C所形成的两个向量之间的夹角
/// 向量AB，BC形成的夹角
/// </summary>
/// <param name="preCoord">第一个点</param>
/// <param name="midCoord">中间点</param>
/// <param name="nextCoord">第三个点</param>
/// <returns></returns>
public static double GetIncludedAngle(RPoint preCoord, RPoint midCoord, RPoint nextCoord)
{
double innerProduct = (midCoord.X - preCoord.X) * (nextCoord.X - midCoord.X) + (midCoord.Y - preCoord.Y) * (nextCoord.Y - midCoord.Y);
double mode1 = Math.Sqrt(Math.Pow((midCoord.X - preCoord.X), 2.0) + Math.Pow((midCoord.Y - preCoord.Y), 2.0));
double mode2 = Math.Sqrt(Math.Pow((nextCoord.X - midCoord.X), 2.0) + Math.Pow((nextCoord.Y - midCoord.Y), 2.0));
return Math.Acos(innerProduct / (mode1 * mode2));
}
/// <summary>
/// 获取由两个点所形成的向量的模(长度)
/// </summary>
/// <param name="preCoord">第一个点</param>
/// <param name="nextCoord">第二个点</param>
/// <returns>由两个点所形成的向量的模(长度)</returns>
public static double GetDistance(RPoint preCoord, RPoint nextCoord)
{
return Math.Sqrt(Math.Pow((nextCoord.X - preCoord.X), ) + Math.Pow((nextCoord.Y - preCoord.Y), ));
}
}
}
RMath