Liu's method is a novel algorithm for efficiently testing whether a point is inside, exactly on, or outside a polyhedron defined by a triangular mesh.  The basic idea is that certain faces (triangles) of the mesh are sufficient to perform the point containment test; once one of these determining faces can be found, the test point's relation to this face is sufficient to determine containment.  The method is especially efficient for large-scale problems in which the polyhedron is composed of very many vertices and faces, and in which thousands of points are to be tested against the same polyhedron.  The method is extensible to N dimensions, and has been implemented in 2D and 3D.  A multi-material version has also been implemented; see the August 7 entry below.

 

The method requires the polyhedron to be non-leaky.  I.e., the polyhedron must not have any holes.

 

Liu's method is available below in two forms: the source code for a stand-alone C++ program (both single- and multi-material versions), and an R package (single-material version).  Note that the C++ code assumes that the indexing of points in the triangular faces starts at 0, not 1.  However, the R package assumes indexing starts at 1, following the R indexing convention.

 

Reference: J. Liu, Y.Q. Chen, J.M. Maisog, G. Luta, A new point containment test algorithm based on preprocessing and determining triangles, Computer-Aided Design, Volume 42, Issue 12, December 2010, Pages 1143-1150

 

R package available on R-Forge and on CRAN; see October 31, 2010 post below.  This is a 3D implementation of the original single-material version.

 

C++ source code for a 3D implementation of the stand-alone single-material version available for Windows and for Linux; see August 10, 2010 post below.

 

C++ source code for a 3D multi-material variant available in a ZIP file, a TAR file, and a RAR file; see August 7, 2010 post below.

 

-- Jose M. Maisog 

 

---

(December 29,2010) From this project, I have learned a lot about how to make R packages based on underlying C++ code, and how to use R-Forge.  I have put together a hands-on tutorial on the things I have learned, available here.  The tutorial comes with example C++ code and steps you through the process of building a working R package from the C++ code.  Chapter 6 in the tutorial (version 1.1) is an extension of the material I wrote earlier on how to use R-Forge.

 

 

---

(October 31, 2010) Version 1.9 of the R package is now available for Windows, Linux, and MacOS builds, on R-Forge and on CRAN.  This version fixes the memory leak and the crashing with leaky polyhedra, for all three builds.

 

To install directly from R-Forge, type (within the R environment):

 

install.packages("ptinpoly", repos="http://R-Forge.R-project.org")

 

To install directly from CRAN, type instead:

 

install.packages("ptinpoly", repos="http://cran.r-project.org")

 

Notes on setting up a repository on R-Forge, and on updating that repository.

 

---

 

(August 10, 2010) Here are updated distributions of the stand-alone single-material code.  The previous version would crash if the polyhedron had a topological defect like a hole.  This updated version uses try-catch when attempting to instantiate the PtInPolyhedron object, and if it fails it gives you a warning message rather than crashing.

 

Windows distribution: ptinpoly_081010.zip

Linux distribution: ptinpoly_081010.tgz

 

The two distributions are nearly identical.  The only difference is in the stopwatch class: under Windows, the macro for converting clock tics to seconds is CLK_TCK, whereas under Linux it's CLOCKS_PER_SEC.

 

To try the program on the simple cube data, run this command:

 

ptinpoly_081010.exe vertices_cube.txt faces_cube.txt queries_cube.txt

 

It should return the numbers 1, 0, 0, 0, -1, indicating that the first query point in queries_cube.txt was inside the cube, the next three were exactly on the cube, and the last was outside the cube.

 

Next, consider the case of a defective polyhedron; the file faces_cube_hole.txt is the same cube data, but with one of the faces removed.  So, this polyhedron is "leaky".  Run this command:

 

ptinpoly_081010.exe vertices_cube.txt faces_cube_hole.txt queries_cube.txt

This time, the program should warn you that the polyhedron has a topological defect and exit gracefully, rather than crashing ignominiously.

 

The distribution also includes a very simple program called checkPoly which performs two very basic checks on the polyhedron: each vertex must be used in the polyhedron at least three times, and each edge must be used exactly two times.  Feed this program the faces text file as input.

 

Under Linux, you can try running the C-shell script BuildAndRun.csh. This script will attempt to run the build, and then it'll try running the program.

 

---

 

(August 8, 2010) Version 1.9 of the R package (single material) is now available.  This fixes a memory leak problem.

 

ptinpoly v. 1.9 Windows build: ptinpoly_1.9.zip

ptinpoly v. 1.9 Linux build (package source): ptinpoly_1.9.tar.gz

 

---

 

(August 7, 2010) Multi-material version C++ code for a stand-alone program available: ZIP file; TAR file; RAR file.

 

When you run this executable, it will prompt you for the name of an input text file.  Try feeding it one of the two sample text files provided with the code (box.txt or twonestedboxes.txt).  Then the program will ask you how many random test points do you want to try; feed it some large integer, e.g. 10000.  Finally, to quit simply feed it any character, e.g. 'a'.

 

 

---

 

(July 13, 2010) Version 1.8 of the R package available!  This fixes the problem of pip3d crashing R if the polyhedron is "leaky".  If the polyhedron is indeed leaky, -2 is returned.  Also, pip3d would hang if more than two vertices were "very close" to one another; this threshold has been increased to three.  (Thanks to Jessica W. for providing motivation for these fixes.)  I also fixed the simple cube data that comes with the package; the verts and faces sample matrices appeared to represent some other polyhedron, not a cube.

 

ptinpoly v. 1.8 Windows build: ptinpoly_1.8.zip

ptinpoly v. 1.8 Linux build (package source): ptinpoly_1.8.tar.gz

 

I will update the official ptinpoly CRAN page by updating the R-Forge SVN respository; this will make a MacOS X build available, too.  And I'll update the R-Forge SVN repository when I figure out how to check files in and out in TortoiseSVN.

 

Here are some references for Dr. Liu's method:

 

J. Liu, Y.Q. Chen, J.M. Maisog, G. Luta, A new point containment test algorithm based on preprocessing and determining triangles, Computer-Aided Design, accepted.

W.P. Horn and D.L. Taylor, A theorem to determine the spatial containment of a point in a planar polygon, Computer Vision, Graphics and Image Processing, vol. 45, pp. 106-116,1989.

 

S. Nordbeck and B. Rysedt, Computer cartography point-in-polygon programs, BIT, vol. 7, pp. 39-64, 1967.

 

---

 

(September 13, 2009)  A ZIP file containing the C++ code for a stand-alone version of Liu's Point-In-Polyhedron code can be found HERE.  When unzipped, you should see the following files:

 

tmain.cpp

pinpolyhedron.h

pinpolyhedron.cpp

kodtree.h

kodtree.cpp 

i.8

JianfeiPtInPoly.dev

Makefile.win

 

tmain.cpp is the main program, that reads in data.  pinpolyhedron.h, pinpolyhedron.cpp, kodtree.h, and kodtree.cpp implement the algorithm.  i.8 is a text file containing the data for a simple polyhedron defining a cube, as well as five test points (data is originally from O'Rourke).  JianfeiPtInPoly.dev is the Dev-C++ project file, analogous to a .sln file in Visual Studio.  Makefile.win is a makefile that was automatically generated by Dev-C++; I haven't tried it, but it may work under Linux.  (Dev-C++ is a free C/C++ compiler; I have written a tutorial on it here.)

 

Let's assume that you have built the stand-alone executable, and it is named JianfeiPtInPoly.exe.  To run it on the data in i.8, simply use file redirection to feed i.8 to the program:

 

JianfeiPtInPoly.exe < i.8

 

It should generate the following output:

 

Polyhedron Vertices:

     i      x      y      z

     1      0      0      0

     2     10      0      0

     3      0     10      0

     4     10     10      0

     5      0      0     10

     6     10      0     10

     7      0     10     10

     8     10     10     10

numvert =   8 vertices read

About to read in  12 faces...

 

query point: [  5,  5,  5 ] : inside

query point: [ 10, 10, 10 ] : exactly on

query point: [  3,  0,  0 ] : exactly on

query point: [  3,  2, 10 ] : exactly on

query point: [ 11, 12, -3 ] : outside