# 雁过无痕

 C++博客 :: 首页 :: 新随笔 :: 联系 :: 聚合  :: 管理 ::

ab × ac = ab × ap + ap × ac + pb × pc

|ab × ac| = |ab × ap| + |ap × ac| + |pb × pc|

#include<cfloat>

template<typename T> class Vec3 {

T x, y, z;

public:

Vec3(T xx, T yy, T zz) : x(xx), y(yy), z(zz) {}

Vec3 operator-(const Vec3& v) const { return Vec3(x - v.x, y - v.y, z - v.z); }

T dot(const Vec3& v) const { return x * v.x + y * v.y + z * v.z; }

Vec3 cross(const Vec3& v) const {

return Vec3(y * v.z - z * v.y, z * v.x - x * v.z,x * v.y - y * v.x);

}

};

typedef Vec3<double> V3d;

//方法一

bool is_in_triangle3a(const V3d& a, const V3d& b, const V3d& c, const V3d& p)

{

V3d ab(b - a), ac(c - a), ap(p -a);

if (fabs(ab.cross(ac).dot(ap)) >= DBL_EPSILON) return false; //四点不共面

V3d abc = ab.cross(ac), abp = ab.cross(ap), apc = ap.cross(ac);

double t0 = abc.dot(abc), t1 = abp.dot(abc), t2 = apc.dot(abc);

//t1 >= 0     t2 >= 0    t1 + t2 <= t0       t0肯定大于0

// return (t1 >= -DBL_EPSILON) & (t2 >= -DBL_EPSILON) & (t0 - t1 - t2 >= -DBL_EPSILON);

double delta = fabs(t1) + fabs(t2) + fabs(t0 - t1 - t2) - t0;

return fabs(delta) < DBL_EPSILON;

}

(a × b)· c = (c × a)· b　　　    （混合积）

a × (b × c) = b（a·c） c（a·b）  (拉格朗日公式)

(a × b)·(a × c) = ((a × c) × a)·b = ((a·a)c – (a·c)a)·b

= (a·a) * (b·c) – (a·c) * (a·b)

//方法二

bool is_in_triangle3b(const V3d& a, const V3d& b, const V3d& c, const V3d& p)

{

V3d ab(b - a), ac(c - a), ap(p -a);

if (fabs(ab.cross(ac).dot(ap)) >= DBL_EPSILON) return false; //四点不共面

V3d abc = ab.cross(ac), abp = ab.cross(ap), apc = ap.cross(ac);

//double t0 = abc.dot(abc), t1 = abp.dot(abc), t2 = apc.dot(abc); //对这三个计算公式进行展开

double v11 = ab.dot(ab), v22 = ac.dot(ac), v12 = ab.dot(ac);

double v13 = ab.dot(ap), v23 = ac.dot(ap);

double t0 = v11 * v22 - v12 * v12;

double t1 = v11 * v23 - v12 * v13;

double t2 = v22 * v13 - v12 * v23;

double delta = fabs(t1) + fabs(t2) + fabs(t0 - t1 - t2) - t0;

return fabs(delta) < DBL_EPSILON;

}

posted on 2011-07-14 23:28 flyinghearts 阅读(1469) 评论(0)  编辑 收藏 引用 所属分类: 算法