Topic: Triangulation (re-visited)

[SOFITO_SOFT]

Hello:
In my pc, I can see. Sorry.
Moving the cursor over the LWPOLYLINES, clicking a point, it is estimated the 2 polylines closest point is interpolated and the original. It displays the Z. .. Now I can do with 2 sets of polylines in the layer "original" and "futuro".  The lsp displays the Z and differences and time calculation
I attached a static image of the video. and a dwg with the new method of double layer/mdt.
Greetings.

[ElpanovEvgeniy]

I 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

A very strange sentence!
Individually written program, will always be significantly more expensive than to sell in large quantities. Free software has its price.

I sometimes are participating in such projects if the projects are necessary to me as well. 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. The project is done under a separate technology, which is uncommon and can not be maintained well-known programs. The problem of triangulation, very far from me...

[SOFITO_SOFT]

  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

[ElpanovEvgeniy]

 
You've worked with solids "ACIS"? They are solid autocad. Are well documented.
very very interesting.
Greetings

Yes, I use a simple command entmakex, entmode, entupd...

[SOFITO_SOFT]

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. 

[ElpanovEvgeniy]

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.  

All this is not the decoding and modeling - Change, for example, add or remove a port complex shape, or simply add / remove a point on one edge. This is the easiest level, and then - the cross section, sections, species in its representation of geometry. The problem of automatically creating species and the cross section with a complete adaptability - changed item, change views, changed views, changed body.

ps. arx not give full access to the geometry. Using arx, you can not add or delete a single point in the body ...

[ymg]

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

[ElpanovEvgeniy]

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

I think you will release the project documentation.
I manage the release of automation software.

If you need to develop small programs or large system, you can contact me with offers. Until now, I work only in the domestic market, but is willing to work with other countries...

[ymg]

Evgenyi,

To go back to your triangulation, I did follow part of the code.

When you setup your first triangle in list s it is too small.

The vertex of that triangle are inside the circumcircle of other points on the outside.
In the litterature some authors advocate infinity for the so called supertriangle.

I have attached a picture showing the interference of your s triangle with some
of the points in the cloud.

If you fix that, I believe we have a Delaunay triangulator. 

ymg

[ymg]

Here I have revised your code:

Code: [Select]

(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
I don't think anything I did will affect the speed.

Attached a gif of the result.

ymg

[SOFITO_SOFT]

Hello:
ymg :
Another way of navigating in 3d over a MDT ... for a XY point which are the points that cut the 3dsolids generates by 3DFACEs.
In this case, draw lines between points found in the vertical of 3 3dsolid.

Time for the example: 0.374973 seg. for 6 points founds....To improve...


Greetings.

[SOFITO_SOFT]

Oops...
after some small problems with sub-functions releases  ....the LSP for navigating in 3d over a MDT ( several 3dsolids )
( c:ll ) ...for initialize and...
( c:pto_en_Z_solid ) for run...
Greetings 

[ymg]

Here is Evgeniy's triangulation commented for analysis.

The problem with first triangle is corrected, but note that in certain case it could be too small.

I have renamed most variables and commented.

Code: [Select]

(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)
      )
   )
)


ymg

[ElpanovEvgeniy]

Hi ymg!

Now, the program still seems complicated?

[ymg]

Now, the program still seems complicated?

No, since it is basically the same algoritmn as Piazza's program.

However I have to agree with Vovka's comment about your variable naming style 

The Divide and conquer algorithm could be faster than this one if written properly.

ymg