;******************************************************************;
; TRIANGULATE - Lisp command to create a TIN from 3D points. ;
; =========== ;
; ;
; Written by Daniele Piazza, ADN member Mechanical Solution s.r.l. ;
; http://pdcode.com/code.htm ;
; ;
; Original C coding "Triangulate" written by PAUL BOURKE ;
; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/ ;
; ;
; This program triangulates an irregular set of points. ;
; You can replace some code (sorting, list manipulation,...) with ;
; VLisp functions to reduce the execution time. ;
; ;
; This code is not seriously tested, if you find a bug...sorry!! ;
; Goodbye, Daniele ;
;*******************************************************************
;
(defun C:TRIANGULATE (/ fuzzy nulllist ss1 ptlst nv supertriangle trianglelst i j k edgelst
circle pt flag perc)
(setq OLDCMD (getvar "CMDECHO"))
(setvar "CMDECHO" 0)
(command ".UNDO" "GROUP")
(setq OLDSNAP (getvar "OSMODE"))
(setvar "OSMODE" 0)
(setq fuzzy 1e-8) ; tolerance in equality test
(setq nulllist nil)
(princ "\nSelect points...")
(setq ss1 (ssget '((0 . "POINT"))))
(setq start (getvar "date") THINK-CNT 0) ; initiate timer & Progress Spinner Counter
(setq ptlst (getptlist ss1)) ; convert selection set to point list
(setq ptlst (xsort ptlst)) ; sort point list by X co-ordinate
(setq nv (length ptlst)) ; number of points
(setq supertriangle (findsupertriangle ptlst)) ; find super triangle
(setq ptlst (append ptlst supertriangle)) ; append coordinates to the end of vertex list
(setq trianglelst (list (list supertriangle nil))) ; add supertriangle to the triangle list
(setq i 0)
(while (< i nv)
(THINKING (strcat "Processing TIN - " (itoa (/ (* i 100) nv)) "% ")) ; update progress spinner
(setq pt (nth i ptlst))
(setq edgelst nil) ; initialize edge buffer
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(setq flag T)
(if (not (cadr (nth j trianglelst)))
(progn
(setq circle (getcircircumcircle triangle)) ; calculate circumcircle
(if (< (+ (caar circle) (cadr circle)) (car pt)) ; test point x and (pt) location
(setq trianglelst (nth_subst j (list (car (nth j trianglelst)) T) trianglelst))
)
(if (isinside pt circle)
(setq edgelst (addtriangleedges triangle edgelst)
trianglelst (nth_del j trianglelst)
flag nil
)
)
) ; end progn
) ; end if
(if flag (setq j (1+ j)) )
) ; end while loop
(setq edgelst (removedoublyedges edgelst fuzzy nulllist)) ; remove all doubly specified edges
(setq trianglelst (addnewtriangles pt edgelst trianglelst)) ; form new triangles for current point
(setq i (1+ i)) ; get next vertex
) ; end while loop
(setq trianglelst (purgetrianglelst trianglelst supertriangle fuzzy)) ; remove triangles with supertriangles edges
(foreach triangle (mapcar 'car trianglelst) ; draw triangles
(drawtriangle triangle)
)
(setvar "OSMODE" OLDSNAP)
(setq OLDSNAP nil)
(command ".UNDO" "END")
(setq stop (getvar "date"))
(princ (strcat "\r TIN Complete - Elapsed time: " (rtos (* 86400.0 (- stop start)) 2 2) " secs."))
(setvar "CMDECHO" OLDCMD)
(princ)
)
; XSORT - Original Shell Sort function replaced with VLISP sort (much quicker ;
; ;
(defun XSORT ( PTLST /)
(vl-sort PTLST (function (lambda (e1 e2) (< (car e1) (car e2)) ) ) )
)
; NTH_DEL ;
; ;
; delete the n item in the list (by position, not by value!!) ;
; ;
; Elimina l'oggetto che si trova nella posizione N della lista LST. L'utilizzo di ;
; funzioni ricorsive,oltre a non assicurare maggiore velocità, può creare problemi;
; di overflow dello stack in caso di liste molto lunghe. ;
(defun NTH_DEL (N LST / l)
(repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l)(cdr lst))
)
; NTH_SUBST ;
; ;
; Replace the index element in the list with new element. This function is ;
; recursive this is not a great solution with a large amount of data. ;
; ;
(defun NTH_SUBST (index new Alist)
(cond
((minusp index) Alist)
((zerop index)(cons new (cdr Alist)))
(T (cons (car Alist)(nth_subst (1- index) new (cdr Alist))))
)
)
; GETPTLIST ;
; ;
; sset -> list (p1 p2 p3 ... pn) ;
; ;
(defun GETPTLIST (ss1 / i pt ptlst)
(if (not (zerop (sslength ss1)))
(progn
(setq i 0)
(while
(setq pt (ssname ss1 i))
(setq ptlst (cons (cdr (assoc 10 (entget pt))) ptlst))
(setq i (1+ i))
)
)
)
ptlst
)
; FINDSUPERTRIANGLE ;
; ;
; Search the supertriangle that contain all points in the data set ;
; ;
(defun FINDSUPERTRIANGLE (ptlst / xmax xmin ymax ymin dx dy dmax xmid ymid
trx1 trx2 trx3 try1 try2 try3 trz1 trz2 trz3
)
(setq xmax (apply 'max (mapcar 'car ptlst))
xmin (apply 'min (mapcar 'car ptlst))
ymax (apply 'max (mapcar 'cadr ptlst))
ymin (apply 'min (mapcar 'cadr ptlst))
dx (- xmax xmin)
dy (- ymax ymin)
dmax (max dx dy)
xmid (* (+ xmax xmin) 0.5)
ymid (* (+ ymax ymin) 0.5)
trx1 (- xmid (* dmax 2.0))
try1 (- ymid dmax)
trz1 0.0
trx2 xmid
try2 (+ ymid dmax)
trz2 0.0
trx3 (+ xmid (* dmax 2.0))
try3 (- ymid dmax)
trz3 0.0
)
(list (list trx1 try1 trz1)
(list trx2 try2 trz2)
(list trx3 try3 trz3)
)
)
; GETCIRCIRCUMCIRCLE ;
; ;
; Calculate the circumcircle (center, radius) of triangle in input ;
; ;
(defun GETCIRCIRCUMCIRCLE (triangle / p1 p2 p3 p1x p2x p3x p1y p2y p3y d xc yc rad)
(setq p1 (car triangle)
p2 (cadr triangle)
p3 (caddr triangle)
p1x (car p1) p1y (cadr p1)
p2x (car p2) p2y (cadr p2)
p3x (car p3) p3y (cadr p3)
d (* 2.0 (+ (* p1y p3x)
(* p2y p1x)
(- (* p2y p3x))
(- (* p1y p2x))
(- (* p3y p1x))
(* p3y p2x)
)
)
xc (/ (+ (* p2y p1x p1x )
(- (* p3y p1x p1x))
(- (* p2y p2y p1y))
(* p3y p3y p1y)
(* p2x p2x p3y)
(* p1y p1y p2y)
(* p3x p3x p1y)
(- (* p3y p3y p2y))
(- (* p3x p3x p2y))
(- (* p2x p2x p1y))
(* p2y p2y p3y)
(- (* p1y p1y p3y))
)
d
)
yc (/ (+ (* p1x p1x p3x)
(* p1y p1y p3x)
(* p2x p2x p1x)
(- (* p2x p2x p3x))
(* p2y p2y p1x)
(- (* p2y p2y p3x))
(- (* p1x p1x p2x))
(- (* p1y p1y p2x))
(- (* p3x p3x p1x))
(* p3x p3x p2x)
(- (* p3y p3y p1x))
(* p3y p3y p2x)
)
d
)
rad (sqrt (+ (* (- p1x xc)(- p1x xc))
(* (- p1y yc)(- p1y yc))
)
)
)
(list (list xc yc) rad)
)
; ISINSIDE ;
; ;
; test if pt is inside a circle ;
; ;
(defun ISINSIDE (pt circle)
(setq ctr (car circle)
rad (cadr circle)
)
(< (distance pt ctr) rad)
)
; ADDTRIANGLEEDGES ;
; ;
; add triangle edges at the edge queue ;
; ;
(defun ADDTRIANGLEEDGES (triangle edgelst)
(append edgelst (list (list (car triangle) (cadr triangle))
(list (cadr triangle) (caddr triangle))
(list (caddr triangle)(car triangle))
)
)
)
; DRAWTRIANGLE ;
; ;
; the fun side of the algorithm. Draw triangulation. ;
; ;
(defun DRAWTRIANGLE (triangle)
(entmake (list (cons 0 "3DFACE") (cons 10 (car triangle)) (cons 11 (caddr triangle))
(cons 12 (cadr triangle)) (cons 13 (cadr triangle))))
)
; EQUALMEMBER ;
; ;
; Check if "item" is in "lista" or not by equality test. With real number the ;
; standard fuction "member" not work correctly. ;
; ;
(defun EQUALMEMBER (item lista fuzzy /)
(apply 'or (mapcar '(lambda (x) (equal x item fuzzy)) lista))
)
; REMOVEDOUBLYEDGES ;
; ;
; Test the edge queue to remove duplicates (warning CW & CCW!) ;
; ;
(defun REMOVEDOUBLYEDGES (edgelst fuzzy nulllist /)
(setq j 0)
(while (< j (length edgelst))
(setq k (1+ j))
(while (< k (length edgelst))
(if
(or (and (equal (car (nth j edgelst)) (car (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
)
(and (equal (car (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (car (nth k edgelst)) fuzzy)
)
)
(setq edgelst (nth_subst j nulllist edgelst)
edgelst (nth_subst k nulllist edgelst)
)
)
(setq k (1+ k))
)
(setq j (1+ j))
)
edgelst
)
; ADDNEWTRIANGLES ;
; ;
; Add new triangle generated by pt to triangle list. ;
; ;
(defun ADDNEWTRIANGLES (pt edgelst trianglelst / j triangle )
(setq j 0)
(while (< j (length edgelst))
(if (nth j edgelst)
(setq triangle (cons pt (nth j edgelst))
trianglelst (cons (list triangle nil) trianglelst)
)
)
(setq j (1+ j))
)
trianglelst
)
; PURGETRIANGLELST ;
; ;
; replace all triangles that share a vertex with supertriangle ;
; ;
(defun PURGETRIANGLELST (trianglelst supertriangle fuzzy /)
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(if (apply 'or
(mapcar '(lambda (x) (equalmember x supertriangle fuzzy))
triangle
)
)
(setq trianglelst (nth_del j trianglelst))
(setq j (1+ j))
)
)
)
; ;
; THINKING - STANDARD PROGRESS SPINNER ;
; ;
(defun THINKING (prmpt)
(setq THINK-CNT (1+ THINK-CNT))
(princ (strcat "\r" (nth (rem THINK-CNT 4) '("\|" "\/" "\-" "\\")) prmpt))
)
; ********************************* END OF CODING *******************************************
(princ "\n'TRIANGULATE' Loaded \n")
(princ)
Command: list
1 found
3D FACE Layer: "GS-ELEV"
Space: Model space
Handle = AB8
first point, X=539837.13728 Y=293829.99543 Z= 2.82000
second point, X=539828.96328 Y=293844.62043 Z= 2.82000
third point, X=539826.73328 Y=293843.57243 Z= 2.87000
fourth point, X=539826.73328 Y=293843.57243 Z= 2.87000
Command:
Command: list
1 found
3D FACE Layer: "GS-ELEV"
Space: Model space
Handle = AB9
first point, X=539837.13728 Y=293829.99543 Z= 2.82000
second point, X=539826.73328 Y=293843.57243 Z= 2.87000
third point, X=539835.71028 Y=293829.00243 Z= 2.87000
fourth point, X=539835.71028 Y=293829.00243 Z= 2.87000;******************************************************************;
; TRIANGULATE - Lisp command to create a TIN from 3D points. ;
; =========== ;
; ;
; Written by Daniele Piazza, ADN member Mechanical Solution s.r.l. ;
; http://pdcode.com/code.htm ;
; ;
; Original C coding "Triangulate" written by PAUL BOURKE ;
; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/ ;
; ;
; This program triangulates an irregular set of points. ;
; You can replace some code (sorting, list manipulation,...) with ;
; VLisp functions to reduce the execution time. ;
; ;
; This code is not seriously tested, if you find a bug...sorry!! ;
; Goodbye, Daniele ;
;*******************************************************************
;
;;
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;
(defun C:TRIANGULATE (/ fuzzy nulllist ss1 ptlst nv supertriangle trianglelst i j k edgelst
circle pt flag perc)
(setq OLDCMD (getvar "CMDECHO"))
(setvar "CMDECHO" 0)
(command ".UNDO" "GROUP")
(setq OLDSNAP (getvar "OSMODE"))
(setvar "OSMODE" 0)
(setq fuzzy 1e-8) ; tolerance in equality test
(setq nulllist nil)
(princ "\nSelect points...")
(setq ss1 (ssget '((0 . "POINT"))))
(setq start (getvar "date") THINK-CNT 0) ; initiate timer & Progress Spinner Counter
(setq ptlst (getptlist ss1)) ; convert selection set to point list
(setq ptlst (xsort ptlst)) ; sort point list by X co-ordinate
(setq nv (length ptlst)) ; number of points
(setq supertriangle (findsupertriangle ptlst)) ; find super triangle
(setq ptlst (append ptlst supertriangle)) ; append coordinates to the end of vertex list
(setq trianglelst (list (list supertriangle nil))) ; add supertriangle to the triangle list
(setq i 0)
(setq cab 0) ; CAB debug
(while (< i nv)
(THINKING (strcat "Processing TIN - " (itoa (/ (* i 100) nv)) "% ")) ; update progress spinner
(setq pt (nth i ptlst))
(setq edgelst nil) ; initialize edge buffer
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(setq flag T)
(if (not (cadr (nth j trianglelst)))
(progn
(setq circle (getcircircumcircle triangle)) ; calculate circumcircle
(if (< (+ (caar circle) (cadr circle)) (car pt)) ; test point x and (pt) location
(setq trianglelst (nth_subst j (list (car (nth j trianglelst)) T) trianglelst))
)
(if (isinside pt circle)
(setq edgelst (addtriangleedges triangle edgelst)
trianglelst (nth_del j trianglelst)
flag nil
)
)
) ; end progn
) ; end if
(if flag (setq j (1+ j)) )
) ; end while loop
(setq edgelst (removedoublyedges edgelst fuzzy nulllist)) ; remove all doubly specified edges
(setq trianglelst (addnewtriangles pt edgelst trianglelst)) ; form new triangles for current point
(setq i (1+ i)) ; get next vertex
) ; end while loop
(setq trianglelst (purgetrianglelst trianglelst supertriangle fuzzy)) ; remove triangles with supertriangles edges
(foreach triangle (mapcar 'car trianglelst) ; draw triangles
(drawtriangle triangle)
)
(setvar "OSMODE" OLDSNAP)
(setq OLDSNAP nil)
(command ".UNDO" "END")
(setq stop (getvar "date"))
(princ (strcat "\r TIN Complete - Elapsed time: " (rtos (* 86400.0 (- stop start)) 2 2) " secs."))
(setvar "CMDECHO" OLDCMD)
(princ)
)
; XSORT - Original Shell Sort function replaced with VLISP sort (much quicker :-) ;
; ;
(defun XSORT ( PTLST /)
(vl-sort PTLST (function (lambda (e1 e2) (< (car e1) (car e2)) ) ) )
)
; NTH_DEL ;
; ;
; delete the n item in the list (by position, not by value!!) ;
; ;
; Elimina l'oggetto che si trova nella posizione N della lista LST. L'utilizzo di ;
; funzioni ricorsive,oltre a non assicurare maggiore velocità, può creare problemi;
; di overflow dello stack in caso di liste molto lunghe. ;
(defun NTH_DEL (N LST / l)
(repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l)(cdr lst))
)
; NTH_SUBST ;
; ;
; Replace the index element in the list with new element. This function is ;
; recursive this is not a great solution with a large amount of data. ;
; ;
(defun NTH_SUBST (index new Alist)
(cond
((minusp index) Alist)
((zerop index)(cons new (cdr Alist)))
(T (cons (car Alist)(nth_subst (1- index) new (cdr Alist))))
)
)
; GETPTLIST ;
; ;
; sset -> list (p1 p2 p3 ... pn) ;
; ;
(defun GETPTLIST (ss1 / i pt ptlst)
(if (not (zerop (sslength ss1)))
(progn
(setq i 0)
(while
(setq pt (ssname ss1 i))
(setq ptlst (cons (cdr (assoc 10 (entget pt))) ptlst))
(setq i (1+ i))
)
)
)
ptlst
)
; FINDSUPERTRIANGLE ;
; ;
; Search the supertriangle that contain all points in the data set ;
; ;
(defun FINDSUPERTRIANGLE (ptlst / xmax xmin ymax ymin dx dy dmax xmid ymid
trx1 trx2 trx3 try1 try2 try3 trz1 trz2 trz3
)
(setq xmax (apply 'max (mapcar 'car ptlst))
xmin (apply 'min (mapcar 'car ptlst))
ymax (apply 'max (mapcar 'cadr ptlst))
ymin (apply 'min (mapcar 'cadr ptlst))
dx (- xmax xmin)
dy (- ymax ymin)
dmax (max dx dy)
xmid (* (+ xmax xmin) 0.5)
ymid (* (+ ymax ymin) 0.5)
trx1 (- xmid (* dmax 2.0))
try1 (- ymid dmax)
trz1 0.0
trx2 xmid
try2 (+ ymid dmax)
trz2 0.0
trx3 (+ xmid (* dmax 2.0))
try3 (- ymid dmax)
trz3 0.0
)
(list (list trx1 try1 trz1)
(list trx2 try2 trz2)
(list trx3 try3 trz3)
)
)
;;=============================================================
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;=============================================================
(defun getcircircumcircle (triangle / p1 p2 p3 pr1 pr2 cen rad bisector)
;; return a pt list for a perpendicular bisector 20 units long
(defun bisector (p1 p2 / perp_ang midpt)
(setq p1 (list (car p1) (cadr p1)) ; make sure 2d point
perp_ang (+ (angle p1 p2) (/ pi 2.0))) ; perpendicular angle
(setq midpt (mapcar '(lambda (pa pb) (+ (/ (- pb pa) 2.0) pa)) p1 p2))
(list (polar midpt perp_ang 10) (polar midpt (+ pi perp_ang) 10))
)
(setq p1 (car triangle)
p2 (cadr triangle)
p3 (caddr triangle)
pr1 (bisector p1 p2)
pr2 (bisector p1 p3)
cen (inters (car pr1) (cadr pr1) (car pr2) (cadr pr2) nil)
rad (distance cen p1)
)
(list cen rad)
)
;;=============================================================
; ISINSIDE ;
; ;
; test if pt is inside a circle ;
; ;
(defun ISINSIDE (pt circle)
(setq ctr (car circle)
rad (cadr circle)
)
(< (distance pt ctr) rad)
)
; ADDTRIANGLEEDGES ;
; ;
; add triangle edges at the edge queue ;
; ;
(defun ADDTRIANGLEEDGES (triangle edgelst)
(append edgelst (list (list (car triangle) (cadr triangle))
(list (cadr triangle) (caddr triangle))
(list (caddr triangle)(car triangle))
)
)
)
; DRAWTRIANGLE ;
; ;
; the fun side of the algorithm. Draw triangulation. ;
; ;
(defun DRAWTRIANGLE (triangle)
(entmake (list (cons 0 "3DFACE") (cons 10 (car triangle)) (cons 11 (caddr triangle))
(cons 12 (cadr triangle)) (cons 13 (cadr triangle))))
)
; EQUALMEMBER ;
; ;
; Check if "item" is in "lista" or not by equality test. With real number the ;
; standard fuction "member" not work correctly. ;
; ;
(defun EQUALMEMBER (item lista fuzzy /)
(apply 'or (mapcar '(lambda (x) (equal x item fuzzy)) lista))
)
; REMOVEDOUBLYEDGES ;
; ;
; Test the edge queue to remove duplicates (warning CW & CCW!) ;
; ;
(defun REMOVEDOUBLYEDGES (edgelst fuzzy nulllist /)
(setq j 0)
(while (< j (length edgelst))
(setq k (1+ j))
(while (< k (length edgelst))
(if
(or (and (equal (car (nth j edgelst)) (car (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
)
(and (equal (car (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (car (nth k edgelst)) fuzzy)
)
)
(setq edgelst (nth_subst j nulllist edgelst)
edgelst (nth_subst k nulllist edgelst)
)
)
(setq k (1+ k))
)
(setq j (1+ j))
)
edgelst
)
; ADDNEWTRIANGLES ;
; ;
; Add new triangle generated by pt to triangle list. ;
; ;
(defun ADDNEWTRIANGLES (pt edgelst trianglelst / j triangle )
(setq j 0)
(while (< j (length edgelst))
(if (nth j edgelst)
(setq triangle (cons pt (nth j edgelst))
trianglelst (cons (list triangle nil) trianglelst)
)
)
(setq j (1+ j))
)
trianglelst
)
; PURGETRIANGLELST ;
; ;
; replace all triangles that share a vertex with supertriangle ;
; ;
(defun PURGETRIANGLELST (trianglelst supertriangle fuzzy /)
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(if (apply 'or
(mapcar '(lambda (x) (equalmember x supertriangle fuzzy))
triangle
)
)
(setq trianglelst (nth_del j trianglelst))
(setq j (1+ j))
)
)
)
; ;
; THINKING - STANDARD PROGRESS SPINNER ;
; ;
(defun THINKING (prmpt)
(setq THINK-CNT (1+ THINK-CNT))
(princ (strcat "\r" (nth (rem THINK-CNT 4) '("\|" "\/" "\-" "\\")) prmpt))
)
; ********************************* END OF CODING *******************************************
(princ "\n'TRIANGULATE' Loaded \n")
(princ)
;;; Purpose: To Generate the triangle by incresing point
;;; Version: 0.1
;;; Credit to Paul Bourke (pbourke@swin.edu.au) for the original C
Program :))
;;; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/
;;; The following codes are translate from C to AutoLisp by QJCHEN
;;; South China University of Technology
;;; Thanks : Eachy at xdcad.net introduce the source code pages
;;; and the STDLIB Function of Reini Urban at http://xarch.tu-
graz.ac.at/autocad/stdlib/archive/
;;; 2006.06.30
(defun c:test (/ tpoints temp howmany ij p1 p2 p3)
(setq tpoints 1
vertex (givever)
triangle (givetri)
edges (giveedg)
)
(while (setq temp (getpoint))
(setq vertex (qj-setnmth (nth 0 temp) tpoints 1 vertex))
(setq vertex (qj-setnmth (nth 1 temp) tpoints 2 vertex))
(if (> tpoints 2)
(progn
(setq howmany (Triangulate tpoints))
)
)
(setq tpoints (1+ tpoints))
(setq ij 0)
(command "redraw")
(if (>= tpoints 4)
(progn
(repeat howmany
(setq ij (1+ ij))
(setq p1 (nth (1- (nth 0 (nth (1- ij) triangle))) vertex))
(setq p2 (nth (1- (nth 1 (nth (1- ij) triangle))) vertex))
(setq p3 (nth (1- (nth 2 (nth (1- ij) triangle))) vertex))
(grdraw p2 p1 1)
(grdraw p1 p3 1)
(grdraw p2 p3 1)
)
)
) ; (grdraw p1 p3 1)
; (grdraw p2 p3 1)
; (grdraw p3 p1 1)
)
)
;|The main function|;
(defun Triangulate (nvert / xmin ymin xmax ymax i dx dy xmid ymid
complete
ntri inc nedge i j Triangulate1
)
(setq xmin (xofv vertex 1))
(setq ymin (yofv vertex 1))
(setq xmax xmin
ymax ymin
)
(setq i 2)
(while (<= i nvert)
(if (< (xofv vertex i) xmin)
(setq xmin (xofv vertex i))
)
(if (> (xofv vertex i) xmax)
(setq xmax (xofv vertex i))
)
(if (< (yofv vertex i) ymin)
(setq ymin (yofv vertex i))
)
(if (> (yofv vertex i) ymax)
(setq ymax (yofv vertex i))
)
(setq i (1+ i))
)
(setq dx (- xmax xmin))
(setq dy (- ymax ymin))
(if (> dx dy)
(setq dmax dx)
(setq dmax dy)
)
(setq xmid (/ (+ xmax xmin) 2))
(setq ymid (/ (+ ymax ymin) 2))
(setq vertex (qj-setnmth (- xmid (* dmax 2)) (1+ nvert) 1 vertex))
(setq vertex (qj-setnmth (- ymid dmax) (1+ nvert) 2 vertex))
(setq vertex (qj-setnmth xmid (+ nvert 2) 1 vertex))
(setq vertex (qj-setnmth (+ ymid (* 2 dmax)) (+ nvert 2) 2 vertex))
(setq vertex (qj-setnmth (+ xmid (* 2 dmax)) (+ nvert 3) 1 vertex))
(setq vertex (qj-setnmth (- ymid dmax) (+ nvert 3) 2 vertex))
(setq triangle (qj-setnmth (+ nvert 1) 1 1 triangle))
(setq triangle (qj-setnmth (+ nvert 2) 1 2 triangle))
(setq triangle (qj-setnmth (+ nvert 3) 1 3 triangle))
(setq complete (append
complete
(list nil)
)
)
(setq ntri 1);;;;;;;;;;;start loop i
(setq i 1)
(while (<= i nvert)
(setq nedge 0);;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
(setq j 0
temp (- 1)
)
(while (< temp ntri)
(setq j (1+ j)
temp j
)
(if (/= (nth (1- j) complete) T)
(progn
(setq inc (InCircle1 (xofv vertex i) (yofv vertex i) (xof vertex
triangle
j 1
)
(yof vertex triangle j 1) (xof vertex
triangle j 2
) (yof vertex
triangle j 2
) (xof vertex
triangle j
3
) (yof vertex triangle
j 3
)
)
)
)
)
(if inc
(progn
(setq edges (qj-setnmth (nth 0 (nth (1- j) triangle)) 1
(+ nedge 1) edges
)
)
(setq edges (qj-setnmth (nth 1 (nth (1- j) triangle)) 2
(+ nedge 1) edges
)
)
(setq edges (qj-setnmth (nth 1 (nth (1- j) triangle)) 1
(+ nedge 2) edges
)
)
(setq edges (qj-setnmth (nth 2 (nth (1- j) triangle)) 2
(+ nedge 2) edges
)
)
(setq edges (qj-setnmth (nth 2 (nth (1- j) triangle)) 1
(+ nedge 3) edges
)
)
(setq edges (qj-setnmth (nth 0 (nth (1- j) triangle)) 2
(+ nedge 3) edges
)
)
(setq Nedge (+ Nedge 3))
(setq triangle (qj-setnmth ([n,m] triangle ntri 1) j 1 triangle))
(setq triangle (qj-setnmth ([n,m] triangle ntri 2) j 2 triangle))
(setq triangle (qj-setnmth ([n,m] triangle ntri 3) j 3 triangle))
(setq complete (std-setnth (nth (1- ntri) complete) (1- j)
complete
)
)
(setq j (1- j)
temp j
)
(setq ntri (1- ntri))
)
)
);;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
(setq j 1)
(while (<= j (1- Nedge))
(if (and
(/= ([n,m] edges 1 j) 0)
(/= ([n,m] edges 2 j) 0)
)
(progn
(setq k (1+ j))
(while (<= k Nedge)
(if (and
(/= ([n,m] edges 1 k) 0)
(/= ([n,m] edges 2 k) 0)
)
(if (= ([n,m] edges 1 j) ([n,m] edges 2 k))
(if (= ([n,m] edges 2 j) ([n,m] edges 1 k))
(progn
(setq edges (qj-setnmth 0 1 j edges))
(setq edges (qj-setnmth 0 2 j edges))
(setq edges (qj-setnmth 0 1 k edges))
(setq edges (qj-setnmth 0 1 k edges))
)
)
)
)
(setq k (1+ k))
)
)
)
(setq j (1+ j))
);;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;;;;;;;;;;;;;;;;;;;;;;;;;;;
(setq j 1)
(while (<= j Nedge)
(if (and
(/= ([n,m] edges 1 j) 0)
(/= ([n,m] edges 2 j) 0)
)
(progn
(setq ntri (1+ ntri))
(setq triangle (qj-setnmth ([n,m] edges 1 j) ntri 1 triangle))
(setq triangle (qj-setnmth ([n,m] edges 2 j) ntri 2 triangle))
(setq triangle (qj-setnmth i ntri 3 triangle))
(setq complete (std-setnth nil (1- ntri) complete))
)
)
(setq j (1+ j))
);;;;;;;;;;;;;;;;;;;;;;;;;;;
(setq i (1+ i))
);;;;;end of loop i
(setq i 0
temp (- 1)
)
(while (< temp ntri)
(setq i (1+ i)
temp i
)
(if (or
(> ([n,m] triangle i 1) nvert)
(> ([n,m] triangle i 2) nvert)
(> ([n,m] triangle i 3) nvert)
)
(progn
(setq triangle (qj-setnmth ([n,m] triangle ntri 1) i 1 triangle))
(setq triangle (qj-setnmth ([n,m] triangle ntri 2) i 2 triangle))
(setq triangle (qj-setnmth ([n,m] triangle ntri 3) i 3 triangle))
(setq i (1- i)
temp i
)
(setq ntri (1- ntri))
)
)
)
(setq Triangulate1 ntri)
Triangulate1
)
;;; stdlib - to substitute the i element of the list to new value
(defun std-%setnth (new i lst / fst len)
(cond
((minusp i)
lst
)
((> i (setq len (length lst)))
lst
)
((> i (/ len 2))
(reverse (std-%setnth new (1- (- len i)) (reverse lst)))
)
(t
(append
(progn
(setq fst nil) ; ; possible vl lsa compiler bug
(repeat (rem i 4)
(setq fst (cons (car lst) fst)
lst (cdr lst)
)
)
(repeat (/ i 4)
(setq fst (cons (cadddr lst) (cons (caddr lst) (cons
(cadr lst)
(cons
(car lst)
fst
)
)
)
)
lst (cddddr lst)
)
)
(reverse fst)
)
(if (listp new)
new
(list new)
) ; v0.4001
(cdr lst)
)
)
)
)
(defun std-setnth (new i lst)
(std-%setnth (list new) i lst)
)
;;; substitute the i row j column of a 2 dimension array
(defun qj-setnmth (new i j lst / listb lista)
(setq listb lst)
(setq i (1- i))
(setq j (1- j))
(setq lista (nth i lst))
(setq lista (std-setnth new j lista))
(setq listb (std-setnth lista i listb))
listb
)
;;; get the n row m column of the two dimension array
(defun [n,m] (a n m / i)
(setq i (nth (1- m) (nth (1- n) a)))
i
)
;;; get the n row of the one dimension array
(defun [n] (a n / i)
(setq i (nth (1- n) a))
i
)
;|Vertex has the form '((x1 y1)(x2 y2)(x3 y3)(x4 y4)())
The function xofv is to get the x value of the i element,i start from 1|;
(defun xofv (vertex i / res)
(setq res (nth 0 (nth (- i 1) vertex)))
res
)
;|Vertex has the form '((x1 y1)(x2 y2)(x3 y3)(x4 y4)())
The function yofv is to get the y value of the i element,i start from 1|;
(defun yofv (vertex i / res)
(setq res (nth 1 (nth (- i 1) vertex)))
res
)
;|Lis has the form '(((x11 y11)(x12 y12)(x13 y13))((x21 y21)(x22 y22)(x23
y23))(()()()))
The function xof is to get the x value of the i,j element,i and j start from
1
and j is the outer sequence, and i is the inter sequence, total 3|;
(defun xof (lisa lisb j v123 / res1 res2 res)
(setq res1 (nth (1- j) lisb))
(setq res2 (nth (1- v123) res1))
(setq res3 (nth (1- res2) lisa))
(setq res (nth 0 res3))
res
)
;|Lis has the form '(((x11 y11)(x12 y12)(x13 y13))((x21 y21)(x22 y22)(x23
y23))(()()()))
The function xof is to get the y value of the i,j element,i and j start from
1
and j is the outer sequence, and i is the inter sequence, total 3|;
(defun yof (lisa lisb j v123 / res1 res2 res)
(setq res1 (nth (1- j) lisb))
(setq res2 (nth (1- v123) res1))
(setq res3 (nth (1- res2) lisa))
(setq res (nth 1 res3))
res
)
;(defun append1 (new n lis / res1 res2 res)
;
; (setq res1 (nth (1- n) lis))
; (setq res2 (append
; res1
; (list new)
; )
; )
; (setq res (std-setnth res2 (1- n) lis))
; res
;)
;
;|Return TRUE if the point (xp,yp) lies inside the circumcircle
made up by points (x1,y1) (x2,y2) (x3,y3)
The circumcircle centre is returned in (xc,yc) and the radius r
NOTE: A point on the edge is inside the circumcircle|;
(defun InCircle1 (xp yp x1 y1 x2 y2 x3 y3 / InCircle eps mx2 my2 xc yc
m1
mx1 my1 m2 mx2 my2 dx dy rsqr r drsqr
)
(setq eps 0.000001)
(setq InCircle nil)
(if (and
(< (abs (- y1 y2)) eps)
(< (abs (- y2 y3)) eps)
)
(alert "INCIRCUM - F - Points are coincident !!")
(progn
(cond
((< (abs (- y2 y1)) eps)
(setq m2 (/ (- x2 x3) (- y3 y2)))
(setq mx2 (/ (+ x2 x3) 2))
(setq my2 (/ (+ y2 y3) 2))
(setq xc (/ (+ x2 x1) 2))
(setq yc (+ my2 (* m2 (- xc mx2))))
)
((< (abs (- y3 y2)) eps)
(setq m1 (/ (- x1 x2) (- y2 y1)))
(setq mx1 (/ (+ x1 x2) 2))
(setq my1 (/ (+ y1 y2) 2))
(setq xc (/ (+ x3 x2) 2))
(setq yc (+ my1 (* m1 (- xc mx1))))
)
(T
(setq m1 (/ (- x1 x2) (- y2 y1)))
(setq m2 (/ (- x2 x3) (- y3 y2)))
(setq mx1 (/ (+ x1 x2) 2))
(setq mx2 (/ (+ x2 x3) 2))
(setq my1 (/ (+ y1 y2) 2))
(setq my2 (/ (+ y2 y3) 2))
(setq xc (/ (- (+ (* m1 mx1) my2) my1 (* m2 mx2)) (- m1
m2)))
(setq yc (+ my1 (* m1 (- xc mx1))))
)
)
(setq dx (- x2 xc))
(setq dy (- y2 yc))
(setq rsqr (+ (* dx dx) (* dy dy)))
(setq r (sqrt rsqr))
(setq dx (- xp xc))
(setq dy (- yp yc))
(setq drsqr (+ (* dx dx) (* dy dy)))
(if (<= drsqr rsqr)
(setq InCircle T)
)
)
)
InCircle
)
;|Determines which side of a line the point (xp,yp) lies.
The line goes from (x1,y1) to (x2,y2)
Returns -1 for a point to the left
0 for a point on the line
+1 for a point to the right|;
(defun whichside (xp yp x1 y1 x2 y2 / equation)
(setq equation (- (* (- yp y1) (- x2 x1)) (* (- y2 y1) (- xp x1))))
(cond
((> equation 0)
(setq whichside (- 0 1))
)
((= equation 0)
(setq whichside 0)
)
(T
(setq whichside 1)
)
)
whichside
)
(defun givetri (/ lis)
(repeat 200
(setq lis (append
lis
(list (list nil nil nil))
)
)
)
lis
)
(defun givever (/ lis)
(repeat 200
(setq lis (append
lis
(list (list nil nil))
)
)
)
lis
)
(defun giveedg (/ lis lis1 lis2)
(repeat 200
(setq lis1 (append
lis1
(list nil)
)
)
)
(setq lis2 lis1)
(setq lis (append
lis
(list lis1)
)
)
(setq lis (append
lis
(list lis2)
)
)
lis
)
effect(defun c:test (/ I L S)
(princ (strcat "\n select points"))
(if (setq i 0
s (ssget '((0 . "POINT")))
) ;_ setq
(progn (repeat (sslength s)
(setq l (cons (cdr (assoc 10 (entget (ssname s i)))) l)
i (1+ i)
) ;_ setq
) ;_ repeat
(eea-delone-triangulate i l)
) ;_ progn
) ;_ if
) ;_ defun
(defun eea-delone-triangulate (i1 L / A A1 A2 A3
I I2 L1 L2 L3 LP MA
MI P S TI TR X1 X2
Y1 Y2
)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; Program triangulate an irregular set of 3d points.
;;
;;*********************************************************
(if l
(progn
(setq ti (car (_VL-TIMES))
i 1
i1 (/ i1 100.)
i2 0
l (vl-sort
(mapcar
(function (lambda (p)
(list (/ (fix (* (car p) 1000)) 1000.)
(/ (fix (* (cadr p) 1000)) 1000.)
(caddr p)
) ;_ list
) ;_ lambda
) ;_ function
l
) ;_ mapcar
(function (lambda (a b) (>= (car a) (car b))))
) ;_ vl-sort
x2 (caar l)
y1 (cadar l)
y2 y1
) ;_ setq
(while l
(setq a (fix (caar l))
a1 (list (car l))
l (cdr l)
) ;_ setq
(while (and l (= (fix (caar l)) a))
(setq a2 (car l))
(if (<= (cadr a2) y1)
(setq y1 (cadr a2))
(if (> (cadr a2) y2)
(setq y2 (cadr a2))
) ;_ if
) ;_ if
(setq a1 (cons (car l) (vl-remove a2 a1))
l (cdr l)
) ;_ setq
) ;_ while
(foreach a a1 (setq lp (cons a lp)))
) ;_ while
(setq x1 (caar lp)
a (list (/ (+ x1 x2) 2) (/ (+ y1 y2) 2))
a1 (distance a (list x1 y1))
ma (+ (car a) a1 a1)
mi (- (car a) a1)
s (list (list ma (cadr a) 0)
(list mi (+ (cadr a) a1 a1) 0)
(list (- (car a) a1) (- (cadr a) a1 a1) 0)
) ;_ list
l (list (cons x2 (cons a (cons (+ a1 a1) s))))
ma (1- ma)
mi (1+ mi)
) ;_ setq
(while lp
(setq p (car lp)
lp (cdr lp)
l1 nil
) ;_ setq
(while l
(setq tr (car l)
l (cdr l)
) ;_ setq
(cond
((< (car tr) (car p)) (setq l2 (cons (cdddr tr) l2)))
((< (distance p (cadr tr)) (caddr tr))
(setq tr (cdddr tr)
a1 (car tr)
a2 (cadr tr)
a3 (caddr tr)
l1 (cons (list (+ (car a1) (car a2))
(+ (cadr a1) (cadr a2))
a1
a2
) ;_ list
(cons (list (+ (car a2) (car a3))
(+ (cadr a2) (cadr a3))
a2
a3
) ;_ list
(cons (list (+ (car a3) (car a1))
(+ (cadr a3) (cadr a1))
a3
a1
) ;_ list
l1
) ;_ cons
) ;_ cons
) ;_ cons
) ;_ setq
)
(t (setq l3 (cons tr l3)))
) ;_ cond
) ;_ while
(setq l l3
l3 nil
l1 (vl-sort l1
(function (lambda (a b)
(if (= (car a) (car b))
(<= (cadr a) (cadr b))
(< (car a) (car b))
) ;_ if
) ;_ lambda
) ;_ function
) ;_ vl-sort
) ;_ setq
(while l1
(if (and (= (caar l1) (caadr l1))
(= (cadar l1) (cadadr l1))
) ;_ and
(setq l1 (cddr l1))
(setq l (cons (eea-data-triangle p (cddar l1)) l)
l1 (cdr l1)
) ;_ setq
) ;_ if
) ;_ while
(if (and (< (setq i (1- i)) 1) (< i2 100))
(progn
(setvar
"MODEMACRO"
(strcat
" "
(itoa (setq i2 (1+ i2)))
" % "
(substr
"||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||"
1
i2
) ;_ substr
(substr
"..."
1
(- 100 i2)
) ;_ substr
) ;_ strcat
) ;_ setvar
(setq i i1)
) ;_ progn
) ;_ if
) ;_ while
(foreach a l (setq l2 (cons (cdddr a) l2)))
(setq
l2 (vl-remove-if-not
(function
(lambda (a)
(and (< mi (caadr a) ma) (< mi (caaddr a) ma))
) ;_ lambda
) ;_ function
l2
) ;_ vl-remove-if
) ;_ setq
(foreach a l2
(entmake (list (cons 0 "3DFACE")
(cons 10 (car a))
(cons 11 (car a))
(cons 12 (cadr a))
(cons 13 (caddr a))
) ;_ list
) ;_ entmake
) ;_ foreach
) ;_ progn
) ;_ if
(setvar "MODEMACRO" "")
(princ (strcat "\n "
(rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4)
" secs."
) ;_ strcat
) ;_ princ
(princ)
) ;_ defun
(defun eea-data-triangle (P1 l / A A1 P2 P3 P4 S)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; Calculation of the centre of a circle and circle radius
;; for program triangulate
;;
;; (eea-data-triangle (getpoint)(list(getpoint)(getpoint)))
;;*********************************************************
(setq p2 (car l)
p3 (cadr l)
p4 (list (car p3) (cadr p3))
) ;_ setq
(if
(not
(zerop
(setq s (sin (setq a (- (angle p2 p4) (angle p2 p1)))))
) ;_ zerop
) ;_ not
(progn (setq a (polar p4
(+ -1.570796326794896 (angle p4 p1) a)
(setq a1 (/ (distance p1 p4) s 2.))
) ;_ polar
a1 (abs a1)
) ;_ setq
(list (+ (car a) a1) a a1 p1 p2 p3)
) ;_ progn
) ;_ if
) ;_ defunCommand:
Command: (LOAD "D:/Work/triangle/test/test.VLX") nil
Command: test
select points
Select objects: Specify opposite corner: 100000 found
Select objects:
39.0630 secs.
Command:Test this:Code: [Select];******************************************************************;
; TRIANGULATE - Lisp command to create a TIN from 3D points. ;
; =========== ;
; ;
; Written by Daniele Piazza, ADN member Mechanical Solution s.r.l. ;
; http://pdcode.com/code.htm ;
; ;
; Original C coding "Triangulate" written by PAUL BOURKE ;
; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/ ;
; ;
; This program triangulates an irregular set of points. ;
; You can replace some code (sorting, list manipulation,...) with ;
; VLisp functions to reduce the execution time. ;
; ;
; This code is not seriously tested, if you find a bug...sorry!! ;
; Goodbye, Daniele ;
;*******************************************************************
;
;;
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;
(defun C:TRIANGULATE (/ fuzzy nulllist ss1 ptlst nv supertriangle trianglelst i j k edgelst
circle pt flag perc)
(setq OLDCMD (getvar "CMDECHO"))
(setvar "CMDECHO" 0)
(command ".UNDO" "GROUP")
(setq OLDSNAP (getvar "OSMODE"))
(setvar "OSMODE" 0)
(setq fuzzy 1e-8) ; tolerance in equality test
(setq nulllist nil)
(princ "\nSelect points...")
(setq ss1 (ssget '((0 . "POINT"))))
(setq start (getvar "date") THINK-CNT 0) ; initiate timer & Progress Spinner Counter
(setq ptlst (getptlist ss1)) ; convert selection set to point list
(setq ptlst (xsort ptlst)) ; sort point list by X co-ordinate
(setq nv (length ptlst)) ; number of points
(setq supertriangle (findsupertriangle ptlst)) ; find super triangle
(setq ptlst (append ptlst supertriangle)) ; append coordinates to the end of vertex list
(setq trianglelst (list (list supertriangle nil))) ; add supertriangle to the triangle list
(setq i 0)
(setq cab 0) ; CAB debug
(while (< i nv)
(THINKING (strcat "Processing TIN - " (itoa (/ (* i 100) nv)) "% ")) ; update progress spinner
(setq pt (nth i ptlst))
(setq edgelst nil) ; initialize edge buffer
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(setq flag T)
(if (not (cadr (nth j trianglelst)))
(progn
(setq circle (getcircircumcircle triangle)) ; calculate circumcircle
(if (< (+ (caar circle) (cadr circle)) (car pt)) ; test point x and (pt) location
(setq trianglelst (nth_subst j (list (car (nth j trianglelst)) T) trianglelst))
)
(if (isinside pt circle)
(setq edgelst (addtriangleedges triangle edgelst)
trianglelst (nth_del j trianglelst)
flag nil
)
)
) ; end progn
) ; end if
(if flag (setq j (1+ j)) )
) ; end while loop
(setq edgelst (removedoublyedges edgelst fuzzy nulllist)) ; remove all doubly specified edges
(setq trianglelst (addnewtriangles pt edgelst trianglelst)) ; form new triangles for current point
(setq i (1+ i)) ; get next vertex
) ; end while loop
(setq trianglelst (purgetrianglelst trianglelst supertriangle fuzzy)) ; remove triangles with supertriangles edges
(foreach triangle (mapcar 'car trianglelst) ; draw triangles
(drawtriangle triangle)
)
(setvar "OSMODE" OLDSNAP)
(setq OLDSNAP nil)
(command ".UNDO" "END")
(setq stop (getvar "date"))
(princ (strcat "\r TIN Complete - Elapsed time: " (rtos (* 86400.0 (- stop start)) 2 2) " secs."))
(setvar "CMDECHO" OLDCMD)
(princ)
)
; XSORT - Original Shell Sort function replaced with VLISP sort (much quicker :-) ;
; ;
(defun XSORT ( PTLST /)
(vl-sort PTLST (function (lambda (e1 e2) (< (car e1) (car e2)) ) ) )
)
; NTH_DEL ;
; ;
; delete the n item in the list (by position, not by value!!) ;
; ;
; Elimina l'oggetto che si trova nella posizione N della lista LST. L'utilizzo di ;
; funzioni ricorsive,oltre a non assicurare maggiore velocità, può creare problemi;
; di overflow dello stack in caso di liste molto lunghe. ;
(defun NTH_DEL (N LST / l)
(repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l)(cdr lst))
)
; NTH_SUBST ;
; ;
; Replace the index element in the list with new element. This function is ;
; recursive this is not a great solution with a large amount of data. ;
; ;
(defun NTH_SUBST (index new Alist)
(cond
((minusp index) Alist)
((zerop index)(cons new (cdr Alist)))
(T (cons (car Alist)(nth_subst (1- index) new (cdr Alist))))
)
)
; GETPTLIST ;
; ;
; sset -> list (p1 p2 p3 ... pn) ;
; ;
(defun GETPTLIST (ss1 / i pt ptlst)
(if (not (zerop (sslength ss1)))
(progn
(setq i 0)
(while
(setq pt (ssname ss1 i))
(setq ptlst (cons (cdr (assoc 10 (entget pt))) ptlst))
(setq i (1+ i))
)
)
)
ptlst
)
; FINDSUPERTRIANGLE ;
; ;
; Search the supertriangle that contain all points in the data set ;
; ;
(defun FINDSUPERTRIANGLE (ptlst / xmax xmin ymax ymin dx dy dmax xmid ymid
trx1 trx2 trx3 try1 try2 try3 trz1 trz2 trz3
)
(setq xmax (apply 'max (mapcar 'car ptlst))
xmin (apply 'min (mapcar 'car ptlst))
ymax (apply 'max (mapcar 'cadr ptlst))
ymin (apply 'min (mapcar 'cadr ptlst))
dx (- xmax xmin)
dy (- ymax ymin)
dmax (max dx dy)
xmid (* (+ xmax xmin) 0.5)
ymid (* (+ ymax ymin) 0.5)
trx1 (- xmid (* dmax 2.0))
try1 (- ymid dmax)
trz1 0.0
trx2 xmid
try2 (+ ymid dmax)
trz2 0.0
trx3 (+ xmid (* dmax 2.0))
try3 (- ymid dmax)
trz3 0.0
)
(list (list trx1 try1 trz1)
(list trx2 try2 trz2)
(list trx3 try3 trz3)
)
)
;;=============================================================
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;=============================================================
(defun getcircircumcircle (triangle / p1 p2 p3 pr1 pr2 cen rad bisector)
;; return a pt list for a perpendicular bisector 20 units long
(defun bisector (p1 p2 / perp_ang midpt)
(setq p1 (list (car p1) (cadr p1)) ; make sure 2d point
perp_ang (+ (angle p1 p2) (/ pi 2.0))) ; perpendicular angle
(setq midpt (mapcar '(lambda (pa pb) (+ (/ (- pb pa) 2.0) pa)) p1 p2))
(list (polar midpt perp_ang 10) (polar midpt (+ pi perp_ang) 10))
)
(setq p1 (car triangle)
p2 (cadr triangle)
p3 (caddr triangle)
pr1 (bisector p1 p2)
pr2 (bisector p1 p3)
cen (inters (car pr1) (cadr pr1) (car pr2) (cadr pr2) nil)
rad (distance cen p1)
)
(list cen rad)
)
;;=============================================================
; ISINSIDE ;
; ;
; test if pt is inside a circle ;
; ;
(defun ISINSIDE (pt circle)
(setq ctr (car circle)
rad (cadr circle)
)
(< (distance pt ctr) rad)
)
; ADDTRIANGLEEDGES ;
; ;
; add triangle edges at the edge queue ;
; ;
(defun ADDTRIANGLEEDGES (triangle edgelst)
(append edgelst (list (list (car triangle) (cadr triangle))
(list (cadr triangle) (caddr triangle))
(list (caddr triangle)(car triangle))
)
)
)
; DRAWTRIANGLE ;
; ;
; the fun side of the algorithm. Draw triangulation. ;
; ;
(defun DRAWTRIANGLE (triangle)
(entmake (list (cons 0 "3DFACE") (cons 10 (car triangle)) (cons 11 (caddr triangle))
(cons 12 (cadr triangle)) (cons 13 (cadr triangle))))
)
; EQUALMEMBER ;
; ;
; Check if "item" is in "lista" or not by equality test. With real number the ;
; standard fuction "member" not work correctly. ;
; ;
(defun EQUALMEMBER (item lista fuzzy /)
(apply 'or (mapcar '(lambda (x) (equal x item fuzzy)) lista))
)
; REMOVEDOUBLYEDGES ;
; ;
; Test the edge queue to remove duplicates (warning CW & CCW!) ;
; ;
(defun REMOVEDOUBLYEDGES (edgelst fuzzy nulllist /)
(setq j 0)
(while (< j (length edgelst))
(setq k (1+ j))
(while (< k (length edgelst))
(if
(or (and (equal (car (nth j edgelst)) (car (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
)
(and (equal (car (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (car (nth k edgelst)) fuzzy)
)
)
(setq edgelst (nth_subst j nulllist edgelst)
edgelst (nth_subst k nulllist edgelst)
)
)
(setq k (1+ k))
)
(setq j (1+ j))
)
edgelst
)
; ADDNEWTRIANGLES ;
; ;
; Add new triangle generated by pt to triangle list. ;
; ;
(defun ADDNEWTRIANGLES (pt edgelst trianglelst / j triangle )
(setq j 0)
(while (< j (length edgelst))
(if (nth j edgelst)
(setq triangle (cons pt (nth j edgelst))
trianglelst (cons (list triangle nil) trianglelst)
)
)
(setq j (1+ j))
)
trianglelst
)
; PURGETRIANGLELST ;
; ;
; replace all triangles that share a vertex with supertriangle ;
; ;
(defun PURGETRIANGLELST (trianglelst supertriangle fuzzy /)
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(if (apply 'or
(mapcar '(lambda (x) (equalmember x supertriangle fuzzy))
triangle
)
)
(setq trianglelst (nth_del j trianglelst))
(setq j (1+ j))
)
)
)
; ;
; THINKING - STANDARD PROGRESS SPINNER ;
; ;
(defun THINKING (prmpt)
(setq THINK-CNT (1+ THINK-CNT))
(princ (strcat "\r" (nth (rem THINK-CNT 4) '("\|" "\/" "\-" "\\")) prmpt))
)
; ********************************* END OF CODING *******************************************
(princ "\n'TRIANGULATE' Loaded \n")
(princ)
;******************************************************************;
; TRIANGULATE - Lisp command to create a TIN from 3D points. ;
; =========== ;
; ;
; Written by Daniele Piazza, ADN member Mechanical Solution s.r.l. ;
; http://pdcode.com/code.htm ;
; ;
; Original C coding "Triangulate" written by PAUL BOURKE ;
; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/ ;
; ;
; This program triangulates an irregular set of points. ;
; You can replace some code (sorting, list manipulation,...) with ;
; VLisp functions to reduce the execution time. ;
; ;
; This code is not seriously tested, if you find a bug...sorry!! ;
; Goodbye, Daniele ;
;*******************************************************************
;;
;;
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;
;; Change by ymg 05/19/2011
;;Modified FINDSUPERTRIANGLE and PURGETRIANGLELST
;;Remove recursion in NTH_SUBST
(defun C:TRIANGULATE (/ fuzzy nulllist ss1
nv ptlst trianglelst
i j k circle
pt flag perc edgelst
)
(setq OLDCMD (getvar "CMDECHO"))
(setvar "CMDECHO" 0)
(command ".UNDO" "GROUP")
(setq OLDSNAP (getvar "OSMODE"))
(setvar "OSMODE" 0)
(setq fuzzy 1e-8) ; tolerance in equality test
(setq nulllist nil)
(princ "\nSelect points...")
(setq ss1 (ssget '((0 . "POINT"))))
(setq start (getvar "date")
THINK-CNT 0
)
(setq ptlst (getptlist ss1))
(setq ptlst (xsort ptlst))
(setq nv (length ptlst))
(setq supertriangle (findsupertriangle ptlst)) ; find super triangle
(setq ptlst (append ptlst supertriangle)) ; append coordinates to the end of vertex list
(setq trianglelst (list (list supertriangle nil))) ; add supertriangle to the triangle list
(setq i 0)
(while (< i nv)
(THINKING (strcat "Processing TIN - " (itoa (/ (* i 100) nv)) "% "))
(setq pt (nth i ptlst))
(setq edgelst nil) ; initialize edge buffer
(setq j 0)
(while
(and trianglelst (setq triangle (car (nth j trianglelst))))
(setq flag T)
(if (not (cadr (nth j trianglelst)))
(progn
(setq circle (getcircircumcircle triangle)) ; calculate circumcircle
(if (< (+ (caar circle) (cadr circle)) (car pt)) ; test point x and (pt) location
(setq trianglelst (nth_subst j (list (car (nth j trianglelst)) T) trianglelst))
)
(if (isinside pt circle)
(setq edgelst (addtriangleedges triangle edgelst)
trianglelst (nth_del j trianglelst)
flag nil
)
)
) ; end progn
) ; end if
(if flag
(setq j (1+ j))
)
)
(setq edgelst (removedoublyedges edgelst fuzzy nulllist))
(setq trianglelst (addnewtriangles pt edgelst trianglelst))
(setq i (1+ i))
)
(setq trianglelst (purgetrianglelst trianglelst supertriangle fuzzy)) ; remove triangles with supertriangles edges
(foreach triangle (mapcar 'car trianglelst)
(drawtriangle triangle)
)
(setvar "OSMODE" OLDSNAP)
(setq OLDSNAP nil)
(command ".UNDO" "END")
(setq stop (getvar "date"))
(princ (strcat "\r TIN Complete - Elapsed time: "
(rtos (* 86400.0 (- stop start)) 2 2)
" secs."
)
)
(setvar "CMDECHO" OLDCMD)
(princ)
)
; XSORT - Original Shell Sort function replaced with VLISP sort (much quicker :-) ;
; ;
(defun XSORT ( PTLST /)
(vl-sort PTLST (function (lambda (a b) (< (car a) (car b)))))
)
; NTH_DEL ;
; ;
; delete the n item in the list (by position, not by value!!) ;
; ;
; Elimina l'oggetto che si trova nella posizione N della lista LST. L'utilizzo di ;
; funzioni ricorsive,oltre a non assicurare maggiore velocità, può creare problemi;
; di overflow dello stack in caso di liste molto lunghe. ;
(defun NTH_DEL (N LST / l)
(repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l)(cdr lst))
)
; NTH_SUBST ;
; ;
; Replace the index element in the list with new element. This function is ;
; recursive this is not a great solution with a large amount of data. ;
; ;
;(defun NTH_SUBST (index new Alist)
; (cond
; ((minusp index) Alist)
; ((zerop index) (cons new (cdr Alist)))
; (T
; (cons (car Alist) (nth_subst (1- index) new (cdr Alist)))
; )
; )
;)
;;;Removed the recursion ymg --------------------------
(defun NTH_SUBST (index new Alist / temp)
(cond
((minusp index) Alist)
((zerop index)(cons new (cdr Alist)))
((> index 0) (while (> index 0)
(setq temp (cons (car alist) temp)
alist (cdr alist)
index (1- index)
)
)
(append (reverse temp) (list new) (cdr alist))
)
)
)
; GETPTLIST ;
; ;
; sset -> list (p1 p2 p3 ... pn) ;
; ;
(defun GETPTLIST (ss1 / i pt ptlst)
(if (not (zerop (sslength ss1)))
(progn
(setq i 0)
(while (setq pt (ssname ss1 i))
(setq ptlst (cons (cdr (assoc 10 (entget pt))) ptlst))
(setq i (1+ i))
)
)
)
ptlst
)
; FINDSUPERTRIANGLE ;
; ;
; Search the supertriangle that contain all points in the data set ;
; ;
(defun FINDSUPERTRIANGLE (ptlst / xmax xmin ymax ymin
dx dy dmax xmid ymid trx1
trx2 trx3 try1 try2 try3 trz1
trz2 trz3
)
(setq xmax (apply 'max (mapcar 'car ptlst))
xmin (apply 'min (mapcar 'car ptlst))
ymax (apply 'max (mapcar 'cadr ptlst))
ymin (apply 'min (mapcar 'cadr ptlst))
dx (- xmax xmin)
dy (- ymax ymin)
dmax (max dx dy)
xmid (* (+ xmax xmin) 0.5)
ymid (* (+ ymax ymin) 0.5)
trx1 (- xmid (* dmax 20.0)) ;modified ymg
try1 (- ymid dmax)
trz1 0.0
trx2 xmid
try2 (+ ymid (* dmax 20.0)) ;modified ymg
trz2 0.0
trx3 (+ xmid (* dmax 20.0)) ;modified ymg
try3 (- ymid dmax)
trz3 0.0
)
(list (list trx1 try1 trz1)
(list trx2 try2 trz2)
(list trx3 try3 trz3)
)
)
;;=============================================================
;; Changes by CAB 03/13/06
;; replaced the GETCIRCIRCUMCIRCLE routine
;;=============================================================
(defun getcircircumcircle (triangle / p1 p2 p3 pr1 pr2 cen rad bisector)
;; return a pt list for a perpendicular bisector 20 units long
(defun bisector (p1 p2 / perp_ang midpt)
(setq p1 (list (car p1) (cadr p1)) ; make sure 2d point
perp_ang (+ (angle p1 p2) (/ pi 2.0))) ; perpendicular angle
(setq midpt (mapcar '(lambda (pa pb) (+ (/ (- pb pa) 2.0) pa)) p1 p2))
(list (polar midpt perp_ang 10) (polar midpt (+ pi perp_ang) 10))
)
(setq p1 (car triangle)
p2 (cadr triangle)
p3 (caddr triangle)
pr1 (bisector p1 p2)
pr2 (bisector p1 p3)
cen (inters (car pr1) (cadr pr1) (car pr2) (cadr pr2) nil)
rad (distance cen p1)
)
(list cen rad)
)
;;=============================================================
; ISINSIDE ;
; ;
; test if pt is inside a circle ;
; ;
(defun ISINSIDE (pt circle)
(setq ctr (car circle)
rad (cadr circle)
)
(< (distance pt ctr) rad)
)
; ADDTRIANGLEEDGES ;
; ;
; add triangle edges at the edge queue ;
; ;
(defun ADDTRIANGLEEDGES (triangle edgelst)
(append edgelst (list (list (car triangle) (cadr triangle))
(list (cadr triangle) (caddr triangle))
(list (caddr triangle)(car triangle))
)
)
)
; DRAWTRIANGLE ;
; ;
; the fun side of the algorithm. Draw triangulation. ;
; ;
(defun DRAWTRIANGLE (triangle)
(entmake (list (cons 0 "3DFACE") (cons 10 (car triangle)) (cons 11 (caddr triangle))
(cons 12 (cadr triangle)) (cons 13 (cadr triangle))
)
)
)
; EQUALMEMBER ;
; ;
; Check if "item" is in "lista" or not by equality test. With real number the ;
; standard fuction "member" not work correctly. ;
; ;
(defun EQUALMEMBER (item lista fuzzy /)
(apply 'or (mapcar '(lambda (x) (equal x item fuzzy)) lista))
)
; REMOVEDOUBLYEDGES ;
; ;
; Test the edge queue to remove duplicates (warning CW & CCW!) ;
; ;
(defun REMOVEDOUBLYEDGES (edgelst fuzzy nulllist /)
(setq j 0)
(while (< j (length edgelst))
(setq k (1+ j))
(while (< k (length edgelst))
(if
(or (and (equal (car (nth j edgelst)) (car (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
)
(and (equal (car (nth j edgelst)) (cadr (nth k edgelst)) fuzzy)
(equal (cadr (nth j edgelst)) (car (nth k edgelst)) fuzzy)
)
)
(setq edgelst (nth_subst j nulllist edgelst)
edgelst (nth_subst k nulllist edgelst)
)
)
(setq k (1+ k))
)
(setq j (1+ j))
)
edgelst
)
; ADDNEWTRIANGLES ;
; ;
; Add new triangle generated by pt to triangle list. ;
; ;
(defun ADDNEWTRIANGLES (pt edgelst trianglelst / j triangle)
(setq j 0)
(while (< j (length edgelst))
(if (nth j edgelst)
(setq triangle (cons pt (nth j edgelst))
trianglelst (cons (list triangle nil) trianglelst)
)
)
(setq j (1+ j))
)
trianglelst
)
; PURGETRIANGLELST ;
; ;
; replace all triangles that share a vertex with supertriangle ;
; ;
(defun PURGETRIANGLELST (trianglelst supertriangle fuzzy /)
(setq j 0)
(while (and trianglelst (setq triangle (car (nth j trianglelst))))
(if
(apply 'or
(mapcar '(lambda (x) (equalmember x supertriangle fuzzy))
triangle
)
)
(setq trianglelst (nth_del j trianglelst))
(setq j (1+ j))
)
)
(setq trianglelst trianglelst); modified to return trianglelst was returning J in certain case ymg
)
; ;
; THINKING - STANDARD PROGRESS SPINNER ;
; ;
(defun THINKING (prmpt)
(setq THINK-CNT (1+ THINK-CNT))
(princ (strcat "\r" (nth (rem THINK-CNT 4) '("\|" "\/" "\-" "\\")) prmpt))
)
; ********************************* END OF CODING *******************************************
(princ "\n'TRIANGULATE' Loaded \n")
(princ)((1652.17 356.759 446.623) (1666.15 431.163 -353.053) (1688.64 379.861 -372.616) (1708.17 888.849 489.959)
(1763.96 799.643 117.206) (1811.9 678.149 -387.295) (1818.56 140.657 -256.432) (1883.13 226.078 -79.1498)
(1888.23 124.665 122.761) (1900.26 864.281 -41.0016) (1950.15 730.671 -164.785) (1979.73 671.496 -19.8523)
(2031.64 260.656 -497.925) (2069.42 69.732 -278.069) (2071.19 123.139 -183.401) (2096.73 383.737 -280.053)
(2173.55 135.927 283.044) (2241.51 1048.67 47.7767) (2298.4 460.399 211.447) (2304.6 871.301 -156.27)
(2441.41 517.957 -411.649) (2455.6 695.636 -390.896) (2462.35 249.15 99.3225) (2585.43 387.857 201.498)
(2591.77 477.032 100.238) (-17456.9 -419.739 0.0) (2121.97 20138.0 0.0) (21700.8 -419.739 0.0))Hello :
In other word we need to create more than just the 3dface.
With such a structure it becomes easy to walk on the TIN, go to the neighbor of a 3dface , contour etc.
So far none of the triangulation that we have do it well.
ymg
;Create a contour slice through a 3 vertex facet pa, pb, pc
;The contour "level" is a horizontal plane perpendicular to the z axis,
;ie: The equation of the contour plane Ax + By + Cz + D = 0
;has A = 0, B = 0, C = 1, D = -level
; Return:
; ((nil nil nil) (nil nil nil) 0) if the contour plane doesn't cut the facet
; ( (x y z) (x y z) 2) if it does cut the facet
; ((nil nil nil) (nil nil nil) -1) for an unexpected occurrence
;If a vertex touches the contour plane nothing need to be drawn
(defun contourfacet (triangle level / pa pb pc sidea sideb sidec p1x p1y p1z p2x p2y p2z)
(setq pa (car triangle)
pb (cadr triangle)
pc (caddr triangle)
)
(setq sidea (-(caddr pa)level)
sideb (-(caddr pb)level)
sidec (-(caddr pc)level)
)
(cond
;---Are all the vertices on one side----------
((and (>= sidea 0) (>= sideb 0) (>= sidec 0))
(setq ret 0)
)
((and (<= sidea 0) (<= sideb 0) (<= sidec 0))
(setq ret 0)
)
;---Is pa the only point on a side by itself---
((and (/=(sign sidea) (sign sideb)) (/=(sign sidea) (sign sidec)))
(setq p1x (- (car pa) (* sidea (/ (- (car pc) (car pa)) (- sidec sidea)))))
(setq p1y (- (cadr pa) (* sidea (/ (- (cadr pc) (cadr pa)) (- sidec sidea)))))
(setq p1z (- (caddr pa) (* sidea (/ (- (caddr pc) (caddr pa)) (- sidec sidea)))))
(setq p2x (- (car pa) (* sidea (/ (- (car pb) (car pa)) (- sideb sidea)))))
(setq p2y (- (cadr pa) (* sidea (/ (- (cadr pb) (cadr pa)) (- sideb sidea)))))
(setq p2z (- (caddr pa) (* sidea (/ (- (caddr pb) (caddr pa)) (- sideb sidea)))))
(setq ret 2)
)
;---Is pb the only point on a side by itself---
((and (/=(sign sideb) (sign sidea)) (/=(sign sideb) (sign sidec)))
(setq p1x (- (car pb) (* sideb (/ (- (car pc) (car pb)) (- sidec sideb)))))
(setq p1y (- (cadr pb) (* sideb (/ (- (cadr pc) (cadr pb)) (- sidec sideb)))))
(setq p1z (- (caddr pb) (* sideb (/ (- (caddr pc) (caddr pb)) (- sidec sideb)))))
(setq p2x (- (car pb) (* sideb (/ (- (car pa) (car pb)) (- sidea sideb)))))
(setq p2y (- (cadr pb) (* sideb (/ (- (cadr pa) (cadr pb)) (- sidea sideb)))))
(setq p2z (- (caddr pb) (* sideb (/ (- (caddr pa) (caddr pb)) (- sidea sideb)))))
(setq ret 2)
)
;---Is pc the only point on a side by itself---
((and (/=(sign sidec) (sign sidea)) (/=(sign sidec) (sign sideb)))
(setq p1x (- (car pc) (* sidec (/ (- (car pa) (car pc)) (- sidea sidec)))))
(setq p1y (- (cadr pc) (* sidec (/ (- (cadr pa) (cadr pc)) (- sidea sidec)))))
(setq p1z (- (caddr pc) (* sidec (/ (- (caddr pa) (caddr pc)) (- sidea sidec)))))
(setq p2x (- (car pc) (* sidec (/ (- (car pb) (car pc)) (- sideb sidec)))))
(setq p2y (- (cadr pc) (* sidec (/ (- (cadr pb) (cadr pc)) (- sideb sidec)))))
(setq p2z (- (caddr pc) (* sidec (/ (- (caddr pb) (caddr pc)) (- sideb sidec)))))
(setq ret 2)
)
(t (setq ret -1))
)
(list (list p1x p1y p1z) (list p2x p2y p2z) ret)
)
;--Define the Signum function-----
(defun sign (x)
( cond
(( minusp x ) -1 )
(( zerop x) 0 )
( t 1 )
)
)Now, especially checked - Playing a random mistake. I accidentally deleted a few lines. The code is similar to the line and the code is understandable, what is missing! If, indeed, no one can understand this program, I promise next week to lay out the full code! But I will be sorry that such a code just copy and do not try to understand...hmmm, your variables naming style scares us off, Evgeniy :)
:-o
(vl-remove-if-not
(function (lambda (a) (and (< mi (caadr a) ma) (< mi (caaddr a) ma))))
l2
)
to this(vl-remove-if
(function
(lambda (a)
(or (vl-some (function (lambda (c) (vl-position c s))) a) (null a))
)
)
l2
)
to make it workSofito_Soft,Hello:
That's an interesting way to do it.
How fast would it be if we wanted to generate a contour map from it ?
I will look at it.
Thanks,
ymg
Sofito_soft,Hello :
That is a lot of time in order to unite al those 3dface.
Whereby if we would preserve the data structure of the triangulation in a double edged list it would take very
little time to contour. Most of the work is already done by the triangulation.
ymg
Hello :
That's true. If you just want to contour, no need to tape "3DFACE" as each 3DFACE have what it takes to calculate its intersection with a plane parallel to XY with Z determined. I think what I have done. I'll look and tell you something.
Sofito_Soft,QuoteHello :
That's true. If you just want to contour, no need to tape "3DFACE" as each 3DFACE have what it takes to calculate its intersection with a plane parallel to XY with Z determined. I think what I have done. I'll look and tell you something.
What is more difficult is to trace that contour from face to face. We do not want to go by brute force
because It's gonna take way too long.
With my prehistoric face to face LSP:
Command: 3DF->CURVA_NIVEL
Selecciona las 3DFACE
Select objects: Specify opposite corner: 3073 found
Select objects:
Cronometro a 00:00
Tiempo para la funcion "c:3df->curva_nivel" : 20.985 seg.
With a intervale of 0.25 m in Z.-> 33000 contour level lines (separates )
Greetings.
Tracing contours between triangles may be an inefficient process, since triangles that span multiple elevation levels will have to be processed multiple times. It may be more efficient to process all the contour segments for each triangle, then process those segments into buckets (for lack of a better term) for later processing. Optimization would likely involve a pre-sort of the triangles based on the lowest (bottom to top) or highest (top to bottom). Or, for even more optimization creation of the contour segments as each triangle is created.
What is more difficult is to trace that contour from face to face. We do not want to go by brute force
because It's gonna take way too long.
With my prehistoric face to face LSP:
If you put the polilines in Z, then you can filter to select all under a certain height. It is easy then weld all.
In my LSP
( defun c:3df->SEG_curva_nivel ( / ) : 3dfaces a segment of polilines in deteminate Z
(defun c:test (/ I L S tmp a b )
(princ (strcat "\n select points"))
(if (setq i 0
s (ssget '((0 . "POINT")))
) ;_ setq
(progn (repeat (sslength s)
( setq tmp (cdr (assoc 10 (entget (ssname s i)))))
( setq l (cons (list (/ (fix (* (car tmp ) 1000)) 1000. )
(/ (fix (* (cadr tmp) 1000)) 1000. )
(caddr tmp)
)
l)
i (1+ i)
) ;_ setq
) ;_ repeat
( setq l ( remove_doubles l ) )
( setq l (vl-sort
l
(function (lambda (a b) (>= (car a) (car b))))
) ;_ vl-sort
)
( eea-delone-triangulate l )
) ;_ progn
) ;_ if
) ;_ defun
(defun eea-delone-triangulate ( L / A A1 A2 A3
I I2 L1 L2 L3 LP MA
MI P S TI TR X1 X2
Y1 Y2
)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; Program triangulate an irregular set of 3d points.
;;
;;*********************************************************
(if l
(progn
(setq ti (car (_VL-TIMES))
i1 ( length l )
i 1
i1 (/ i1 100.)
i2 0
x2 (caar l)
y1 (cadar l)
y2 y1
) ;_ setq
(while l
(setq a (fix (caar l))
a1 (list (car l))
; l (cdr l) ; <<<<<<<<<<<<<<<<<<<<<<< !!!!
) ;_ setq
.......
(defun c:test (/ I L S)
(princ (strcat "\n select points"))
(if (setq i 0
s (ssget '((0 . "POINT")))
) ;_ setq
(progn (repeat (sslength s)
(setq l (cons (cdr (assoc 10 (entget (ssname s i)))) l)
i (1+ i)
) ;_ setq
) ;_ repeat
(eea-delone-triangulate i l)
) ;_ progn
) ;_ if
) ;_ defun
(defun eea-delone-triangulate
(i1 L / A A1 A2 A3 I I2 L1 L2 L3 LP MA MI P S TI TR X1 X2 Y1 Y2)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; edit 20.05.2011
;; Program triangulate an irregular set of 3d points.
;;
;;*********************************************************
(if l
(progn
(setq ti (car (_VL-TIMES))
i 1
i1 (/ i1 100.)
i2 0
l (vl-sort (mapcar (function (lambda (p)
(list (/ (fix (* (car p) 1000)) 1000.)
(/ (fix (* (cadr p) 1000)) 1000.)
(caddr p)
) ;_ list
) ;_ lambda
) ;_ function
l
) ;_ mapcar
(function (lambda (a b) (>= (car a) (car b))))
) ;_ vl-sort
x2 (caar l)
y1 (cadar l)
y2 y1
) ;_ setq
(while l
(setq a2 (car l))
(if (<= (cadr a2) y1)
(setq y1 (cadr a2))
(if (> (cadr a2) y2)
(setq y2 (cadr a2))
)
)
(setq a (fix (caar l))
a1 (list (car l))
l (cdr l)
) ;_ setq
(while (and l (= (fix (caar l)) a))
(setq a2 (car l))
(if (<= (cadr a2) y1)
(setq y1 (cadr a2))
(if (> (cadr a2) y2)
(setq y2 (cadr a2))
) ;_ if
) ;_ if
(setq a1 (cons (car l) (vl-remove a2 a1))
l (cdr l)
) ;_ setq
) ;_ while
(foreach a a1 (setq lp (cons a lp)))
) ;_ while
(setq x1 (caar lp)
a (list (/ (+ x1 x2) 2) (/ (+ y1 y2) 2))
a1 (distance a (list x1 y1))
ma (+ (car a) a1 a1)
mi (- (car a) a1)
s (list (list ma (cadr a) 0)
(list mi (+ (cadr a) a1 a1) 0)
(list (- (car a) a1) (- (cadr a) a1 a1) 0)
) ;_ list
l (list (cons x2 (cons a (cons (+ a1 a1) s))))
ma (1- ma)
mi (1+ mi)
) ;_ setq
(while lp
(setq p (car lp)
lp (cdr lp)
l1 nil
) ;_ setq
(while l
(setq tr (car l)
l (cdr l)
) ;_ setq
(cond ((< (car tr) (car p)) (setq l2 (cons (cdddr tr) l2)))
((< (distance p (cadr tr)) (caddr tr))
(setq tr (cdddr tr)
a1 (car tr)
a2 (cadr tr)
a3 (caddr tr)
l1 (cons (list (+ (car a1) (car a2)) (+ (cadr a1) (cadr a2)) a1 a2)
(cons (list (+ (car a2) (car a3)) (+ (cadr a2) (cadr a3)) a2 a3)
(cons (list (+ (car a3) (car a1)) (+ (cadr a3) (cadr a1)) a3 a1) l1)
) ;_ cons
) ;_ cons
) ;_ setq
)
(t (setq l3 (cons tr l3)))
) ;_ cond
) ;_ while
(setq l l3
l3 nil
l1 (vl-sort l1
(function (lambda (a b)
(if (= (car a) (car b))
(<= (cadr a) (cadr b))
(< (car a) (car b))
) ;_ if
) ;_ lambda
) ;_ function
) ;_ vl-sort
) ;_ setq
(while l1
(if (and (= (caar l1) (caadr l1)) (= (cadar l1) (cadadr l1)))
(setq l1 (cddr l1))
(setq l (cons (eea-data-triangle p (cddar l1)) l)
l1 (cdr l1)
) ;_ setq
) ;_ if
) ;_ while
(if (and (< (setq i (1- i)) 1) (< i2 100))
(progn
(setvar
"MODEMACRO"
(strcat
" "
(itoa (setq i2 (1+ i2)))
" % "
(substr
"||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||"
1
i2
) ;_ substr
(substr "..." 1 (- 100 i2))
) ;_ strcat
) ;_ setvar
(setq i i1)
) ;_ progn
) ;_ if
) ;_ while
(foreach a l (setq l2 (cons (cdddr a) l2)))
(setq l2 (vl-remove-if-not
(function (lambda (a) (and (< mi (caadr a) ma) (< mi (caaddr a) ma))))
l2
) ;_ vl-remove-if
) ;_ setq
(foreach a l2
(entmake (list (cons 0 "3DFACE")
(cons 10 (car a))
(cons 11 (car a))
(cons 12 (cadr a))
(cons 13 (caddr a))
) ;_ list
) ;_ entmake
) ;_ foreach
) ;_ progn
) ;_ if
(setvar "MODEMACRO" "")
(princ (strcat "\n " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) " secs."))
(princ)
) ;_ defun
(defun eea-data-triangle (P1 l / A A1 P2 P3 P4 S)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; Calculation of the centre of a circle and circle radius
;; for program triangulate
;;
;; (eea-data-triangle (getpoint)(list(getpoint)(getpoint)))
;;*********************************************************
(setq p2 (car l)
p3 (cadr l)
p4 (list (car p3) (cadr p3))
) ;_ setq
(if (not (zerop (setq s (sin (setq a (- (angle p2 p4) (angle p2 p1)))))))
(progn (setq a (polar p4
(+ -1.570796326794896 (angle p4 p1) a)
(setq a1 (/ (distance p1 p4) s 2.))
) ;_ polar
a1 (abs a1)
) ;_ setq
(list (+ (car a) a1) a a1 p1 p2 p3)
) ;_ progn
) ;_ if
) ;_ defunprogram with the corrections:
You are still missing triangles on the convex hull Evgenyi.
ymg
((1652.17 356.759 446.623) (1666.15 431.163 -353.053) (1688.64 379.861 -372.616) (1708.17 888.849 489.959)
(1763.96 799.643 117.206) (1811.9 678.149 -387.295) (1818.56 140.657 -256.432) (1883.13 226.078 -79.1498)
(1888.23 124.665 122.761) (1900.26 864.281 -41.0016) (1950.15 730.671 -164.785) (1979.73 671.496 -19.8523)
(2031.64 260.656 -497.925) (2069.42 69.732 -278.069) (2071.19 123.139 -183.401) (2096.73 383.737 -280.053)
(2173.55 135.927 283.044) (2241.51 1048.67 47.7767) (2298.4 460.399 211.447) (2304.6 871.301 -156.27)
(2441.41 517.957 -411.649) (2455.6 695.636 -390.896) (2462.35 249.15 99.3225) (2585.43 387.857 201.498)
(2591.77 477.032 100.238) (-17456.9 -419.739 0.0) (2121.97 20138.0 0.0) (21700.8 -419.739 0.0))program with the corrections:can you provide a test points list which will show the difference between old and new versions?
can you provide a test points list which will show the difference between old and new versions?
(foreach p '((2516.06 458.26 0.0)
(2520.33 508.832 0.0)
(2563.56 516.802 0.0)
(2562.18 470.832 0.0)
(2594.14 435.458 0.0)
(2614.14 445.705 0.0)
(2611.41 488.645 0.0)
)
(entmakex (list '(0 . "point") (cons 10 p)))
)
;*****************************************************************************;
; TRIANGULATE - Lisp command to create a TIN from 3D points. ;
; =========== ;
; ;
; Written by Daniele Piazza, ADN member Mechanical Solution s.r.l. ;
; http://pdcode.com/code.htm ;
; ;
; Original C coding "Triangulate" written by PAUL BOURKE ;
; http://astronomy.swin.edu.au/~pbourke/modelling/triangulate/ ;
; ;
; This program triangulates an irregular set of points. ;
; You can replace some code (sorting, list manipulation,...) with ;
; VLisp functions to reduce the execution time. ;
; ;
; This code is not seriously tested, if you find a bug...sorry!! ;
; Goodbye, Daniele ;
;*****************************************************************************;
;; ;
;; ;
;; Changes by CAB 03/13/06 ;
;; Replaced the GETCIRCIRCUMCIRCLE routine ;
;; ;
;; Changes by ymg : ;
;; ;
;;Modified FINDSUPERTRIANGLE 19/05/2011 ;
;;and PURGETRIANGLELST 19/05/2011 ;
;; ;
;;Removed recursion in NTH_SUBST 19/05/2011 ;
;; ;
;;Reverted to original GETCIRCIRCUMCIRCLE routine and changed it's ;
;;name to GETCIRCUMCIRCLE 22/05/2011 ;
;; ;
;;Modified so that trianglelst and edgelst are now list of ;
;;indices into ptlst. 25/05/2011 ;
;; ;
;;Removed EQUALMEMBER 25/05/2011 ;
;; ;
;;For Contour generation: ;
;; ;
;;Added GETEDGELST 28/05/2011 ;
;;Added GETCROSSEDEGDGE 28/05/2011 ;
;;Added MAKE_CONTOUR 28/05/2011 ;
;; ;
(defun C:TRIANGULATE (/ ss1 nv
i j k
circle pt flag
ptlst edgelst trianglelst
oldcmd oldsnap supertriangle
triangle start think-cnt
intv zmax zmin
)
(setq OLDCMD (getvar "CMDECHO"))
(setvar "CMDECHO" 0)
(command ".UNDO" "GROUP")
(setq OLDSNAP (getvar "OSMODE"))
(setvar "OSMODE" 0)
(princ "\nSelect points...")
(setq ss1 (ssget '((0 . "POINT"))))
(setq start (car (_VL-TIMES))
THINK-CNT 0
)
(setq ptlst (getptlist ss1))
(setq ptlst (xsort ptlst))
(setq nv (length ptlst))
(setq supertriangle (findsupertriangle ptlst))
(setq ptlst (append ptlst supertriangle))
(setq supertriangle (list nv (1+ nv) (+ 2 nv)))
(setq trianglelst (list(list supertriangle nil)))
(setq i 0)
(while (< i nv)
(THINKING (strcat "Processing TIN - " (itoa (/ (* i 100) nv)) "% ") think-cnt)
(setq pt (nth i ptlst))
(setq edgelst nil)
(setq j 0)
(while (and trianglelst (setq triangle (car(nth j trianglelst))))
(setq flag T)
(if (not (cadr (nth j trianglelst)))
(progn
(setq circle (GETCIRCUMCIRCLE triangle ptlst))
(cond
((< (+ (caar circle) (cadr circle)) (car pt))
(setq trianglelst (nth_subst j (list (car (nth j trianglelst)) T) trianglelst))
)
((isinside pt circle)
(setq edgelst (addtriangleedges triangle edgelst)
trianglelst (nth_del j trianglelst)
flag nil
)
)
)
)
)
(if flag (setq j (1+ j)))
)
(setq edgelst (removedoublyedges edgelst))
(setq trianglelst (addnewtriangles i edgelst trianglelst))
(setq i (1+ i))
)
(setq trianglelst (purgetrianglelst trianglelst supertriangle))
(repeat 3
(setq ptlst (vl-remove (last ptlst) ptlst))
)
(foreach triangle (mapcar 'car trianglelst)
(drawtriangle triangle ptlst)
)
;; Now we set-up to trace the contour ;
(setq edgelst (GETEDGELST trianglelst)
zmax (apply 'max (mapcar 'caddr ptlst))
zmin (apply 'min (mapcar 'caddr ptlst))
intv 1
zmin (+ (fix zmin) intv)
zmax (fix zmax)
)
(MAKE_CONTOUR zmin zmax intv edgelst ptlst)
(setvar "OSMODE" OLDSNAP)
(setq OLDSNAP nil)
(command ".UNDO" "END")
(princ (strcat "\r TIN Complete - Elapsed time: " (rtos (/ (- (car (_VL-TIMES)) start) 1000.) 2 4) " secs."))
(princ)
(setvar "CMDECHO" OLDCMD)
(princ)
)
;; XSORT ;
;; ;
;; Original Shell Sort function replaced with VLISP sort (much quicker :-) ;
;; ;
(defun XSORT ( PTLST /)
(vl-sort PTLST (function (lambda (a b) (< (car a) (car b)))))
)
; NTH_DEL ;
; ;
; delete the n item in the list (by position, not by value!!) ;
; ;
(defun NTH_DEL (N LST / l)
(repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l) (cdr lst))
)
;(defun NTH_DEL (index lista)
; (setq i -1)
; (apply 'append (mapcar '(lambda (x)
; (if (= (setq i (+ 1 i)) index)
; nil
; (list x)
; )
; )
; lista
; )
; )
;)
; NTH_SUBST ;
; ;
; Replace the index element in the list with new element. ;
; ;
; ;
(defun NTH_SUBST (n new lst / l)
(cond
((minusp n) lst)
((zerop n)(cons new (cdr lst)))
(t (repeat n
(setq l (cons (car lst) l)
lst (cdr lst)
)
)
(append (reverse l) (list new) (cdr lst))
)
)
)
;Need to check that one for speed
;(defun NTH_SUBST (n new lst / i)
; (setq i -1)
; (mapcar '(lambda (x)
; (if (= (setq i (+ 1 i)) n)
; new
; x
; )
; )
; lst
; )
;)
; GETPTLIST ;
; ;
; sset -> list (p1 p2 p3 ... pn) ;
; ;
(defun GETPTLIST (ss1 / i pt ptlst)
(setq i 0)
(if (not (zerop (sslength ss1)))
(progn
(repeat (sslength ss1)
(setq ptlst (cons (cdr (assoc 10 (entget (ssname ss1 i)))) ptlst)
i (1+ i)
)
)
)
)
ptlst
)
; FINDSUPERTRIANGLE ;
; ;
; Search the supertriangle that contain all points in the data set ;
; ;
(defun FINDSUPERTRIANGLE (ptlst / xmax xmin ymax ymin
dmax xmid ymid x1 x2 x3
y1 y2 y3
)
(setq xmax (apply 'max (mapcar 'car ptlst))
xmin (apply 'min (mapcar 'car ptlst))
ymax (apply 'max (mapcar 'cadr ptlst))
ymin (apply 'min (mapcar 'cadr ptlst))
dmax (max (- xmax xmin) (- ymax ymin))
xmid (* (+ xmax xmin) 0.5)
ymid (* (+ ymax ymin) 0.5)
x1 (- xmid (* dmax 20.0)) ;modified ymg
y1 (- ymid dmax)
x2 xmid
y2 (+ ymid (* dmax 20.0)) ;modified ymg
x3 (+ xmid (* dmax 20.0)) ;modified ymg
y3 (- ymid dmax)
)
(list (list x1 y1 0.0) (list x2 y2 0.0) (list x3 y3 0.0))
)
;; GETCIRCUMCIRCLE ;
;; ;
;; Find circle passing through three points ;
;; ;
(defun GETCIRCUMCIRCLE (triangle ptlst / p1 p2 p3 p1x p2x p3x
p1y p2y p3y d xc yc rad)
(setq p1 (nth (car triangle) ptlst)
p2 (nth (cadr triangle) ptlst)
p3 (nth (caddr triangle) ptlst)
p1x (car p1) p1y (cadr p1)
p2x (car p2) p2y (cadr p2)
p3x (car p3) p3y (cadr p3)
d (* 2.0 (+ (* p1y p3x)
(* p2y p1x)
(- (* p2y p3x))
(- (* p1y p2x))
(- (* p3y p1x))
(* p3y p2x)
)
)
xc (/ (+ (* p2y p1x p1x )
(- (* p3y p1x p1x))
(- (* p2y p2y p1y))
(* p3y p3y p1y)
(* p2x p2x p3y)
(* p1y p1y p2y)
(* p3x p3x p1y)
(- (* p3y p3y p2y))
(- (* p3x p3x p2y))
(- (* p2x p2x p1y))
(* p2y p2y p3y)
(- (* p1y p1y p3y))
)
d
)
yc (/ (+ (* p1x p1x p3x)
(* p1y p1y p3x)
(* p2x p2x p1x)
(- (* p2x p2x p3x))
(* p2y p2y p1x)
(- (* p2y p2y p3x))
(- (* p1x p1x p2x))
(- (* p1y p1y p2x))
(- (* p3x p3x p1x))
(* p3x p3x p2x)
(- (* p3y p3y p1x))
(* p3y p3y p2x)
)
d
)
rad (sqrt (+ (* (- p1x xc)(- p1x xc))
(* (- p1y yc)(- p1y yc))
)
)
)
(list (list xc yc) rad)
)
; ISINSIDE ;
; ;
; test if pt is inside a circle ;
; ;
(defun ISINSIDE (pt circle)
(< (distance pt (car circle)) (cadr circle))
)
; ADDTRIANGLEEDGES ;
; ;
; add triangle edges at the edge queue ;
; ;
(defun ADDTRIANGLEEDGES (triangle edgelst)
(append edgelst (list (list (car triangle) (cadr triangle))
(list (cadr triangle) (caddr triangle))
(list (caddr triangle) (car triangle))
)
)
)
; DRAWTRIANGLE ;
; ;
; the fun side of the algorithm. Draw triangulation. ;
; ;
(defun DRAWTRIANGLE (triangle ptlst /)
(entmake (list (cons 0 "3DFACE")
(cons 8 "TIN")
(cons 10 (nth (car triangle) ptlst))
(cons 11 (nth (caddr triangle) ptlst))
(cons 12 (nth (cadr triangle) ptlst))
(cons 13 (nth (cadr triangle) ptlst))
)
)
)
; REMOVEDOUBLYEDGES ;
; ;
; Test the edge queue to remove duplicates (warning CW & CCW!) ;
; ;
(defun REMOVEDOUBLYEDGES (edgelst / j k)
(setq j 0)
(while (< j (length edgelst))
(setq k (1+ j))
(while (< k (length edgelst))
(if
(or (and (= (car (nth j edgelst)) (car (nth k edgelst)))
(= (cadr (nth j edgelst)) (cadr (nth k edgelst)))
)
(and (= (car (nth j edgelst)) (cadr (nth k edgelst)))
(= (cadr (nth j edgelst)) (car (nth k edgelst)))
)
)
(setq edgelst (nth_subst j nil edgelst)
edgelst (nth_subst k nil edgelst)
)
)
(setq k (1+ k))
)
(setq j (1+ j))
)
edgelst
)
; ADDNEWTRIANGLES ;
; ;
; Add new triangle generated by pt to triangle list. ;
; ;
(defun ADDNEWTRIANGLES (pt edgelst trianglelst / j triangle)
(setq j 0)
(repeat (length edgelst)
(if (nth j edgelst)
(setq triangle (append (cons pt (nth j edgelst)))
trianglelst (cons (list triangle nil) trianglelst)
)
)
(setq j (1+ j))
)
trianglelst
)
; PURGETRIANGLELST ;
; ;
; replace all triangles that share a vertex with supertriangle ;
; ;
(defun PURGETRIANGLELST (trianglelst supertriangle / a b j triangle)
(setq j 0)
(repeat (length trianglelst)
(setq triangle (car (nth j trianglelst)))
(if
(apply 'or (mapcar '(lambda (a) (apply 'or (mapcar '(lambda (b) (= b a )) supertriangle)))
triangle
)
)
(setq trianglelst (nth_del j trianglelst))
(setq j (1+ j))
)
)
trianglelst
)
;; GETEDGELST ;
;; ;
;; Create list of EDGES of all triangles in trianglelst. ;
;; ;
(defun GETEDGELST (trianglelst / edgelst neighlst)
(setq edgelst nil)
(foreach tr trianglelst
(setq edgelst (cons (list (caar tr)(cadar tr)) edgelst)
edgelst (cons (list (cadar tr)(caddar tr)) edgelst)
edgelst (cons (list (caddar tr)(caar tr)) edgelst)
)
)
)
;; GETCROSSEDEGDGE ;
;; ;
;;Traverses the edges list and creates a list of edges that are crossed by a ;
;;contour level. It then follows contlst from neigbour to neighbor until it ;
;;reaches the exterior of the TIN or the starting point. ;
;; ;
;;Returns crossedlst which contains list of edges. ;
;; ;
(defun GETCROSSEDEDGE (lev edgelst ptlst / e pos nxt contlst lwplst
z1 z2 crossedlst)
(setq contlst nil)
(foreach e edgelst
(setq z1 (caddr (nth (car e) ptlst))
z2 (caddr (nth (cadr e) ptlst))
)
(if (= z1 lev) (setq z1 (+ z1 1e-8)))
(if (= z2 lev) (setq z2 (+ z2 1e-8)))
(if (or (< z1 lev z2)(> z1 lev z2))
(setq contlst (cons e contlst))
)
)
(setq crossedlst nil)
(while contlst
(setq pos 0)
;;Find an edge on the convex hull and start from it if none start at last pos in contlst
(while (and (member (reverse (nth pos contlst)) contlst) (< pos (1- (length contlst))))
(setq pos (1+ pos))
)
(setq lwplst (list (nth pos contlst))
contlst (nth_del pos contlst)
pos (- pos (rem pos 2))
lwplst (cons (nth pos contlst) lwplst)
nxt (reverse (car lwplst))
contlst (nth_del pos contlst)
)
(while (and (setq pos (vl-position nxt contlst)) (not (member nxt lwplst)))
(cond
((> (length contlst) 1)
(setq contlst (nth_del pos contlst)
pos (- pos (rem pos 2))
lwplst (cons (nth pos contlst) lwplst)
nxt (reverse (car lwplst))
contlst (nth_del pos contlst)
)
)
(t (setq lwplst (cons (nth pos contlst) lwplst)
contlst (nth_del pos contlst)
)
)
)
)
(setq crossedlst (cons lwplst crossedlst))
)
crossedlst
)
;; MAKE_CONTOUR ;
;; ;
;; Creates LWPolylines from the edgelst ;
;; ;
(defun MAKE_CONTOUR (zmin zmax intv edgelst ptlst /
lwplst lev edge p1 p2 d pt tmp)
(setq lev zmin)
(repeat (fix (/ (- zmax zmin) intv))
(setq lwplst (GETCROSSEDEDGE lev edgelst ptlst))
(foreach pline lwplst
(setq tmp nil)
(foreach edge pline
(setq p1 (nth (car edge) ptlst)
p2 (nth (cadr edge) ptlst)
r (/ (- l (caddr p1)) (- (caddr p2) (caddr p1)));Modified june 2013
p1 (list (car p1)(cadr p1)) ; Calculation of contour was incorrect
p2 (list (car p2)(cadr p2))
d (* (distance p1 p2) r)
pt (polar p1 (angle p1 p2) d)
pt (list (car pt)(cadr pt))
tmp (cons (cons 10 pt) tmp)
)
)
(setq tmp (cons (cons 38 lev) tmp)
tmp (cons (cons 43 0.0) tmp)
tmp (cons (cons 70 0) tmp)
tmp (cons (cons 90 (length pline)) tmp)
tmp (cons (cons 100 "AcDbPolyline") tmp)
tmp (cons (cons 100 "AcDbEntity") tmp)
tmp (cons (cons 8 "Contour") tmp)
tmp (cons (cons 0 "LWPOLYLINE") tmp)
)
(entmake tmp)
)
(setq lev (+ lev intv))
)
)
;; THINKING ;
;; ;
;; Standard progress spinner ;
;; ;
(defun THINKING (prmpt think-cnt)
(setq think-cnt (1+ think-cnt))
(princ (strcat "\r" (nth (rem think-cnt 4) '("\|" "\/" "\-" "\\")) prmpt))
)
; ********************************* END OF CODING *****************************
(princ "\n'TRIANGULATE' Loaded \n")
(princ) ( (1652.17 356.759 16.4)
(1666.15 431.163 16.1)
(1688.64 379.861 16.2)
(1708.17 888.849 13.8)
(1763.96 799.643 14.1)
(1818.56 140.657 18.3)
(1883.13 226.078 18.1)
(1888.23 124.665 18.8)
(1900.26 864.281 12.8)
(1950.15 730.671 12.5)
(1979.73 671.496 12.1)
(2031.64 260.656 15.8)
(2069.42 69.732 22.2)
(2071.19 123.139 19.8)
(2096.73 383.737 15.4)
(2173.55 135.927 20.7)
(2241.51 1048.67 15.1)
(2298.4 460.399 12.3)
(2304.6 871.301 12.0)
(2441.41 517.957 11.5)
(2455.6 695.636 10.0)
(2462.35 249.15 18.3)
(2585.43 387.857 9.8)
(2591.77 477.032 9.5)
)Holas:
Very small differences, but actually looks better.
The best test is made with points that are in one of the sought contour. So, if the new curve (polyline) boundary passes over the points, is that the process can be taken for good.
Eugeny: another "extravaganza" of your algorithm ... please watch the yellow circle in extreme of wellow arrow......
I hope my words do not destroy the great atmosphere of communication and mutual learning on the forum! :?
You've contributed a great piece of code ElpanovEvgeniy, Спасибо.
You've contributed a great piece of code ElpanovEvgeniy, Спасибо.
You wrote this code in another language and your code is much faster...
Also I am sure that this could be optimized <<<< There are 3 or 4 methods ( p.e.geometric general calculation ) that are not highly purified. If you have any interest in something really fast ( productive? ), please write me to e-mail private.
Very good idea! Bravo.
I am manage a process of comprehensive automation of CAD, CAM and CAE.
I have done some testing ...with one MDT based on LWPOLYLINE. I post a micro-video... I fail to generate animated GIF, so herewith a AVI
Greetings
QuoteI am manage a process of comprehensive automation of CAD, CAM and CAE.
Evgenyi,
I got back into triangulation and contouring as part of a hobby of mine which is writing GCode
for router machine.
There exist very good program like Vectric Aspire but they are so pricey!!
If you want to do 3d routing, thr problem is essantially contouring.
Cheers!
ymg
Now I write 3D kernel for the direct management of vertices, edges and faces 3dsolid. All this is necessary in a project of precast ferroconcrete.
You've worked with solids "ACIS"? They are solid autocad. Are well documented.
very very interesting.
Greetings
Eugeny:
I thought rather directly access to 3DSOLID with the "SAT" files .
Giles opened the door, with a wonderful ACISDECODE subfunction.
From there, you can get the points, edges, vertices, faces, normal, holes, etc, of the ACIS solids.
I have read that it is developed also for NET.
It is always a good thing that new developments are compatible with something already established.
I do not know if your head warm.
A greeting. :-)
A very strange sentence!
Individually written program, will always be significantly more expensive than to sell in large quantities. Free software has its price.
Of course developing custom solutions is always more expensive.
However If you are doing it as a hobby, buying every nice program out there can break the bank.
Plus programming is another hobby.
ymg
(defun c:test (/ I L S)
(princ (strcat "\n select points"))
(if (setq i 0
s (ssget '((0 . "POINT")))
) ;_ setq
(progn (repeat (sslength s)
(setq l (cons (cdr (assoc 10 (entget (ssname s i)))) l)
i (1+ i)
) ;_ setq
) ;_ repeat
(eea-delone-triangulate i l)
) ;_ progn
) ;_ if
) ;_ defun
(defun eea-delone-triangulate
(i1 L / A A1 A2 A3 I I2 L1
L2 L3 LP MA MI P S TI TR X1
X2 Y1 Y2
)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; edit 20.05.2011
;; Program triangulate an irregular set of 3d points.
;;
;;*********************************************************
(if l
(progn
(setq ti (car (_VL-TIMES))
i 1
i1 (/ i1 100.)
i2 0
l (vl-sort
(mapcar
(function (lambda (p)
(list (/ (fix (* (car p) 1000)) 1000.)
(/ (fix (* (cadr p) 1000)) 1000.)
(caddr p)
) ;_ list
) ;_ lambda
) ;_ function
l
) ;_ mapcar
(function (lambda (a b) (>= (car a) (car b))))
) ;_ vl-sort
x2 (caar l)
y1 (cadar l)
y2 y1
) ;_ setq
(while l
(setq a2 (car l))
(if (<= (cadr a2) y1)
(setq y1 (cadr a2))
(if (> (cadr a2) y2)
(setq y2 (cadr a2))
)
)
(setq a (fix (caar l))
a1 (list (car l))
l (cdr l)
) ;_ setq
(while (and l (= (fix (caar l)) a))
(setq a2 (car l))
(if (<= (cadr a2) y1)
(setq y1 (cadr a2))
(if (> (cadr a2) y2)
(setq y2 (cadr a2))
) ;_ if
) ;_ if
(setq a1 (cons (car l) (vl-remove a2 a1))
l (cdr l)
) ;_ setq
) ;_ while
(foreach a a1 (setq lp (cons a lp)))
) ;_ while
(setq x1 (caar lp)
a (list (/ (+ x1 x2) 2) (/ (+ y1 y2) 2))
a1 (distance a (list x1 y1))
ma (+ (car a)(* 20 a1));changed ymg
mi (- (car a)(* 20 a1));changed ymg
s (list (list ma (cadr a) 0)
(list mi (+ (cadr a)(* 20 a1)) 0);changed ymg
(list mi (- (cadr a)(* 20 a1)) 0);changed ymg
) ;_ supertriangle
l (list (cons x2 (cons a (cons (+ a1 a1) s))))
;mi (- (car a)(* 20 a1));added ymg
ma (1- ma)
mi (1+ mi)
) ;_ setq
(while lp
(setq p (car lp)
lp (cdr lp)
l1 nil
) ;_ setq
(while l
(setq tr (car l)
l (cdr l)
) ;_ setq
(cond
((< (car tr) (car p)) (setq l2 (cons (cdddr tr) l2)))
((< (distance p (cadr tr)) (caddr tr))
(setq tr (cdddr tr)
a1 (car tr)
a2 (cadr tr)
a3 (caddr tr)
l1 (cons (list (+ (car a1) (car a2))
(+ (cadr a1) (cadr a2))
a1
a2
)
(cons (list (+ (car a2) (car a3))
(+ (cadr a2) (cadr a3))
a2
a3
)
(cons (list (+ (car a3) (car a1))
(+ (cadr a3) (cadr a1))
a3
a1
)
l1
)
) ;_ cons
) ;_ cons
) ;_ setq
)
(t (setq l3 (cons tr l3)))
) ;_ cond
) ;_ while
(setq l l3
l3 nil
l1 (vl-sort l1
(function (lambda (a b)
(if (= (car a) (car b))
(<= (cadr a) (cadr b))
(< (car a) (car b))
) ;_ if
) ;_ lambda
) ;_ function
) ;_ vl-sort
) ;_ setq
(while l1
(if (and (= (caar l1) (caadr l1))
(= (cadar l1) (cadadr l1))
)
(setq l1 (cddr l1))
(setq l (cons (eea-data-triangle p (cddar l1)) l)
l1 (cdr l1)
) ;_ setq
) ;_ if
) ;_ while
(if (and (< (setq i (1- i)) 1) (< i2 100))
(progn
(setvar
"MODEMACRO"
(strcat
" "
(itoa (setq i2 (1+ i2)))
" % "
(substr
"||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||"
1
i2
) ;_ substr
(substr "..." 1 (- 100 i2))
) ;_ strcat
) ;_ setvar
(setq i i1)
) ;_ progn
) ;_ if
) ;_ while
(foreach a l (setq l2 (cons (cdddr a) l2)))
;purge triangle list
(setq
l2 (vl-remove-if-not
(function
(lambda (a)
(and (< mi (caadr a) ma) (< mi (caaddr a) ma))
)
)
l2
) ;_ vl-remove-if
) ;_ setq
(foreach a l2
(entmake (list (cons 0 "3DFACE")
(cons 10 (car a))
(cons 11 (car a))
(cons 12 (cadr a))
(cons 13 (caddr a))
) ;_ list
) ;_ entmake
) ;_ foreach
) ;_ progn
) ;_ if
(setvar "MODEMACRO" "")
(princ (strcat "\n "
(rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4)
" secs."
)
)
(princ)
) ;_ defun
(defun eea-data-triangle (P1 l / A A1 P2 P3 P4 S)
;;*********************************************************
;;
;; Written by ElpanovEvgeniy
;; 17.10.2008
;; Calculation of the centre of a circle and circle radius
;; for program triangulate
;;
;; (eea-data-triangle (getpoint)(list(getpoint)(getpoint)))
;;*********************************************************
(setq p2 (car l)
p3 (cadr l)
p4 (list (car p3) (cadr p3))
) ;_ setq
(if (not
(zerop
(setq s (sin (setq a (- (angle p2 p4) (angle p2 p1)))))
)
)
(progn (setq a (polar p4
(+ -1.570796326794896 (angle p4 p1) a)
(setq a1 (/ (distance p1 p4) s 2.))
) ;_ polar
a1 (abs a1)
) ;_ setq
(list (+ (car a) a1) a a1 p1 p2 p3)
) ;_ progn
) ;_ if
) ;_ defun(defun c:test (/ i pl s)
(princ (strcat "\n select points"))
(if (setq i 0
s (ssget '((0 . "POINT")))
)
(progn (repeat (sslength s)
(setq pl (cons (cdr (assoc 10 (entget (ssname s i)))) pl)
i (1+ i)
)
)
(triangulate pl)
)
)
)
;;********************************************************;
;; ;
;; Written by ElpanovEvgeniy ;
;; 17.10.2008 ;
;; edit 20.05.2011 ;
;; Program triangulate an irregular set of 3d points. ;
;; Modified by ymg June 2011 ;
;;********************************************************;
(defun triangulate (pl / a b c i i1 i2
bb sl al el tl l ma mi
ti tr x1 x2 y1 y2 p r cp
)
(if pl
(progn
(setq ti (car (_VL-TIMES))
i 1
i1 (/ (length pl) 100.)
i2 0
pl (vl-sort pl
(function (lambda (a b) (< (car a) (car b))))
)
bb (list (apply 'mapcar (cons 'min pl))
(apply 'mapcar (cons 'max pl))
)
;Replaced code to get minY and maxY with 3d Bounding Box Routine;
;A bit slower but clearer. minZ and maxZ kept for contouring ;
x1 (caar bb) ;minX
x2 (caadr bb) ;maxX
y1 (cadar bb) ;minY
y2 (cadadr bb) ;maxY
)
(setq cp (list (/ (+ x1 x2) 2.0) (/ (+ y1 y2) 2.0)); Midpoint of points cloud and center point of circumcircle through supertriangle.
r (* (distance cp (list x1 y1)) 20) ;This could still be too small in certain case. No harm if we make it bigger.
ma (+ (car cp) r);ma is maxX of supertriangle
mi (- (car cp) r);mi is minX of supertriangle
sl (list (list ma (cadr cp) 0);list of 3 points defining the supertriangle
(list mi (+ (cadr cp) r) 0)
(list mi (- (cadr cp) r) 0)
)
al (list (cons x2 (cons cp (cons (* 20 r) sl))))
;al is a work list that contains active triangles each item is a list that contains: ;
; item 0: Xmax of points in triangle. ;
; item 1: List 2d coordinates of center of circle circumscribing triangle. ;
; item 2: Radius of above circle. ;
; item 3: List of 3 points defining the triangle ;
ma (1- ma);Reducing ma to prepare for removal of triangles having a vertex
mi (1+ mi);common with supertriangle. Done once triangulation is completed.
) ;Increasing mi for same reason.
;Begin insertion of points
(repeat (length pl)
(setq p (car pl);Get one point from point list
pl (cdr pl);Remove point from point list
el nil ;Initialize edge list
)
(while al ;Active triangle list
(setq tr (car al);Get one triangle from active triangle list.
al (cdr al);Remove the triangle from the active list.
)
(cond
((< (car tr) (car p)) (setq tl (cons (cdddr tr) tl)));This triangle inactive. We store it's 3 vertex in tl (Final triangle list).
((< (distance p (cadr tr)) (caddr tr));p is inside the triangle.
(setq tr (cdddr tr) ;Trim tr to vertex of triangle only.
a (car tr) ; First point.
b (cadr tr) ; Second point.
c (caddr tr) ; Third point.
el (cons (list (+ (car a) (car b)) ;We create 3 new edges ab, bc and ca,
(+ (cadr a) (cadr b)) ;and store them in edge list.
a
b
)
(cons (list (+ (car b) (car c))
(+ (cadr b) (cadr c))
b
c
)
(cons (list (+ (car c) (car a))
(+ (cadr c) (cadr a))
c
a
)
el
)
)
)
)
)
(t (setq l (cons tr l))) ;tr did not meet any cond so it is still active.
) ;we store it a temporary list.
);Go to next triangle of active list.
(setq al l ;Restore active triangle list from the temporary list.
l nil ;Re-initialize the temporary list to prepare for next insertion.
el (vl-sort el ;Sort the edges list on X and Y
(function (lambda (a b)
(if (= (car a) (car b))
(<= (cadr a) (cadr b))
(< (car a) (car b))
)
)
)
)
)
;Removes doubled edges, form new triangles, calculates circumcircles and add them to active list.
(while el
(if (and (= (caar el) (caadr el))
(= (cadar el) (cadadr el))
)
(setq el (cddr el))
(setq al (cons (getcircumcircle p (cddar el)) al)
el (cdr el)
)
)
)
;Spinner to show progress
(if (and (< (setq i (1- i)) 1) (< i2 100))
(progn
(setvar
"MODEMACRO"
(strcat
" "
(itoa (setq i2 (1+ i2)))
" % "
(substr
"||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||"
1
i2
)
(substr "..." 1 (- 100 i2))
)
)
(setq i i1)
)
)
) ;Go to insert next point.
;We are done with the triangulation
(foreach tr al (setq tl (cons (cdddr tr) tl)));What's left in active list is added to triangle list
;Purge triangle list of any triangle that has a common vertex with supertriangle.
(setq
tl (vl-remove-if-not
(function
(lambda (a)
(and (< mi (caadr a) ma) (< mi (caaddr a) ma))
)
)
tl
)
)
;Create a layer and Draw the triangulation
(or (tblsearch "LAYER" "TIN")
(entmake (list
'(0 . "LAYER")
'(100 . "AcDbSymbolTableRecord")
'(100 . "AcDbLayerTableRecord")
'(2 . "TIN")
'(70 . 0)
'(62 . 8)
'(6 . "Continuous")
'(290 . 1)
'(370 . -3)
)
)
)
(setvar "CLAYER" "TIN")
(foreach tr tl
(entmake (list (cons 0 "3DFACE")
(cons 10 (car tr))
(cons 11 (car tr))
(cons 12 (cadr tr))
(cons 13 (caddr tr))
)
)
)
)
)
(setvar "MODEMACRO" "")
(princ (strcat "\n "
(rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4)
" secs."
)
)
(princ)
)
;;*********************************************************;
;; ;
;; Written by ElpanovEvgeniy ;
;; 17.10.2008 ;
;; Calculation of the centre of a circle and circle radius ;
;; for program triangulate ;
;; ;
;; Modified ymg june 2011 (renamed variables) ;
;;*********************************************************;
(defun getcircumcircle (a el / b c c2 cp r ang)
(setq b (car el)
c (cadr el)
c2 (list (car c) (cadr c)) ;c2 is point c but in 2d
)
(if (not
(zerop
(setq ang (- (angle b c) (angle b a)))
)
)
(progn (setq cp (polar c2
(+ -1.570796326794896 (angle c a) ang)
(setq r (/ (distance a c2) (sin ang) 2.0))
)
r (abs r)
)
(list (+ (car cp) r) cp r a b c)
)
)
)Now, the program still seems complicated?
No, since it is basically the same algoritmn as Piazza's program.
The first is GA for TSP .
Can be a little more detail?
I first heard about this algorithm ...
For this program, I developed an algorithm by yourself!]
The Divide and conquer algorithm could be faster than this one if written properly.
do you propose to continue the contest in Lisp?
It is a challenge?
(grread t); Dummy just to force a screen update and avoid that Autocad hangs while calculating
(princ (strcat "\n TIN Completed - Elapsed time: " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) " secs."))
(princ (strcat "\n Generated " (itoa (length tl)) " 3DFACES"))
(princ "\n")
; Begin Calculation of Contour line
;
; The routine follows each contour from edge to edge resulting in Connected LWPOLYLINE.
;
; 1. Sort triangle list tl on the maximum z value of vertices of each triangles.
; 2. Create el, a list containing all edges of all triangles. ((a b)(b c) (c a).......)
; 3. For each desired contour level l, we traverse el and create list cl containing all
; edges crossed by level l. At this point cl contains 2 edges per triangle.
; As we go through el, we can destroy any triplets of edges whose max z value is below
; the current contour level, thus limiting traversal for next level.
; 4. We now process cl, following from edge to edge until either we've closed the polyline
; or reached the convex hull. If we reached the convex hull, we insure that the beginning
; of the polyline is also on the convex hull. Otherwise we reverse the polyline and
; continue processing. Process continues until cl is empty.
; As we complete each polyline list xl is formed.
; 5. We entmake each element of list xl.
; 6. We go back to step 3 until all levels are completed and el is empty.
;
; Step 4 and 5 could be combined but it is easier to follow contour in index form
;
;
; An alternate way to do this would be compute all al segment between two edges joining with
; with a line for all contour level and at end Join everything together.
;
; Setup for Contouring
(setq ti (car (_VL-TIMES)) ; Re-initialize timer for Contouring
i 1
i1 (/ (length pl) 100.)
i2 0
intv 1 ; Interval between contour
zmin (+ (fix zmin) intv) ; zmin was calculated during triangulation
zmax (fix zmax) ; z2 also calculated at beginning
l zmin ; Initialize Contour level
el nil ; el will be list of all edges
vc 0 ; Vertices Count
pc 0 ; LWPOLYLINE Count
)
(setq tl (vl-sort tl
(function
(lambda (a b)
(< (max (caddr (nth (car a) pl))
(caddr (nth (cadr a) pl)) ;Gotta be a more concise way
(caddr (nth (caddr a) pl)) ;to write this probably with apply.
) ;Help!!!
(max (caddr (nth (car b) pl))
(caddr (nth (cadr b) pl))
(caddr (nth (caddr b) pl))
)
)
)
)
)
)
; Extract all triangle edges from tl and form list el
(foreach tr tl
(setq el (cons (list (car tr)(cadr tr)) el)
el (cons (list (cadr tr)(caddr tr)) el)
el (cons (list (caddr tr)(car tr)) el)
)
)
(repeat (+(fix (/ (- zmax zmin) intv)) 1) ;Main Loop through each contour level
(setq cl nil ; cl will be list of all edges crossed at current level l
j 0 ; Index into edge list el
zm 1e-8 ; zmax value for a given triplets of edges
)
(repeat (length el)
(setq e (nth j el)
z1 (caddr (nth (car e) pl)) ; Get elevation of edges from the point list.
z2 (caddr (nth (cadr e) pl))
zm (max zm z1 z2)
j (1+ j)
)
(if (and (= (rem j 3) 0) (< zm (+ l intv))) ; Reduce size of el on zmax criteria.
(setq j (- j 3)
el (vl-remove (nth j el) el)
el (vl-remove (nth j el) el)
el (vl-remove (nth j el) el)
zm 1e-8
)
)
(if (= z1 l) (setq z1 (- z1 1e-8))); If vertex is equal to l we disturb
(if (= z2 l) (setq z2 (- z2 1e-8))); the z value a little.
(if (or (< z1 l z2)(> z1 l z2))
(setq cl (cons e cl)) ; Edge is added to Crossed List
)
);end foreach e
; cl now contains all edges where all contours at level l passes.
(setq xl nil)
(while cl
;We Initialize a Polyline
(setq pol (list (cadr cl)(car cl)) ; We go reverse as we will cons the polyline
nxt (reverse (car pol)) ; nxt will be our next edge
cl (cddr cl) ; Remove first two edges from cl
)
(if (not (member nxt cl)) ;The previous edge was on convex hull
(setq pol (reverse pol) ;We reverse our Polyline
nxt (reverse(car pol)) ;and adjust our next edge
)
)
(while (setq n (vl-position nxt cl))
(setq cl (vl-remove (nth n cl) cl)
n (- n (rem n 2))
pol (cons (nth n cl) pol)
cl (vl-remove (nth n cl) cl)
)
(if (member nxt pol)
(setq nxt nil)
(setq nxt (reverse (car pol)))
)
(if (not (vl-position nxt cl))
(setq pol (reverse pol)
nxt (reverse (car pol))
)
)
);end while
(setq xl (cons pol xl))
);end while cl
(setq pc (+ pc (length xl)))
(foreach p xl
(setq ent nil
vc (+ vc (length p))
)
(if (equal (car p) (reverse (last p)))
(setq isclosed 1
p (cddr p)
)
(setq isclosed 0)
)
(foreach e p
(setq p1 (nth (car e) pl)
p2 (nth (cadr e) pl)
r (/ (- l (caddr p1)) (- (caddr p2) (caddr p1)))
p1 (list (car p1)(cadr p1))
p2 (list (car p2)(cadr p2))
d (* (distance p1 p2) r)
pt (polar p1 (angle p1 p2) d)
ent (cons (cons 10 pt) ent)
)
)
(setq ent (cons (cons 38 l) ent)
ent (cons (cons 43 0.0) ent)
ent (cons (cons 70 isclosed) ent)
ent (cons (cons 90 (length p)) ent)
ent (cons (cons 100 "AcDbPolyline") ent)
ent (cons (cons 100 "AcDbEntity") ent)
ent (cons (cons 8 "Contour") ent)
ent (cons (cons 0 "LWPOLYLINE") ent)
)
(entmake ent)
);end foreach p
(setq l (+ l intv))
);end repeat
(setvar "MODEMACRO" "")
(princ (strcat "\n CONTOUR Completed - Elapsed time: " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) " secs."))
(princ (strcat "\n Generated " (itoa pc) " LWPOLYLINE."))
(princ (strcat "\n Total " (itoa vc) " Vertices."))
(princ)
);end defun triangulate
(foreach p xl
(setq ent nil
vc (+ vc (length p))
)
(if (equal (car p) (reverse (last p)))
(setq isclosed 1
p (cdr p) ;Was wrong (cddr p)
)
(setq isclosed 0)
)
I inadvertently introduce a bug in the last posting of tin.lsp, :ugly:
; Contour by ymg (July 2013)
;
; The routine follows each contour from edge to edge resulting in Connected LWPOLYLINE.
; Basic method is per:
; http://www.originlab.com/www/helponline/Origin/en/UserGuide/Creating_Contour_Graphs.html
;
; Contours are then smoothed using Christensen's Eclectic Method.
; For details see:
; http://www.asprs.org/a/publications/pers/2001journal/april/2001_apr_511-517.pdf
; http://www.google.com.mx/patents/US5333248
;
; General Flow:
;
; 1. Sort triangle list tl on the maximum z value of vertices of each triangles.
;
; 2. Create el, a list containing all edges of all triangles. ((a b)(b c) (c a).......)
;
; 3. For each desired contour level l, we traverse el and create list cl containing all
; edges crossed by level l. At this point cl contains 2 edges per triangle.
; As we go through el, we can destroy any triplets of edges whose max z value is below
; the current contour level, thus limiting traversal for next level.
;
; 4. We now process cl, following from edge to edge until either, we've closed the polyline
; or reached the convex hull. If we reached the convex hull, we insure that the beginning
; of the polyline is also on the convex hull. Otherwise we reverse the polyline and
; continue processing until cl is empty.
; As we complete each polyline list xl is formed and contains a series of list tracing all the
; polylines for a given level.
;
; 5. We entmake each element of list xl, smoothing vertices of polylines when conditions are right.
; Variable mang controls if an edge can be smoothed or not preventing crossing of contours.
; Smoothing is controlled by two variables lang and lres which limits the number of segments
; tracing the parabolas of an edge.
;
; 6. We go back to step 3 until all levels are completed and el is empty.
;
;
;
(defun contour (tl pl / a a1 a2 a3 b bb c cl cn1 cn2 e el ent hfac i i1 i2 intv isclosed
j l lang lres mang n nxt p1 p2 p3 p4 pc pl pol prv seg sm ti
v1 v2 v3 vcr vcs xl z1 z2 zm zmax zmin
)
; Setup for Contouring
(setq ti (car (_VL-TIMES)) ; Re-initialize timer for Contouring
i 1
i1 (/ (length pl) 100.)
i2 0
bb (list
(apply 'mapcar (cons 'min pl))
(apply 'mapcar (cons 'max pl))
) ; bb, bounding box of Point List.
zmin (caddar bb) ; zmin, minimum z of Point List.
zmax (caddr(cadr bb)) ; zmax, maximum z of Point List.
intv 1 ; Interval between contour.
zmin (+ (fix zmin) intv) ; zmin, first level for contour.
zmax (fix zmax) ; zmax; last level for contour.
l zmin ; Initialize Current Contour level.
intv 1 ; Interval between contour.
el nil ; el, will be list of all edges in triangle list.
vcr 0 ; Vertices Count for Raw Contours.
vcs 0 ; Vertices Count for Smoothed Contours.
pc 0 ; LWPOLYLINE Count.
lres 5 ; Linear resolution for smoothing.
lang (dtr 5) ; Angular resolution for smoothing.
mang (* lang 1.5) ; Min angle controls if an edge is suitable for smoothing.
hfac 0.5 ; Degree of smoothing, Max is 0.5, Min would be 1.
)
(setq tl (vl-sort tl
(function
(lambda (a b)
(< (max (caddr (nth (car a) pl))
(caddr (nth (cadr a) pl)) ;Gotta be a more concise way
(caddr (nth (caddr a) pl)) ;to write this probably with apply.
) ;Help!!!
(max (caddr (nth (car b) pl))
(caddr (nth (cadr b) pl))
(caddr (nth (caddr b) pl))
)
)
)
)
)
)
; Extract all triangle edges from tl and form edges list el
(foreach tr tl
(setq el (cons (list (car tr)(cadr tr)) el)
el (cons (list (cadr tr)(caddr tr)) el)
el (cons (list (caddr tr)(car tr)) el)
)
)
(repeat (+(fix (/ (- zmax zmin) intv)) 1) ;Main Loop through each contour level
(setq cl nil ; cl will be list of all edges crossed at current level l
j 0 ; Index into edge list el
zm 1e-8 ; zmax value for a given triplets of edges
)
(repeat (length el)
(setq e (nth j el)
z1 (caddr (nth (car e) pl)) ; Get elevation of edges from the point list.
z2 (caddr (nth (cadr e) pl))
zm (max zm z1 z2)
j (1+ j)
)
(if (and (= (rem j 3) 0) (< zm (+ l intv))) ; Reduce size of el on zmax criteria.
(setq j (- j 3)
el (vl-remove (nth j el) el)
el (vl-remove (nth j el) el)
el (vl-remove (nth j el) el)
zm 1e-8
)
)
(if (= z1 l) (setq z1 (- z1 1e-8))); If vertex is equal to l we disturb
(if (= z2 l) (setq z2 (- z2 1e-8))); the z value a little.
(if (or (< z1 l z2)(> z1 l z2))
(setq cl (cons e cl)) ; Edge is added to Crossed List
)
);end foreach e
; cl now contains all edges where all contours at level l passes.
(setq xl nil)
(while cl
;We Initialize a Polyline
(setq pol (list (cadr cl)(car cl)) ; We go reverse as we will cons the polyline
nxt (reverse (car pol)) ; nxt will be our next edge
cl (cddr cl) ; Remove first two edges from cl
)
(if (not (member nxt cl)) ;The previous edge was on convex hull
(setq pol (reverse pol) ;We reverse our Polyline
nxt (reverse(car pol)) ;and adjust our next edge
)
)
(while (setq n (vl-position nxt cl))
(setq cl (vl-remove (nth n cl) cl)
n (- n (rem n 2))
pol (cons (nth n cl) pol)
cl (vl-remove (nth n cl) cl)
)
(if (member nxt pol)
(setq nxt nil)
(setq nxt (reverse (car pol)))
)
(if (not (vl-position nxt cl))
(setq pol (reverse pol)
nxt (reverse (car pol))
)
)
);end while
(setq xl (cons pol xl))
);end while cl
(setq pc (+ pc (length xl)))
(foreach p xl
(setq ent nil)
(if (equal (car p) (reverse (last p)))
(setq isclosed 1
p (append p (list(cadr p)))
)
(setq isclosed 0
ent (list (cons 10 (clv l (car p) pl)))
)
)
(while (> (length p) 2)
(setq v1 (clv l (car p) pl)
v2 (clv l (cadr p) pl)
)
(setq v3 (clv l (caddr p) pl)
prv (car p)
nxt (caddr p)
p (cdr p)
p1 (nth (caar p) pl)
p3 (nth (cadar p) pl)
p2 (nth (car prv) pl)
p4 (nth (car nxt) pl)
)
(if (or (equal p2 p1) (equal p2 p3))
(setq p2 (nth (cadr prv) pl))
)
(if (or (equal p4 p1) (equal p4 p3))
(setq p4 (nth (cadr nxt) pl))
)
; Need to test here if edge p1 p3 is suitable for smoothing
; Unclear to me, what mang should if I set it to same value as lang, we still have
; Contours crossing each other. With (* lang 1.5) no crossing!!
(cond
((and (or (< (defl p1 p2 p3) mang) (< (defl p1 p4 p3) mang))
(> (/ (abs (- (caddr p1) (caddr p3))) intv) 1)
)
(setq sm (list v2 )); Unsuitable we set v2 as the vertices on this edge.
)
(t (setq cn1 (mapcar '(lambda (a b c) (/ (+ a b c) 3.0)) p1 p2 p3)
cn2 (mapcar '(lambda (a b c) (/ (+ a b c) 3.0)) p1 p3 p4)
a1 (cond
((inters cn1 p1 v1 v2))
((inters cn1 p3 v1 v2))
;((inters cn1 p2 v1 v2))
)
a2 v2
a3 (cond
((inters cn2 p1 v2 v3))
((inters cn2 p3 v2 v3))
;((inters cn2 p4 v2 v3))
)
sm (smooth a1 a2 a3); Probably would be better not to call a function here
); end setq
); end t
); end cond
(foreach e sm
(setq ent (cons (cons 10 e) ent))
)
); end while p
(if (= isclosed 0) (setq ent (cons(cons 10 v3) ent)))
(setq seg (length ent)
vcs (+ vcs seg)
ent (cons (cons 38 l) ent)
ent (cons (cons 43 0.0) ent)
ent (cons (cons 70 isclosed) ent)
ent (cons (cons 90 seg) ent)
ent (cons (cons 100 "AcDbPolyline") ent)
ent (cons (cons 100 "AcDbEntity") ent)
ent (cons (cons 8 "Contour") ent)
ent (cons (cons 0 "LWPOLYLINE") ent)
)
(entmake ent)
);end foreach p
(setq l (+ l intv))
);end repeat
(setvar "MODEMACRO" "")
(princ (strcat "\n CONTOUR - Elapsed time: " (rtos (/ (- (car (_VL-TIMES)) ti) 1000.) 2 4) " secs, "
(itoa pc) " LWPOLYLINES, "
(itoa vcs) " Vertices."
)
)
(princ)
);end defun contour
; Degree to Radian Conversion
(defun dtr (a) (* pi (/ a 180.0)))
; Given: lev, (level of contour being traced).
; ed, (edge defined by a list of two index into the point list.,
; pl, Point List
;
; Returns: Point where contour cross the edge.
(defun clv (lev ed pl / r d p1 p2)
(setq p1 (nth (car ed) pl)
p2 (nth (cadr ed) pl)
r (/ (- lev (caddr p1))
(- (caddr p2) (caddr p1))
)
p1 (list (car p1)(cadr p1))
p2 (list (car p2)(cadr p2))
d (* (distance p1 p2) r)
)
(polar p1 (angle p1 p2) d)
)
; Smooth an edge defined by 3 points
; Returns a list of points starting at a1 and ending at a3
; which traces two parabolas joining a1 to a3
; Global variables lang and lres limits number of segments on parabolas.
; Global variable hfac controls how much smoothing we apply (Not implemented yet)
(defun smooth (a1 a2 a3 / cp f ps tmp)
(setq ps (subdiv a1 a2 a3)
f nil
)
(while ps
(setq cp (car ps)
ps (cdr ps)
)
(while (and (> (min (distance (car cp) (cadr cp))
(distance (cadr cp) (caddr cp))
)
lres
)
(< (defl (car cp) (cadr cp) (caddr cp))
(- pi lang)
)
)
(setq tmp (subdiv (caddr cp)(cadr cp)(car cp))
cp (car tmp)
ps (cons (cadr tmp) ps)
)
) ; end while lres lang
(setq f (cons (car cp) f))
) ; end while ps
(cons a1 f)
)
; Subdivide a pair defined by 3 points (a1 b1 c1)
; Return a list of 6 points (a1 d1 e1 e1 d2 c1)
(defun subdiv (a1 b1 c1 / m1 m2 m3)
(setq m1 (mid a1 b1)
m3 (mid b1 c1)
m2 (mid m1 m3)
)
(list (list c1 m3 m2)(list m2 m1 a1))
)
; Midpoint of two points
(defun mid (a b)
(mapcar '(lambda (a b) (/ (+ a b) 2.0)) a b)
)
; Given 3 points, return the angle smaller than pi between a b and c b
(defun defl (a b c / d)
(abs (- (setq d (- (angle b a)(angle b c))) (* (fix (/ d pi)) (+ pi pi))))
)
;Generate a bunch of points for testing on layer named points.
(defun c:gen ( / n rangex rangey rangez)
(setq n (getint "\nNumber of points: ")
rangex 5000 ; Extent in X for the points
rangey (/ rangex 1.6); Extent in Y * Golden Ratio
rangez 10 ; Extent in Z
)
(while (> n 0)
(entmake
(list
'(0 . "POINT")
'(8 . "Points")
(cons 10 (list (* rangex (ran)) (* rangey (ran)) (* rangez (ran))))
)
)
(setq n (1- n))
)
)
; Random number generator
(defun ran ()
(setq seed (if seed (rem (+ (* seed 15625.7) 0.21137152) 1) 0.3171943))
)
I have also found a bug, the contour iteration is not completing all contours
Just as a thought, if you where to adjust the lres variable (line resolution) at the same time the smoothing factor is
Unacceptable, 3 lashes with a wet hypotenuse.
Thank you ymg . Where i can find the new version ?
intvL '(" 0.02" " 0.1" " 0.2" " 1" " 2" " 5" " 10")
majcntL '(" 0.1" " 0.5" " 1" " 5" " 10" " 50" "100")
hfacL '(" Max" " -20%" " -40%" " -60%" " -80%" " Min" "None")intvL '(" 0.02" " 0.1" " 0.2"......................)
majcntL '(" 0.1" " 0.5" " 1" ..........................")
hfacL '(" Max" " -20%" " -40%" " -60%" " -80%" " Min" "None")"Here is an update to the TIN program.Great Program , Yang !
It is still very much a work in progress, so I would appreciate
bug reports and enhancement suggestion.
The Contouring part has been accelerated a little, and now work with
fractionnal contours.
You can also call CONT to do contours on a Selection of existing 3DFACE.
Express Tools is required, as I make use of the Progress Bar.
Still missing Breaklines handling for the triangulation,
and Highlighting of Depression Contours for the contouring section.
ymg
"You can also call CONT to do contours on a Selection of existing 3DFACE." , I hope It can
chlh_jd,I'm sorry for Just finished reading your code in *gress5.lsp and foud add cont command .Quote"You can also call CONT to do contours on a Selection of existing 3DFACE." , I hope It can
Thanks for the feedback.
Agree that my wording is not too clear. What was meant by that is in the case where you already have
a selection of 3dfaces created by another program, you can use command "CONT" and do the contouring.
Normal way would be a set of points that you triangulate and then contours.Sorry for my poor English , I mean that , the zone I mark , smoothed contours less than ideal .
I do not understand exactly what you mean from the attached image. Could you clarify?
I've been playing around a bit with it and manage to accelerate it some. Shall post it soon.Not long ago, I also wrote a contour program is based on existing 3DFACES, and not for smoothing; and I'v lost in smooth annoying.
We still need to handle the marking of "Depression Contour" currently not available.
I highly recommend Christensen's paper, it is well worth the reading.
ymg
I believe that this zone would get crossing contours unless you do it this way.
In other word this is inherent to the method.
ymg
chlh_jd,Thanks Ymg , it did fix the bug in the case "intv=0.25 majcnt=1.0"; but error in the case "intv=0.1 majcnt=0.5".
Here is a revision with the bugs in contour generation removed.
Bear in mind, that this is still a work in progress, so the routine does not clean themselves nicely on completion.
ymg
Thanks Ymg , it did fix the bug in the case "intv=0.25 majcnt=1.0"; but error in the case "intv=0.1 majcnt=0.5".
...
(setq ti (car (_VL-TIMES)) ; Re-initialize timer for Contouring
zl (mapcar (function caddr) pl) ;_only for zlist
zmin (apply (function min) zl) ;_improve speed by GSLS(SS) 2014.3.2
zmax (apply (function max) zl)
zmin (* (fix (/ zmin intv)) intv) ;_zmin, first level for contour. (* (+ (fix (/ zmin intv)) 1) intv)
zmax (* (fix (/ zmax intv)) intv) ;_ zmax; last level for contour.
z zmin ; Initialize Current Contour level.
...
)
(setq tl (vl-sort tl
(function
(lambda (a b)
(< (max (nth (car a) zl)
(nth (cadr a) zl)
(nth (caddr a) zl))
(max (nth (car b) zl)
(nth (cadr b) zl)
(nth (caddr b) zl)))))))
...
In my head , Longer code still requires a larger computer memory
Command: prof
Select a Linear Entity:
Error: bad argument type: fixnump: nil
(setq pl (vl-sort pl (function (lambda (a b) (< (cadr a) (cadr b))))) ;; Sort on Y Coordinate
pl (vl-sort pl (function (lambda (a b) (< (car a) (car b))))) ;; Sort on X Coordinate
pl (vl-remove nil ;; Remove Duplicate Points
(mapcar '(lambda (a b) (if (not (equal (cdr (reverse a)) (cdr (reverse b)) 1e-03)) a))
pl (append (cdr pl) (list (car pl)))
)
)
)
(contour pl tl intv majcnt majcolor mincolor hfac norm)2. In the "get_ss" function , will be error , if the 3dface is manual change ( because of the osnap accuracy ), it can be solve by two way : (setq pl (cons (mapcar (function (lambda (a) (/ (fix (* a 1e6)) 1e6))) p) pl)3. In the contour function , some position may get error :?getContour@AcGeCurveBoundary@@QBEXAAHPAPAPAVAcGeEntity3d@@PAPAPAVAcGeCurve2d@@PAPAH3@Z
?getContour@AcGeImpCurveBoundary@@QBEXAAHPAPAPAVAcGeImpEntity3d@@PAPAPAVAcGeImpCurve3d@@PAPAH3@Z
?getContours@AcGeExternalBoundedSurface@@QBEXAAHAAPAVAcGeCurveBoundary@@@Z
?getContours@AcGeImpExternalBoundedSurface@@QBEXAAHAAPAVAcGeImpCurveBoundary@@@Z
ymg,
is a fantastic tool.
However, PROF function is not working. Always returns the same error.QuoteCommand: prof
Select a Linear Entity:
Error: bad argument type: fixnump: nil
I am using a lwpolyline as a linear entity. Not sure what am I doing wrong.
Attached the sample file
(if C:cal (princ "\nGeomcal loaded..")(arxload "geomcal"));** <-- this is optional
(if C:cal ;if loaded
(defun getz (p t1 t2 t3 / ptt ) ;defun opted by hanhphuc
(setq
ptt (list (car p)(cadr p)(1+(caddr p))))
(cal "ilp(p,ptt,t1,t2,t3)")
); getz
(defun getz .... code unchanged ) ;defun by ymg
) ;c:cal