Author Topic: Triangulation (re-visited)  (Read 312784 times)

0 Members and 1 Guest are viewing this topic.

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #105 on: July 04, 2013, 10:59:36 AM »
Hello Evgeniy,

I've look into your algorythm.  Sure did appreciate the subject of the study.
Gives a whole new meaning to Model Space  :-)

Your process as far as I can tell is patch based (As opposed to line smoothing)
and applicable to matrix. 

If you look into Christensen's paper,

   http://www.asprs.org/a/publications/pers/2001journal/april/2001_apr_511-517.pdf

you will get some of the pluses and minuses of this method.

He then goes on proposing a so called eclectic method which actually is kind of an hybrid (line and patch)
solutions.

This is what I propose to implement.

Hope everything is well with you.


ymg


ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #106 on: July 04, 2013, 12:54:10 PM »
I inadvertently introduce a bug in the last posting of tin.lsp,    :ugly:
as I checked for isclosed polyline the list should be reduced by one element
instead of two.

Apologies !!

ymg

Code: [Select]
(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)
   )

MP

  • Seagull
  • Posts: 17750
  • Have thousands of dwgs to process? Contact me.
Re: Triangulation (re-visited)
« Reply #107 on: July 04, 2013, 12:57:46 PM »
I inadvertently introduce a bug in the last posting of tin.lsp,    :ugly:

Unacceptable, 3 lashes with a wet hypotenuse.
Engineering Technologist • CAD Automation Practitioner
Automation ▸ Design ▸ Drafting ▸ Document Control ▸ Client
cadanalyst@gmail.comhttp://cadanalyst.slack.comhttp://linkedin.com/in/cadanalyst

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #108 on: July 20, 2013, 07:07:31 PM »
Here is Contour with Christensen's Eclectic Method for smoothing implemented.

I still have not implemented the variable smoothing that he describe in his paper.

Would appreciate comments, bug(s) report and suggestion for speeding up.

Code: [Select]
; 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))
)




Still to do:
  • Code Cean-up
  • Variable smoothing
  • Adding Constraints to triangulation
  • Putting Major and Minor contour on separate layers
  • Coloring with a gradient (ramp) the contours
  • Labelling of contours
  • Calculating Volume from the tin
  • etc... etc...


To test, downlad tins.lsp attached to this post.
Appload it to Autocad.
Issue command gen and creates a bunch of point.
Issue command tin and select all points created on previous step.

Triangulation and Smoothed Contours will be created,
on layer "TIN" and layer "Contour".

ymg

snownut2

  • Swamp Rat
  • Posts: 971
  • Bricscad 22 Ultimate
Re: Triangulation (re-visited)
« Reply #109 on: July 21, 2013, 03:42:20 PM »

Added Dialog Interface for Selecting Major & Minor Colors, Major & Minor Contour Intervals & Smoothing Factor
(all contained within LISP so there is only one file)

This also works in Bricscad 13......

See attached file for code.

I have also found a bug, the contour iteration is not completing all contours.  I have attached a dwg file that will not complete all contours.


« Last Edit: July 21, 2013, 08:54:54 PM by snownut2 »

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #110 on: July 22, 2013, 02:12:12 PM »
Snowmut,

Thanks for the dcl for contour.

Quote
I have also found a bug, the contour iteration is not completing all contours

Bug is confirmed, need to do a test or disturb position of point when a point is exactly on Contour level.
In your example be bug at point # 3 which is exactly at elevation 100.

I'll add the necessary and repost.

Thanks for the testing.

ymg

snownut2

  • Swamp Rat
  • Posts: 971
  • Bricscad 22 Ultimate
Re: Triangulation (re-visited)
« Reply #111 on: July 22, 2013, 05:34:14 PM »
Found another bug of sorts, the contours where not starting out at an elevation that was a multiple of the interval.  I have fixed that along with the contour layer placement (related) that added yesterday. 

See attached file (seems the entire code is to long to place in this box)

Revised code, added some error trapping in DCL and activated smoothing in DCL, also added some layer control to make sure appropriate layers are on.

Revised to;
       Accept either Blocks or Points
       Added TIN Name option to Dialog

Updated Dialog Interface to;
      Selection for Blocks or Point Objects to Contour (radio buttons)
      Added Error checking to ensure all required items have been done
     
I have also incorporated all the latest updates from ymg's latest file.
     

« Last Edit: July 24, 2013, 12:49:06 PM by snownut2 »

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #112 on: July 23, 2013, 03:57:42 PM »
Find attached tinsd.lsp which fixes the bugs found so far.

Also implemented variable smoothing as per Christensen's

A dialogue box compliment of Snownut2 is included for setting
various parameters for contouring.  Some of this seems to be
from LeeMac.

I need to get going on adding constraints to the triangulation.

Again would appreciate suggestion for speeding up, bug report
as well as better way to write some section.

ymg

P.S. Take the latest file below
« Last Edit: July 23, 2013, 10:39:37 PM by ymg »

snownut2

  • Swamp Rat
  • Posts: 971
  • Bricscad 22 Ultimate
Re: Triangulation (re-visited)
« Reply #113 on: July 23, 2013, 04:06:10 PM »
Indeed the LM:Popup function is included in the code, it has not been modified from Lee's code.

Bruce

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #114 on: July 23, 2013, 11:04:13 PM »
Here is two images of the output from tinsd.

Only a 100 points and 20 contours.

On these images I have superimposed output with max smoothing,
min smoothing and no smoothing.

Not sure the variable smoothing is worth the trouble.

ymg




snownut2

  • Swamp Rat
  • Posts: 971
  • Bricscad 22 Ultimate
Re: Triangulation (re-visited)
« Reply #115 on: July 24, 2013, 12:27:37 AM »
I am not surprised at what you are calling minor differences, if the lines where shifted/relocated much more then the contour accuracy would suffer.  Typically as a rule in the civil world you cannot be more then 1/2 the contour interval off at any location in the contoured area.  That is on the outside you would like to see the contours closer than that, if the lines where adjusted anymore with the smoothing then you could be outside of that accuracy range.

Just as a thought, if you where to adjust the lres variable (line resolution) at the same time the smoothing factor is adjusted, making the resolution smaller as the smoothing factor is increased, the results are definitely more noticeable.

Regarding triangulation constraints, it appears that the very large angle >160deg triangles lay at the outer edges of the TIN.  These triangles cause the contour lines to take a sharp turn and head off into an un-natural direction.  How about not allowing triangles with any angle greater than 160deg. (or some other limit that makes sense)
« Last Edit: July 24, 2013, 01:08:08 AM by snownut2 »

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #116 on: July 24, 2013, 01:05:59 PM »
Quote
Just as a thought, if you where to adjust the lres variable (line resolution) at the same time the smoothing factor is

That suggestion makes perfect sense.

The second one about triangle could prove to be more involved.
However if we implement constraints with an outer boundary
that problem would go away.

ymg

ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #117 on: July 24, 2013, 02:40:13 PM »

Here is a case that tinsd does nor handle properly.

Has to do when a contour has only one edge.

Quote
Unacceptable, 3 lashes with a wet hypotenuse.

Will test and correct.  By now the wet hypotenuse
of MP is used up. So I will use the opposite side.   :pissed:

ymg


ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #118 on: July 24, 2013, 04:22:53 PM »
Bug is corrected.

Find below revised file.


Code - Auto/Visual Lisp: [Select]
  1. (foreach p xl
  2.         (setq ent nil)
  3.         (if (equal (car p) (reverse (last p)))
  4.                  (setq isclosed 1
  5.                               p (append p (list(cadr p)))
  6.                  )
  7.                  (setq isclosed 0
  8.                             ent (list (cons 10 (clv l (car p) pl)))
  9.                  )
  10.         )
  11.        
  12.        (if (< (length p) 3); Added to handle the case where polyline is of lenght 2.              
  13.            (setq  v3 (clv l (cadr p) pl))
  14.         )
  15.         (while (> (length p) 2)
  16.                
  17.            (setq  v1 (clv l (car p) pl)
  18.                   v2 (clv l (cadr p) pl)
  19.                   v3 (clv l (caddr p) pl)
  20.                  prv (car p)
  21.  
  22. etc. etc.
  23.  


ymg

  • Guest
Re: Triangulation (re-visited)
« Reply #119 on: August 30, 2013, 06:37:31 PM »
I have updated the triangulation so that an adjacency list or neighbour list of triangle is generated.

Also modified so that we generate Voronoi Diagram plus the Delaunay triangulation.

The neigbour list now permits us to locate ourself efficiently in the tin via a new routine called Triloc.
This is actually Lawson's walk who invented the procedure.  Since then many different scheme have been devised
to accelerate it when you have huge triangulation.

Here's a paper showing some of them:  Walking Location Algorithms  by Roman Soukal.

Another paper Walking in a triangulation by Olivier Devillers analyzes the speed of the  various schemes.

Here are the added function to the triangulation:

Code - Auto/Visual Lisp: [Select]
  1. ;;****************************************************************************;
  2. ;; (getneighbour pl tl)                                                                                                                   ;
  3. ;;                                                                                                                                                 ;
  4. ;; Parameters: pl, Points List                                                                                                        ;
  5. ;;             tl, Triangle List (index into pl)                                                                                       ;
  6. ;;    Returns: nl, List of neigbours for each edges of each triangles.                                               ;
  7. ;;                                                                                                                                                 ;
  8. ;; Should be done by the triangulation.                                                                                         ;
  9. ;;****************************************************************************;
  10.  
  11. (defun getneighbour (pl tl / nl pos tmp)
  12.  
  13.              (setq el nil)
  14.              (foreach tr tl
  15.                 (setq el (cons (list (car tr)(cadr tr)) el)
  16.                       el (cons (list (cadr tr)(caddr tr)) el)
  17.                       el (cons (list (caddr tr)(car tr)) el)
  18.                 )
  19.              )
  20.              (setq el (reverse el))
  21.  
  22.             (setq tmp nil)
  23.             (foreach e el
  24.                 (setq pos (vl-position (reverse e) el))
  25.                 (if pos
  26.                    (setq tmp (cons (/ pos 3) tmp)); Integer Division          ;
  27.                    (setq tmp (cons nil tmp))
  28.                 )
  29.             )
  30.            
  31.             (setq nl nil)
  32.             (while tmp
  33.                 (setq nl (cons (list (caddr tmp)(cadr tmp)(car tmp)) nl)
  34.                      tmp (cdddr tmp)
  35.                 )
  36.             )
  37.             nl
  38. )
  39.  
  40.  
  41. ;;****************************************************************************;
  42. ;; (triloc p)                                                                                                                                   ;
  43. ;;                                                                                                                                                ;
  44. ;; Locates triangle which encloses point p using Lawson's Walk.                                                    ;
  45. ;;                                                                                                                                                ;
  46. ;; Given p a point, Returns Index in tl of triangle containing the point.                                          ;
  47. ;; If outside the triangulation Return is nil.                                                                                   ;
  48. ;;                                                                                                                                                ;
  49. ;; Point List pl and Neigbour List nl are defined outside this routine.                                              ;
  50. ;; by ymg  August 2013                                                                                                               ;
  51. ;;****************************************************************************;
  52.  
  53. (defun triloc (p / v v1 v2 v3)
  54.  
  55.     ; Negative return, point is on right of v1->v2                                                                            ;
  56.     ; Positive return, point is on left of v1->v2                                                                                ;
  57.     ; 0 return, point is smack on the vector.                                                                                    ;
  58.     (defun onside (p v1 v2 /)
  59.         (setq v1 (list (car v1)(cadr v1))); We want 2d points                                                             ;
  60.         (apply '- (mapcar '* (mapcar '- v2 v1) (reverse (mapcar '-  p v1))))
  61.     )
  62.  
  63.     (setq found nil)
  64.     (if (not tn) (setq tn (/ (length tl) 2)))
  65.     (while (and (not found) tl)
  66.         (setq  v (nth tn tl)
  67.               v1 (nth (car v) pl)
  68.               v2 (nth (cadr v) pl)
  69.               v3 (nth (caddr v) pl)
  70.         )      
  71.         (cond
  72.            ((minusp (onside p v1 v2))(setq tn (car (nth tn nl))))
  73.            ((minusp (onside p v2 v3))(setq tn (cadr (nth tn nl))))
  74.            ((minusp (onside p v3 v1))(setq tn (caddr (nth tn nl))))
  75.            (t (setq found t))      ; We are inside triangle tn                ;
  76.         )
  77.         (if (not tn)(setq found t)); We are outside the triangulation         ;
  78.     )
  79.     tn ;
  80. )
  81.  
  82.  

The body of the triangulation has been modified to make use of the getneighbour function.
Note that contour has been removed from the attachment as I am currently doing a major
rewrite to make use of the neighbour list and hopefully gain some speed.

I have also put together a demoz routine to demonstrate what can be done with triloc and the neighbour list.
It is heavily inspired (leeched  :?)  by LeeMac's grtext.lsp

Here is the code to that one followed by an animated gif showing it in action.

Code - Auto/Visual Lisp: [Select]
  1. ;;****************************************************************************;
  2. ;; Demoz                                                                                                                                     ;
  3. ;; Demonstrate Lawson's Walk, to locate a triangle in a triangulation                                             ;
  4. ;; As you hover around the screen, the triangle under the screen is                                               ;
  5. ;; highlighted and the coordinates of the cursor position including Z, plus                                      ;
  6. ;; the triangle number are written to the screen.                                                                           ;
  7. ;;                                                                                                                                                 ;
  8. ;; Inspired by LeeMac's Grtext.lsp                                                                                                 ;
  9. ;;                                                                                                                                                 ;
  10. ;; by ymg        August 2013                                                                                                          ;
  11. ;;****************************************************************************;
  12.  
  13. (defun c:demoZ ( / *error* pnt str col ent etxt p prevtxt prevhatch
  14.                    prevtn scl t1 t2 t3 tr txt txtp xof yof)
  15.  
  16.     (defun *error* ( m ) (princ m) (redraw) (princ))
  17.  
  18.     (toglay '("Points" "Voronoi"))
  19.      
  20.     (while (= 5 (car (setq pnt (grread nil 13 0))))
  21.         (redraw)
  22.         (setq p (cadr pnt))
  23.         (setq tn (triloc p))
  24.         (if tn
  25.           (progn
  26.             (setq tr (nth tn tl)
  27.                   t1 (nth (car tr) pl)
  28.                   t2 (nth (cadr tr) pl)
  29.                   t3 (nth (caddr tr) pl)
  30.                    p (getz p t1 t2 t3)
  31.             )    
  32.             (if (not (= tn prevtn))
  33.                (progn
  34.                   ; Built with MakeEntmake.lsp by Cab at TheSwamp             ;
  35.                   (setq ent (list (cons 0 "HATCH")
  36.                                   (cons 100 "AcDbEntity")  
  37.                                   (cons 8 "Tin")
  38.                                   (cons 62 1)
  39.                                   (cons 440 33554534)
  40.                                   (cons 100 "AcDbHatch")
  41.                                   (cons 10 (list 0.0 0.0 0.0))
  42.                                   (cons 210 (list 0.0 0.0 1.0))      
  43.                                   (cons 2 "SOLID")
  44.                                   (cons 70 1)
  45.                                   (cons 71 0)
  46.                                   (cons 91 1)
  47.                                   (cons 92 1)
  48.                                   (cons 93 3)
  49.                                   (cons 72 1)
  50.                                   (cons 10 (list (car t1)(cadr t1) 0.0))
  51.                                   (cons 11 (list (car t2)(cadr t2) 0.0))
  52.                                   (cons 72 1)
  53.                                   (cons 10 (list (car t2)(cadr t2) 0.0))
  54.                                   (cons 11 (list (car t3)(cadr t3) 0.0))
  55.                                   (cons 72 1)
  56.                                   (cons 10 (list (car t3)(cadr t3) 0.0))
  57.                                   (cons 11 (list (car t1)(cadr t1) 0.0))
  58.                                   (cons 97 0)  
  59.                                   (cons 75 1)
  60.                                   (cons 76 1)
  61.                                   (cons 98 1)
  62.                                   (cons 10 (list 0.0 0.0 0.0))
  63.                                   (cons 450 0)
  64.                                   (cons 451 1)
  65.                                   (cons 460 0.0)
  66.                                   (cons 461 0.0)
  67.                                   (cons 452 0)
  68.                                   (cons 462 1.0)
  69.                                   (cons 453 2)
  70.                                   (cons 463 0.0)
  71.                                   (cons 63 5)
  72.                                   (cons 421 255)
  73.                                   (cons 463 1.0)
  74.                                   (cons 63 2)
  75.                                   (cons 421 16776960)
  76.                                   (cons 470 "LINEAR")    
  77.                             )
  78.                   )
  79.                   (if prevhatch (entdel prevhatch))
  80.                   (entmake ent)
  81.                   (setq prevhatch (entlast)
  82.                         prevtn tn
  83.                   )
  84.                   (redraw)
  85.                )
  86.             )
  87.           )
  88.           (if prevhatch (progn (entdel prevhatch)(setq prevhatch nil)))
  89.         )
  90.         ; Leeched from LeeMac (LM:GrText)                                     ;
  91.         ; http://www.lee-mac.com/lisp/html/GrTextDemo.html                    ;
  92.         (setq  col 2
  93.                xof 30
  94.                yof -30
  95.                scl (/ (getvar 'viewsize) (cadr (getvar 'screensize)))
  96.                  p (trans p 1 2)
  97.                str (mapcar 'rtos (trans p 1 0))
  98.               txtp (list (+ (car  p) (* xof scl)) (+ (cadr p) (* yof scl)))
  99.                txt (strcat  "T#: " (if tn (itoa tn) "Outside") "\\PX= " (car str) "\\PY= " (cadr str) "\\PZ= " (caddr str))
  100.         )
  101.         (setq etxt (list (cons 0 "MTEXT")
  102.                          (cons 100 "AcDbEntity")
  103.                          (cons 67 0)
  104.                          (cons 8 "Tin")
  105.                          (cons 62 col)
  106.                          (cons 100 "AcDbMText")
  107.                          (cons 10 txtp)
  108.                          (cons 40 30)
  109.                          (cons 41 360.0)
  110.                          (cons 46 0.0)
  111.                          (cons 72 5)
  112.                          (cons 1 txt)
  113.                          (cons 7 "Standard")
  114.                          (cons 11 '(1.0 0.0 0.0))
  115.                          (cons 42 237.0122)
  116.                          (cons 43 181.39154)
  117.                          (cons 50 0)
  118.                          (cons 73 1)
  119.                          (cons 44 1.0)
  120.                   )
  121.         )
  122.         (if prevtxt (entdel prevtxt))
  123.         (entmake etxt)
  124.         (setq prevtxt (entlast))
  125.     )
  126.     (if prevtxt (entdel prevtxt))
  127.     (if prevhatch (progn (entdel prevhatch)(setq prevhatch nil)))
  128.     (toglay '("Points" "Voronoi"))
  129.     (redraw)
  130.    
  131.    
  132.   (princ)
  133. )
  134.  
  135.  
  136.  
  137.  
  138. ;;****************************************************************************;
  139. ;; (getz p t1 t2 t3)                                                                                                                        ;
  140. ;; Given point p and triangle defined by points t1, t2, t3                                                                ;
  141. ;; Returns: (x y z) where z is on face of triangle.                                                                           ;
  142. ;;                                                                                                                                                  ;
  143. ;; By ymg  August 2013                                                                                                                ;
  144. ;;****************************************************************************;
  145.  
  146. (defun getz (p t1 t2 t3 /  v1 v2)
  147.      
  148.    ; Calculating a normal vector with gile's functions in line.               ;
  149.    ; Note that I do not calculate the unit vector for the normal n.           ;
  150.    (setq v1 (mapcar '- t2 t1)
  151.          v2 (mapcar '- t3 t1)
  152.           n (list (- (* (cadr v1) (caddr v2)) (* (caddr v1) (cadr v2)))
  153.                   (- (* (caddr v1) (car v2)) (* (car v1) (caddr v2)))  
  154.                   (- (* (car v1) (cadr v2)) (* (cadr v1) (car v2)))    
  155.             )
  156.    )
  157.    (list (car p)(cadr p)(/ (apply '+ (mapcar '* n (mapcar '- t1 p)))(caddr n)))
  158.  )
  159.  
  160. ;; (toglay lst)                                                                                               ;
  161. ;; Toggles (on/off) a list of Layer names.                                                      ;
  162. ;; From a Fishing Story by Cab at TheSwamp >-=((((°>                              ;
  163. ;;                                                                                                                ;
  164. ;; by ymg               August 2013                                                                 ;
  165.  
  166. (defun toglay (lst / entlist l)
  167.     (foreach l lst
  168.          (and (setq entlist (tblobjname "LAYER" l))
  169.               (setq entlist  (entget entlist))
  170.               (entmod (subst (cons 62 (- (cdr (assoc 62 entlist))))
  171.                              (assoc 62 entlist)
  172.                              entlist
  173.                       )
  174.               )
  175.          )
  176.     )
  177. )
  178.  
  179.  
  180.  
  181.  

To run the demo, load the attachment demoz.lsp, Issue command GEN and enter number of points desired
in the triangulation.

Isue command VOR, the Voronoi Diagram and the Delaunay Triangulation will generate.

Issue command DEMOZ, as you hover in the drawing, the triangle number plus coordinates X,Y and Z of point under cursor
is displayed and the currently occupied triangle is hatched in solid Red.

Sorry!, about the long post.

ymg



 
« Last Edit: September 01, 2013, 10:15:19 AM by ymg »