椭圆论

Highflybird

收集了一些有关椭圆的程序和链接：

一、椭圆的基本要素与其创建和更新
椭圆的几个基本参数：      
DXF 组码中的几个数值      
10  椭圆中心        
11  长轴矢量        
40  长短比率        
41  起点角度        
42  终点角度        
210  法线矢量        
这些参数可以由entget得到，也可以由activeX方法得到。   
注: 起点角度和终点角度的计算，如果是不是椭圆弧，是全椭圆，则这个数值为：
    O 和2Pi，如果是椭圆弧，则如图1；     


    需要注意的是vla-get-StartAngle 和 vla-get-StartParameter 在一般情况
    下，是两个不同数值，它们所代表的角度的意义如上图。同样去理解:
    vla-get-EngAngle 和 vla-get-EndParameter. 但是:                     
    vla-get-StartParameter = vlax-curve-getStartParam               
    vla-get-EndParameter = vlax-curve-getEndParam                    
    vlax-curve-getParamAtPoint 其参数可按照上述方式去理解。   

[pcode=lisp,true]
;;;=====================================================================
;;; 功能: 画一个椭圆                                      
;;; 输入: 中心，长轴，短轴和旋转角度                                    
;;; 输出: 成功返回椭圆的图元名，否则返回nil                             
;;;=====================================================================
(defun Ent:Make_Ellipse (cen a b ang  / maj rat)
  (if (> b a)
    (setq maj (polar '(0 0 0) (+ ang (/ pi 2)) b)
   rat (/ a b 1.0)
    )
    (setq maj (polar '(0 0 0) ang a)
   rat (/ b a 1.0)
    )
  )
  (entmakeX
    (list
      '(0 . "ELLIPSE")
      '(100 . "AcDbEntity")
      '(100 . "AcDbEllipse")
      (cons 10 cen)
      (cons 11 maj)
      (cons 40 rat)
    )
  )
)
;;;=====================================================================
;;; 功能: 画一个空三维空间的椭圆                          
;;; 输入: 中心，长轴，短轴、旋转角度和法线矢量            
;;; 输出: 成功返回椭圆的图元名，否则返回nil               
;;;=====================================================================
(defun Ent:Make_Ellipse_3d (cen a b ang an1 an2 Normal / maj rat)
  (if (> b a)
    (setq maj (polar '(0 0 0) (+ ang (/ pi 2)) b)
   rat (/ a b 1.0)
    )
    (setq maj (polar '(0 0 0) ang a)
   rat (/ b a 1.0)
    )
  )
  (setq Normal (mat:unit Normal))
  (setq maj (trans maj Normal 0 T))
  (entmakeX
    (list
      '(0 . "ELLIPSE")
      '(100 . "AcDbEntity")
      '(100 . "AcDbEllipse")
      (cons 10 cen)
      (cons 11 maj)
      (cons 40 rat)
      (cons 41 an1)
      (cons 42 an2)
      (cons 210 Normal)
    )
  )
)
;;;=====================================================================
;;; 功能: 画一个椭圆弧                                
;;; 输入: 中心，长轴，短轴和旋转角度，起始角度，终点角度  
;;; 输出: 成功返回椭圆弧的图元名，否则返回nil            
;;;=====================================================================
(defun Ent:Make_EllipseArc (cen a b ang an1 an2 / maj rat)
  (if (> b a)
    (setq maj (polar '(0 0 0) (+ ang (/ pi 2)) b)
   rat (/ a b 1.0)
    )
    (setq maj (polar '(0 0 0) ang a)
   rat (/ b a 1.0)
    )
  )
  (entmakeX
    (list
      '(0 . "ELLIPSE")
      '(100 . "AcDbEntity")
      '(100 . "AcDbEllipse")
      (cons 10 cen)
      (cons 11 maj)
      (cons 40 rat)
      (cons 41 an1)
      (cons 42 an2)
    )
  )
)
;;;=====================================================================
;;; 功能: 从某一点获取椭圆的参数(纯数学法)                              
;;; 输入: 椭圆上的一点和构成椭圆的参数：中心，半长轴，半短轴和旋转角   
;;; 输出: 此点的参数值                                             
;;;=====================================================================
(defun ELL:GetParam (pt Cen a b ang / x y p)            
  (setq p (mapcar '- pt cen))
  (setq x (car p))
  (setq y (cadr p))
  (atan
    (* a (- (* y (cos ang))(* x (sin ang))))
    (* b (+ (* x (cos ang))(* y (sin ang))))
  )
)
;;;=====================================================================
;;; 功能: 更新平面椭圆                        
;;; 输入: 要更新的实体名，更新后的中心，半长轴，半短轴和旋转角         
;;; 输出: 成功返回该实体名，否则返回nil                  
;;;=====================================================================
(defun ELL:Update (ent cen a b ang / D C X R m n an)
  (if (> b a)
    (setq m b n a an (+ ang (/ pi 2)))
    (setq m a n b an ang)
  )
  (setq D (entget ent))
  (setq C (cons 10 cen))
  (setq X (cons 11 (polar '(0. 0. 0.) an m)))
  (setq R (cons 40 (/ n m 1.0)))
  (setq D (subst C (assoc 10 D) D))
  (setq D (subst X (assoc 11 D) D))
  (setq D (subst R (assoc 40 D) D))
  (entmod D)
  (entupd ent)
)
[/pcode]

Highflybird

二、椭圆的曲率、弧长、面积和离心率   

[pcode=lisp,true]
;;;*********************************************************************
;;; 二、椭圆的曲率、弧长、面积和离心率                  
;;;*********************************************************************
;;;=====================================================================
;;;功能: 获取椭圆上一点处的曲率和转弯半径                 
;;;参数: 椭圆实体和椭圆上的一点                           
;;;返回: 此处离心圆圆心、转弯半径及其曲率（离心率）      
;;;说明: 如果要在CAD中几何作图，可以参考此贴：   
;;;      http://bbs.mjtd.com/thread-62980-1-1.html   
;;;=====================================================================
(defun ELL:GetCurvature (en pt / dxf maj rat a b p par x y k v1 v2 rad cen)
  (setq dxf (entget en))
  (setq maj (cdr (assoc 11 dxf)))
  (setq rat (cdr (assoc 40 dxf)))
  (setq a   (distance '(0 0) maj))
  (setq b   (* a rat))
  (setq p   (vlax-curve-getclosestpointto en pt))
  (setq par (vlax-curve-getParamAtPoint en p))
  (setq v1  (vlax-curve-getFirstDeriv en par))
  (setq v2  (list (- (cadr v1)) (car v1) (caddr v1)))
  (setq x   (* a (cos par)))
  (setq y   (* b (sin par)))
  (setq k   (expt rat 4))
  (setq rad (/ (expt (+ (* y y) (* k x x)) 1.5) rat rat b b))
  (setq cen (polar p (angle '(0 0 0) v2) rad))
  (list cen rad (/ 1 rad))
)
;;;测试椭圆的曲率:
(defun c:GetCurvature (/ sel ent dxf pnt ret)
  (setq sel (nentselp "\n选取椭圆："))
  (if (and (setq ent (car sel))
    (setq dxf (entget ent))
    (= "ELLIPSE" (cdr (assoc 0 dxf)))
      )
    (progn
      (setq pnt (cadr sel))
      (setq ret (ELL:GetCurvature ent (trans pnt 1 0)))
      (princ ret)
      (entmakeX
(list
   '(0 . "CIRCLE")
   (cons 10 (car ret))
   (cons 40 (cadr ret))
)
      )
    )
  )
  (princ)
)
;;;=====================================================================
;;; 椭圆的离心率由椭圆的长短比率得出。(sqrt  (- 1 (* ratio ratio)))     
;;; 离心率越大圆就越扁，越小则越接近于圆   
;;;---------------------------------------------------------------------
;;; 功能: 获取椭圆的离心率，值越大形状越扁，越小越接近圆  
;;; 输入: 椭圆的图元名                                    
;;; 输出: 椭圆的离心率                                    
;;;=====================================================================
(defun ELL:Eccentricity (e / ratio)
  (setq ratio (cdr (assoc 40 (entget e))))
  (sqrt (- 1 (* ratio ratio)))
)
;;;=====================================================================
;;; 说明: 此函数为纯数学计算，意味可不用图元参与，采用了龙贝塔积分法。  
;;;       适用于某些特殊情况。一般可用vlax-curve函数获得。见样例。
;;; 功能: 获取椭圆弧的长度。               
;;; 输入: 椭圆的长半轴，短半轴，参数1，参数2和精度                     
;;; 输出: 椭圆弧的长度                                          
;;;=====================================================================
(defun ELL:Length (la lb p1 p2 eps / Func)
  (defun Func (x / cx ee)
    (setq cx (cos x))
    (setq ee (* (+ la lb) (- la lb)))
    (sqrt (- (* la la) (* ee cx cx)))
  )
  (Math:Romberg p1 p2 eps)
)
;;;=====================================================================
;;; 说明: 此函数为纯数学计算，意味可不用图元参与，速度居然不慢。        
;;;       适用于某些特殊情况。一般可用vlax-curve-getarea函数获得，也可以
;;;       用vla-get-area获得。测试见样例。                              
;;; 功能: 获取椭圆弧的面积。               
;;; 输入: 椭圆的长半轴，短半轴，参数1，参数2                           
;;; 输出: 椭圆弧的面积                                          
;;;=====================================================================
(defun ELL:Area (la lb A1 A2 / a k)
  (if (> A1 A2)
    (setq A (- (+ pi pi A2) A1))
    (setq A (- A2 A1))
  )
  (setq k (sin (* 0.5 A)))
  (setq k (* 1.333333333333333333333 la k k k))
  (* 0.5 la lb (- A (sin A)))
)
;;;测试椭圆弧长和面积样例
(defun c:TestLA (/ e o la lb par1 par2 l1 l2 s1 s2 s3)
  (setq e  (car (entsel "\n请选择椭圆弧: ")))
  (setq o  (vlax-ename->vla-object e))
  (setq la (vla-get-MajorRadius o))
  (setq lb (vla-get-MinorRadius o))
  (setq par1 (vlax-curve-getStartParam e))
  (setq par2 (vlax-curve-getEndParam e))
  (setq l1 (ell:length la lb par1 par2 1e-14))    ; 用数学计算获得椭圆弧长
  (setq l2 (vlax-curve-getDistAtParam e par2))      ; 用vlax-curve的函数获得椭圆弧长
  (princ "\n用vlax-curve函数计算的结果和积分法计算相差:")
  (princ (rtos (- l1 l2) 2 20))
  (misc:test 1001
    '((ell:length la lb par1 par2 1e-8)                                 ; 因为椭圆弧长没有数学公式，所以慢了
      (vlax-curve-getDistAtParam e par2)
    )
  )
  (setq s1 (vla-get-area o))
  (setq s2 (ELL:Area la lb par1 par2))
  (setq s3 (vlax-curve-getarea e))
  (princ "\n用vla-get-area函数计算的结果和数学法计算相差:")
  (princ (rtos (- s1 s2) 2 20))
  (misc:test 1001
    '((ELL:Area la lb par1 par2)     ; 居然用这个方法最快??
      (vla-get-area o)
      (vlax-curve-getarea e)
    )
  )
  (princ)
)
[/pcode]

Highflybird

三、空间椭圆的变换矩阵及其法线和平面投影

[pcode=lisp,true]
;;;*********************************************************************
;;; 三、空间椭圆的变换矩阵及其法线和平面投影                                
;;;*********************************************************************

;;;=====================================================================
;;; 功能: 获取椭圆的变换矩阵                                             
;;; 输入: 椭圆的图元名                                                  
;;; 输出: 椭圆的变换矩阵及逆矩阵                                       
;;; 说明: 由椭圆的长轴矢量，和短轴矢量以及法线矢量可以构成椭圆自身的变换
;;;       矩阵。如果椭圆位于图块内或者嵌套块内，关于其变换矩阵可以参考这
;;;       个帖子的讨论：http://bbs.mjtd.com/thread-93828-1-1.html.      
;;;       此帖中对椭圆的变换讨论的比较深入。                                   
;;;=====================================================================
(defun ELL:TransMatrix (e / Obj Cen DX DY DZ mat)
  (setq obj (vlax-ename->vla-object e))
  (setq Cen (vlax-get obj 'Center))                                        ; 中心点
  (setq DX  (Mat:unit (vlax-get obj 'MajorAxis)))                       ; OCS的X轴方向 长轴方向
  (setq DY  (Mat:unit (vlax-get obj 'MinorAxis)))                        ; OCS的Y轴方向 短轴方向
  (setq DZ  (Mat:unit (vlax-get obj 'Normal)))                                ; OCS的Z轴方向 法线方向
  (setq mat (list DX DY DZ))                                                    ; 椭圆的旋转矩阵
  (list
    (Mat:DispToMatrix (Mat:trp mat) cen)                                ; 椭圆的变换矩阵 = (trans Pt En 0)
    (Mat:DispToMatrix mat (mapcar '- (Mat:mxv mat cen)))                ; 椭圆的逆变换矩阵 = (trans Pt 0 En)
  )
)

;;;=====================================================================
;;; 功能: 获取椭圆的旋转矩阵                                             
;;; 输入: 椭圆的长轴矢量和椭圆的法线                                         
;;; 输出: 椭圆的旋转变换矩阵及逆矩阵                                    
;;;=====================================================================
(defun ELL:RotationMatrix (Maj Nrm / DX DY DZ mat)        
  (setq DX (Mat:unit Maj))                                                     ; OCS的X轴方向
  (setq DZ (Mat:unit Nrm))                                                ; OCS的Z轴方向
  (setq DY (Mat:Unit (Mat:vxv DZ DX)))                                        ; OCS的Y轴方向
  (setq mat (list DX DY DZ))                                                    ; 椭圆的旋转矩阵的逆矩阵
  (list (Mat:trp mat) mat)                                                 ; 椭圆的旋转矩阵和逆矩阵
)

;;;以下为测试程序
(defun c:testMatrix (/ sel ent mat obj doc)
  (if (setq sel (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq doc (vla-get-ActiveDocument (vlax-get-acad-object)))
      (vla-StartUndoMark doc)
      (setq ent (ssname sel 0))
      (setq mat (ELL:TransMatrix ent))
      (setq obj (vlax-ename->vla-object ent))
      (vla-transformby obj (vlax-tmatrix (cadr mat)))
      (command "Select" ent pause)
      (vla-transformby obj (vlax-tmatrix (car mat)))
      (vla-EndUndoMark doc)
      (princ)
    )
  )
)

;;;=====================================================================
;;; 功能: 修改一个椭圆的法线方向(不能直接用ActiveX方法)                  
;;; 输入: 椭圆的图元名，目标法线方向                                    
;;; 输出: 成功返回椭圆的DXF表，否则返回nil                              
;;;=====================================================================
(defun ELL:Put_Normal (ent Nrm / DXF Maj)
  (setq Nrm (mat:unit Nrm))
  (setq dxf (entget ent))
  (setq maj (cdr (assoc 11 dxf)))                                        ; 首先取得主轴方向
  (setq maj (trans maj 0 (cdr (assoc 210 dxf)) T))                        ; 主轴方向在OCS中的方向
  (setq maj (trans Maj Nrm 0 T))                                        ; 把OCS中的方向转化为目标方向
  (setq dxf (subst (cons 11 maj) (assoc 11 dxf) dxf))                        ; 椭圆的主轴方向
  (setq dxf (subst (cons 210 Nrm) (assoc 210 dxf) dxf))                        ; 椭圆的法线方向
  (entmod dxf)                                                                ; 更新图元
)

;;;=====================================================================
;;; highflybird  2012.11.1 第一版  2013.5.25 修改于深圳                 
;;; 功能: 给定一个椭圆或者圆，画出在平面上的投影（一般为椭圆）                 
;;; 输入: (椭圆、椭圆弧、圆弧、圆)实体E，平面的法线Normal和其一点Origin
;;; 输出: 投影椭圆(创建3d椭圆的几个参数)                                    
;;; 参考: http://bbs.mjtd.com/forum.php?mod=viewthread&tid=84527        
;;;       http://www.theswamp.org/index.php?topic=43031.0               
;;;=====================================================================
(defun ELL:Projection (E Normal Origin /
                       A1  AN1 AN2 ANG B1  C0  C1  C2  ENT OBJ PA0
                       PA1 PB0 PB1 PE0 PE1 PR1 PR2 PS0 PS1 RET)
  (setq C0  (cdr (assoc 10 (entget E))))                                ; 中心点(用vla-get-center可能出错？)
  (setq C0  (trans C0 e 0))                                                ; 此处很重要!!!
  (setq PS0 (vlax-curve-getStartPoint E))                                 ; 这个可以通过角度计算出来StartAngle
  (setq PE0 (vlax-curve-getEndPoint E))                                        ; 这个可以通过角度计算出来EndAngle
  (setq Pr1 (vlax-curve-getStartParam E))
  (setq Pr2 (vlax-curve-getEndParam E))
  (if (LINE:Colinearity C0 PS0 PE0)                                        ; 防止中心，起点，终点共线情况出现
    (setq PA0 (vlax-curve-getPointAtParam E (1+ pr1))                   ; 起点参数加一
          pB0 (vlax-curve-getpointatparam E (1- pr2))                        ; 终点参数加一
    )
    (setq pA0 (vlax-curve-getpointatparam E (* 0.5 (+ pr1 pr2)))        ; 取中点为第二点
          pB0 PE0                                                        ; 终点为第三点
    )
  )
  
  (mapcar
    (function (lambda (x y) (set x (trans y 0 Normal))))                ; 把这些点投影到平面上
    (quote (PA1 PB1 PS1 PE1 C1))
    (list PA0 PB0 PS0 PE0 C0)                                                ; 椭圆上取的三点、起始点和端点
  )
  (if (setq RET (ELL:C3P C1 PS1 PA1 PB1))                                ; 中心和三点画椭圆法
    (progn
      (setq A1  (cadr ret))                                                 ; 投影椭圆的长轴
      (setq B1  (caddr ret))                                                ; 投影椭圆的短轴
      (setq ang (cadddr ret))                                                ; 投影椭圆的旋转角
      (setq C2  (PLANE:Perpendicular_Foot c0 Origin Normal))                ; 椭圆的实际中心
      (if (vlax-curve-isClosed E)
        (list C2 A1 B1 ang 0 (+ pi pi) Normal)                                ; 返回创建3d椭圆的几个参数
        (progn
          (setq an1 (ELL:GetParam PS1 C1 A1 B1 ang))                        ; 投影椭圆的起始参数
          (setq an2 (ELL:GetParam PE1 C1 A1 B1 ang))                        ; 投影椭圆的终点参数
          (if (> (TRI:Det3P pS1 pA1 pB1) 0)                                ; 此处的确要作判断
            (list C2 a1 b1 ang an1 an2 Normal)                            ; 此时是逆时针
            (list C2 a1 b1 ang an2 an1 Normal)                                 ; 顺时针要交换起始角和终止角
          )
        )
      )
    )
  )
)

;;;用于投影椭圆的测试:
(defun C:Projection (/ sel ent normal origin i ret)
  (setq sel (ssget '((0 . "ELLIPSE,ARC,CIRCLE"))))
  (setq normal (trans '(0 0 1) 1 0 T))
  (setq origin (getvar 'ucsorg))
  (if sel
    (repeat (setq i (sslength sel))
      (setq ent (ssname sel (setq i (1- i))))
      (if (setq ret (ELL:Projection ent normal origin))
        (apply 'Ent:Make_ELLIPSE_3D ret)
      )
    )
  )
)
[/pcode]

Highflybird

四、椭圆的展点和椭圆锥的展开曲线

[pcode=lisp,true]
;;;*********************************************************************
;;; 四、椭圆的展点和椭圆锥的展开曲线                                       
;;;*********************************************************************

;;;=====================================================================
;;; 功能: 椭圆展点。                                                         
;;; 输入: 无                                                            
;;; 输出: 无                                                               
;;;=====================================================================
(defun C:ExpandEllipse(/ A B D ENT F I L LEN N OBJ P PAR S SEG)
  (initget 15)
  (setq a (getdist "\n请输入椭圆的长轴:"))        
  (initget 15)
  (setq b (getdist "\n请输入椭圆的短轴:"))
  (initget 7)
  (setq n (getint "\n请输入等分数目:"))
  (setq ent (ENT:MAKE_ELLIPSE '(0 0 0) a b 0))
  (setq obj (vlax-ename->vla-object ent))
  (setq par (vlax-curve-getendparam obj))
  (setq len (vlax-curve-getdistatparam obj par))
  (setq seg (/ len n))
  (setq l nil)
  (setq d 0)
  (setq i 0)
  (repeat n
    (setq d (* i seg))
    (setq p (vlax-curve-getPointAtDist obj d))
    ;(Ent:Make_Point p)
    (setq l (cons p L))
    (setq i (1+ i))
  )
  (setq l (reverse l))
  (vla-erase obj)
  (if (setq f (getfiled "输入坐标文件名" "D:/" "TXT" 1))
    (progn
      (setq f (open f "W"))
      (foreach p l
        (setq s (strcat (rtos (car p) 2 8) "\t" (rtos (cadr p) 2 8)))
        (write-line s f)
      )
      (close f)
    )
  )
)

;;;=====================================================================
;;; 功能: 展开椭圆锥体。                                                
;;; 输入: 椭圆锥底椭圆的长轴a，短轴b，椭圆锥的高h和分弧精度n            
;;; 输出: 展开图的曲线的坐标点(可据此形成展开面)                        
;;;=====================================================================
(defun ELL:ExpandCone (a b h n / Pts Ang D Param DivAng P0 P1 R0 R1 X Y)
  (setq DivAng (/ pi n 0.5))                                                ;等分角度
  (setq Param 0)                                                        ;开始参数为0
  (setq p0 (list a 0))                                                        ;从椭圆最右边的点开始
  (setq r0 (distance '(0 0 0) (list a 0 h)))                                ;开始的锥顶最右边点的距离
  (setq Pts (list (list r0 0) '(0 0)))                                        ;把最开始的两点加入到数据表
  (setq Ang 0)                                                                ;开始角度也为0
  (repeat n
    (setq Param  (+ Param DivAng))
    (setq x   (* a (cos Param)))                                        ;椭圆上的x坐标
    (setq y   (* b (sin Param)))                                        ;椭圆上的y坐标
    (setq p1  (list x y))                                                ;点的坐标
    (setq r1  (distance '(0 0 0) (list x y h)))                                ;展开面的锥线长度
    (setq d   (distance p0 p1))                                                ;椭圆上的上一点到这点的距离
    (setq Ang (+ Ang (car (TRI:CosinesLaw d r0 r1))))                        ;坐标角度
    (setq Pts (cons (polar '(0 0) Ang r1) Pts))                                ;得到了坐标，并把它加入到表中
    (setq p0 p1)                                                        ;用新的坐标位置替换旧位置
    (setq r0 r1)                                                        ;用新的锥线长度替换旧长度
  )
  (reverse Pts)                                                                ;逆转表，得到正序的数据
)

;;; 测试椭圆锥的展开
(defun c:ExpandCone(/ h N e d r c m a b l e o)
  (initget 1)
  (setq h (getdist "\n输入椭圆锥高: "))
  (initget 7)
  (setq N (getint "\n等分数:"))
  (prompt "\n选择椭圆:")
  (if (setq s (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq e (ssname s 0))
      (setq d (entget e))
      (setq r (cdr (assoc 40 d)))
      (setq c (cdr (assoc 10 d)))
      (setq m (cdr (assoc 11 d)))
      (setq a (distance '(0 0 0) m))
      (setq b (* r a))
      (setq l (ELL:ExpandCone a b h n))
      (setq o (vlax-ename->vla-object (Ent:Make_LWPoly l T)))
      (vla-move o (vlax-3d-point '(0 0 0)) (vlax-3d-point c))
      (princ)
    )
  )
)[/pcode]

Highflybird

五、椭圆的相切和相交

[pcode=lisp,true]
;;;*********************************************************************
;;; 五、椭圆的相切和相交                                                
;;;*********************************************************************

;;;=====================================================================
;;; 功能: 数学法求椭圆上一点的切线矢量(相当于vlax-curve-GetFirstDeriv)  
;;; 输入: 这点的参数Param,椭圆中心C,长轴方向M，椭圆比率R和椭圆的法线N   
;;; 输出: 这点的切线矢量                                                
;;;=====================================================================
(defun ELL:PointTangentOn (Param C M R N / a b v)
  (setq a (distance '(0 0 0) M))
  (setq b (* R a))
  (setq v (list (* (- a) (sin Param)) (* b (cos Param)) 0))
  (Mat:mxv (car (ELL:RotationMatrix M N)) v)
)

;;; 切线矢量函数测试
(defun C:PointTangentOn ( / CEN DXF ENT MAJ NRM PAR PNT RAT RET SEL)
  (if (setq sel (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq ent (ssname sel 0))
      (setq dxf (entget ent))
      (setq cen (cdr (assoc 10 dxf)))
      (setq maj (cdr (assoc 11 dxf)))
      (setq rat (cdr (assoc 40 dxf)))
      (setq Nrm (cdr (assoc 210 DXF)))
      (while (setq Pnt (getpoint "\n选取点: "))
        (setq Pnt (trans Pnt 1 0))
        (setq pnt (vlax-curve-getclosestpointto ent pnt))
        (setq par (vlax-curve-getParamAtPoint ent pnt))
        (setq ret (ELL:PointTangentOn par cen maj rat Nrm))
        (ent:make_line Pnt (mapcar '+ Pnt ret))
      )
      (princ)
    )
  )
)

;;; 上面是也可以由某个点的参数通过数学方法算出。
;;; 一般来说可以由vlax-curve-getFirstDeriv 算出，面是测试如下:
(defun C:TestDeriv (/ e o p param f1 f2)
  (if (setq e (car (entsel)))
    (progn
      (setq o (vlax-ename->vla-object e))
      (while (setq p (getpoint "\n测试点:"))
        (setq p (trans p 1 0))
        (setq p (vlax-curve-getclosestpointto e p))
        (setq param (vlax-curve-getParamAtPoint e p))
        (setq f1 (vlax-curve-getFirstDeriv e param))
        (setq f2 (vlax-curve-getSecondDeriv e param))
        (Ent:make_line p (mapcar '+ p f1))
        (Ent:make_line p (mapcar '+ p f2))
      )
    )
  )
)

;;;=====================================================================
;;; 功能: 求椭圆(平面)外一点到椭圆的切线                                
;;; 输入: 椭圆外一点Point,椭圆中心C,椭圆半长轴a,半短轴b,和椭圆的旋转角  
;;; 输出: 此点到与椭圆相切的切点(一个，两个或者nil)                        
;;;=====================================================================
(defun ELL:PointTangentTo (Point C a b ang / AA BB K1 K2 KK P X Y Z XX YY)
  (setq Point (mapcar '- Point C))
  (setq z (polar '(0 0 0) ang a))
  (setq P (trans point 0 z T))
  (setq x (caddr p))
  (setq y (car p))
  (if (equal y 0 1e-8)
    (if (> (abs x) a)
      (progn
        (setq kk (/ a x))
        (setq xx (* a kk))
        (setq yy (* b (sqrt (- 1 (* kk kk)))))
        (list
          (mapcar '+ C (trans (list yy 0 xx) z 0 T))
          (mapcar '+ C (trans (list (- yy) 0 xx) z 0 T))
        )
      )
    )
    (progn
      (setq aa (* a a))
      (setq bb (* b b))
      (setq xx (* x x))
      (setq yy (* y y))
      (setq k1 (+ (* aa yy) (* bb xx)))
      (setq k2 (- k1 (* aa bb)))
      (if (> k2 0)
        (progn
          (defun GetTan (A B C bb K1 X Y Z KK / d1 k4 k5 s1 x1 y1 p1)
            (setq d1 (+ (* bb x) kk))
            (setq k4 (/ (* d1 x) k1))
            (setq k5 (/ (* d1 a) k1))
            (setq s1 (atan (/ (* b (- 1 k4)) y) k5))
            (setq x1 (* a (cos s1)))
            (setq y1 (* b (sin s1)))
            (setq P1 (trans (list y1 0 x1) z 0 T))
            (setq p1 (mapcar '+ P1 C))
          )
          (setq kk (* (sqrt k2) (abs y)))
          (list (GetTan A B C bb K1 X Y Z kk)
                (GetTan A B C bb K1 X Y Z (- kk))
          )
        )
      )
    )
  )
)

;;; 椭圆的切线函数测试
(defun C:PointTangentTo (/ A ANG B C DXF ENT MAJ PNT RET SEL)
  (if (setq sel (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq ent (ssname sel 0))
      (setq dxf (entget ent))
      (setq maj (cdr (assoc 11 dxf)))
      (setq a   (distance '(0 0 0) maj))
      (setq b   (* a (cdr (assoc 40 dxf))))
      (setq c   (cdr (assoc 10 dxf)))
      (setq ang (angle '(0 0 0) maj))
      (while (setq Pnt (getpoint "\n选取点: "))
        (setq Pnt (trans Pnt 1 0))
        (setq ret (ELL:PointTangentTo Pnt C a b ang))
        (foreach p ret
          (Ent:Make_line p Pnt)
        )
      )
    )
  )
)

;;;=====================================================================
;;; 关于椭圆的共切线问题请参考：                                       
;;; http://bbs.mjtd.com/forum.php?mod=viewthread&tid=82900               
;;; 另外可参考我的附件: ellipse-tan-solve1.LSP                                 
;;;=====================================================================

;;;=====================================================================
;;; 功能: 求椭圆与直线的交点(计算几何方式1)                             
;;; 输入: 椭圆的长轴a,短轴b，两点P1,P2.                                 
;;; 输出: 直线P1P2与椭圆的交点。                                       
;;; 参考: http://bbs.mjtd.com/thread-62003-2-1.html                        
;;;       RootOf((M^2*a^2+N^2*b^2)*_Z^2-M^2*a^2+L^2+2*L*N*b*_Z)*b      
;;;=====================================================================
(Defun ELL:Inters_Ellipse_Line (a b P1 P2 / D E K L M N P PA PB X1 X2 Y1 Y2)
  (setq k  (/ b (float a)))
  (setq X1 (car P1))
  (setq X2 (car P2))
  (setq Y1 (/ (cadr P1) k))
  (setq Y2 (/ (cadr P2) k))
  (Setq M  (- Y1 Y2))                                                   ;直线方程系数M
  (Setq N  (- X2 X1))                                                    ;直线方程系数N
  (setq D  (/ (- (* Y2 X1) (* Y1 X2)) (sqrt (+ (* M M) (* N N)))))        ;垂距
  (setq e  (angle (list x1 y1) (list x2 y2)))                                 ;直线的与X轴线的交角
  (setq p  (polar '(0 0) (- e (* pi 0.5)) d))                                ;圆心到直线的垂足
  (setq d  (abs d))
  (if (equal d a 1e-8)                                                         ;如果垂距等于半径
    (list (list (car p) (* k (cadr p))))                                ;相切
    (if (< d a)                                                               ;如果垂距小于半径
      (progn
        (setq L  (sqrt (* (+ a d)(- a d))))                                  ;半弦长
        (setq Pa (polar p e (- L)))
        (setq Pb (polar p e L))
        (setq pa (list (car Pa) (* k (cadr Pa))))
        (setq pb (list (car Pb) (* k (cadr Pb))))
        (list pa pb)                                                          ;有两个交点
      )
    )
  )
)

;;;=====================================================================
;;; 功能: 求椭圆与直线的交点(计算几何方式2,似乎比方式1慢)               
;;; 输入: 椭圆的长轴a,短轴b，两点P1,P2.                                 
;;; 输出: 直线P1P2与椭圆的交点。                                       
;;;=====================================================================
(Defun ELL:Inters_Ellipse_Line_1 (a b P1 P2 / K1 K2 K3 M N L MA NB PS X X1 X2 Y Y1 Y2)
  (Setq X1 (Car  P1))
  (Setq Y1 (Cadr P1))
  (Setq X2 (Car  P2))
  (Setq Y2 (Cadr P2))
  (Setq M  (float (- Y1 Y2)))                                            ;直线方程系数M
  (Setq N  (float (- X2 X1)))                                            ;直线方程系数N
  (Setq L  (float (- (* Y2 X1) (* Y1 X2))))                              ;直线方程系数L
  (setq Nb (* N b))
  (if (equal M 0 1e-8)                                                        ;水平直线
    (cond
      ( (equal (setq K1 (* (+ Nb L) (- Nb L))) 0 1e-8)
        (list (list 0 (- (/ L N))))
      )
      ( (> K1 0)
        (setq x (/ (* a (sqrt k1)) Nb))
        (setq y (- (/ L N)))
        (list (list x y) (list (- x) y))
      )
    )
    (progn
      (setq Ma (* M a))
      (setq k1 (+ (* Ma Ma) (* Nb Nb)))
      (setq k2 (* 2 L Nb))
      (setq k3 (- (* L L) (* Ma Ma)))
      (foreach e (Math:Quadratic_Equation_1 k1 k2 k3)
        (setq x (/ (+ L (* N e b)) (- M)))
        (setq y (* e b))
        (setq ps (cons (list x y) ps))
      )
    )
  )
)

;;;=====================================================================
;;; 功能: 求平面椭圆与直线的交点                                       
;;; 输入: 椭圆的中心Center,长轴方向major,椭率ratio，两点P1,P2.         
;;; 输出: 直线P1P2与椭圆的交点。                                       
;;;=====================================================================
(Defun ELL:Inters_Ellipse_Line_2D (Center major Ratio P1 P2 / v1 v2 a b s)
  (setq v1 (mapcar '- p1 center))
  (setq v2 (mapcar '- p2 center))
  (setq v1 (trans v1 0 major T))
  (setq v2 (trans v2 0 major T))
  (setq v1 (list (caddr v1) (car v1)))
  (setq v2 (list (caddr v2) (car v2)))
  (setq a  (distance '(0 0 0) major))
  (setq b  (* ratio a))
  (foreach p (ELL:Inters_Ellipse_Line_1 a b v1 v2)
    (setq p (trans (list (cadr p) 0 (car p)) major 0 T))
    (setq p (mapcar '+ center p))
    (setq s (cons p s))
  )
)

;;;=====================================================================
;;; 功能: 求空间椭圆与直线的交点（注意，这个交点为在椭圆平面上的交点）  
;;; 输入: 椭圆的中心Center,长轴方向major,椭率ratio,法线Normal,两点P1,P2.
;;; 输出: 直线P1P2与椭圆的交点。                                       
;;;=====================================================================
(Defun ELL:Inters_Ellipse_Line_3D (Center major Ratio Normal P1 P2 / A B M1 M2 S V1 V2)
  (setq v1 (mapcar '- p1 center))
  (setq v2 (mapcar '- p2 center))
  (setq m1 (ELL:RotationMatrix major Normal))
  (setq m2 (cadr m1))
  (setq m1 (car m1))
  (setq v1 (mat:mxv m2 v1))
  (setq v2 (mat:mxv m2 v2))
  (setq a  (distance '(0 0 0) major))
  (setq b  (* ratio a))
  (foreach p (ELL:Inters_Ellipse_Line_1 a b v1 v2)
    (setq p (mat:mxv m1 (append p '(0))))
    (setq p (mapcar '+ center p))
    (setq s (cons p s))
  )
)

;;; 椭圆的与直线相交函数的测试
(defun C:iel (/ p1 p2 sel ent dxf cen maj rat nrm a b ret)
  (initget 1)
  (setq P1 (getpoint "\n选取点1: "))
  (initget 1)
  (setq P2 (getpoint P1 "\n选取点2: "))
  (setq p1 (trans p1 1 0))
  (setq p2 (trans p2 1 0))
  (Ent:Make_Line p1 p2)
  (if (setq sel (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq ent (ssname sel 0))
      (setq dxf (entget ent))
      (setq cen (cdr (assoc 10 dxf)))
      (setq maj (cdr (assoc 11 dxf)))
      (setq rat (cdr (assoc 40 dxf)))
      (setq Nrm (cdr (assoc 210 DXF)))
      (setq a   (distance '(0 0 0) maj))
      (setq b   (* rat a))
      (MISC:test 100
        '((ELL:Inters_Ellipse_Line_2d cen maj rat P1 P2)
          (ELL:Inters_Ellipse_Line_3d cen maj rat Nrm P1 P2)))
      ;(setq ret (ELL:Inters_Ellipse_Line_2d cen maj rat P1 P2))
      (setq ret (ELL:Inters_Ellipse_Line_3d cen maj rat Nrm P1 P2))
      (and ret (mapcar 'Ent:Make_Point ret))
    )
  )
)
[/pcode]
            

Highflybird

六、与椭圆的最近点和最远点
[pcode=lisp,true]
;;;*********************************************************************
;;; 六、与椭圆的最近点和最远点                                               
;;;*********************************************************************

;;;=====================================================================
;;; 求平面上一点到椭圆的距离的极值（包括最小值和最大值）方式1           
;;; 输入：椭圆的长轴a，短轴b和平面的一点P                              
;;; 输出：椭圆的跟平面的P点距离为极值的点                              
;;;=====================================================================
(defun ELL:GetExtremumDist (a b P / am cc K1 K2 K3 m n v x y ps)
  (setq m (car p))
  (setq n (cadr p))
  (setq cc (* (+ a b) (- a b)))
  (setq k1 (* n b))
  (setq am (* a m))
  (setq k2 (+ am cc))
  (setq k2 (+ k2 k2))
  (setq k3 (- am cc))
  (setq k3 (+ k3 k3))
  (foreach s (Math:Quartic_Equation k1 k2 0 k3 (- k1))
    (if        (equal (cadr s) 0 1e-8)
      (progn
        (setq s (car s))
        (setq v (* s s))
        (setq x (/ (* a (- 1 v)) (1+ v)))
        (setq y (/ (* b (+ s s)) (1+ v)))
        (setq ps (cons (list x y) ps))
      )
    )
  )
  (if (equal k1 0 1e-8)
    (setq ps (cons (list (- a) 0) ps))
    (reverse ps)
  )
)

;;;=====================================================================
;;; 求平面上一点到椭圆的距离的极值（包括最小值和最大值）方式2           
;;; 输入：椭圆的长轴a，短轴b和平面的一点P                              
;;; 输出：椭圆的跟平面的P点距离为极值的点                              
;;;=====================================================================
(defun ELL:GetExtremumDist1 (a b P / c cc e ee f K1 K2 K3 K4 K5 m n x y ps)
  (setq m (car p))
  (setq n (cadr p))
  (if (equal m 0 1e-10)
    (list (list 0 (- b)) (list 0 b))
    (progn
      (setq cc (* (+ a b) (- a b)))
      (setq c  (sqrt cc))                               
      (setq e  (/ c a))
      (setq f  (/ b a))
      (setq ee (* e e))
      (setq k1 (* ee ee))
      (setq k2 (* -2 m ee))
      (setq k3 (- (+ (* n n f f) (* m m)) (* ee cc)))
      (setq k4 (* 2 m cc))
      (setq k5 (- (* a a m m)))
      (foreach s (Math:Quartic_Equation k1 k2 k3 k4 k5)
        (if (equal (cadr s) 0 1e-6)
          (progn
            (setq x (car s))
            (setq y (/ (* n b b) (- (/ (* m a a) x) cc)))
            (setq ps (cons (list x y) ps))
          )
        )
      )
      ps
    )
  )
)

;;;=====================================================================
;;; 功能: 求空间上一点到空间椭圆的距离的极值（包括最小值和最大值）      
;;; 输入: 椭圆的中心Center,长轴矢量major,椭率ratio,法线Normal,空间一点P.
;;; 输出: 椭圆的跟空间的P点距离为极值的点                              
;;;=====================================================================
(defun ELL:GetExtremumDist_3D (Center major Ratio Normal P / A B M1 M2 S V1)
  (setq v1 (mapcar '- p center))
  (setq m1 (ELL:RotationMatrix major Normal))
  (setq m2 (cadr m1))
  (setq m1 (car m1))
  (setq v1 (mat:mxv m2 v1))
  (setq a  (distance '(0 0 0) major))
  (setq b  (* ratio a))
  (foreach n (ELL:GetExtremumDist a b v1)
    (setq n (mat:mxv m1 (append n '(0))))
    (setq n (mapcar '+ center n))
    (setq s (cons n s))
  )
)

;;; 测试求极值距离程序
(defun C:GetExtremumDist (/ p sel ent dxf cen maj rat nrm a b)
  (initget 1)
  (setq P (getpoint "\n选取点1: "))
  (setq p (trans p 1 0))
  (Ent:Make_Point p)
  (if (setq sel (ssget ":S" '((0 . "ELLIPSE"))))
    (progn
      (setq ent (ssname sel 0))
      (setq dxf (entget ent))
      (setq cen (cdr (assoc 10 dxf)))
      (setq maj (cdr (assoc 11 dxf)))
      (setq rat (cdr (assoc 40 dxf)))
      (setq Nrm (cdr (assoc 210 DXF)))
      (setq a   (distance '(0 0 0) maj))
      (setq b   (* rat a))
;;;      (misc:test
;;;        1001
;;;        '((ELL:GetExtremumDist1 a b P)
;;;          (ELL:GetExtremumDist a b P)
;;;          (vlax-curve-getclosestpointto ent p)
;;;        )
;;;      )
      ;(mapcar 'Ent:Make_Point (ELL:GetExtremumDist1 a b P))
      (mapcar 'Ent:Make_Point (ELL:GetExtremumDist_3D cen maj rat Nrm P))
    )
  )
)[/pcode]

Highflybird

以后把完整的代码和测试样例附上来。
不要走开，后续更精彩。。。

Lispboy

也算沙发吧，大作又横空出世，拜读中...

牢固

精品文章！认真学习ing...

zyhandw

见了楼主的帖子，总是忍不住要收藏:lol

liuyun242

很给力，要收藏，多多学习

Lisphk

拜读，恶补知识去了。

ysq101

前排支持！！！！！！！这个得多高的学问啊！！学习ING

sicky111

楼主的贴子，一定要学习。

GTJ116600

HighFlyBird版主的精品帖子，收藏