[简介] 重心坐标的应用
本帖最后由 qjchen 于 2011-1-27 23:41 编辑最近，为了求一些三角形的内心、外心坐标学习了一下重心坐标，随手记下一些心得~
比如，假如知道了 三角形三个点 (x1,y1,z1) (x2,y2,z2) (x3 y3 z3) 请问，其内心坐标是多少，外心坐标是多少。
乍一看这问题很简单，但真正求起来还颇麻烦。
于是google之，发现了 Barycentric coordinates 这个名词
其定义可参见
中文wiki
http://zh.wikipedia.org/zh/%E9%87%8D%E5%BF%83%E5%9D%90%E6%A0%87
英文wiki（比较详细）
http://en.wikipedia.org/wiki/Barycentric_coordinates_(mathematics)
cut-the-knot（比较简单，有些链接）
http://www.cut-the-knot.org/triangle/glasses.shtml
Mathworld （对各心有一个比较详细的公式，好用）
http://mathworld.wolfram.com/BarycentricCoordinates.html
那么，也就是说，只要按 http://mathworld.wolfram.com/BarycentricCoordinates.html 这里的公式，计算出各点所对应的权重，那么 就可以根据 http://en.wikipedia.org/wiki/Barycentric_coordinates_(mathematics) 的公式
x=λ1x1+λ2x2+λ3x3.....（其中λ是之前的重心坐标的归一化，也就是λ1:λ2:λ3=重心坐标比，和=1）
就可以得到各种心的坐标了。
公式多数是轮换的，如下图，就很方便用计算软件来进行一下推导验证某些其他性质了。
附上一个
原来发在此处http://www.theswamp.org/index.php?topic=36786.0的LISP代码
;;;;vec plus
(defun q:vec:+(v1 v2) (mapcar '+ v1 v2))
;;;;vec substract
(defun q:vec:-(v1 v2) (mapcar '- v1 v2))
;;;;vec plus constant
(defun q:vec:*c(v a)(mapcar '(lambda(x) (* x a)) v))
;;;;vec dot product
(defun q:vec:dot*(v1 v2) (apply '+ (mapcar '* v1 v2)))
;;;;vec cross product
(defun q:vec:cross*(v1 v2)
(list (q:det:2 (cadr v1) (caddr v1) (cadr v2) (caddr v2))
(q:det:2 (caddr v1) (car v1) (caddr v2) (car v2))
(q:det:2 (car v1) (cadr v1) (car v2) (cadr v2)))
)
;;;;Normalize a vec
(defun q:vec:Norm(v / l)
(if (not (zerop (setq l (distance '(0 0 0) v))))
(mapcar '(lambda(x) (/ x l)) v))
)
;;;; Vector Length
(defun q:vec:Length(v) (sqrt (apply '+ (mapcar '(lambda(x) (expt x 2)) v))))
;;;;determinant library
;;;;cal determinant
;;;;|a1 a2|
;;;;|b1 b2|
(defun q:det:2(a1 a2 b1 b2)
(- (* a1 b2) (* a2 b1))
)
;;;; Geometry code
;;;; qjchen@gmail.com
;;;;http://en.wikipedia.org/wiki/Distance_from_a_point_to_a_line
;;;; distance of point per to line
(defun q:geo:point-dis-to-line(p p1 p2 / n)
(setq n (q:vec:Norm (q:vec:- p2 p1)))
(q:vec:Length (q:vec:- (q:vec:+ p1 (q:vec:*c n (q:vec:dot* (q:vec:- p p1) n))) p))
)
;;;;;;;;;;;;;;;;;;;;;;;;;;;;Use Barycentric Coordinates to calculate multiple center
;;;http://mathworld.wolfram.com/BarycentricCoordinates.html
;;;http://en.wikipedia.org/wiki/Barycentric_coordinates_(mathematics)#Converting_to_barycentric_coordinates
(defun q:geo:barycentric-normalize(a b c)
(mapcar '(lambda(x) (/ x (+ a b c))) (list a b c))
)
;;;;transform from barycentric-to-Cartesian coordinate
(defun q:geo:barycentric-to-Cartesian(tlst P1 P2 P3)
(q:vec:+(q:vec:+ (q:vec:*c P1 (car tlst)) (q:vec:*c P2 (cadr tlst))) (q:vec:*c P3 (caddr tlst)))
)
;;Find out the incircle of three point (three distinct, noncollinear WCS points)
;;;; qjchen@gmail.com
(defun q:geo:incircle3d(p1 p2 p3 / a b c p)
(setq a (distance p2 p3) b (distance p1 p3) c (distance p1 p2))
(setq p(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize a b c) P1 P2 P3))
(list p (q:geo:point-dis-to-line p p1 p2))
)
;;Find out the circumcircle of three point (three distinct, noncollinear WCS points)
(defun q:geo:circumcircle3d(p1 p2 p3 / a b c p temp)
(setq a (distance p2 p3) b (distance p1 p3) c (distance p1 p2))
(defun temp(a b c) (* a a (+ (* b b) (* c c) (* a (- a)))))
(setq p(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize (temp a b c) (temp b c a) (temp c a b)) P1 P2 P3))
(list p (distance p1 p))
)
;;Find out the triangle centroid of three point (three distinct, noncollinear WCS points)
(defun q:geo:centroid3d(p1 p2 p3 / p)
(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize 1.0 1.0 1.0) P1 P2 P3)
)
;;Find out the symmedian point of three point (three distinct, noncollinear WCS points)
(defun q:geo:symmedianpoint3d(p1 p2 p3 / a b c)
(setq a (distance p2 p3) b (distance p1 p3) c (distance p1 p2))
(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize (* a a) (* b b) (* c c)) P1 P2 P3)
)
;;Find out the orthocenter of three point (three distinct, noncollinear WCS points)
(defun q:geo:orthocenterpoint3d(p1 p2 p3 / a b c temp)
(setq a (distance p2 p3) b (distance p1 p3) c (distance p1 p2))
(defun temp(a b c) (* (+ (* a a) (* b b) (* c (- c)))(+ (* c c) (* a a) (* b (- b)))))
(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize (temp a b c) (temp b c a) (temp c a b)) P1 P2 P3)
)
;;Find out the Nagel point of three point (three distinct, noncollinear WCS points)
(defun q:geo:Nagelpoint3d(p1 p2 p3 / a b c temp)
(setq a (distance p2 p3) b (distance p1 p3) c (distance p1 p2))
(defun temp(x) (- (/ (+ a b c) 2.0) x))
(q:geo:barycentric-to-Cartesian (q:geo:barycentric-normalize (temp a) (temp b) (temp c)) P1 P2 P3)
)
(defun c:test1 (/ p1 p2 p3 res1 res2 res3)
(command "UCS" "W")
(setvar "OSMODE" 1)
(setvar "PDMODE" 2)
(setq p1 (getpoint)
p2 (getpoint)
p3 (getpoint)
)
(setvar "OSMODE" 0)
(setq res1 (q:geo:incircle3d p1 p2 p3))
(setq res2 (q:geo:circumcircle3d p1 p2 p3))
(setq res3 (q:geo:orthocenterpoint3d p1 p2 p3))
(command "UCS" "3" p1 p2 p3)
(command "POINT" (trans (car res1) 0 1))
(command "CIRCLE" (trans (car res1) 0 1) (cadr res1))
(command "POINT" (trans (car res2) 0 1))
(command "CIRCLE" (trans (car res2) 0 1) (cadr res2))
(command "POINT" (trans res3 0 1))
(command "UCS" "W")
(grdraw (list 0 0 0) res3 1)
)
支持了！ 本帖最后由 Ea 于 2011-2-20 22:45 编辑
内切圆，几何作图法(角平分线)(defun incir (p1 p2 p3 / an12 an13 an21 an23 d p22 p11 p01 pcen r)
(setq an12 (angle p1 p2)
an13 (angle p1 p3)
an21 (angle p2 p1)
an23 (angle p2 p3)
d (distance p1 p2)
p22(polar p1 an13 d)
p11(polar p2 an23 d)
p00(list (/ (+ (car p2) (car p22)) 2)
(/ (+ (cadr p2) (cadr p22)) 2)
)
p01(list (/ (+ (car p1) (car p11)) 2)
(/ (+ (cadr p1) (cadr p11)) 2)
)
pcen (inters p1 p00 p01 p2 nil)
r (distance pcen
(inters p1 p2 pcen (polar pcen (+ an21 (/ pi 2)) d) nil)
)
)
(list pcen r)
)
(defun c:test (/ p1 p2 p3 p)
(if (and (setq p1 (getpoint "\nFirst Point: "))
(setq p2 (getpoint p1 "\nSecond Point: "))
(setq p3 (getpoint p1 "\nThird Point: "))
)
(progn
(vl-cmdf ".pline" "_non" p1 "_non" p2 "_non" p3 "c")
(setq p (incir p1 p2 p3))
(vl-cmdf ".circle" (car p) (cadr p))
)
)
(princ)
) 支持了！ 本帖最后由 qjchen 于 2011-2-20 23:06 编辑
回复 Ea 的帖子
谢谢Eachy版主~
当时也写过一段角平分线法的，既写了，也发一下，和eachy版主的算法应该是差不多的。那个重心坐标的好处就在于可以求出空间状态下的圆心。
;;2d angle bisector
(defun q:geo:angle-bisector(p1 p2 p3)
(list p1 (polar p1 (/ (+ (angle p1 p2) (angle p1 p3)) 2.0) 1.0))
)
;;2d distance from a point to a line
(defun q:geo:dis-point-to-2dline(p p1 p2)
(distance p (inters p (polar p (+ (angle p1 p2) (/ pi 2.0)) 1.0) p1 p2 nil))
)
;;Find out the incircle of three point (three distinct, noncollinear WCS points)
(defun q:geo:incircle(p1 p2 p3 / l1 l2 p)
(setq l1 (q:geo:angle-bisector p1 p2 p3) l2 (q:geo:angle-bisector p2 p3 p1))
(setq p (inters (car l1) (cadr l1) (car l2) (cadr l2) nil))
(list p (q:geo:dis-point-to-2dline p p1 p2))
)
;;;测试函数
;;;;;Test function
(defun c:test(/ p1 p2 p3)
(setq p1 (getpoint) p2 (getpoint) p3 (getpoint))
(q:entmake:circle (car (q:geo:incircle p1 p2 p3)) (cadr (q:geo:incircle p1 p2 p3)) (getvar "clayer"))
)
;;entmake circle
(defun q:entmake:circle (center rad layer)
(entmake (list (cons 0 "CIRCLE") (cons 8 layer) (cons 10 center) (cons 40 rad)))
)
页:
[1]