Python Bindings for Cork

Makefile

@@ -62,7 +62,7 @@ INC       := $(INC) $(GMPINC)
 # use the second line to disable profiling instrumentation
 # PROFILING := -pg
 PROFILING :=
-CCFLAGS   := -Wall $(INC) $(CONFIG) -O2 -DNDEBUG $(PROFILING)
+CCFLAGS   := -Wall $(INC) $(CONFIG) -O2 -DNDEBUG $(PROFILING) -fPIC
 CXXFLAGS  := $(CCFLAGS) $(CPP11_FLAGS)
 CCDFLAGS  := -Wall $(INC) $(CONFIG) -ggdb
 CXXDFLAGS := $(CCDFLAGS)
@@ -195,6 +195,7 @@ obj/isct/triangle.o: src/isct/triangle.c
 	@echo "Compiling the Triangle library"
 	@$(CC) -O2 -DNO_TIMER \
                -DREDUCED \
+               -fPIC \
                -DCDT_ONLY -DTRILIBRARY \
                -Wall -DANSI_DECLARATORS \
                -o obj/isct/triangle.o -c src/isct/triangle.c

python/cork.i

@@ -0,0 +1,196 @@
+/* This file is part of the Cork library.
+
+ * Cork is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as
+ * published by the Free Software Foundation, either version 3 of
+ * the License, or (at your option) any later version.
+
+ * Cork is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU Lesser General Public License for more details.
+
+ * You should have received a copy 
+ * of the GNU Lesser General Public License
+ * along with Cork.  If not, see <http://www.gnu.org/licenses/>.
+ */
+/*
+ * @author Stephen Dawson-Haggerty <[email protected]>
+ */
+
+%module cork
+%{
+#include <numpy/arrayobject.h>
+#include "cork.h"
+
+
+/* create a "base object" that lets us use the memory returned by the
+ *  cork library.  this is DECREFed when the numpy array we return is
+ *  deleted, and we can use that to free the c++ array using the
+ *  appropriate deallocator.
+ */
+typedef struct {
+  PyObject_HEAD
+  uint *triangles;
+  float *vertices;
+} _DeleteObject;
+
+static void _DeleteObject_dealloc(PyObject *self) {
+  _DeleteObject *real_self = (_DeleteObject *)self;
+  if (real_self->triangles) {
+    delete real_self->triangles;
+  }
+  if (real_self->vertices) {
+    delete real_self->vertices;
+  }
+  self->ob_type->tp_free((PyObject *)self);
+}
+
+static PyTypeObject _DeleteObjectType = {
+  PyObject_HEAD_INIT(NULL)
+  0, 
+  "delete_deallocator",
+  sizeof(_DeleteObject),
+  0,
+  _DeleteObject_dealloc,         /* tp_dealloc */
+  0, /*tp_print*/
+  0, /*tp_getattr*/
+  0, /*tp_setattr*/
+  0, /*tp_compare*/
+  0, /*tp_repr*/
+  0, /*tp_as_number*/
+  0, /*tp_as_sequence*/
+  0, /*tp_as_mapping*/
+  0, /*tp_hash */
+  0, /*tp_call*/
+  0, /*tp_str*/
+  0, /*tp_getattro*/
+  0, /*tp_setattro*/
+  0, /*tp_as_buffer*/
+  Py_TPFLAGS_DEFAULT, /*tp_flags*/
+  "c++ delete deallocator", /* tp_doc */
+};
+
+%}
+
+%init %{
+  /* do this, or segfault */
+  import_array();
+
+  /* add the new deallocator type to the system */
+  _DeleteObjectType.tp_new = PyType_GenericNew;
+  if (PyType_Ready(&_DeleteObjectType) < 0) {
+    return;
+  }
+
+  Py_INCREF(&_DeleteObjectType);
+  PyModule_AddObject(m, "_cork", (PyObject *)&_DeleteObjectType);
+%}
+
+/* Create typemaps for mapping python types to/from a CorkTriMesh 
+ *
+ * We allocate temporary PyObjects so that we can use the
+ * PyArray_FromAny to accept a lot of different potential kinds of
+ * matrixes.  These get cleaned up by the freearg typemap.
+ */
+%typemap(in) CorkTriMesh (PyObject *a = NULL, PyObject *b = NULL) {
+  npy_intp *d1, *d2;
+
+  if (!PyTuple_Check($input) || PyTuple_Size($input) != 2) {
+    PyErr_SetString(PyExc_TypeError, "argument must be a tuple of length 2");
+    return NULL;
+  }
+  a = PyTuple_GetItem($input, 0);
+  b = PyTuple_GetItem($input, 1);
+  if (a == NULL || b == NULL) {
+    PyErr_SetString(PyExc_ValueError, "argument must not be none");
+    a = b = NULL;
+    return NULL;
+  }
+
+  /* create the PyArrays */
+  a = PyArray_FromAny(a, PyArray_DescrFromType(NPY_UINT), 2, 2, NPY_ARRAY_CARRAY, NULL);
+  if (a == NULL) {
+    return NULL;
+  }
+  b = PyArray_FromAny(b, PyArray_DescrFromType(NPY_FLOAT), 2, 2, NPY_ARRAY_CARRAY, NULL);
+  if (b == NULL) {
+    Py_DECREF(a);
+    a = b = NULL;
+    return NULL;
+  }
+
+  /* check the dimensions and get the heights... */
+  d1 = PyArray_DIMS(a);
+  d2 = PyArray_DIMS(b);
+  if (d1[1] != 3 || d2[1] != 3) {
+    PyErr_SetString(PyExc_ValueError, "arrays must be Nx3");
+    return NULL;
+  }
+
+  /* how to track allocations.. we can't free this since  */
+  $1.n_triangles = d1[0];
+  $1.n_vertices = d2[0];
+  $1.triangles = (uint *)PyArray_DATA(a);
+  $1.vertices = (float *)PyArray_DATA(b);
+}
+
+%typemap(freearg) CorkTriMesh {
+  if (a$argnum) {
+    Py_DECREF(a$argnum);
+  }
+  if (b$argnum) {
+    Py_DECREF(b$argnum);
+  }
+}
+
+/* allocate a temporary CorkTriMesh on the way in for use as an output
+   parameter */
+%typemap(in, numinputs=0) CorkTriMesh *OUTPUT (CorkTriMesh temp) {
+  memset(&temp, 0, sizeof(temp));
+  $1 = &temp;
+}
+
+/* build numpy matrices for the return type reusing the memory
+   allocated by Cork. */
+%typemap(argout) CorkTriMesh *OUTPUT {
+  /* return some numpy arrays with the coordinates in there */
+  PyObject *a, *b;
+  _DeleteObject *base_a, *base_b;
+  npy_intp dims[2];
+  dims[1] = 3;
+
+
+
+  dims[0] = temp$argnum.n_triangles;
+  a = PyArray_SimpleNewFromData(2, dims, NPY_UINT, temp$argnum.triangles);
+
+  dims[0] = temp$argnum.n_vertices;
+  b = PyArray_SimpleNewFromData(2, dims, NPY_FLOAT, temp$argnum.vertices);
+
+  /* set a base object on the numpy array that will deallocate our c++ buffer */
+  base_a = PyObject_New(_DeleteObject, &_DeleteObjectType);
+  base_a->triangles = temp$argnum.triangles;
+  base_a->vertices = NULL;
+  PyArray_SetBaseObject((PyArrayObject *)a, (PyObject *)base_a);
+
+  base_b = PyObject_New(_DeleteObject, &_DeleteObjectType);
+  base_b->triangles = NULL;
+  base_b->vertices = temp$argnum.vertices;
+  PyArray_SetBaseObject((PyArrayObject *)b, (PyObject *)base_b);
+
+
+  $result = Py_BuildValue("(OO)", a, b);
+  /* BuildValue calls Py_INCREF on its arguments */
+  Py_DECREF(a);
+  Py_DECREF(b);
+}
+
+/* Export the Cork API */
+extern bool isSolid(CorkTriMesh mesh);
+extern void computeUnion(CorkTriMesh in0, CorkTriMesh in1, CorkTriMesh *OUTPUT);
+extern void computeDifference(CorkTriMesh in0, CorkTriMesh in1, CorkTriMesh *OUTPUT);
+extern void computeIntersection(CorkTriMesh in0, CorkTriMesh in1, CorkTriMesh *OUTPUT);
+extern void computeSymmetricDifference(
+                        CorkTriMesh in0, CorkTriMesh in1, CorkTriMesh *OUTPUT);
+extern void resolveIntersections(CorkTriMesh in0, CorkTriMesh in1, CorkTriMesh *OUTPUT);

python/test.py

@@ -0,0 +1,112 @@
+"""
+This file is part of the Cork library.
+
+Cork is free software: you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as
+published by the Free Software Foundation, either version 3 of
+the License, or (at your option) any later version.
+
+Cork is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy 
+of the GNU Lesser General Public License
+along with Cork.  If not, see <http://www.gnu.org/licenses/>.
+"""
+"""
+Test that the python bindings produce the same answer as the c program
+on the sample input.
+
+For this to work, you must edit src/util/prelude.h and change initRand
+to have an srand(0) in it, so we can get the same results from both
+runs.
+
+@author Stephen Dawson-Haggerty <[email protected]>
+"""
+
+import os
+import unittest
+import subprocess
+import numpy as np
+import ctypes
+import ctypes.util
+import _cork
+
+class TestPythonBindings(unittest.TestCase):
+    """Compare the output of the cork binary to the python library
+    """
+    BALL_A = os.path.join(os.path.dirname(__file__), 
+                          "../samples/ballA.off")
+    BALL_B = os.path.join(os.path.dirname(__file__), 
+                          "../samples/ballB.off")
+
+    @staticmethod
+    def read_mesh(f):
+        fp = open(f, "r")
+        fp.readline()
+        dims = map(int, fp.readline().split(" "))
+        triangles, vertices = [], []
+        for i in xrange(0, dims[0]):
+            vertices.append(map(float, fp.readline().split(" ") ))
+
+        for i in xrange(0, dims[1]):
+            triangles.append(map(float, fp.readline().split(" ")[1:] ))
+        return triangles, vertices
+
+    def call_cork(self, method):
+        args = [os.path.join(os.path.dirname(__file__), '../bin/cork'),
+                '-' + method,
+                self.BALL_A, self.BALL_B, '/tmp/_temp_mesh.off']
+        subprocess.check_call(args)
+        return self.read_mesh('/tmp/_temp_mesh.off')
+
+    def assertTriMeshEquals(self, a, b):
+        t1, v1 = a
+        t2, v2 = b
+
+        # assert these are close
+        self.assertTrue((np.sum(np.abs(v1 - v2)) / np.sum(np.abs(v1))) < 1e-6)
+        # these should be equal since they're integers to begin with
+        self.assertEqual(np.sum(np.abs(t1 - t2)), 0)
+
+    def setUp(self):
+        self.ballA = self.read_mesh(self.BALL_A)
+        self.ballB = self.read_mesh(self.BALL_B)
+
+        # are we having fun yet?  cork uses srand to seed the crappy
+        # PRNG, but doesn't have a way to reset it.  This way we can
+        # call it between runs so we get the same answer every time.
+        libc = ctypes.cdll.LoadLibrary(ctypes.util.find_library('c'))
+        libc.srand(0)
+
+    def testIsSolid(self):
+        """Both demo objects are solid
+        """
+        self.assertTrue(_cork.isSolid(self.ballA))
+        self.assertTrue(_cork.isSolid(self.ballB))
+
+    def testUnion(self):
+        self.assertTriMeshEquals(_cork.computeUnion(self.ballA, self.ballB),
+                                 self.call_cork('union'))
+
+    def testDifference(self):
+        self.assertTriMeshEquals(_cork.computeDifference(self.ballA, self.ballB),
+                                 self.call_cork('diff'))
+
+    def testIntersection(self):
+        self.assertTriMeshEquals(_cork.computeIntersection(self.ballA, self.ballB),
+                                 self.call_cork('isct'))
+
+    def testSymmetricDifference(self):
+        self.assertTriMeshEquals(_cork.computeSymmetricDifference(self.ballA, self.ballB),
+                                 self.call_cork('xor'))
+
+    def testResolve(self):
+        self.assertTriMeshEquals(_cork.resolveIntersections(self.ballA, self.ballB),
+                                 self.call_cork('resolve'))
+
+
+if __name__ == '__main__':
+    unittest.main()

setup.py

@@ -0,0 +1,68 @@
+"""
+This file is part of the Cork library.
+
+Cork is free software: you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as
+published by the Free Software Foundation, either version 3 of
+the License, or (at your option) any later version.
+
+Cork is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU Lesser General Public License for more details.
+
+You should have received a copy 
+of the GNU Lesser General Public License
+along with Cork.  If not, see <http://www.gnu.org/licenses/>.
+
+"""
+"""
+The setup.py file, and the associated python/cork.i file provide the
+typemaps necessary to call cork from Python.
+
+Before building the Python module, you must have compiled the cork
+library; it should be placed in lib/ (the default makefile does this.
+
+The swig interface contains typemaps which compile to a single python
+module, named _cork.  It exports the definitions from cork.h, except
+for freeCorkTriMesh, which is not necessary.
+
+All of these function use the CorkTriMesh type for both input and
+output, which is very convenient.  In Python, this type is represented
+as a tuple of of (triangles, vertices).  For input, these can be
+anything implementing the python iterator protocol.  If you input
+numpy matrices, with triangles represented as a uint32 Nx3 matrix, and
+the vertices as a float32 Nx3 matrix, it should be possible to work
+without a copy.  
+
+The return values are also represented as numpy matrices.  
+
+Currently, input/output from .off files is not implemented.
+
+@author Stephen Dawson-Haggerty <[email protected]>
+"""
+
+from distutils.core import setup, Extension
+import numpy as np
+
+cork_module = Extension('_cork',
+  sources=['python/cork.i'],
+    language="c++",
+
+  # build extension wrapper with c++11 support
+  swig_opts=['-c++', '-threads'],
+  extra_compile_args=['-std=c++11'],
+
+  libraries=['cork', 'gmp'],
+  library_dirs=['lib/'],
+  include_dirs=['src', np.get_include()])
+
+setup(
+    name='cork',
+    version='0.1',
+    author='Stephen Dawson-Haggerty',
+    description=('Python interface to the cork library'),
+    license='GPLv3',
+    requires=['numpy'],
+    ext_modules=[cork_module],
+)