Author Topic: intersection  (Read 11623 times)

0 Members and 1 Guest are viewing this topic.

Fuccaro

  • Guest
intersection
« on: April 21, 2004, 08:41:34 AM »
Hello Lispers!

I use only AutoLisp for the moment. I try to write a program to calculate the illumination on surfaces.
I search for the intersection point between a 3dface and a ray (line) defined by two points. This involves a lot of mathematics: first I have to calculate the intersection between the line and the plane passing trough the 3 corners of the 3Dface, than I must check if that point is inside of the triangle (face).
It can be done, but probable you, the VL wizards, just take out from the hat a magic formula, something like vla_get_intersection_between_that_and_this...

hendie

  • Guest
intersection
« Reply #1 on: April 21, 2004, 09:17:09 AM »
hi Fuccaro
do you mean something like
Code: [Select]
(setq P1 (vla-intersectwith Object1 Object2 acextendnone))

 ~ Object1 & Object2 need to be valid VLA-objects, so you use something like
Code: [Select]
(setq Object1 (vlax-ename->vla-object (entlast)))

but I'm not sure how that would apply in your case since you are talking about the intersection of a line and a planar region. According to the help file it works with all objects except Except Pviewport and PolygonMesh.
If need be, I suppose you could always create a temporary "region" to test on.
If the two objects do not intersect, no data is returned.

I'm sure some of our more esteemed colleagues will have an answer for you

hendie

  • Guest
intersection
« Reply #2 on: April 21, 2004, 10:03:04 AM »
okay, I could be wrong.
I can't get it to work. I have my VLA object line and I have my VLA object region and I know for a fact that they do intersect.
I get a variant array returned but the value is always Nil.

bloody VL-alchemy  :genie:

anyone tell me where I'm going wrong ?


*edited* ~ tried it with a 3DFace also but still get the same error

hendie

  • Guest
intersection
« Reply #3 on: April 21, 2004, 10:35:39 AM »
Code: [Select]

_$ (setq MyLine (car (entsel "\nPlease select LINE : ")))
(setq Vline (vlax-ename->vla-object MyLine))

<Entity name: 4002cf78>
#<VLA-OBJECT IAcadLine 028db734>

_$ (setq MyRegion (car (entsel "\nPlease select REGION : ")))
(setq Vreg (vlax-ename->vla-object MyRegion))

<Entity name: 4002cf60>
#<VLA-OBJECT IAcadRegion 028d9664>

_$ (setq P1 (vla-intersectwith Vline Vreg acextendnone))

#<variant 8197 ...>



and when I do an inspect on p1 I get:
Code: [Select]

<type> Array of double
<value> # <safearray>


and the value is NIL


I even tried extruding the region but got an ; error: Automation Error. Not implemented yet message



and it all seemed so simple  :crazy:

hendie

  • Guest
intersection
« Reply #4 on: April 21, 2004, 10:51:16 AM »
I did a google search and found a few other Vlisps for the intersectwith function. All of them returned a variant array of Nil.
perhaps the object types are incompatible for this method.

Mark

  • Custom Title
  • Seagull
  • Posts: 28753
intersection
« Reply #5 on: April 21, 2004, 10:58:33 AM »
What do you get when you do this:
Code: [Select]
(vlax-safearray->list  (vlax-variant-value P1))
TheSwamp.org  (serving the CAD community since 2003)

Mark

  • Custom Title
  • Seagull
  • Posts: 28753
intersection
« Reply #6 on: April 21, 2004, 11:02:33 AM »
how about a dwg we can work with. upload to
http://theswamp.org/lilly.pond/
TheSwamp.org  (serving the CAD community since 2003)

hendie

  • Guest
intersection
« Reply #7 on: April 21, 2004, 11:04:49 AM »
Quote from: Mark Thomas
What do you get when you do this:
Code: [Select]
(vlax-safearray->list  (vlax-variant-value P1))


Code: [Select]

$ (vlax-safearray->list  (vlax-variant-value P1))
; error: ActiveX Server returned an error: Invalid index



 :horror:

DEVITG

  • Bull Frog
  • Posts: 479
finding intersections even with out lisp
« Reply #8 on: April 21, 2004, 10:08:52 PM »
On ACAD there is a marvelous tool

geomcal .
Quote


Command line: cal (or 'cal for transparent use)

CAL is an online geometry calculator that evaluates point (vector), real, or integer expressions. The expressions can access existing geometry using the object snap functions such as CEN, END, and INS. You can insert AutoLISP® variables into the arithmetic expression and assign the value of the expression back to an AutoLISP variable. You can use these arithmetic and vector expressions in any AutoCAD® command that expects points, vectors, or numbers.

CAL topics:

Understanding Syntax of Expressions
Formatting Feet and Inches
Formatting Angles
Using Points and Vectors
Using AutoLISP Variables
Using AutoCAD System Variables
Converting Units of Measurement
Using Standard Numeric Functions
Calculating a Vector from Two Points
Calculating the Length of a Vector
Obtaining a Point by Cursor
Obtaining the Last-Specified Point
Using AutoCAD Snap Modes in Arithmetic Expressions
Converting Points between UCS and WCS
Calculating a Point on a Line
Rotating a Point About an Axis
Obtaining an Intersection Point
Calculating a Distance
Obtaining a Radius
Obtaining an Angle
Calculating a Normal Vector
Using Shortcut Functions

It has a lot of tools and one just fited to yuor need
Quote
The ill and ilp functions determine intersection points.

ill(p1,p2,p3,p4)

Determines the intersection point between two lines (p1,p2) and (p3,p4). AutoCAD considers all points three-dimensional.

ilp(p1,p2,p3,p4,p5)

Determines the intersection point between a line (p1,p2) and a plane passing through three points (p3,p4,p5).

 


It can be used on Lisp

Code: [Select]


 (Setq ach (c:cal "pld(pa1,pa2,dst)")) ; locate a point on line [pa1 pa2] at dist

(SETQ ANPab (C:CAL "ang(INAB,pA2,pb2)")) ; calc a angle between two line's END point ant it's apex or vertice, Does not matter wher it is located






Quote
dist(p1,p2)

Determines the distance between two points, p1 and p2. This is the same as the vector expression abs(p1-p2).

dpl(p,p1,p2)

Determines the shortest distance between point p and the line passing through points p1 and p2.

 

dpp(p,p1,p2,p3)

Determines the distance from a point p to a plane defined by three points (p1,p2,p3).

dist(p1,p2)

Determines the distance between two points p1 and p2. This is the same as the vector expression abs(p1-p2).

 

The following example returns half the distance between the centers of two selected objects:

dist(cen,cen)/2

The following example finds the distance between the point 3,2,4 and a plane you define by selecting three endpoints:

dpp([3,2,4],end, end, end)





Really I do not know why it is not used more often.
Location @ Córdoba Argentina Using ACAD 2019  at Window 10

Fuccaro

  • Guest
intersection
« Reply #9 on: April 22, 2004, 02:12:40 AM »
GREAT!
Thank you DEVITG.

The ILP function fit exactly to my needs. It will return the intersection between the plane and the ray (let’s call X that point). Having the intersection point it is a 2D problem to check if it is inside of the triangle ABC. With the ILL I will have the intersection point (A1) between AX and BC. If the distance AX is greater than AA1, the point is outside. Repeating 3 times this test (for each point A, B and C) I can decide if the point is inside.

They are so smart people in Argentina... Thank you again!
And thanks for the other answers too!

hendie

  • Guest
intersection
« Reply #10 on: April 22, 2004, 03:24:36 AM »
brilliant ! glad you got it sorted Fuccaro, and thanks to Devitg for the solution.

it's still buggin me why I couldn't get there with the vla-intersectwith method though

DEVITG

  • Bull Frog
  • Posts: 479
the invisible it is allways at our reach
« Reply #11 on: April 22, 2004, 07:46:01 AM »
I discover the CAL funtion just browsing , I can not remember when I find
it but it is a simple tool , it can be used as a direct acad command from the command line and it can be used transparently.

just 'cal

and remember the famous acronym

KISS : KEEP IT SIMPLE SIR



 :D  :)
Location @ Córdoba Argentina Using ACAD 2019  at Window 10

hendie

  • Guest
intersection
« Reply #12 on: April 22, 2004, 08:04:04 AM »
I knew about the 'CAL tool, but I was not aware of the ILL & ILP options within it.

I must remember to store that in my brain somewhere in case I ever need it again

CAB

  • Global Moderator
  • Seagull
  • Posts: 10401
intersection
« Reply #13 on: April 22, 2004, 08:21:48 AM »
Fuccaro
You may find some use for these functions that were in my collection from the forums.

Code: [Select]
;|Here are a set of LISP functions that do the job for old-style polylines
(I did it a few years ago and it may not be perfect LISP). You will have
to modify the function PLINEVERTICES to be able to use new-style plines.
INPOLY does not calculate bulges, but gives correct results even for
self-intersecting plines. I used the fact that the sum of angles from
the checkpoint to all vertices is equal 0 (360) degrees if the
checkpoint is inside and any othre value if it is outside. The FUZZY
value is to avoid rounding errors.
Code is free to use by anybody

Tom
|;

(setq imr_fuzzy 0.0001)

(defun inpoly (pt en /)
  (if (equal 0.0 (anglesum pt (plinevertices en)) imr_fuzzy)
    nil
    t
  )
)


(defun c:inpoly (/ pt en)
  (setvar "cmdecho" 0)
  (if (= "POLYLINE" (cdr (assoc 0 (entget (setq en (car (entsel)))))))
    (while (setq pt (getpoint "\npoint to check: "))
      (if (inpoly pt en)
        (princ "INSIDE  polyline!")
        (princ "OUTSIDE polyline!")
      )
    )
  )
  (prin1)
)

(defun angd (scheitel p1 p2) (r2d (angr scheitel p1 p2)))

(defun r2d (winkel) (* (/ (float winkel) pi) 180.0))

(defun angr (scheitel p1 p2 / alf bet)
  (setq bet (angle scheitel p1)
        alf (angle scheitel p2)
        alf (- alf bet)
  )
  (if (< alf 0)
    (setq alf (+ (* 2 pi) alf))
  )
  alf
)

(defun plistisplanar (plist / isplanar)
  (setq isplanar t
        pl (cddr plist)
  )
  (while (and isplanar (< 3 (length plist)))
    (while (< 1 (length pl))
      (if (inters (car plist) (cadr plist) (car pl) (cadr pl))
        (setq isplanar nil
              pl nil
              plist nil
        )
      )
      (setq pl (cdr pl))
    )
    (if isplanar
      (setq plist (cdr plist)
            pl    (cddr plist)
      )
      (setq plist nil)
    )
  )
  isplanar
)

(defun anglesum (pt plist / as p1 p2 an1 an2 an)
  (setq as    0.0
        stp   (car plist)
        count 1
  )
  (while (< 1 (length plist))
    (setq p1    (car plist)
          p2    (cadr plist)
          plist (cdr plist)
          an    (angd pt p1 p2)
          an    (if (< 180.0 an)
                  (- an 360.0)
                  an
                )
          as    (+ as an)
    )
    (setq count (1+ count))
  )
  (setq
    an (angd pt p2 stp)
    an (if (< 180.0 an)
         (- an 360.0)
         an
       )
  )
  (+ as an)
)

(defun plinevertices (en / el vl)
  (while (= "VERTEX"
            (cdr (assoc 0
                        (setq el (entget (setq en (entnext
                                                    en
                                                  )
                                         )
                                 )
                        )
                 )
            )
         )
    (setq vl (cons (cdr (assoc 10 el)) vl))
  )
  (if (= (car vl) (car (reverse vl)))
    (reverse (cdr (reverse vl)))
    vl
  )
)




Code: [Select]
;|Re: Point inside the perimeter of a closed polyline ?
Subject: Re: Point inside the perimeter of a closed polyline ?
From: "Michael Doerr" <mdoerr@cebacus.de>
Date: Tue, 14 Dec 1999 06:47:25 +0000
Newsgroups: autodesk.autocad.customization
Hi Dominique,

the following code calculates the sum of the angles from a given point to
the points of the polyline and decides wether the point is inside or outside
in that way. 'PointInQuestion' means the point which has to be tested and
'Point_List' is the list of Points which represent the polyline. If the
function returns 'T' the given point is inside your polyline.

Greetings, Michael|;

(defun punktinpolylinie (pointinquestion
                         point_list
                         /
                        )
  (if (equal 0.0 (pipwinkelsumme pointinquestion point_list) 0.0001)
    nil
    t
  )
)

(defun pipwinkelsumme (pointinquestion           point_list
                       /            count        p1
                       p2           scheitel     winkeleins
                       winkelzwei
                      )
  (setq winkeleins 0.0
        scheitel   (car point_list)
        count      1
  )
  (while (< 1 (length point_list))
    (setq p1         (car point_list)
          p2         (cadr point_list)
          point_list (cdr point_list)
          winkelzwei (pipwinkelhilfe pointinquestion p1 p2)
          winkelzwei (if (< 180.0 winkelzwei)
                       (- winkelzwei 360.0)
                       winkelzwei
                     )
          winkeleins (+ winkeleins winkelzwei)
    )
    (setq count (1+ count))
  )
  (setq winkelzwei (pipwinkelhilfe pointinquestion p2 scheitel)
        winkelzwei (if (< 180.0 winkelzwei)
                     (- winkelzwei 360.0)
                     winkelzwei
                   )
  )
  (+ winkeleins winkelzwei)
)

(defun pipwinkelhilfe (pointinquestion p1 p2 / alpha beta)
  (setq beta  (angle pointinquestion p1)
        alpha (angle pointinquestion p2)
        alpha (- alpha beta)
  )
  (if (< alpha 0)
    (setq alpha (+ (* 2 pi) alpha))
  )
  (* (/ (float alpha) pi) 180.0)
)




Code: [Select]
;; Updated (03-07-03)
;; Given:
;;-------------------------------------------------------
;; Function to find or create an invisible ray VLA-Object
;; on Layer "0" an to create a global symbol $cv_ray.
;; It will leave the Ray object in the drawing to save time.
(defun @cv_ray ( / e ss)
  (if (setq ss (ssget "X" '((0 . "RAY")(8 . "0")(60 . 1))))
    (setq $cv_ray (vlax-ename->vla-object (ssname ss 0)))
    (setq $cv_ray
      (entmakex
        (list
         '(0 . "RAY")'(100 . "AcDbEntity")'(100 . "AcDbRay")
         '(8 . "0")'(60 . 1)'(10 0.0 0.0 0.0)'(11 1.0 0.0 0.0)
        )
      )
      $cv_ray (vlax-ename->vla-object $cv_ray)
    )
  )
)
;; And:
;;-------------------------------------------------
;; Function originated by Ken Alexander (03-05-03)
;; that is 10X faster than my @cv_parse_list,
;; to group data into triplets.
;; Thanks, Ken!
;;
(defun @cv_triple_up (old / new)
  (while
    (setq new (cons (list (car old)(cadr old)(caddr old)) new)
          old (cdddr old)
    )
  )
  (reverse new)
)
;; And:
;;-------------------------------------------------------------------
;; Function to determine if a point <PIQ> is inside a closed polyline
;; based on the number of intersections found between a ray whose
;; basepoint is the PIQ, and that PIQ is not one of the intersection
;; points.
;; Arguments:
;;   PIQ = 3D Point in WCS
;;   Outer = Outer Polyline VLA-Object
;; Returns:
;;   either T (inside) or nil (on or outside)
(defun @cv_inside (PIQ Outer / Points)
  (vl-load-com)
  (and
    (> (vl-list-length PIQ) 1)
    (vl-every 'numberp PIQ)
    (= (type Outer) 'VLA-Object)
    (vlax-property-available-p Outer 'Closed)
    (= (vla-get-closed Outer) :vlax-true)
    (or $cv_ray (@cv_ray))
    (or (vlax-put $cv_ray "Basepoint" PIQ) T)
    (setq Points (vlax-invoke Outer "IntersectWith" $cv_ray acExtendNone))
    (setq Points (@cv_triple_up Points))
    (= (rem (length Points) 2) 1)
    (not (equal PIQ (vlax-curve-getclosestpointto Object PIQ) 1e-11))
  )
)

;; Then:
;; ... define your point and Outer object ...
(setq Inside (@cv_inside PIQ Outer))



Code: [Select]
Re: Point inside the perimeter of a closed polyline ?
Subject: Re: Point inside the perimeter of a closed polyline ?
From: kboutora@francemel.com (Kamal Boutora)
Date: Thu, 16 Dec 1999 15:40:45 +0000
Newsgroups: autodesk.autocad.customization
In article <838ik5$kje17@adesknews2.autodesk.com>,
not.robertb@mwengineers.com says...
> Could you post code, or a link to some?
>

  From the faq of comp.graphics.algorithms (a very interesting faq that
  every programmer should have IMHO).

    you can look at the latest HTML version at either of these two sites:

http://www.exaflop.org/docs/cgafaq
http://www.cis.ohio-state.edu/hypertext/faq/usenet/graphics/algorithms-
faq/faq.html

    The code below is from Wm. Randolph Franklin <wrf@ecse.rpi.edu>

    int pnpoly(int npol, float *xp, float *yp, float x, float y)
    {
      int i, j, c = 0;
      for (i = 0, j = npol-1; i < npol; j = i++) {
        if ((((yp[i]<=y) && (y<yp[j])) ||
             ((yp[j]<=y) && (y<yp[i]))) &&
            (x < (xp[j] - xp[i]) * (y - yp[i]) / (yp[j] - yp[i]) + xp[i]))

          c = !c;
      }
      return c;
    }




Subject: Re: About Polyline
Author: Horst Kraemer <horst.kraemer@snafu.de>
Date Posted: Oct 2 1999 9:15:21:000AM


On Sat, 02 Oct 1999 09:34:42 GMT, philipma@ms5.hinet.net (Philip Ma)
wrote:

> Hi,
>
> Maybe the questions are stupid, but I need the answers. Thank you in
> advance!
>
> There are some questions about 2-D closed polyline as following:
>   1. How to know a polyline is crisscross (intersection itself) ?
>   2. If the polyline is not crisscross, how to know the polyline is
> clockwise or counterclockwise?
>   3. If the polyline is not crisscross, how to know a point is inside,
> outside or on the polyline?
>   4. If the polyline is not crisscross, how to know the area of the
> polyline?

In the order of ascending computational cost:

If

          P0=(x0,y0),P1=(x1,y1),....,PN=(xN,yN)=P0

is your closed polygon consisting of a sequence of N distinct points
P1,..,PN then the _oriented_ area is 1/2 of the sum of

        x_i*y_{i+1}-y_i*x_{i+1}

for i from 0 to N-1, i.e. for all adjacent pairs of vertices of the
closed polygon. This area may be positive or negative, according to
the orientation of the sequence relative to the coordinate system. In
any case its magnitude is the surface area.

If the coordinate system is oriented in a way that turning the
positive x-axis counterclockwise by 90deg will turn it into the
positive y-axis, then a positive area means that the polygon is
oriented counterclockwise. In fact "positive" doesn't mean
"counterclockwise". It just means "the same orientation as the
coordinate system". In fact in 2D geometry "clockwise" has no meaning.
If you look at a plane from the other side "clockwise" will turn into
"counterclockwise". You have to define additionally which side of the
plane is the "upside" and this definition can only be done by "looking
at it" from 3D. (Sorry for the side-step...)

There are many methods to test if a point of the plane is inside a
simple polygon. None of them is technical trivial.

The method I prefer is this one. It has been published by the german
mathematician Ahlemeyer, although he may not be the first who invented
it.  Let's assume for simplicity that the point to be tested is the
origin.

Divide the plane into 3 regions.

S0 consists only of the origin (0,0)

S+ is the "upper sector", consisting of all points with y>0 + all
points (x,0) with x>0

S- is the "lower sector", consisting of all points with y<0 + all
points (x,0) with x<0.

Every point of the plane belongs to exactly one sector.

Initialize a counter with 0.

Now make a round trip as before through all pairs of vertex points. If
any vertex is in S0 then (0,0) is a a vertex of the polygon and you
make your decision depending on the detail if an edge/vertex is
"inside" or "outside" the polygon.

Otherwise
  if two adjacent vertices (x_i,y_i),(x_{i+1},y_{i+1}) are in
  the same sector, do nothing.
  else if they are in two distinct sectors then calculate

      d =  x_i*y_{i+1) - y_i*x_{i+1)

  if d=0 then the origin is situated on the edge and you
  make your decision
  else if d>0 then add 1 to the counter
  else (d<0) subtract 1 from the counter

If no final decision could be taken in the middle then the point is
outside if the counter is ZERO, otherwise the point is outside. (Then
the counter is either 2 or -2).


Testing if a polygon is simple (not crossing itself) is probably the
most difficult task in terms of computation time. The crude and
straightforward approach is to test for each pair of edges if they
intersect or not. This will only be a problem for _very_ many
vertices.

Some time ago we hat a long thread about this problem in
sci.math.research, but nobody could come up with a _working_ algorithm
which was not O(n^2), i.e. which did need significantly less than
n^2/2 distinct tests. There were a lot of "convincing" proposals but
none of them was waterproof.


Regards
Horst
I've reached the age where the happy hour is a nap. (°¿°)
Windows 10 core i7 4790k 4Ghz 32GB GTX 970
Please support this web site.

Fuccaro

  • Guest
intersection
« Reply #14 on: April 22, 2004, 10:01:53 AM »
CAB
Receiving Devitg’s message I thought that I am on the right way and in short time I will finish my lisp. Now I have a lot of information to analyze... Just joking!
Before to ask help in this forum, I thought at two different solutions to find if a point X is inside of a triangle ABC. One was the method you sent me now: by summing the angles AXB+BXC+CXA. If the sum is =2*PI the point is inside. I anticipated some complications because in fact we are working in 3D and I prefer not to change the UCS.
Look at me: I write you long explanations, when in fact I wish to say you a single thing: thank you CAB!