Lee Mac 发表于 2022-7-5 18:04:54

 
谢谢Stefan!
您的代码产生了良好的结果。

Tyke 发表于 2022-7-5 18:08:04

@李
 
您的代码适用于三维点还是仅适用于二维点?我在一批3D点上试用了它,结果很奇怪,而Stefan的结果似乎很管用。
 
经过进一步检查,Stefan的代码生成了一条穿过三维点的连接线,但它只是一个二维实体。
 
有没有可能仅仅通过点簇创建一条与XLINE或光线平行的线?

Lee Mac 发表于 2022-7-5 18:11:38

 
嗨,小家伙,
 
由于我的代码使用简单线性回归,它将只适用于2变量系统,即2D点;当然,更高维度的线性回归是可能的,但要复杂一些。

Tyke 发表于 2022-7-5 18:13:47

 
谢谢李,我想可能是这样的。我有来自3D点云的点子集,大约200个点,需要获得点的最佳3D拟合线。有什么想法可以实现吗?

Tyke 发表于 2022-7-5 18:17:54

我的LISP确实很生疏,所以我在VBA中实现了它,并在点簇中画了一条线。我现在打算做进一步的线性回归,但使用点的X和Z坐标来获得通过点的线的梯度。三维线性回归所需的数学知识超出了我的能力。

ymg3 发表于 2022-7-5 18:21:03

在原始问题的背景下,我认为Stefan的解决方案优于
对直线应用正交回归,并最小化点到直线的距离。
 
Lee的解决方案(杰出的实现)是另一个beast,它将沿Y坐标到直线的距离缩小,当X被视为独立变量时,它是完全正确的。
 
该问题的解决方案首次发表在R.J.ADCOCK出版的《分析员》第5卷第2期(1878年3月)第53-54页。
这是本文的链接:最小二乘法中的一个问题
 
它也称为戴明回归或正交回归。
 
ymg公司
 
这是一张来自维基百科的图片,展示了斯特凡的解决方案:
http://upload.wikimedia.org/wikipedia/commons/thumb/8/81/Total_least_squares.svg/220px-Total_least_squares.svg.png
 
以下是李正在最小化的内容:
 

Lee Mac 发表于 2022-7-5 18:24:11

非常有趣且信息丰富的比较ymg-感谢您提供的信息!

ymg3 发表于 2022-7-5 18:26:18

泰克,
 
这是3D中的正交回归。
 
可能已经太晚了,但我正忙于其他事情。
 
代码还可以,但我仍然需要添加一个残差表来分析结果。
 
公式推导见:Jean Jacquelin的三维线性回归
 
ymg公司
 

;; fl3d, Orthogonal Regression in 3D   by ymg                                  ;
;;                                                                           ;
;; For Derivation of Formulas See Following Paper:                           ;
;;      3-D LINEAR REGRESSIONby Jean Jacquelin                           ;
;; http://www.scribd.com/doc/31477970/Regressions-et-trajectoires-3D         ;
;;                                                                           ;

(defun c:fl3d (/ *acaddoc* a avg b c0 c1 c2 cos2a cosa errl go k00 k01 k10 k11 k12 k22
               n p phi pl q r rho s2m sin2a sina ss sum sxx sxy sxz syy syz szz u v w xm ym zm)
(vl-load-com)
(defun *error* (msg)
           (mapcar 'eval errl)
(if (and msg (not (wcmatch (strcase msg) "*BREAK*,*CANCEL*,*EXIT*")))
          (princ (strcat "\nError: " msg))
       )
(and *AcadDoc* (vla-endundomark *AcadDoc*))
       (princ)
)

(defun acos (z)(atan (sqrt (- 1.0 (* z z))) z))

(setq errl '("CMDECHO" "DIMZIN")
      errl (mapcar (function (lambda (a) (list 'setvar a (getvar a)))) errl)
)
   
   
(or *AcadDoc*
   (setq *AcadDoc* (vla-get-activedocument (vlax-get-acad-object)))
)

(princ (strcat "\n Select Points"))
(if (setqn 0
   pl nil
   ss (ssget '((0 . "POINT")))
            go (> (sslength ss) 2)
      )
   (progn
      (vla-startundomark *AcadDoc*)          
      (setvar 'CMDECHO 0)
      (setvar 'DIMZIN 0)
(repeat (sslength ss)
         (setq pl (cons (cdr (assoc 10 (entget (ssname ss n)))) pl)
                n (1+ n)
         )
)
      (setqavg(mapcar '/(apply 'mapcar (cons '+ pl)) (list n n n))
               ;dif(mapcar '(lambda (a) (mapcar '- a avg)) pl)
               sum(apply
                     'mapcar
                        (cons '+
                           (mapcar
                              '(lambda (a)
                                     (list (* (car a) (car a))       ;XX
                                           (* (car a) (cadr a))      ;XY
                                           (* (cadr a) (cadr a))   ;YY
                                           (* (car a) (caddr a))   ;XZ
                                           (* (caddr a) (caddr a))   ;ZZ
                                           (* (cadr a) (caddr a))    ;YZ
                                     )
                                 )
                                 pl
                           )
                        )
                  )
               sum(mapcar '/ sum (list n n n n n n))
               ;; avg contains:(Xm Ym Zm)                                  ;
               ;; sum contains:(sXX sXY sYY sXZ sZZ sYZ) Divided by n      ;
                  Xm (car avg) Ym (cadr avg) Zm (caddr avg)
               Sxx (- (car sum)(* Xm Xm))
               Sxy (- (cadr sum)(* Xm Ym))
               Syy (- (caddr sum)(* Ym Ym))
               Sxz (- (cadddr sum)(* Xm Zm))
               Szz (- (car (cddddr sum))(* Zm Zm))
               Syz (- (cadr (cddddr sum))(* Ym Zm))
                   a (/ (atan (/ (* Sxy 2.)(- Sxx Syy))) 2.)         
                Cosa (cos a)      Sina (sin a)
               Cos2a (* Cosa Cosa) Sin2a (* Sina Sina)            
               K11 (+ (* (+ Syy Szz) Cos2a) (* (+ Sxx Szz) Sin2a) (* -2. Sxy Cosa Sina))
               K22 (+ (* (+ Syy Szz) Sin2a) (* (+ Sxx Szz) Cos2a) (*2. Sxy Cosa Sina))      
               K12 (+ (* (* Sxy -1.)(- Cos2a Sin2a))(* (- Sxx Syy) Cosa Sina))
               K10 (+ (* Sxz Cosa)(* Syz Sina))
               K01 (+ (* Sxz Sina -1.)(* Syz Cosa))
               K00 (+ Sxx Syy)         
                  c2 (* (+ K00 K11 K22) -1.)
                  c1 (+ (* K00 K11)(* K00 K22)(* K11 K22)(* K01 K01 -1.)(* K10 K10 -1.))
                  c0 (+ (* K01 K01 K11)(* K10 K10 K22)(* -1. K00 K11 K22))
                   p (- c1 (/ (* c2 c2) 3.))
                   q (+ (/ (* 2.0 c2 c2 c2) 27.)(/ (* -1. c1 c2) 3.) c0)
                   R (+ (/ (* q q) 4.)(/ (* p p p) 27.))         
      )
      (if (minusp R)
         (progn
               (setq rho (sqrt (/ (* -1. p p p) 27.))
                     phi (acos (/ q (* -2. rho)))
                     s2m (min (+ (/ c2 -3.)(* 2. (expt rho (/ 1. 3.))(cos (/ phi 3.))))
                              (+ (/ c2 -3.)(* 2. (expt rho (/ 1. 3.))(cos (/ (+ phi pi pi) 3.))))
                              (+ (/ c2 -3.)(* 2. (expt rho (/ 1. 3.))(cos (/ (+ phi pi pi pi pi) 3.))))
                         )
               )      
         )
         (setq s2m (+ (/ c2 -3.)(expt (+ (/ q -2.)(sqrt R)) (/ 1. 3.))(expt (- (/ q -2.)(sqrt R)) (/ 1. 3.))))
      )

      (setq a (+ (* (/ (* -1. K10)(- K11 s2m)) Cosa) (* (/ K01 (- K22 s2m)) Sina))
            b (+ (* (/ (* -1. K10)(- K11 s2m)) Sina) (* (/ (* -1. K01) (- K22 s2m)) Cosa))
            u (* (/ 1. (+ 1. (* a a) (* b b))) (+ (* (+ 1. (* b b)) Xm) (* -1. a b Ym) (* a Zm)))
            v (* (/ 1. (+ 1. (* a a) (* b b))) (+ (* -1 a b Xm) (* (+ 1 (* a a)) Ym) (* b Zm)))
            w (* (/ 1. (+ 1. (* a a) (* b b))) (+ (* a Xm) (* b Ym) (* (+ (* a a) (* b b)) Zm)))
      )      
      (entmakex (list '(0 . "LINE") (cons 10 avg) (cons 11 (list u v w))))
      
      
   )
)
(*error* nil)
)

(princ "\nFit a 3D Line to Selection Set of Points, Type FL3D to run")


prodromosm 发表于 2022-7-5 18:28:26

我使用这个lisp代码
 
5

ymg3 发表于 2022-7-5 18:32:33

Prodromosm,
 
我知道这一个,它来自Cadalyst,做普通最小二乘法。
 
这并不理想,因为您希望将一条线拟合到您测量的点。
斯特凡的那部更好。
 
新的用于三维拟合,也适用于二维线,
然而,可能对舍入更加敏感。我没有检查。
 
ymg公司
页: 1 [2]
查看完整版本: 从f建立一条平均线