;;;************************************
;;;求最小包圍圓的lisp程序--------------
;;;其算法為參見了有關文獻--------------
;;;這種算法在退化很嚴重的情況結果也正確
;;;其中程序主段是核心算法，其他的附加程
;;;序為取點，畫點，畫圓和半徑，用來測試
;;;************************************
(defun C:test (/ olderr en errmsg oce oldmode ss ptlist x cen radius ptmax)
  ;;定義錯誤函數和預處理---------------
  (setvar "errno" 0)
  (setq olderr *error*)
  (defun *error* (msg / en errmsg)
    (setq en (getvar "errno"))
    (setq errmsg (strcat "errno=" (itoa en) "\nError:" msg))
    (alert errmsg)
    (setq *error* olderr)
  )
  (graphscr)
  (setq oldmode (getvar "osmode"))
  (setq oce (getvar "cmdecho"))
  (setvar "cmdecho" 0)
  (command ".ucs" "W")
  ;;取點，畫點，并對函數用時計算-------
  (setq sl ' ( (-4 . "<OR" )
	      (0 . "POINT")
	      (0 . "LINE")
	      (0 . "POLYLINE")
	      (0 . "LWPOLYLINE")
	      (-4 . "OR>" )))
  (setq ss (ssget sl))
  (setq ptlist (ssgetpoint1 ss))
  (setq t1 (getvar "CDATE"))
  (setq x (mincir ptlist))
  (setq t2 (getvar "CDATE"))
  (setq usetime (* (- t2 t1) 1e6))
  (princ (strcat "\n用時=" (rtos usetime 2 6) "秒"))
  (if (= nil x)
    (progn 
      (alert "點的有效數目太小，請重新輸入!")
      (command ".ucs" "p")
      (setvar "osmode" oldmode)
      (setvar "cmdecho" oce)
      (princ "\n")
      (princ)
    )
    (progn 
      (setq cen (car x) radius (cadr x) ptmax (caddr x))
      ;;;畫圓及半徑，列出圓的圓心半徑值
      (entmake
        (append
          '((0 . "circle") (100 . "AcDbEntity") (100 . "AcDbCircle")) 
          (list (cons 10 cen))(list (cons 40 radius))(list (cons 62 1))
        )
      )
      (entmake
        (append
          '((0 . "line") (100 . "AcDbEntity") (100 . "AcDbLine")) 
          (list (cons 10 cen))(list (cons 11 ptmax))(list (cons 62 1))
        )
      )
      (command ".ucs" "p")
      (setvar "osmode" oldmode)
      (setvar "cmdecho" oce)
      (princ "\n")
      (list cen radius)
    )
  )
)
;;;************************************
;;;求最小包圍圓的函數，空集返回空集，否
;;;則返回最小圓的圓心，半徑和圓上的一點
;;;這是程序的主段----------------------
;;;************************************
(defun mincir (ptlist / p1 p2 p3 ptmax cen_r cen radius)
  ;;定義中點函數，本來R2004版中無須定義
  ;;但不知道為什么到R2006版沒有定義了。
  (defun mid (p1 p2)
    (polar p1 (angle p1 p2) (* (distance p1 p2) 0.5)) 
  )
  ;;判斷有效點個數---------------------
  (cond
    ((= (length ptlist) 0)
     nil
    )
    ((= (length ptlist) 1)
     (progn
       (alert "點集合為一點,最小圓半徑為0")
       (list (car ptlist) 0 (car ptlist))
     )
    )
    ((= (length ptlist) 2)
     (progn
       (alert "點集合為兩點，最小圓直徑為其兩點距離，\n圓心為其連線中點")
       (setq cen (mid (car ptlist) (cadr ptlist)) radius (/ (distance (car ptlist) (cadr ptlist)) 2))
       (list cen radius (car ptlist))
     )
    )
    (t
     (progn
       ;;上面啰嗦的一大段在實際情況中一般不會出現
       ;;判斷點是否在圓內------------------------
       (defun in1 (pt cen r)
	 (if (< (- (distance pt cen) r) 1e-8)
	   t
	   nil
	 )
       )
       ;;判斷點集是否在圓內----------------------
       (defun in2 (ptl cen r)
	 (if (apply 'and (mapcar '(lambda (x) (in1 x cen r))  ptl))
	   t
	   nil
	 )
       )
       ;;定義三點最小圓圓心及其半徑，若是銳角三角
       ;;形，則是其三點圓，否則是其最大邊的直徑圓
       (defun 3pc (pa pb pc / a b c l p ja jb jc ppa ppb ppc cen radius)
	 (setq a (list (distance pb pc) pa))
	 (setq b (list (distance pc pa) pb))
	 (setq c (list (distance pa pb) pc))
	 (setq l (list a b c))
	 (setq p (/ (+ (car a) (car b) (car c)) 2))
	 (setq a (nth (car (vl-sort-i (mapcar 'car l) '>)) l))
	 (setq b (nth (cadr (vl-sort-i (mapcar 'car l) '>)) l))
	 (setq c (nth (caddr (vl-sort-i (mapcar 'car l) '>)) l))
	 (setq l (+ (* (car b) (car b)) (* (car c) (car c)) (* (car a) (car a) -1.0)))
	 ;;上面l利用了余弦定理作為判斷-----------
	 (if (< l 1e-8)
	   (list (mid (cadr b) (cadr c))(/ (car a) 2)(list (cadr b) (cadr c) (cadr a)))
	   (progn 
	     (setq ja (angle pb pc))
	     (setq jb (angle pc pa))
	     (setq jc (angle pa pb))
	     (setq ppc (polar (mid pa pb) (+ (/ pi 2) jc) p))
	     (setq ppa (polar (mid pb pc) (+ (/ pi 2) ja) p))
	     (setq ppb (polar (mid pc pa) (+ (/ pi 2) jb) p))
	     (setq cen (inters ppc (mid pa pb) ppa (mid pb pc) nil))
	     (setq radius (distance cen pa))
	     (list cen radius (list pa pb pc))
	   )
	 )
       )
       ;;定義四點的最小圓圓心半徑，并返回三點坐標
       (defun 4pc (p1 p2 p3 ptmax / pts 3pt)
	 (setq pts (list (3pc p1 p2 p3) (3pc p1 p2 ptmax) (3pc p1 p3 ptmax) (3pc p2 p3 ptmax)))
	 (setq 3pt (vl-sort-i (mapcar 'cadr pts) '<))
	 (setq pts (list (nth (car 3pt) pts)  (nth (cadr 3pt) pts)
		         (nth (caddr 3pt) pts)(nth (cadddr 3pt) pts)))
	 (nth (vl-position t (mapcar '(lambda (x) (in2 (list p1 p2 p3 ptmax) (car x) (cadr x))) pts)) pts)
       )
       ;;定義求點集中離圓心最遠的點的函數--------
       (defun maxd-cir (ptl cen / distl)
	 (setq distl (mapcar '(lambda (x) (distance x cen)) ptl))
	 (nth (car (vl-sort-i distl '>)) ptl)
       )
       ;;開始遞歸運算----------------------------
       (setq p1	(car ptlist) p2 (cadr ptlist) p3 (caddr ptlist))
       (setq cen_r (3pc p1 p2 p3))
       (setq ptmax (maxd-cir ptlist (car cen_r)))
       (while (= nil (in1 ptmax (car cen_r) (cadr cen_r)))
         (setq cen_r (4pc p1 p2 p3 ptmax)) 
         (setq p1 (car (caddr cen_r)) p2 (cadr (caddr cen_r)) p3 (caddr (caddr cen_r)))
	 (setq ptmax (maxd-cir ptlist (car cen_r)))
       ) 
       (list (car cen_r) (cadr cen_r) ptmax)
      );;for progn
    );;  for t    
  );;    for cond 
);;      for defun
;;以下代碼來自曉東
;;定義取點函數----
(defun ssgetpoint (ss / i listpp a b c) 
  (setq i 0 listpp nil ) 
  (if ss 
    (repeat (sslength ss) 
      (setq a (ssname ss i)) 
      (setq b (entget a)) 
      (setq c (cdr (assoc 10 b))) 
      (setq listpp (cons c listpp)) 
      (setq i (1+ i))  
    ) 
  ) 
  listpp 
)
(defun ssgetpoint1 (ss / i listpp a b c) 
  (setq i 0 listpp nil ) 
  (if ss 
    (repeat (sslength ss) 
      (setq a (ssname ss i)) 
      (setq b (entget a))
      (setq ename (cdr (assoc 0 b)))
      (cond
	( (or nil (= ename  "POLYLINE") (= ename "LWPOLYLINE"))
	  (progn
	    (setq c (GetListOfPline a))
	    (setq listpp (append c listpp))
	  )
        )
	( (= ename "LINE")
	  (progn 
            (setq c (cdr (assoc 10 b)))
            (setq d (cdr (assoc 11 b)))
            (setq listpp (cons c listpp))
            (setq listpp (cons d listpp))
	  )
	)
	( (= ename "POINT")
	  (progn
	    (setq c (cdr (assoc 10 b)))
	    (setq listpp (cons c listpp))
	  )
	)
      )
      (setq i (1+ i))
    )	    
  )  
  listpp 
)
;; Get all nodes of the LWPolyline, Polyline.
(defun GetListOfPline (EntityName / SSE_Pline N newEntityName)
  (setq SSE_Pline (entget EntityName))
  (setq LastList nil)
  (if (= (cdr (assoc 0 SSE_Pline)) "LWPOLYLINE")
    (progn
      (setq LastList (LIST (LIST 0 0)))
      (setq N 0)
      (while (/= (nth N SSE_Pline) nil)
	(if (= (car (nth N SSE_Pline)) 10)
	  (setq	LastList (append LastList
				 (list (list (cadr (nth N SSE_Pline))
					     (caddr (nth N SSE_Pline))
					     0
				       )
				 )
			 )
	  )
	)
	(setq N (+ N 1))
      )
      (setq LastList (cdr LastList))
    )
  )
  (if (= (cdr (ASSOC 0 SSE_Pline)) "POLYLINE")
    (PROGN
      (setq LastList (list (list 0 0)))
      (setq newEntityName (entnext EntityName))
      (while (= (cdr (assoc 0 (entget newEntityName))) "VERTEX")
	(setq
	  LastList (append
		     LastList
		     (list
		       (list (cadr (assoc 10 (entget newEntityName)))
			     (caddr (assoc 10 (entget newEntityName)))
			     0
		       )
		     )
		   )
	)
	(setq newEntityName (entnext newEntityName))
      )
      (setq LastList (cdr LastList))
    )
  )
  (setq LastList LastList)
) ;_defun