三角网插值

ndgxy · 发表于 2017-11-1 10:41:39

那太慢了呀，如果三角形上万个，甚至几十万个，那效率太低了呀
有没其他办法呀？

求指点

newer · 发表于 2017-11-1 10:52:58

本帖最后由 newer 于 2017-11-1 10:54 编辑

ndgxy 发表于 2017-11-1 10:41

一点也不慢，你试下 xdrx_points_delaunay创建三角网，几十万的，也不到1秒就创建出来了，这还是创建实体不是搜索，搜索就更快了。
如果你只想通过点就得到，没有别的办法，除非你交互的时候不只给点，还点选下三角网的3D FACE，一点也不慢。遍历几十万实体的数据库，几毫秒。ACAD很多操作都是遍历数据库。如果你想只给一个点，还想改进想快点，那么就把屏幕内看见的三角网选出来，只遍历屏幕内的。

fools · 发表于 2017-11-1 12:56:11

本帖最后由 fools 于 2017-11-1 13:11 编辑

加了注释，帮助你理解

;;                                                                         ;
;; triloc       by ymg  August 2013                                        ;
;;                                                                         ;
;; Locates triangle which encloses point p using Lawson's Walk.             ;
;;                                                                         ;
;; Given p a point, Returns Index in tl of triangle containing the point.    ;
;; If outside the triangulation Return is nil.                               ;
;;                                                                         ;
;; Point list pl, Neigbour list nl and Triangle list tl are defined          ;
;; outside this routine by calling program.                                  ;
;;                                                                         ;
;; Optimized Speed and re-organized code January 2014                      ;
;; Nice but get lost when triangulation is disjointed.                      ;
;; Modified to start on a random edge and operate in a CDT                   ;
;;                                                                         ;
;;pl-端点集，每条数据储存一个三维点
;;tl-三角形集，每条数据储存三角形三个顶点在端点集中的索引号
;;nl-相邻三角形集，每条数据是每个三角形的三条边与哪个三角形相邻
(defun triloc (p / e i notfound p1 p2 p3 x x1 x2 x3 y y1 y2 y3)
  (if (not tn)
(setq tn (/ (length tl) 2)) ;_对半查找，提高效率
  )
  (setq x (car p) ;_横坐标
y (cadr p) ;_纵坐标
notfound t ;_查找标记
  )
  (while (and notfound tn)
(setq i  (nth tn tl) ;_定位一个三角形
   p1 (nth (car i) pl) ;_提取该三角形的三个顶点
   p2 (nth (cadr i) pl)
   p3 (nth (caddr i) pl)
   x1 (- (car p1) x) ;_三个端点到插值点的X,Y分量
   y1 (- (cadr p1) y)
   x2 (- (car p2) x)
   y2 (- (cadr p2) y)
   x3 (- (car p3) x)
   y3 (- (cadr p3) y)
   e  (fix (mrand 3)) ;_调用随机数提高查找效率，减少发生最不利查找的概率
)
;;下述语句原理一致，即不断对半查找
;;原理：两点确定一条直线，一条无限延伸的直线相当于将平面一分为二
;;通过判断在哪条边的外侧后提取相邻三角形，则另一侧的三角形均不需要遍历
;;相当于算法复杂度从O(n)降为O(LogN)
(cond ((= e 0)
   (cond ((minusp (- (* x1 y2) (* y1 x2))) (setq tn (car (nth tn nl))))
   ((minusp (- (* x2 y3) (* y2 x3))) (setq tn (cadr (nth tn nl))))
   ((minusp (- (* x3 y1) (* y3 x1))) (setq tn (caddr (nth tn nl))))
   ((setq notfound nil))
   )
   )
   ((= e 1)
   (cond ((minusp (- (* x2 y3) (* y2 x3))) (setq tn (cadr (nth tn nl))))
   ((minusp (- (* x3 y1) (* y3 x1))) (setq tn (caddr (nth tn nl))))
   ((minusp (- (* x1 y2) (* y1 x2))) (setq tn (car (nth tn nl))))
   ((setq notfound nil))
   )
   )
   ((= e 2)
   (cond ((minusp (- (* x3 y1) (* y3 x1))) (setq tn (caddr (nth tn nl))))
   ((minusp (- (* x1 y2) (* y1 x2))) (setq tn (car (nth tn nl))))
   ((minusp (- (* x2 y3) (* y2 x3))) (setq tn (cadr (nth tn nl))))
   ((setq notfound nil))
   )
   )
)
  )
  tn
)

ndgxy · 发表于 2017-11-2 08:33:47

感谢回答
请问p,tl,nl从哪里得到呀？TRIANGULATE吗？没看到呀

求指点，谢谢了

fools · 发表于 2017-11-2 12:27:14

p c:demoZ -> (setq pnt (grread nil 13 0)) -> (setq p (cadr pnt))
pl c:tin -> (setq pl (remduppoint (sortxy pl) fuzz))
tl triangulate -> (setq tl (cons (cadddr tr) tl))
nl (get_neighbour tl)

ndgxy · 发表于 2017-11-2 16:59:31

谢谢了，我试试看

fools · 发表于 2017-11-2 17:28:52

本帖最后由 fools 于 2017-11-3 21:33 编辑

不客气，triloc选用的算法很好。我测试过10000点，20000个三角形，1000个插值测试点，找寻三角形的平均次数为124.7次，结果秒出，远远快于遍历20000个三角形。

ndgxy · 发表于 2017-11-3 09:31:29

本帖最后由 ndgxy 于 2017-11-3 10:00 编辑

我加载了“测试里的文件”，用gen命令算出100点测试，用tin命令绘制出三角网和等高线，
并计算出nl (setq nl(getneighbour pl tl))
但是插值好像不对
里边的命令cz是我自己写的
这是为什么呀？求指点
我在triangulate函数里，把tl改为全局变量，估计pl tl nl的参数有问题，不知道问题出在哪里
求指点，多谢了

fools · 发表于 2017-11-3 21:12:03

以下是基于原程序写的简单插值示例，如果想改写原程序尽量少用全局变量，建议用XRECORD保存

;;gen->tin后
(defun c:tt (/ PTS RANGEX RANGEY TN TR XYZ )
(setq rangex 5000 ; Extent in X for the points ;
rangey (/ rangex 1.6) ; Extent in Y / Golden Ratio ;
)
(and (not nl) (setq nl (get_neighbour tl))) ;_生成相邻三角形集
(repeat 1000 (setq pts (cons (list (mrand rangex) (mrand rangey)) pts))) ;_产生1000个随机插值点
(foreach pt pts
(and (setq tn (triloc pt)) ;_找寻所在三角形
(setq tr (mapcar '(lambda (a) (nth a pl)) (nth tn tl)) ;_三角形三顶点
xyz (getz pt (car tr) (cadr tr) (caddr tr)) ;_计算插值点高程
)
(entmakex (list '(0 . "TEXT")
(cons 1 (rtos (caddr xyz) 2 3))
(cons 10 xyz)
(CONS 72 4)
(cons 11 xyz)
(cons 40 10)
)
) ;_文本值为中心点高程
)
)
)

;;gen->tin后 
(defun c:tt (/  PTS RANGEX RANGEY TN TR XYZ ) 
  (setq        rangex 5000                        ; Extent in X for the points ; 
        rangey (/ rangex 1.6)                ; Extent in Y / Golden Ratio ; 
  ) 
  (and (not nl) (setq nl (get_neighbour tl))) ;_生成相邻三角形集 
  (repeat 1000 (setq pts (cons (list (mrand rangex) (mrand rangey)) pts))) ;_产生1000个随机插值点 
  (foreach pt pts 
    (and (setq tn (triloc pt)) ;_找寻所在三角形 
         (setq tr (mapcar '(lambda (a) (nth a pl)) (nth tn tl)) ;_三角形三顶点 
               xyz  (getz pt (car tr) (cadr tr) (caddr tr)) ;_计算插值点高程 
         ) 
         (entmakex (list '(0 . "TEXT") 
                         (cons 1 (rtos (caddr xyz) 2 3)) 
                         (cons 10 xyz) 
                         (CONS 72 4) 
                         (cons 11 xyz) 
                         (cons 40 10) 
                   ) 
         ) ;_文本值为中心点高程 
    ) 
  ) 
)

lrd1861 · 发表于 2017-11-5 13:08:21

ndgxy · 发表于 2017-11-7 20:16:31

fools 发表于 2017-11-3 21:12

好的，谢谢了
我试试看。

xyp1964 · 发表于 2017-11-8 22:04:56

账号		自动登录	找回密码
密码			立即注册

[求助] 三角网插值

点评

浏览过的版块

[求助] 三角网 插值

点评

浏览过的版块

[求助] 三角网插值