【越飞越高讲堂12】最小包围盒和最大距离点对 - AutoLISP/Visual LISP 编程
最小包围盒和最大距离点对在实际生活中有很多运用。譬如，最小包围盒指的是能包围一个多边形或者一些多边形的最小面积（可能周长）的这样的一个矩形，这个帖子我只是讨论二维的，所以说是矩形，
求得最小包围和就能合理利用材料，选择合适的截面。
在这个算法中，采用了“游标卡尺”的算法，所以运行速度很快，在求出一个形状的凸包之后，就能迅速得到最小矩形。O(n *log(n))的时间复杂度，
大量运算成为可能。
对于CAD本身来说，用getboundingBox得出的结果仅仅是平行与WCS的坐标系统的矩形，而且有时候很不准确，特别对于spline来说，有时候相差很大。
然而，WCS下的boundingBox一般来说不是最小矩形，如果要得到最小的，必须要旋转这个图形很多次才能得到结果，但这时一个费事费力的事情。
所以，我们有必要考虑一种算法。
这些天通过几个比较晚的晚上，终于得出了算法。为了增加讨论效果，我在这里先贴出arx文件，适用于2004-2006平台。至于LISP，以后在贴出。
速度之快，几乎是0秒。
用法： 先输入命令 Hull,
然后输入分弧精度---指的是对样条曲线或者弧形的分段数目,取值在100-2000比较合理,太大容易引起问题且也不能有效提高精度。
选择图形中的多边形或者样条曲线，然后你就可以看到结果了。最远距离用红线标出。
;;;The procedure for Test
(defun C:test (/ PP PTLIST SEL T0 n)
(princ "\n命令是test(The command is Test).")
(princ "\nPlease select LWPOLYLINE,LINE,SPLINE,POINT")
(setq sel (ssget (list '(0 . "POINT,LWPOLYLINE,LINE,SPLINE")))) ;select curve or point
(initget 7)
(setq n (getint "\nPlease Enter the precision of dividing curve(The resonable value is 100 to 2000):"))
(if sel
(progn
(setq ptlist (getpt sel n)) ;construct the set of points
(setq ptlist (Graham-scan ptlist)) ;construct the CCW Hull of this set.
(if (<= (det (car ptlist) (cadr ptlist) (caddr ptlist)) 0.0) ;ensure the hull is CCW.
(setq ptlist (reverse ptlist)) ;if it isn't CCW,then reverse it
)
(setq t0 (getvar "TDUSRTIMER")) ;The start time of this algorithm
(setq pp (car (MinAreaRectangle ptlist))) ;start calculating
(princ "\nIt takes :")
(princ (* (- (getvar "TDUSRTIMER") t0) 86400)) ;The End time
(princ "seconds")
(if pp
(make-poly pp) ;draw rectangle.
)
)
)
(princ)
)
;;;=======================================================
;;;Function : Find the minimum area of encasing rectangle.
;;;Arguments : A CCW HULL
;;;Return: The Four points of Rectangle and its Area
;;;=======================================================
(defun MinAreaRectangle (ptlist / AA AI BB D1
D2 EDGEI I1X I1Y I2X I2Y
IL INF IX IY J1 J2 MINA
MINHMINWNORHNORM PI1 PI2 PTI0
PTI1PTI2PTJ1PTK1 PTM1PTS1PTS2
PTS3PTS4REC1REC2 REC3REC4RECT
VECHVECLVJ12VM12
)
(setq INF 1e309)
(setq minA INF) ;Initiating the Minimum area is infinite
(setq pti0 (car ptlist)) ;the first point of Hull.
(setq pts1 (append ptlist (list pti0))) ;add the first point at back of Hull
(setq pts2 (cdr (append ptlist ptlist (list pti0)))) ;Construct a loop for the hull
(setq i 0)
;;Find area of encasing rectangle anchored on each edge.
(repeat (length ptlist)
(setq pi1 (car pts1)
pi2 (cadrpts1)
i1x (car pi1)
i1y (cadrpi1)
i2x (car pi2)
i2y (cadrpi2)
ix(- i2x i1x)
iy(- i2y i1y)
il(distance (list ix iy) '(0.0 0.0))
)
;;寻找最左点
;;Find a vertex on on first perpendicular line of support
(while (> (DOTPR ix iy pts2) 0.0)
(setq pts2 (cdr pts2))
)
;;寻找最上点
;;Find a vertex on second perpendicular line of suppoer
(if (= i 0)
(setq pts3 pts2)
)
(while (> (CROSSPR ix iy pts3) 0.0)
(setq pts3 (cdr pts3))
)
;;寻找最右点
;;Find a vertex on second perpendicular line of suppoer
(if (= i 0)
(setq pts4 pts3)
)
(while (< (DOTPR ix iy pts4) 0.0)
(setq pts4 (cdr pts4))
)
;;得出了每边的矩形
;;Find distances between parallel and perpendicular lines of support
(cond
((equal i1x i2x 1e-4) ;如果边两点的X值相同
(setq d1 (- (caar pts3) i1x) ;那么矩形的高就是最上点与边的X的差值
d2 (- (cadar pts4) (cadar pts2)) ;矩形的宽就是最左和最右的Y的差值
)
)
((equal i1y i2y 1e-4) ;如果边两点的Y值相同
(setq d1 (- (cadar pts3) i1y) ;那么矩形的高就是最上点与边的Y的差值
d2 (- (caar pts4) (caar pts2)) ;矩形的宽就是最左和最右的X的差值
)
)
(T
(setq aa (det pi1 pi2 (car pts3))) ;否则计算边和最上点构成的面积的二倍(det)
(setq d1 (/ aa il)) ;高就是det值除以边长
(setq j1 (car pts2)) ;最右边点
(setq j2 (list (- (car j1) iy) (+ (cadr j1) ix))) ;通过最右边点的垂直边的点
(setq bb (det j1 j2 (car pts4))) ;最右边点，上面的点和最左边的点
(setq d2 (/ bb il)) ;这三点的det除以边长就是宽
)
)
;;计算矩形的面积，必要时更新最小面积
;;Compute area of encasing rectangle anchored on current edge.
;;if the area is smaller than the old Minimum area,then update,and record the width,height and five points.
(setq Ai (abs (* d1 d2))) ;面积就是高和宽的积
(if (< Ai MinA) ;如果面积小于先前的最小面积，则记录：
(setq MinA Ai ;更新最小面积
MinH d1 ;最小面积的高
MinW d2 ;最小面积的宽
pti1 pi1 ;最小面积的边的第一个端点
pti2 pi2 ;最小面积的边的第二个端点
ptj1 (car pts2) ;最右边的点
ptk1 (car pts3) ;最上面的点
ptm1 (car pts4) ;最左边的点
)
)
(setq pts1 (cdr pts1)) ;检测下一条边
(setq i (1+ i)) ;计数器加一
)
;;according to the result ,draw the Minimum Area Rectangle
(setq edge (mapcar '- pti2 pti1)) ;最小面积的边对应的向量
(setq VecL (distance edge '(0.0 0.0))) ;最小面积的边的长度
(setq NorH (abs (/ MinH VecL))) ;这边的法线
(setq Norm (list (- (cadr edge)) (car edge))) ;这边的垂直向量
(setq vj12 (mapcar '+ ptj1 Norm)) ;通过最右点的垂直向量
(setq vm12 (mapcar '+ ptm1 Norm)) ;通过最左点的垂直向量
(setq vecH (mapcar '* (list NorH NorH) Norm))
(setq rec1 (inters pti1 pti2 ptj1 vj12 nil)) ;矩形的第一点
(setq rec4 (inters pti1 pti2 ptm1 vm12 nil)) ;矩形的第四点
(setq rec2 (mapcar '+ rec1 vecH)) ;矩形的第二点
(setq rec3 (mapcar '+ rec4 vecH)) ;矩形的第三点
(setq rect (list Rec1 rec2 rec3 rec4)) ;矩形的点表
(cons rect MinA) ;返回这个矩形的点表和最大距离
)
;;;========================================
;;;求凸壳的直径的程序
;;;参数：逆时针的凸壳 H-------注意逆时针!!!
;;;返回值: 直径的两个端点和直径 Pair . MaxD
;;;========================================
(defun Max-distance (H / D M MAXD P PAIR Q U V W)
(setq Q (cdr (append H H (list (car H))))) ;构造一个首尾循环的凸集,且起始点为凸壳的第二点
(setq MaxD 0.0) ;初始化最小距离为0
(foreach U H ;依次检查凸壳的边
(setq V (car Q)) ;循环集的第一点
(setq W (cadr Q)) ;循环集的第二点
(setq M (mid-pt V W)) ;这两点的中点
(while (> (dot M U V) 0.0) ;如果夹角小于90度（即点积大于0）
(setq Q (cdr Q)) ;循环集推进
(setq V (car Q)) ;取下一点
(setq W (cadr Q)) ;下下一点
(setq M (mid-pt V W)) ;这两点的中点
)
(setq D (distance U V)) ;计算这时的最大距离
(if (> D MaxD) ;如果大于前面的最大距离
(setq MaxD D ;就替换前面的最大距离
Pair (list U V) ;并记录这对点
)
)
)
(cons Pair MaxD) ;返回这对点和最大距离
)
;;;==================
;;;Graham扫描法求凸包
;;;==================
(defun Graham-scan (ptlist / hullpt maxXpt sortPt P Q)
(if (< (length ptlist) 3) ;3点以下
ptlist ;是本集合
(progn
(setq maxXpt (assoc (apply 'max (mapcar 'car ptlist)) ptlist)) ;最右边的点
(setq sortPt (sort-by-angle-distance ptlist maxXpt)) ;分类点集
(setq hullPt (list (cadr sortPt) maxXpt)) ;开始的两点
(foreach n (cddr sortPt) ;从第3点开始
(setq hullPt (cons n HullPt)) ;把Pi加入到凸集
(setq P (cadr hullPt)) ;Pi-1
(setq Q (caddr hullPt)) ;Pi-2
(while (and q (> (det n P Q) -1e-6)) ;如果左转
(setq hullPt (cons n (cddr hullPt))) ;删除Pi-1点
(setq P (cadr hullPt)) ;得到新的Pi-1点
(setq Q (caddr hullPt)) ;得到新的Pi-2点
)
)
(reverse hullpt) ;返回凸集
)
)
)
;;;中点函数
(defun mid-pt (p1 p2)
(list
(* (+ (car p1) (car p2)) 0.5)
(* (+ (cadr p1) (cadr p2)) 0.5)
)
)
;;;以某点为基点，按照角度和距离分类点集
(defun sort-by-angle-distance (ptlist pt / )
(vl-sort ptlist
(function
(lambda (e1 e2 / ang1 ang2 )
(setq ang1 (angle pt e1))
(setq ang2 (angle pt e2))
(if (= ang1 ang2)
(< (distance pt e1) (distance pt e2))
(< ang1 ang2)
)
)
)
)
)
;;;点积= x1*x2 + y1*y2
(defun DOTPR (ix iy pts / pt1 pt2)
(setq pt1 (car pts))
(setq pt2 (cadr pts))
(+ (* ix (- (car pt2) (car pt1)))
(* iy (- (cadr pt2) (cadr pt1)))
)
)
;;;叉积= x1*y2 - x2*y1
(defun CROSSPR (ix iy pts / pt1 pt2)
(setq pt1 (car pts))
(setq pt2 (cadr pts))
(- (* ix (- (cadr pt2) (cadr pt1)))
(* iy (- (car pt2) (car pt1)))
)
)
;;;中点函数
(defun mid-pt (p1 p2)
(list
(* (+ (car p1) (car p2)) 0.5)
(* (+ (cadr p1) (cadr p2)) 0.5)
)
)
;;;定义三点的行列式,即三点之倍面积
(defun det (p1 p2 p3 / x2 y2)
(setq x2 (car p2)
y2 (cadr p2)
)
(- (* (- x2 (car p3)) (- y2 (cadr p1)))
(* (- x2 (car p1)) (- y2 (cadr p3)))
)
)
;;;定义向量的点积函数
(defun dot (p1 p2 p3 / x1 y1)
(setq x1 (car p1)
y1 (cadr p1)
)
(+ (* (- (car p2) x1) (- (car p3) x1))
(* (- (cadr p2) y1) (- (cadr p3) y1))
)
)
;;;取点函数2
(defun getpt (ss n / i s a b c d e)
(setq i 0)
(if ss
(repeat (sslength ss)
(setq a (ssname ss i))
(setq b (entget a))
(setq e (cdr (assoc 0 b)))
(cond
((= e "LWPOLYLINE")
(setq c (get-pline-vertexs a n))
(setq s (append c s))
)
((= e "SPLINE")
(setq c (get-spline-vertexs a n))
(setq s (append c s))
)
((= e "LINE")
(setq c (cdr (assoc 10 b)))
(setq d (cdr (assoc 11 b)))
(setq c (list (car c) (cadr c)))
(setq d (list (car d) (cadr d)))
(setq s (cons c s))
(setq s (cons d s))
)
((= e "POINT")
(setq c (cdr (assoc 10 b)))
(setq c (list (car c) (cadr c)))
(setq s (cons c s))
)
)
(setq i (1+ i))
)
)
s
)
;;取得多边形顶点
(defun get-LWpolyline-vertexs (DXF / lst)
(foreach n DXF
(if (= (car n) 10)
(setq lst (cons (cdr n) lst))
)
)
(reverse lst)
)
(defun get-3dpolyline-vertexs ( ent / p )
(if (and (setq ent (entnext ent)) (setq p (cdr (assoc 10 (entget ent)))))
(cons p (get-3dpolyline-vertexs ent))
)
p
)
;;;取得样条曲线的点
(defun get-spline-vertexs (ent n / DIST ENDPAR I LEN OBJ PT PTS SEG)
(setq obj (vlax-ename->vla-object ent))
(setq endpar(vlax-curve-getEndParam obj))
(setq len (vlax-curve-getDistAtParam obj endpar))
(setq seg (/ len n))
(setq dist 0)
(while (< dist len)
(setq pt (vlax-curve-getPointAtDist obj dist))
(setq pts (cons pt pts))
(setq dist (+ seg dist))
)
(if (= (vla-get-closed obj) :vlax-false)
(setq pt (vlax-curve-getEndPoint obj)
pts (cons pt pts)
)
)
(reverse pts)
)
;;;取得含有圆弧的多段线的点
(defun get-pline-vertexs (ent n / BLG DIST ENDPAR I L1 L2 L3 LI OBJ PT PTS VEXNUM)
(setq obj (vlax-ename->vla-object ent))
(setq endpar (vlax-curve-getEndParam obj))
(setq vexNum (fix endPar))
(setq pts nil)
(setq i 0)
(repeat vexNum
(setq pt (vlax-curve-getPointAtParam obj i))
(setq pts (cons pt pts))
(setq blg (vla-getbulge obj i))
(if (/= blg 0.0)
(progn
(setq l1 (vlax-curve-getDistAtParam obj i))
(setq l2 (vlax-curve-getDistAtParam obj (1+ i)))
(setq l3 (- l2 l1))
(setq li (/ l3 n))
(setq dist l1)
(repeat (1- n)
(setq dist (+ dist li))
(setq pt (vlax-curve-getPointAtDist obj dist))
(setq pts (cons pt pts))
)
)
)
(setq i (1+ i))
)
(if (= (vla-get-closed obj) :vlax-false)
(setq pt (vlax-curve-getEndPoint obj)
pts (cons pt pts)
)
)
pts
)
;;;绘制多段线
(defun Make-Poly (pp / x)
(entmake ;画凸包
(append
'((0 . "LWPOLYLINE")
(100 . "AcDbEntity")
(100 . "AcDbPolyline")
)
(list (cons 90 (length pp))) ;顶点个数
(mapcar
(function (lambda (x) (cons 10 x)))
pp
) ;多段线顶点
(list (cons 70 1)) ;闭合的
(list (cons 62 1)) ;红色的
)
)
)
;;;========================================
;;;求凸壳的直径的程序
;;;参数：逆时针的凸壳 H-------注意逆时针!!!
;;;返回值: 直径的两个端点和直径 Pair . MaxD
;;;========================================
(defun Max-distance (H / D M MAXD P PAIR Q U V W)
(setq Q (cdr (append H H (list (car H))))) ;构造一个首尾循环的凸集,且起始点为凸壳的第二点
(setq MaxD 0.0) ;初始化最小距离为0
(foreach U H ;依次检查凸壳的边
(setq V (car Q)) ;循环集的第一点
(setq W (cadr Q)) ;循环集的第二点
(setq M (mid-pt V W)) ;这两点的中点
(while (> (dot M U V) 0.0) ;如果夹角小于90度（即点积大于0）
(setq Q (cdr Q)) ;循环集推进
(setq V (car Q)) ;取下一点
(setq W (cadr Q)) ;下下一点
(setq M (mid-pt V W)) ;这两点的中点
)
(setq D (distance U V)) ;计算这时的最大距离
(if (> D MaxD) ;如果大于前面的最大距离
(setq MaxD D ;就替换前面的最大距离
Pair (list U V) ;并记录这对点
)
)
)
(cons Pair MaxD) ;返回这对点和最大距离
)
页:
[1]