MathGroup Archive 2000

[Date Index] [Thread Index] [Author Index]

Search the Archive

Point inside a plygon (again)

  • To: mathgroup at smc.vnet.net
  • Subject: [mg25402] Point inside a plygon (again)
  • From: Adriano Pascoletti <pascolet at dimi.uniud.it>
  • Date: Fri, 29 Sep 2000 01:06:42 -0400 (EDT)
  • Sender: owner-wri-mathgroup at wolfram.com

Dear MathGroup,

some days ago ([mg25323] Re: [mg25239] Point inside a polygon?) I suggested
Adriano Moreira
>See the constructive proof of Theorem 2.1, Sect. 2.2 of
>F.P.Preparata, M.I.Shamos: Computational Geometry, Springer 1985

Here is an implementation (with some comments) whose accuracy depends on
the accuracy of two tests
i. x==y (x,y reals)
ii. Det[m]>0 (m 3 by 3).


Adriano Pascoletti



(***********************************************************************

                    Mathematica-Compatible Notebook

This notebook can be used on any computer system with Mathematica 4.0,
MathReader 4.0, or any compatible application. The data for the notebook
starts with the line containing stars above.

To get the notebook into a Mathematica-compatible application, do one of
the following:

* Save the data starting with the line of stars above into a file
  with a name ending in .nb, then open the file inside the application;

* Copy the data starting with the line of stars above to the
  clipboard, then use the Paste menu command inside the application.

Data for notebooks contains only printable 7-bit ASCII and can be
sent directly in email or through ftp in text mode.  Newlines can be
CR, LF or CRLF (Unix, Macintosh or MS-DOS style).

NOTE: If you modify the data for this notebook not in a Mathematica-
compatible application, you must delete the line below containing the
word CacheID, otherwise Mathematica-compatible applications may try to
use invalid cache data.

For more information on notebooks and Mathematica-compatible
applications, contact Wolfram Research:
  web: http://www.wolfram.com
  email: info at wolfram.com
  phone: +1-217-398-0700 (U.S.)

Notebook reader applications are available free of charge from
Wolfram Research.
***********************************************************************)

(*CacheID: 232*)


(*NotebookFileLineBreakTest
NotebookFileLineBreakTest*)
(*NotebookOptionsPosition[      7063,        207]*)
(*NotebookOutlinePosition[      7905,        234]*)
(*  CellTagsIndexPosition[      7861,        230]*)
(*WindowFrame->Normal*)



Notebook[{

Cell[CellGroupData[{
Cell["Point location w.r. to a simple polygon", "Section"],

Cell[TextData[{
  "Adriano Pascoletti\n",
  ButtonBox["pascolet at dimi.uniud.it",
    ButtonData:>{
      URL[ "mailto:pascolet at dimi.uniud.it"], None},
    ButtonStyle->"Hyperlink"]
}], "Text"],

Cell[TextData[{
  "PLSPolygon (Point Location w.r.to a simple polygon) takes a list of 2D \
points (poly), a query point (q) and returns ",
  StyleBox["inside",
    FontSlant->"Italic"],
  " if q is in the interior of poly, ",
  StyleBox["outside",
    FontSlant->"Italic"],
  " if q is on the boundary or in the exterior."
}], "Text"],

Cell[TextData[{
  "The approach follows the ",
  Cell[BoxData[
      \(TraditionalForm\`O(n)\)]],
  " method suggested constructive proof of Theorem 2.1, Sect. 2.2 of\n\
F.P.Preparata, M.I.Shamos: Computational Geometry, Springer 1985. \n\nIf a \
horizontal line goes from -\[Infinity] to q then count the edges crossing the \
halfline strictly before q: if there is an odd number of intersections then q \
belongs to the interior of poly.\nA symbolic sufficiently small clockwise \
rotation of the halfline around q  lets us assume that  the halfline doesn't \
pass through any vertex.\nIn view of this it is enough considering \n(1) non \
horizontal edges ",
  Cell[BoxData[
      \(TraditionalForm\`\((a, b)\)\)]],
  " s.t. \n(2) one vertex has the y-coordinate \[GreaterEqual] ",
  Cell[BoxData[
      \(TraditionalForm\`q\_y\)]],
  " and the other ",
  Cell[BoxData[
      \(TraditionalForm\`\(\(<\)\(q\_y\)\)\)]],
  ".\n If ",
  Cell[BoxData[
      \(TraditionalForm\`a\ b\)]],
  " is such an edge one can assume ",
  Cell[BoxData[
      \(TraditionalForm\`a\_y > b\_y\)]],
  "(see 3) and ",
  Cell[BoxData[
      \(TraditionalForm\`a\ b\)]],
  " intersects the halfline through q before q iff the signed area of the \
triangle ",
  Cell[BoxData[
      \(TraditionalForm\`a\ b\ q\)]],
  " is positive. \n\nThe signed double area is computed by means of  ",
  Cell[BoxData[
      FormBox[
        RowBox[{"det", "(", GridBox[{
              {\(a\_x\), \(a\_y\), "1"},
              {\(b\_x\), \(b\_y\), "1"},
              {\(q\_x\), \(q\_y\), "1"}
              }], ")"}], TraditionalForm]]],
  "\nThe algorithm doesn't require floating point divisions.\n\n"
}], "Text"],

Cell[BoxData[
    \(\(PLSPolygon[poly : {\(({_, _})\) .. }, q : {x_, y_}] /;
          Length[poly] \[GreaterEqual] 3 :=
        Module[{edges = \((Append[poly, First[poly]] //
                  Partition[#, 2, 1] &)\), temp},  (*\
            1 : \ drop\ horizontal\ edges\ *) \
          temp = DeleteCases[edges, {{x1_, yy_}, {x2_, yy_}}];  (*\
            2 : \ drop\ edges\ with\ both\ vertex\ ordinates\  \
\[GreaterEqual] y\ or\ with\ both\  < \ y\ *) \ \[IndentingNewLine]temp =
            DeleteCases[
              temp, {{x1_, y1_}, {x2_, y2_}} /;
                Min[y1, y2] \[GreaterEqual] y \[Or] Max[y1, y2] < y];  (*\
            3 : \ edges\ in\ temp\ intersect\ a\ horizontal\ line\ through\ \
q; \ direct\ them\ downwards\ *) \[IndentingNewLine]temp =
            Map[If[#\[LeftDoubleBracket]1,
                      2\[RightDoubleBracket] > #\[LeftDoubleBracket]2,
                      2\[RightDoubleBracket], #, Reverse[#]] &,
              temp]; \[IndentingNewLine] (*\
            4 : \ an\ edge\ \((a, b)\)\ with\ y \((a)\) > \
                y \((b)\)\ intersects\ the\ horizontal\ halfline\ running\ \
from\  - \[Infinity]\ to\ q\ at\ a\ point\  \[NotEqual]
                q\ iff\ the\ signed\ area\ of\ the\ triangle\ \((a, b,
                    q)\)\ is\ positive\ *) \ \ \[IndentingNewLine]If[
            OddQ[Count[temp,
                e_ /; darea[Append[e, q]] >
                    0]], "\<inside\>", "\<outside\>"]];\)\)], "Input"],

Cell[BoxData[
    \(darea[triangle : {p1_, p2_, p3_}] :=
      Det[\(Append[#, 1] &\) /@ triangle]\)], "Input"],

Cell[CellGroupData[{

Cell["test", "Subsection"],

Cell[BoxData[{
    \(\(p = {{0, 0}, {1, \(-1\)}, {2, 0}, {2.5, \(-1.1\)}, {3, 0}, {4,
            0}, {4, 1}, {2.5, 2.5}, {0.5, 1}};\)\), "\[IndentingNewLine]",
    \(\(Show[
        Graphics[{{Thickness[0.007], Line[Append[p, First[p]]]},
            Line[{{\(-1\), 0}, {5, 0}}],
            Line[{{0, \(-1.5\)}, {0, 3}}]}]];\)\)}], "Input"],

Cell[CellGroupData[{

Cell[BoxData[
    \(PLSPolygon[p, {1, 0}]\)], "Input"],

Cell[BoxData[
    \("inside"\)], "Output"]
}, Open  ]],

Cell[BoxData[
    \(\(Show[
        Graphics[{Line[Append[p, First[p]]], PointSize[0.015], Point /@ p,
            RGBColor[1, 0, 0], PointSize[0.02],
            Point /@ Table[{i, 0}, {i, \(-1\)/2, 4 + 1/2,
                  1}]}]];\)\)], "Input"],

Cell[CellGroupData[{

Cell[BoxData[
    \(\(PLSPolygon[p, #] &\) /@
      Table[{i, 0}, {i, \(-1\)/2, 4 + 1/2, 1}]\)], "Input"],

Cell[BoxData[
    \({"outside", "inside", "inside", "inside", "outside",
      "outside"}\)], "Output"]
}, Open  ]],

Cell[CellGroupData[{

Cell[" ", "Subsubsection"],

Cell["\<\
Counterexample by William Self to David Park's solution ([mg25350] \
RE: [mg25239] Point inside a plygon? and [mg24992] )\
\>", "Text"],

Cell[CellGroupData[{

Cell[BoxData[
    \(PLSPolygon[{{0, 1}, {0, \(-1\)}, \ {2, 0}}, {1, 0}]\)], "Input"],

Cell[BoxData[
    \("inside"\)], "Output"]
}, Open  ]]
}, Open  ]]
}, Open  ]]
}, Open  ]]
},
FrontEndVersion->"4.0 for Macintosh",
ScreenRectangle->{{0, 832}, {0, 604}},
WindowSize->{520, 482},
WindowMargins->{{1, Automatic}, {Automatic, 1}},
MacintoshSystemPageSetup->"\<\
00<0001804P000000`d26_oQon at 3:`8g0dL5N`?P0080001804P000000]P2:001
0000I00000400`<30?l00BL?00400 at 0000000000000006P801T1T00000000000
00000000004000000000000000000000\>"
]


(***********************************************************************
Cached data follows.  If you edit this Notebook file directly, not using
Mathematica, you must remove the line containing CacheID at the top of
the file.  The cache data will then be recreated when you save this file
from within Mathematica.
***********************************************************************)

(*CellTagsOutline
CellTagsIndex->{}
*)

(*CellTagsIndex
CellTagsIndex->{}
*)

(*NotebookFileOutline
Notebook[{

Cell[CellGroupData[{
Cell[1739, 51, 58, 0, 56, "Section"],
Cell[1800, 53, 191, 6, 46, "Text"],
Cell[1994, 61, 335, 9, 64, "Text"],
Cell[2332, 72, 1676, 42, 319, "Text"],
Cell[4011, 116, 1491, 25, 379, "Input"],
Cell[5505, 143, 112, 2, 43, "Input"],

Cell[CellGroupData[{
Cell[5642, 149, 26, 0, 46, "Subsection"],
Cell[5671, 151, 346, 6, 123, "Input"],

Cell[CellGroupData[{
Cell[6042, 161, 54, 1, 27, "Input"],
Cell[6099, 164, 42, 1, 26, "Output"]
}, Open  ]],
Cell[6156, 168, 252, 5, 75, "Input"],

Cell[CellGroupData[{
Cell[6433, 177, 106, 2, 27, "Input"],
Cell[6542, 181, 104, 2, 26, "Output"]
}, Open  ]],

Cell[CellGroupData[{
Cell[6683, 188, 26, 0, 42, "Subsubsection"],
Cell[6712, 190, 145, 3, 46, "Text"],

Cell[CellGroupData[{
Cell[6882, 197, 84, 1, 27, "Input"],
Cell[6969, 200, 42, 1, 26, "Output"]
}, Open  ]]
}, Open  ]]
}, Open  ]]
}, Open  ]]
}
]
*)




(***********************************************************************
End of Mathematica Notebook file.
***********************************************************************)



----------------------------------------
Adriano Pascoletti
Dipartimento di Matematica e Informatica
Universita' di Udine
Via delle Scienze 206
I-33100 UDINE
Italy

Tel. +39-0432-558441
Fax. +39-0432-558499

http://www.dimi.uniud.it/~pascolet




  • Prev by Date: troubles with 3D plot
  • Next by Date: Why is Mathematica so slow ?
  • Previous by thread: troubles with 3D plot
  • Next by thread: Why is Mathematica so slow ?