Browse files

Merge pull request #377 from pv/qhull-imp

ENH: spatial/qhull: add convex hull + voronoi qhull wrappers, improve docs, add tutorial

This PR completes the Qhull wrappers in scipy.spatial, adding support
for convex hulls and Voronoi diagrams, and the incremental construction
and furthest-site modes.

A tutorial plus convenience 2D plotting routines are also added.
  • Loading branch information...
2 parents a257657 + 7b500f9 commit d0ddf549ada070e8b6fc83b9d5b62cdea971b2c9 @pv pv committed Dec 17, 2012
Showing with 24,180 additions and 10,293 deletions.
  1. +13 −0 doc/release/0.12.0-notes.rst
  2. +1 −0 doc/source/tutorial/index.rst
  3. +275 −0 doc/source/tutorial/spatial.rst
  4. +224 −192 scipy/interpolate/interpnd.c
  5. +30 −30 scipy/interpolate/interpnd.pyx.in
  6. +2 −1 scipy/spatial/SConscript
  7. +52 −3 scipy/spatial/__init__.py
  8. +160 −0 scipy/spatial/_plotutils.py
  9. +2 −1 scipy/spatial/bscript
  10. +0 −29 scipy/spatial/generate_qhull.py
  11. +21,015 −7,904 scipy/spatial/qhull.c
  12. +1 −1 scipy/spatial/qhull.pxd
  13. +1,165 −163 scipy/spatial/qhull.pyx
  14. +3 −4 scipy/spatial/qhull/Announce.txt
  15. +1 −1 scipy/spatial/qhull/COPYING.txt
  16. +328 −321 scipy/spatial/qhull/README.txt
  17. +0 −1,343 scipy/spatial/qhull/src/Changes.txt
  18. +3 −3 scipy/spatial/qhull/src/geom.c
  19. +3 −3 scipy/spatial/qhull/src/geom.h
  20. +3 −3 scipy/spatial/qhull/src/geom2.c
  21. +11 −8 scipy/spatial/qhull/src/global.c
  22. +8 −10 scipy/spatial/qhull/src/io.c
  23. +3 −3 scipy/spatial/qhull/src/io.h
  24. +5 −5 scipy/spatial/qhull/src/libqhull.c
  25. +15 −7 scipy/spatial/qhull/src/libqhull.h
  26. +8 −4 scipy/spatial/qhull/src/mem.c
  27. +9 −9 scipy/spatial/qhull/src/mem.h
  28. +4 −4 scipy/spatial/qhull/src/merge.c
  29. +3 −3 scipy/spatial/qhull/src/merge.h
  30. +5 −4 scipy/spatial/qhull/src/poly.c
  31. +3 −3 scipy/spatial/qhull/src/poly.h
  32. +17 −6 scipy/spatial/qhull/src/poly2.c
  33. +0 −18 scipy/spatial/qhull/src/qhull.h
  34. +3 −3 scipy/spatial/qhull/src/qhull_a.h
  35. +132 −92 scipy/spatial/qhull/src/qset.c
  36. +10 −6 scipy/spatial/qhull/src/qset.h
  37. +2 −2 scipy/spatial/qhull/src/random.c
  38. +3 −3 scipy/spatial/qhull/src/random.h
  39. +2 −1 scipy/spatial/qhull/src/rboxlib.c
  40. +3 −3 scipy/spatial/qhull/src/stat.c
  41. +11 −4 scipy/spatial/qhull/src/stat.h
  42. +18 −10 scipy/spatial/qhull/src/user.c
  43. +15 −5 scipy/spatial/qhull/src/user.h
  44. +2 −2 scipy/spatial/qhull/src/usermem.c
  45. +1 −17 scipy/spatial/qhull/src/userprintf.c
  46. +53 −0 scipy/spatial/qhull/src/userprintf_rbox.c
  47. +7 −16 scipy/spatial/setup.py
  48. +45 −0 scipy/spatial/tests/test__plotutils.py
  49. +501 −43 scipy/spatial/tests/test_qhull.py
View
13 doc/release/0.12.0-notes.rst
@@ -39,6 +39,19 @@ count_neighbors and sparse_distance_matrix) are between 200 and 1000 times
faster in cKDTree than in KDTree. With very minor caveats, cKDTree has
exactly the same interface as KDTree, and can be used as a drop-in replacement.
+Voronoi diagrams and convex hulls
+---------------------------------
+`scipy.spatial` now contains functionality for computing Voronoi
+diagrams and convex hulls using the Qhull library. (Delaunay
+triangulation was available since Scipy 0.9.0.)
+
+Delaunay improvements
+---------------------
+It's now possible to pass in custom Qhull options in Delaunay
+triangulation. Coplanar points are now also recorded, if present.
+Incremental construction of Delaunay triangulations is now also
+possible.
+
``scipy.optimize`` improvements
-------------------------------
A callback mechanism was added to L-BFGS-B and TNC minimization solvers.
View
1 doc/source/tutorial/index.rst
@@ -18,6 +18,7 @@ SciPy Tutorial
linalg
arpack
csgraph
+ spatial
stats
ndimage
io
View
275 doc/source/tutorial/spatial.rst
@@ -0,0 +1,275 @@
+.. _qhulltutorial:
+
+Spatial data structures and algorithms (`scipy.spatial`)
+========================================================
+
+.. currentmodule:: scipy.spatial
+
+`scipy.spatial` can compute triangulations, Voronoi diagrams, and
+convex hulls of a set of points, by leveraging the `Qhull
+<http://qhull.org/>`__ library.
+
+Moreover, it contains `KDTree` implementations for nearest-neighbor point
+queries, and utilities for distance computations in various metrics.
+
+Delaunay triangulations
+-----------------------
+
+The Delaunay triangulation is a subdivision of a set of points into a
+non-overlapping set of triangles, such that no point is inside the
+circumcircle of any triangle. In practice, such triangulations tend to
+avoid triangles with small angles.
+
+Delaunay triangulation can be computed using `scipy.spatial` as follows:
+
+.. plot::
+
+ >>> from scipy.spatial import Delaunay
+ >>> points = np.array([[0, 0], [0, 1.1], [1, 0], [1, 1]])
+ >>> tri = Delaunay(points)
+
+ We can visualize it:
+
+ >>> import matplotlib.pyplot as plt
+ >>> plt.triplot(points[:,0], points[:,1], tri.simplices.copy())
+ >>> plt.plot(points[:,0], points[:,1], 'o')
+
+ And add some further decorations:
+
+ >>> for j, p in enumerate(points):
+ ... plt.text(p[0]-0.03, p[1]+0.03, j, ha='right') # label the points
+ >>> for j, s in enumerate(tri.simplices):
+ ... p = points[s].mean(axis=0)
+ ... plt.text(p[0], p[1], '#%d' % j, ha='center') # label triangles
+ >>> plt.xlim(-0.5, 1.5); plt.ylim(-0.5, 1.5)
+ >>> plt.show()
+
+The structure of the triangulation is encoded in the following way:
+the ``simplices`` attribute contains the indices of the points in the
+``points`` array that make up the triangle. For instance:
+
+>>> i = 1
+>>> tri.simplices[i,:]
+array([3, 1, 0], dtype=int32)
+>>> points[tri.simplices[i,:]]
+array([[ 1. , 1. ],
+ [ 0. , 1.1],
+ [ 0. , 0. ]])
+
+Moreover, neighboring triangles can also be found out:
+
+>>> tri.neighbors[i]
+array([-1, 0, -1], dtype=int32)
+
+What this tells us is that this triangle has triangle #0 as a neighbor,
+but no other neighbors. Moreover, it tells us that neighbor 0 is
+opposite the vertex 1 of the triangle:
+
+>>> points[tri.simplices[i, 1]]
+array([ 0. , 1.1])
+
+Indeed, from the figure we see that this is the case.
+
+Qhull can also perform tesselations to simplices also for
+higher-dimensional point sets (for instance, subdivision into
+tetrahedra in 3-D).
+
+
+Coplanar points
+^^^^^^^^^^^^^^^
+
+It is important to note that not *all* points necessarily appear as
+vertices of the triangulation, due to numerical precision issues in
+forming the triangulation. Consider the above with a duplicated
+point:
+
+>>> points = np.array([[0, 0], [0, 1], [1, 0], [1, 1], [1, 1]])
+>>> tri = Delaunay(points)
+>>> np.unique(tri.simplices.ravel())
+array([0, 1, 2, 3], dtype=int32)
+
+Observe that point #4, which is a duplicate, does not occur as a
+vertex of the triangulation. That this happened is recorded:
+
+>>> tri.coplanar
+array([[4, 0, 3]], dtype=int32)
+
+This means that point 4 resides near triangle 0 and vertex 3, but is
+not included in the triangulation.
+
+Note that such degeneracies can occur not only because of duplicated
+points, but also for more complicated geometrical reasons, even in
+point sets that at first sight seem well-behaved.
+
+However, Qhull has the "QJ" option, which instructs it to perturb the
+input data randomly until degeneracies are resolved:
+
+>>> tri = Delaunay(points, qhull_options="QJ Pp")
+>>> points[tri.simplices]
+array([[[1, 1],
+ [1, 0],
+ [0, 0]],
+ [[1, 1],
+ [1, 1],
+ [1, 0]],
+ [[0, 1],
+ [1, 1],
+ [0, 0]],
+ [[0, 1],
+ [1, 1],
+ [1, 1]]])
+
+Two new triangles appeared. However, we see that they are degenerate
+and have zero area.
+
+
+Convex hulls
+------------
+
+Convex hull is the smallest convex object containing all points in a
+given point set.
+
+These can be computed via the Qhull wrappers in `scipy.spatial` as
+follows:
+
+.. plot::
+
+ >>> from scipy.spatial import ConvexHull
+ >>> points = np.random.rand(30, 2) # 30 random points in 2-D
+ >>> hull = ConvexHull(points)
+
+ The convex hull is represented as a set of N-1 dimensional simplices,
+ which in 2-D means line segments. The storage scheme is exactly the
+ same as for the simplices in the Delaunay triangulation discussed
+ above.
+
+ We can illustrate the above result:
+
+ >>> import matplotlib.pyplot as plt
+ >>> plt.plot(points[:,0], points[:,1], 'o')
+ >>> for simplex in hull.simplices:
+ >>> plt.plot(points[simplex,0], points[simplex,1], 'k-')
+ >>> plt.show()
+
+The same can be achieved with `scipy.spatial.convex_hull_plot_2d`.
+
+
+Voronoi diagrams
+----------------
+
+A Voronoi diagram is a subdivision of the space into the nearest
+neighborhoods of a given set of points.
+
+There are two ways to approach this object using `scipy.spatial`.
+First, one can use the `KDTree` to answer the question "which of the
+points is closest to this one", and define the regions that way:
+
+.. plot::
+
+ >>> from scipy.spatial import KDTree
+ >>> points = np.array([[0, 0], [0, 1], [0, 2], [1, 0], [1, 1], [1, 2],
+ ... [2, 0], [2, 1], [2, 2]])
+ >>> tree = KDTree(points)
+ >>> tree.query([0.1, 0.1])
+ (0.14142135623730953, 0)
+
+ So the point ``(0.1, 0.1)`` belongs to region ``0``. In color:
+
+ >>> x = np.linspace(-0.5, 2.5, 31)
+ >>> y = np.linspace(-0.5, 2.5, 33)
+ >>> xx, yy = np.meshgrid(x, y)
+ >>> xy = np.c_[xx.ravel(), yy.ravel()]
+ >>> import matplotlib.pyplot as plt
+ >>> plt.pcolor(x, y, tree.query(xy)[1].reshape(33, 31))
+ >>> plt.plot(points[:,0], points[:,1], 'ko')
+ >>> plt.show()
+
+ This does not, however, give the Voronoi diagram as a geometrical
+ object.
+
+ The representation in terms of lines and points can be again
+ obtained via the Qhull wrappers in `scipy.spatial`:
+
+ >>> from scipy.spatial import Voronoi
+ >>> vor = Voronoi(points)
+ >>> vor.vertices
+ array([[ 0.5, 0.5],
+ [ 1.5, 0.5],
+ [ 0.5, 1.5],
+ [ 1.5, 1.5]])
+
+ The Voronoi vertices denote the set of points forming the polygonal
+ edges of the Voronoi regions. In this case, there are 9 different
+ regions:
+
+ >>> vor.regions
+ [[-1, 0], [-1, 1], [1, -1, 0], [3, -1, 2], [-1, 3], [-1, 2], [3, 1, 0, 2], [2, -1, 0], [3, -1, 1]]
+
+ Negative value ``-1`` again indicates a point at infinity. Indeed,
+ only one of the regions, ``[3, 1, 0, 2]``, is bounded. Note here that
+ due to similar numerical precision issues as in Delaunay triangulation
+ above, there may be fewer Voronoi regions than input points.
+
+ The ridges (lines in 2-D) separating the regions are described as a
+ similar collection of simplices as the convex hull pieces:
+
+ >>> vor.ridge_vertices
+ [[-1, 0], [-1, 0], [-1, 1], [-1, 1], [0, 1], [-1, 3], [-1, 2], [2, 3], [-1, 3], [-1, 2], [0, 2], [1, 3]]
+
+ These numbers indicate indices of the Voronoi vertices making up the
+ line segments. ``-1`` is again a point at infinity --- only four of
+ the 12 lines is a bounded line segment while the others extend to
+ infinity.
+
+ The Voronoi ridges are perpendicular to lines drawn between the
+ input points. Which two points each ridge corresponds to is also
+ recorded:
+
+ >>> vor.ridge_points
+ array([[0, 3],
+ [0, 1],
+ [6, 3],
+ [6, 7],
+ [3, 4],
+ [5, 8],
+ [5, 2],
+ [5, 4],
+ [8, 7],
+ [2, 1],
+ [4, 1],
+ [4, 7]], dtype=int32)
+
+ This information, taken together, is enough to construct the full
+ diagram.
+
+ We can plot it as follows. First the points and the Voronoi vertices:
+
+ >>> plt.plot(points[:,0], points[:,1], 'o')
+ >>> plt.plot(vor.vertices[:,0], vor.vertices[:,1], '*')
+ >>> plt.xlim(-1, 3); plt.ylim(-1, 3)
+
+ Plotting the finite line segments goes as for the convex hull,
+ but now we have to guard for the infinite edges:
+
+ >>> for simplex in vor.ridge_vertices:
+ >>> simplex = np.asarray(simplex)
+ >>> if np.all(simplex >= 0):
+ >>> plt.plot(vor.vertices[simplex,0], vor.vertices[simplex,1], 'k-')
+
+ The ridges extending to infinity require a bit more care:
+
+ >>> center = points.mean(axis=0)
+ >>> for pointidx, simplex in zip(vor.ridge_points, vor.ridge_vertices):
+ >>> simplex = np.asarray(simplex)
+ >>> if np.any(simplex < 0):
+ >>> i = simplex[simplex >= 0][0] # finite end Voronoi vertex
+ >>> t = points[pointidx[1]] - points[pointidx[0]] # tangent
+ >>> t /= np.linalg.norm(t)
+ >>> n = np.array([-t[1], t[0]]) # normal
+ >>> midpoint = points[pointidx].mean(axis=0)
+ >>> far_point = vor.vertices[i] + np.sign(np.dot(midpoint - center, n)) * n * 100
+ >>> plt.plot([vor.vertices[i,0], far_point[0]],
+ ... [vor.vertices[i,1], far_point[1]], 'k--')
+ >>> plt.show()
+
+This plot can also be created using `scipy.spatial.voronoi_plot_2d`.
View
416 scipy/interpolate/interpnd.c
@@ -53,12 +53,15 @@
(PyErr_Format(PyExc_TypeError, \
"expected index value, got %.200s", Py_TYPE(o)->tp_name), \
(PyObject*)0))
- #define PyIndex_Check(o) (PyNumber_Check(o) && !PyFloat_Check(o) && !PyComplex_Check(o))
+ #define __Pyx_PyIndex_Check(o) (PyNumber_Check(o) && !PyFloat_Check(o) && \
+ !PyComplex_Check(o))
+ #define PyIndex_Check __Pyx_PyIndex_Check
#define PyErr_WarnEx(category, message, stacklevel) PyErr_Warn(category, message)
#define __PYX_BUILD_PY_SSIZE_T "i"
#else
#define __PYX_BUILD_PY_SSIZE_T "n"
#define CYTHON_FORMAT_SSIZE_T "z"
+ #define __Pyx_PyIndex_Check PyIndex_Check
#endif
#if PY_VERSION_HEX < 0x02060000
#define Py_REFCNT(ob) (((PyObject*)(ob))->ob_refcnt)
@@ -508,7 +511,7 @@ struct __pyx_t_5scipy_7spatial_5qhull_DelaunayInfo_t {
int npoints;
int nsimplex;
double *points;
- int *vertices;
+ int *simplices;
int *neighbors;
double *equations;
double *transform;
@@ -1080,7 +1083,7 @@ static char __pyx_k_35[] = "Format string allocated too short.";
static char __pyx_k_37[] = "\nSimple N-D interpolation\n\n.. versionadded:: 0.9\n\n";
static char __pyx_k_38[] = "scipy.spatial.qhull";
static char __pyx_k_39[] = "*";
-static char __pyx_k_42[] = "/home/rgommers/Code/scipy/scipy/interpolate/interpnd.pyx";
+static char __pyx_k_42[] = "/home/pauli/prj/scipy/scipy/scipy/interpolate/interpnd.pyx";
static char __pyx_k_43[] = "scipy.interpolate.interpnd";
static char __pyx_k_51[] = "\n Common routines for interpolators.\n\n .. versionadded:: 0.9\n\n ";
static char __pyx_k_60[] = "\n LinearNDInterpolator(points, values)\n\n Piecewise linear interpolant in N dimensions.\n\n .. versionadded:: 0.9\n\n Parameters\n ----------\n points : ndarray of floats, shape (npoints, ndims)\n Data point coordinates.\n values : ndarray of float or complex, shape (npoints, ...)\n Data values.\n fill_value : float, optional\n Value used to fill in for requested points outside of the\n convex hull of the input points. If not provided, then\n the default is ``nan``.\n\n Notes\n -----\n The interpolant is constructed by triangulating the input data\n with Qhull [1]_, and on each triangle performing linear\n barycentric interpolation.\n\n References\n ----------\n .. [1] http://www.qhull.org/\n\n ";
@@ -1165,10 +1168,10 @@ static char __pyx_k____init__[] = "__init__";
static char __pyx_k____main__[] = "__main__";
static char __pyx_k____test__[] = "__test__";
static char __pyx_k__isimplex[] = "isimplex";
-static char __pyx_k__vertices[] = "vertices";
static char __pyx_k__warnings[] = "warnings";
static char __pyx_k__enumerate[] = "enumerate";
static char __pyx_k__eps_broad[] = "eps_broad";
+static char __pyx_k__simplices[] = "simplices";
static char __pyx_k__transpose[] = "transpose";
static char __pyx_k__ValueError[] = "ValueError";
static char __pyx_k__asanyarray[] = "asanyarray";
@@ -1278,13 +1281,13 @@ static PyObject *__pyx_n_s__ret;
static PyObject *__pyx_n_s__rg;
static PyObject *__pyx_n_s__self;
static PyObject *__pyx_n_s__shape;
+static PyObject *__pyx_n_s__simplices;
static PyObject *__pyx_n_s__start;
static PyObject *__pyx_n_s__tol;
static PyObject *__pyx_n_s__transpose;
static PyObject *__pyx_n_s__tri;
static PyObject *__pyx_n_s__values;
static PyObject *__pyx_n_s__values_shape;
-static PyObject *__pyx_n_s__vertices;
static PyObject *__pyx_n_s__w;
static PyObject *__pyx_n_s__warn;
static PyObject *__pyx_n_s__warnings;
@@ -3023,7 +3026,7 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
PyArrayObject *__pyx_v_values = 0;
PyArrayObject *__pyx_v_out = 0;
CYTHON_UNUSED PyArrayObject *__pyx_v_points = 0;
- PyArrayObject *__pyx_v_vertices = 0;
+ PyArrayObject *__pyx_v_simplices = 0;
double __pyx_v_c[NPY_MAXDIMS];
double __pyx_v_fill_value;
int __pyx_v_i;
@@ -3041,10 +3044,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__Pyx_Buffer __pyx_pybuffer_out;
__Pyx_LocalBuf_ND __pyx_pybuffernd_points;
__Pyx_Buffer __pyx_pybuffer_points;
+ __Pyx_LocalBuf_ND __pyx_pybuffernd_simplices;
+ __Pyx_Buffer __pyx_pybuffer_simplices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_values;
__Pyx_Buffer __pyx_pybuffer_values;
- __Pyx_LocalBuf_ND __pyx_pybuffernd_vertices;
- __Pyx_Buffer __pyx_pybuffer_vertices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_xi;
__Pyx_Buffer __pyx_pybuffer_xi;
PyObject *__pyx_r = NULL;
@@ -3095,10 +3098,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_pybuffer_points.refcount = 0;
__pyx_pybuffernd_points.data = NULL;
__pyx_pybuffernd_points.rcbuffer = &__pyx_pybuffer_points;
- __pyx_pybuffer_vertices.pybuffer.buf = NULL;
- __pyx_pybuffer_vertices.refcount = 0;
- __pyx_pybuffernd_vertices.data = NULL;
- __pyx_pybuffernd_vertices.rcbuffer = &__pyx_pybuffer_vertices;
+ __pyx_pybuffer_simplices.pybuffer.buf = NULL;
+ __pyx_pybuffer_simplices.refcount = 0;
+ __pyx_pybuffernd_simplices.data = NULL;
+ __pyx_pybuffernd_simplices.rcbuffer = &__pyx_pybuffer_simplices;
__pyx_pybuffer_xi.pybuffer.buf = NULL;
__pyx_pybuffer_xi.refcount = 0;
__pyx_pybuffernd_xi.data = NULL;
@@ -3146,21 +3149,21 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_t_1 = PyObject_GetAttr(__pyx_v_self, __pyx_n_s__tri); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 199; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
- __pyx_t_4 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__vertices); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 199; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ __pyx_t_4 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__simplices); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 199; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
if (!(likely(((__pyx_t_4) == Py_None) || likely(__Pyx_TypeTest(__pyx_t_4, __pyx_ptype_5numpy_ndarray))))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 199; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__pyx_t_5 = ((PyArrayObject *)__pyx_t_4);
{
__Pyx_BufFmt_StackElem __pyx_stack[1];
- if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer, (PyObject*)__pyx_t_5, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
- __pyx_v_vertices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf = NULL;
+ if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer, (PyObject*)__pyx_t_5, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
+ __pyx_v_simplices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf = NULL;
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 199; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
- } else {__pyx_pybuffernd_vertices.diminfo[0].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_vertices.diminfo[0].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_vertices.diminfo[1].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_vertices.diminfo[1].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[1];
+ } else {__pyx_pybuffernd_simplices.diminfo[0].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_simplices.diminfo[0].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_simplices.diminfo[1].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_simplices.diminfo[1].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[1];
}
}
__pyx_t_5 = 0;
- __pyx_v_vertices = ((PyArrayObject *)__pyx_t_4);
+ __pyx_v_simplices = ((PyArrayObject *)__pyx_t_4);
__pyx_t_4 = 0;
@@ -3372,9 +3375,9 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_t_25 = __pyx_v_isimplex;
__pyx_t_26 = __pyx_v_j;
- if (__pyx_t_25 < 0) __pyx_t_25 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_26 < 0) __pyx_t_26 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_v_m = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_25, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_26, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_25 < 0) __pyx_t_25 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_26 < 0) __pyx_t_26 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_v_m = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_25, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_26, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_27 = __pyx_v_m;
@@ -3416,8 +3419,8 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__Pyx_ErrFetch(&__pyx_type, &__pyx_value, &__pyx_tb);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);}
__Pyx_AddTraceback("scipy.interpolate.interpnd.LinearNDInterpolator._evaluate_double", __pyx_clineno, __pyx_lineno, __pyx_filename);
@@ -3426,14 +3429,14 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_L0:;
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__pyx_L2:;
__Pyx_XDECREF((PyObject *)__pyx_v_values);
__Pyx_XDECREF((PyObject *)__pyx_v_out);
__Pyx_XDECREF((PyObject *)__pyx_v_points);
- __Pyx_XDECREF((PyObject *)__pyx_v_vertices);
+ __Pyx_XDECREF((PyObject *)__pyx_v_simplices);
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
@@ -3507,7 +3510,7 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
PyArrayObject *__pyx_v_values = 0;
PyArrayObject *__pyx_v_out = 0;
CYTHON_UNUSED PyArrayObject *__pyx_v_points = 0;
- PyArrayObject *__pyx_v_vertices = 0;
+ PyArrayObject *__pyx_v_simplices = 0;
double __pyx_v_c[NPY_MAXDIMS];
__pyx_t_double_complex __pyx_v_fill_value;
int __pyx_v_i;
@@ -3525,10 +3528,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__Pyx_Buffer __pyx_pybuffer_out;
__Pyx_LocalBuf_ND __pyx_pybuffernd_points;
__Pyx_Buffer __pyx_pybuffer_points;
+ __Pyx_LocalBuf_ND __pyx_pybuffernd_simplices;
+ __Pyx_Buffer __pyx_pybuffer_simplices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_values;
__Pyx_Buffer __pyx_pybuffer_values;
- __Pyx_LocalBuf_ND __pyx_pybuffernd_vertices;
- __Pyx_Buffer __pyx_pybuffer_vertices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_xi;
__Pyx_Buffer __pyx_pybuffer_xi;
PyObject *__pyx_r = NULL;
@@ -3588,10 +3591,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_pybuffer_points.refcount = 0;
__pyx_pybuffernd_points.data = NULL;
__pyx_pybuffernd_points.rcbuffer = &__pyx_pybuffer_points;
- __pyx_pybuffer_vertices.pybuffer.buf = NULL;
- __pyx_pybuffer_vertices.refcount = 0;
- __pyx_pybuffernd_vertices.data = NULL;
- __pyx_pybuffernd_vertices.rcbuffer = &__pyx_pybuffer_vertices;
+ __pyx_pybuffer_simplices.pybuffer.buf = NULL;
+ __pyx_pybuffer_simplices.refcount = 0;
+ __pyx_pybuffernd_simplices.data = NULL;
+ __pyx_pybuffernd_simplices.rcbuffer = &__pyx_pybuffer_simplices;
__pyx_pybuffer_xi.pybuffer.buf = NULL;
__pyx_pybuffer_xi.refcount = 0;
__pyx_pybuffernd_xi.data = NULL;
@@ -3639,21 +3642,21 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_t_1 = PyObject_GetAttr(__pyx_v_self, __pyx_n_s__tri); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 249; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
- __pyx_t_4 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__vertices); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 249; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ __pyx_t_4 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__simplices); if (unlikely(!__pyx_t_4)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 249; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_4);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
if (!(likely(((__pyx_t_4) == Py_None) || likely(__Pyx_TypeTest(__pyx_t_4, __pyx_ptype_5numpy_ndarray))))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 249; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__pyx_t_5 = ((PyArrayObject *)__pyx_t_4);
{
__Pyx_BufFmt_StackElem __pyx_stack[1];
- if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer, (PyObject*)__pyx_t_5, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
- __pyx_v_vertices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf = NULL;
+ if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer, (PyObject*)__pyx_t_5, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
+ __pyx_v_simplices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf = NULL;
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 249; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
- } else {__pyx_pybuffernd_vertices.diminfo[0].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_vertices.diminfo[0].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_vertices.diminfo[1].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_vertices.diminfo[1].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[1];
+ } else {__pyx_pybuffernd_simplices.diminfo[0].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_simplices.diminfo[0].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_simplices.diminfo[1].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_simplices.diminfo[1].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[1];
}
}
__pyx_t_5 = 0;
- __pyx_v_vertices = ((PyArrayObject *)__pyx_t_4);
+ __pyx_v_simplices = ((PyArrayObject *)__pyx_t_4);
__pyx_t_4 = 0;
@@ -3826,18 +3829,20 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_v_k = __pyx_t_19;
+ __pyx_t_15 = __Pyx_CREAL(__pyx_v_fill_value);
__pyx_t_20 = __pyx_v_i;
__pyx_t_21 = __pyx_v_k;
if (__pyx_t_20 < 0) __pyx_t_20 += __pyx_pybuffernd_out.diminfo[0].shape;
if (__pyx_t_21 < 0) __pyx_t_21 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_20, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_21, __pyx_pybuffernd_out.diminfo[1].strides)).real = __Pyx_CREAL(__pyx_v_fill_value);
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_20, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_21, __pyx_pybuffernd_out.diminfo[1].strides)).real = __pyx_t_15;
+ __pyx_t_15 = __Pyx_CIMAG(__pyx_v_fill_value);
__pyx_t_22 = __pyx_v_i;
__pyx_t_23 = __pyx_v_k;
if (__pyx_t_22 < 0) __pyx_t_22 += __pyx_pybuffernd_out.diminfo[0].shape;
if (__pyx_t_23 < 0) __pyx_t_23 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_22, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_23, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __Pyx_CIMAG(__pyx_v_fill_value);
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_22, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_23, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __pyx_t_15;
}
@@ -3879,9 +3884,9 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_t_30 = __pyx_v_isimplex;
__pyx_t_31 = __pyx_v_j;
- if (__pyx_t_30 < 0) __pyx_t_30 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_31 < 0) __pyx_t_31 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_v_m = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_30, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_31, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_30 < 0) __pyx_t_30 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_31 < 0) __pyx_t_31 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_v_m = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_30, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_31, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_32 = __pyx_v_m;
@@ -3934,8 +3939,8 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__Pyx_ErrFetch(&__pyx_type, &__pyx_value, &__pyx_tb);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);}
__Pyx_AddTraceback("scipy.interpolate.interpnd.LinearNDInterpolator._evaluate_complex", __pyx_clineno, __pyx_lineno, __pyx_filename);
@@ -3944,14 +3949,14 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_20LinearNDInterpolator_
__pyx_L0:;
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__pyx_L2:;
__Pyx_XDECREF((PyObject *)__pyx_v_values);
__Pyx_XDECREF((PyObject *)__pyx_v_out);
__Pyx_XDECREF((PyObject *)__pyx_v_points);
- __Pyx_XDECREF((PyObject *)__pyx_v_vertices);
+ __Pyx_XDECREF((PyObject *)__pyx_v_simplices);
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
@@ -4815,22 +4820,22 @@ static double __pyx_f_5scipy_11interpolate_8interpnd__clough_tocher_2d_single_do
int __pyx_t_2;
- __pyx_v_e12x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]));
+ __pyx_v_e12x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]));
- __pyx_v_e12y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]));
+ __pyx_v_e12y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]));
- __pyx_v_e23x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]));
+ __pyx_v_e23x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]));
- __pyx_v_e23y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]));
+ __pyx_v_e23y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]));
- __pyx_v_e31x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]));
+ __pyx_v_e31x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]));
- __pyx_v_e31y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]));
+ __pyx_v_e31y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]));
__pyx_v_e14x = ((__pyx_v_e12x - __pyx_v_e31x) / 3.0);
@@ -4956,10 +4961,10 @@ static double __pyx_f_5scipy_11interpolate_8interpnd__clough_tocher_2d_single_do
__pyx_L5:;
- (__pyx_v_y[0]) = ((((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
+ (__pyx_v_y[0]) = ((((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
- (__pyx_v_y[1]) = ((((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
+ (__pyx_v_y[1]) = ((((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
__pyx_f_5scipy_7spatial_5qhull__barycentric_coordinates(2, (__pyx_v_d->transform + ((__pyx_v_isimplex * 2) * 3)), __pyx_v_y, __pyx_v_c);
@@ -5115,22 +5120,22 @@ static __pyx_t_double_complex __pyx_f_5scipy_11interpolate_8interpnd__clough_toc
int __pyx_t_2;
- __pyx_v_e12x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]));
+ __pyx_v_e12x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]));
- __pyx_v_e12y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]));
+ __pyx_v_e12y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]));
- __pyx_v_e23x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]));
+ __pyx_v_e23x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]));
- __pyx_v_e23y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 1)])))]));
+ __pyx_v_e23y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 1)])))]));
- __pyx_v_e31x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]));
+ __pyx_v_e31x = ((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]));
- __pyx_v_e31y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_isimplex) + 2)])))]));
+ __pyx_v_e31y = ((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 0)])))]) - (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_isimplex) + 2)])))]));
__pyx_v_e14x = ((__pyx_v_e12x - __pyx_v_e31x) / 3.0);
@@ -5256,10 +5261,10 @@ static __pyx_t_double_complex __pyx_f_5scipy_11interpolate_8interpnd__clough_toc
__pyx_L5:;
- (__pyx_v_y[0]) = ((((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
+ (__pyx_v_y[0]) = ((((__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(0 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
- (__pyx_v_y[1]) = ((((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->vertices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
+ (__pyx_v_y[1]) = ((((__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 0)])))]) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 1)])))])) + (__pyx_v_d->points[(1 + (2 * (__pyx_v_d->simplices[((3 * __pyx_v_itri) + 2)])))])) / 3.0);
__pyx_f_5scipy_7spatial_5qhull__barycentric_coordinates(2, (__pyx_v_d->transform + ((__pyx_v_isimplex * 2) * 3)), __pyx_v_y, __pyx_v_c);
@@ -5668,7 +5673,7 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
PyArrayObject *__pyx_v_grad = 0;
PyArrayObject *__pyx_v_out = 0;
CYTHON_UNUSED PyArrayObject *__pyx_v_points = 0;
- PyArrayObject *__pyx_v_vertices = 0;
+ PyArrayObject *__pyx_v_simplices = 0;
double __pyx_v_c[NPY_MAXDIMS];
double __pyx_v_f[(NPY_MAXDIMS + 1)];
double __pyx_v_df[((2 * NPY_MAXDIMS) + 2)];
@@ -5690,10 +5695,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_Buffer __pyx_pybuffer_out;
__Pyx_LocalBuf_ND __pyx_pybuffernd_points;
__Pyx_Buffer __pyx_pybuffer_points;
+ __Pyx_LocalBuf_ND __pyx_pybuffernd_simplices;
+ __Pyx_Buffer __pyx_pybuffer_simplices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_values;
__Pyx_Buffer __pyx_pybuffer_values;
- __Pyx_LocalBuf_ND __pyx_pybuffernd_vertices;
- __Pyx_Buffer __pyx_pybuffer_vertices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_xi;
__Pyx_Buffer __pyx_pybuffer_xi;
PyObject *__pyx_r = NULL;
@@ -5756,10 +5761,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_pybuffer_points.refcount = 0;
__pyx_pybuffernd_points.data = NULL;
__pyx_pybuffernd_points.rcbuffer = &__pyx_pybuffer_points;
- __pyx_pybuffer_vertices.pybuffer.buf = NULL;
- __pyx_pybuffer_vertices.refcount = 0;
- __pyx_pybuffernd_vertices.data = NULL;
- __pyx_pybuffernd_vertices.rcbuffer = &__pyx_pybuffer_vertices;
+ __pyx_pybuffer_simplices.pybuffer.buf = NULL;
+ __pyx_pybuffer_simplices.refcount = 0;
+ __pyx_pybuffernd_simplices.data = NULL;
+ __pyx_pybuffernd_simplices.rcbuffer = &__pyx_pybuffer_simplices;
__pyx_pybuffer_xi.pybuffer.buf = NULL;
__pyx_pybuffer_xi.refcount = 0;
__pyx_pybuffernd_xi.data = NULL;
@@ -5824,21 +5829,21 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_1 = PyObject_GetAttr(__pyx_v_self, __pyx_n_s__tri); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1062; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
- __pyx_t_5 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__vertices); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1062; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ __pyx_t_5 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__simplices); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1062; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_5);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
if (!(likely(((__pyx_t_5) == Py_None) || likely(__Pyx_TypeTest(__pyx_t_5, __pyx_ptype_5numpy_ndarray))))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1062; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__pyx_t_6 = ((PyArrayObject *)__pyx_t_5);
{
__Pyx_BufFmt_StackElem __pyx_stack[1];
- if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer, (PyObject*)__pyx_t_6, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
- __pyx_v_vertices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf = NULL;
+ if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer, (PyObject*)__pyx_t_6, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
+ __pyx_v_simplices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf = NULL;
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 1062; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
- } else {__pyx_pybuffernd_vertices.diminfo[0].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_vertices.diminfo[0].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_vertices.diminfo[1].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_vertices.diminfo[1].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[1];
+ } else {__pyx_pybuffernd_simplices.diminfo[0].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_simplices.diminfo[0].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_simplices.diminfo[1].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_simplices.diminfo[1].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[1];
}
}
__pyx_t_6 = 0;
- __pyx_v_vertices = ((PyArrayObject *)__pyx_t_5);
+ __pyx_v_simplices = ((PyArrayObject *)__pyx_t_5);
__pyx_t_5 = 0;
@@ -6013,9 +6018,9 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_24 = __pyx_v_isimplex;
__pyx_t_25 = __pyx_v_j;
- if (__pyx_t_24 < 0) __pyx_t_24 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_25 < 0) __pyx_t_25 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_26 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_24, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_25, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_24 < 0) __pyx_t_24 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_25 < 0) __pyx_t_25 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_26 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_24, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_25, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_27 = __pyx_v_k;
if (__pyx_t_26 < 0) __pyx_t_26 += __pyx_pybuffernd_values.diminfo[0].shape;
if (__pyx_t_27 < 0) __pyx_t_27 += __pyx_pybuffernd_values.diminfo[1].shape;
@@ -6024,9 +6029,9 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_28 = __pyx_v_isimplex;
__pyx_t_29 = __pyx_v_j;
- if (__pyx_t_28 < 0) __pyx_t_28 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_29 < 0) __pyx_t_29 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_30 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_28, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_29, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_28 < 0) __pyx_t_28 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_29 < 0) __pyx_t_29 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_30 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_28, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_29, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_31 = __pyx_v_k;
__pyx_t_32 = 0;
if (__pyx_t_30 < 0) __pyx_t_30 += __pyx_pybuffernd_grad.diminfo[0].shape;
@@ -6037,9 +6042,9 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_33 = __pyx_v_isimplex;
__pyx_t_34 = __pyx_v_j;
- if (__pyx_t_33 < 0) __pyx_t_33 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_34 < 0) __pyx_t_34 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_35 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_33, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_34, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_33 < 0) __pyx_t_33 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_34 < 0) __pyx_t_34 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_35 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_33, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_34, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_36 = __pyx_v_k;
__pyx_t_37 = 1;
if (__pyx_t_35 < 0) __pyx_t_35 += __pyx_pybuffernd_grad.diminfo[0].shape;
@@ -6087,8 +6092,8 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_grad.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);}
__Pyx_AddTraceback("scipy.interpolate.interpnd.CloughTocher2DInterpolator._evaluate_double", __pyx_clineno, __pyx_lineno, __pyx_filename);
@@ -6098,15 +6103,15 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_grad.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__pyx_L2:;
__Pyx_XDECREF((PyObject *)__pyx_v_values);
__Pyx_XDECREF((PyObject *)__pyx_v_grad);
__Pyx_XDECREF((PyObject *)__pyx_v_out);
__Pyx_XDECREF((PyObject *)__pyx_v_points);
- __Pyx_XDECREF((PyObject *)__pyx_v_vertices);
+ __Pyx_XDECREF((PyObject *)__pyx_v_simplices);
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
@@ -6181,7 +6186,7 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
PyArrayObject *__pyx_v_grad = 0;
PyArrayObject *__pyx_v_out = 0;
CYTHON_UNUSED PyArrayObject *__pyx_v_points = 0;
- PyArrayObject *__pyx_v_vertices = 0;
+ PyArrayObject *__pyx_v_simplices = 0;
double __pyx_v_c[NPY_MAXDIMS];
__pyx_t_double_complex __pyx_v_f[(NPY_MAXDIMS + 1)];
__pyx_t_double_complex __pyx_v_df[((2 * NPY_MAXDIMS) + 2)];
@@ -6203,10 +6208,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_Buffer __pyx_pybuffer_out;
__Pyx_LocalBuf_ND __pyx_pybuffernd_points;
__Pyx_Buffer __pyx_pybuffer_points;
+ __Pyx_LocalBuf_ND __pyx_pybuffernd_simplices;
+ __Pyx_Buffer __pyx_pybuffer_simplices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_values;
__Pyx_Buffer __pyx_pybuffer_values;
- __Pyx_LocalBuf_ND __pyx_pybuffernd_vertices;
- __Pyx_Buffer __pyx_pybuffer_vertices;
__Pyx_LocalBuf_ND __pyx_pybuffernd_xi;
__Pyx_Buffer __pyx_pybuffer_xi;
PyObject *__pyx_r = NULL;
@@ -6245,29 +6250,35 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
int __pyx_t_32;
npy_int __pyx_t_33;
int __pyx_t_34;
- int __pyx_t_35;
+ double __pyx_t_35;
int __pyx_t_36;
- npy_int __pyx_t_37;
- int __pyx_t_38;
- long __pyx_t_39;
- int __pyx_t_40;
- int __pyx_t_41;
- npy_int __pyx_t_42;
+ int __pyx_t_37;
+ npy_int __pyx_t_38;
+ int __pyx_t_39;
+ long __pyx_t_40;
+ double __pyx_t_41;
+ int __pyx_t_42;
int __pyx_t_43;
- long __pyx_t_44;
+ npy_int __pyx_t_44;
int __pyx_t_45;
- int __pyx_t_46;
- npy_int __pyx_t_47;
+ long __pyx_t_46;
+ double __pyx_t_47;
int __pyx_t_48;
- long __pyx_t_49;
- int __pyx_t_50;
+ int __pyx_t_49;
+ npy_int __pyx_t_50;
int __pyx_t_51;
- npy_int __pyx_t_52;
- int __pyx_t_53;
- long __pyx_t_54;
+ long __pyx_t_52;
+ double __pyx_t_53;
+ int __pyx_t_54;
int __pyx_t_55;
- int __pyx_t_56;
+ npy_int __pyx_t_56;
int __pyx_t_57;
+ long __pyx_t_58;
+ double __pyx_t_59;
+ double __pyx_t_60;
+ int __pyx_t_61;
+ int __pyx_t_62;
+ int __pyx_t_63;
int __pyx_lineno = 0;
const char *__pyx_filename = NULL;
int __pyx_clineno = 0;
@@ -6288,10 +6299,10 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_pybuffer_points.refcount = 0;
__pyx_pybuffernd_points.data = NULL;
__pyx_pybuffernd_points.rcbuffer = &__pyx_pybuffer_points;
- __pyx_pybuffer_vertices.pybuffer.buf = NULL;
- __pyx_pybuffer_vertices.refcount = 0;
- __pyx_pybuffernd_vertices.data = NULL;
- __pyx_pybuffernd_vertices.rcbuffer = &__pyx_pybuffer_vertices;
+ __pyx_pybuffer_simplices.pybuffer.buf = NULL;
+ __pyx_pybuffer_simplices.refcount = 0;
+ __pyx_pybuffernd_simplices.data = NULL;
+ __pyx_pybuffernd_simplices.rcbuffer = &__pyx_pybuffer_simplices;
__pyx_pybuffer_xi.pybuffer.buf = NULL;
__pyx_pybuffer_xi.refcount = 0;
__pyx_pybuffernd_xi.data = NULL;
@@ -6356,21 +6367,21 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_1 = PyObject_GetAttr(__pyx_v_self, __pyx_n_s__tri); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1119; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_1);
- __pyx_t_5 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__vertices); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1119; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ __pyx_t_5 = PyObject_GetAttr(__pyx_t_1, __pyx_n_s__simplices); if (unlikely(!__pyx_t_5)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1119; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__Pyx_GOTREF(__pyx_t_5);
__Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
if (!(likely(((__pyx_t_5) == Py_None) || likely(__Pyx_TypeTest(__pyx_t_5, __pyx_ptype_5numpy_ndarray))))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1119; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
__pyx_t_6 = ((PyArrayObject *)__pyx_t_5);
{
__Pyx_BufFmt_StackElem __pyx_stack[1];
- if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer, (PyObject*)__pyx_t_6, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
- __pyx_v_vertices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf = NULL;
+ if (unlikely(__Pyx_GetBufferAndValidate(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer, (PyObject*)__pyx_t_6, &__Pyx_TypeInfo_nn_npy_int, PyBUF_FORMAT| PyBUF_STRIDES, 2, 0, __pyx_stack) == -1)) {
+ __pyx_v_simplices = ((PyArrayObject *)Py_None); __Pyx_INCREF(Py_None); __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf = NULL;
{__pyx_filename = __pyx_f[0]; __pyx_lineno = 1119; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
- } else {__pyx_pybuffernd_vertices.diminfo[0].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_vertices.diminfo[0].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_vertices.diminfo[1].strides = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_vertices.diminfo[1].shape = __pyx_pybuffernd_vertices.rcbuffer->pybuffer.shape[1];
+ } else {__pyx_pybuffernd_simplices.diminfo[0].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[0]; __pyx_pybuffernd_simplices.diminfo[0].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[0]; __pyx_pybuffernd_simplices.diminfo[1].strides = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.strides[1]; __pyx_pybuffernd_simplices.diminfo[1].shape = __pyx_pybuffernd_simplices.rcbuffer->pybuffer.shape[1];
}
}
__pyx_t_6 = 0;
- __pyx_v_vertices = ((PyArrayObject *)__pyx_t_5);
+ __pyx_v_simplices = ((PyArrayObject *)__pyx_t_5);
__pyx_t_5 = 0;
@@ -6519,18 +6530,20 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_v_k = __pyx_t_20;
+ __pyx_t_16 = __Pyx_CREAL(__pyx_v_fill_value);
__pyx_t_21 = __pyx_v_i;
__pyx_t_22 = __pyx_v_k;
if (__pyx_t_21 < 0) __pyx_t_21 += __pyx_pybuffernd_out.diminfo[0].shape;
if (__pyx_t_22 < 0) __pyx_t_22 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_21, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_22, __pyx_pybuffernd_out.diminfo[1].strides)).real = __Pyx_CREAL(__pyx_v_fill_value);
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_21, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_22, __pyx_pybuffernd_out.diminfo[1].strides)).real = __pyx_t_16;
+ __pyx_t_16 = __Pyx_CIMAG(__pyx_v_fill_value);
__pyx_t_23 = __pyx_v_i;
__pyx_t_24 = __pyx_v_k;
if (__pyx_t_23 < 0) __pyx_t_23 += __pyx_pybuffernd_out.diminfo[0].shape;
if (__pyx_t_24 < 0) __pyx_t_24 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_23, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_24, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __Pyx_CIMAG(__pyx_v_fill_value);
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_23, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_24, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __pyx_t_16;
}
@@ -6552,94 +6565,102 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__pyx_t_27 = __pyx_v_isimplex;
__pyx_t_28 = __pyx_v_j;
- if (__pyx_t_27 < 0) __pyx_t_27 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_28 < 0) __pyx_t_28 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_29 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_27, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_28, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_27 < 0) __pyx_t_27 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_28 < 0) __pyx_t_28 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_29 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_27, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_28, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_30 = __pyx_v_k;
if (__pyx_t_29 < 0) __pyx_t_29 += __pyx_pybuffernd_values.diminfo[0].shape;
if (__pyx_t_30 < 0) __pyx_t_30 += __pyx_pybuffernd_values.diminfo[1].shape;
- __Pyx_SET_CREAL((__pyx_v_f[__pyx_v_j]), (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_values.rcbuffer->pybuffer.buf, __pyx_t_29, __pyx_pybuffernd_values.diminfo[0].strides, __pyx_t_30, __pyx_pybuffernd_values.diminfo[1].strides)).real);
+ __pyx_t_16 = (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_values.rcbuffer->pybuffer.buf, __pyx_t_29, __pyx_pybuffernd_values.diminfo[0].strides, __pyx_t_30, __pyx_pybuffernd_values.diminfo[1].strides)).real;
+ __Pyx_SET_CREAL((__pyx_v_f[__pyx_v_j]), __pyx_t_16);
__pyx_t_31 = __pyx_v_isimplex;
__pyx_t_32 = __pyx_v_j;
- if (__pyx_t_31 < 0) __pyx_t_31 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_32 < 0) __pyx_t_32 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_33 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_31, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_32, __pyx_pybuffernd_vertices.diminfo[1].strides));
+ if (__pyx_t_31 < 0) __pyx_t_31 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_32 < 0) __pyx_t_32 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_33 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_31, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_32, __pyx_pybuffernd_simplices.diminfo[1].strides));
__pyx_t_34 = __pyx_v_k;
if (__pyx_t_33 < 0) __pyx_t_33 += __pyx_pybuffernd_values.diminfo[0].shape;
if (__pyx_t_34 < 0) __pyx_t_34 += __pyx_pybuffernd_values.diminfo[1].shape;
- __Pyx_SET_CIMAG((__pyx_v_f[__pyx_v_j]), (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_values.rcbuffer->pybuffer.buf, __pyx_t_33, __pyx_pybuffernd_values.diminfo[0].strides, __pyx_t_34, __pyx_pybuffernd_values.diminfo[1].strides)).imag);
+ __pyx_t_35 = (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_values.rcbuffer->pybuffer.buf, __pyx_t_33, __pyx_pybuffernd_values.diminfo[0].strides, __pyx_t_34, __pyx_pybuffernd_values.diminfo[1].strides)).imag;
+ __Pyx_SET_CIMAG((__pyx_v_f[__pyx_v_j]), __pyx_t_35);
- __pyx_t_35 = __pyx_v_isimplex;
- __pyx_t_36 = __pyx_v_j;
- if (__pyx_t_35 < 0) __pyx_t_35 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_36 < 0) __pyx_t_36 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_37 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_35, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_36, __pyx_pybuffernd_vertices.diminfo[1].strides));
- __pyx_t_38 = __pyx_v_k;
- __pyx_t_39 = 0;
- if (__pyx_t_37 < 0) __pyx_t_37 += __pyx_pybuffernd_grad.diminfo[0].shape;
- if (__pyx_t_38 < 0) __pyx_t_38 += __pyx_pybuffernd_grad.diminfo[1].shape;
- if (__pyx_t_39 < 0) __pyx_t_39 += __pyx_pybuffernd_grad.diminfo[2].shape;
- __Pyx_SET_CREAL((__pyx_v_df[(2 * __pyx_v_j)]), (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_37, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_38, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_39, __pyx_pybuffernd_grad.diminfo[2].strides)).real);
+ __pyx_t_36 = __pyx_v_isimplex;
+ __pyx_t_37 = __pyx_v_j;
+ if (__pyx_t_36 < 0) __pyx_t_36 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_37 < 0) __pyx_t_37 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_38 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_36, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_37, __pyx_pybuffernd_simplices.diminfo[1].strides));
+ __pyx_t_39 = __pyx_v_k;
+ __pyx_t_40 = 0;
+ if (__pyx_t_38 < 0) __pyx_t_38 += __pyx_pybuffernd_grad.diminfo[0].shape;
+ if (__pyx_t_39 < 0) __pyx_t_39 += __pyx_pybuffernd_grad.diminfo[1].shape;
+ if (__pyx_t_40 < 0) __pyx_t_40 += __pyx_pybuffernd_grad.diminfo[2].shape;
+ __pyx_t_41 = (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_38, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_39, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_40, __pyx_pybuffernd_grad.diminfo[2].strides)).real;
+ __Pyx_SET_CREAL((__pyx_v_df[(2 * __pyx_v_j)]), __pyx_t_41);
- __pyx_t_40 = __pyx_v_isimplex;
- __pyx_t_41 = __pyx_v_j;
- if (__pyx_t_40 < 0) __pyx_t_40 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_41 < 0) __pyx_t_41 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_42 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_40, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_41, __pyx_pybuffernd_vertices.diminfo[1].strides));
- __pyx_t_43 = __pyx_v_k;
- __pyx_t_44 = 0;
- if (__pyx_t_42 < 0) __pyx_t_42 += __pyx_pybuffernd_grad.diminfo[0].shape;
- if (__pyx_t_43 < 0) __pyx_t_43 += __pyx_pybuffernd_grad.diminfo[1].shape;
- if (__pyx_t_44 < 0) __pyx_t_44 += __pyx_pybuffernd_grad.diminfo[2].shape;
- __Pyx_SET_CIMAG((__pyx_v_df[(2 * __pyx_v_j)]), (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_42, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_43, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_44, __pyx_pybuffernd_grad.diminfo[2].strides)).imag);
+ __pyx_t_42 = __pyx_v_isimplex;
+ __pyx_t_43 = __pyx_v_j;
+ if (__pyx_t_42 < 0) __pyx_t_42 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_43 < 0) __pyx_t_43 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_44 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_42, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_43, __pyx_pybuffernd_simplices.diminfo[1].strides));
+ __pyx_t_45 = __pyx_v_k;
+ __pyx_t_46 = 0;
+ if (__pyx_t_44 < 0) __pyx_t_44 += __pyx_pybuffernd_grad.diminfo[0].shape;
+ if (__pyx_t_45 < 0) __pyx_t_45 += __pyx_pybuffernd_grad.diminfo[1].shape;
+ if (__pyx_t_46 < 0) __pyx_t_46 += __pyx_pybuffernd_grad.diminfo[2].shape;
+ __pyx_t_47 = (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_44, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_45, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_46, __pyx_pybuffernd_grad.diminfo[2].strides)).imag;
+ __Pyx_SET_CIMAG((__pyx_v_df[(2 * __pyx_v_j)]), __pyx_t_47);
- __pyx_t_45 = __pyx_v_isimplex;
- __pyx_t_46 = __pyx_v_j;
- if (__pyx_t_45 < 0) __pyx_t_45 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_46 < 0) __pyx_t_46 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_47 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_45, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_46, __pyx_pybuffernd_vertices.diminfo[1].strides));
- __pyx_t_48 = __pyx_v_k;
- __pyx_t_49 = 1;
- if (__pyx_t_47 < 0) __pyx_t_47 += __pyx_pybuffernd_grad.diminfo[0].shape;
- if (__pyx_t_48 < 0) __pyx_t_48 += __pyx_pybuffernd_grad.diminfo[1].shape;
- if (__pyx_t_49 < 0) __pyx_t_49 += __pyx_pybuffernd_grad.diminfo[2].shape;
- __Pyx_SET_CREAL((__pyx_v_df[((2 * __pyx_v_j) + 1)]), (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_47, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_48, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_49, __pyx_pybuffernd_grad.diminfo[2].strides)).real);
+ __pyx_t_48 = __pyx_v_isimplex;
+ __pyx_t_49 = __pyx_v_j;
+ if (__pyx_t_48 < 0) __pyx_t_48 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_49 < 0) __pyx_t_49 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_50 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_48, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_49, __pyx_pybuffernd_simplices.diminfo[1].strides));
+ __pyx_t_51 = __pyx_v_k;
+ __pyx_t_52 = 1;
+ if (__pyx_t_50 < 0) __pyx_t_50 += __pyx_pybuffernd_grad.diminfo[0].shape;
+ if (__pyx_t_51 < 0) __pyx_t_51 += __pyx_pybuffernd_grad.diminfo[1].shape;
+ if (__pyx_t_52 < 0) __pyx_t_52 += __pyx_pybuffernd_grad.diminfo[2].shape;
+ __pyx_t_53 = (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_50, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_51, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_52, __pyx_pybuffernd_grad.diminfo[2].strides)).real;
+ __Pyx_SET_CREAL((__pyx_v_df[((2 * __pyx_v_j) + 1)]), __pyx_t_53);
- __pyx_t_50 = __pyx_v_isimplex;
- __pyx_t_51 = __pyx_v_j;
- if (__pyx_t_50 < 0) __pyx_t_50 += __pyx_pybuffernd_vertices.diminfo[0].shape;
- if (__pyx_t_51 < 0) __pyx_t_51 += __pyx_pybuffernd_vertices.diminfo[1].shape;
- __pyx_t_52 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_vertices.rcbuffer->pybuffer.buf, __pyx_t_50, __pyx_pybuffernd_vertices.diminfo[0].strides, __pyx_t_51, __pyx_pybuffernd_vertices.diminfo[1].strides));
- __pyx_t_53 = __pyx_v_k;
- __pyx_t_54 = 1;
- if (__pyx_t_52 < 0) __pyx_t_52 += __pyx_pybuffernd_grad.diminfo[0].shape;
- if (__pyx_t_53 < 0) __pyx_t_53 += __pyx_pybuffernd_grad.diminfo[1].shape;
- if (__pyx_t_54 < 0) __pyx_t_54 += __pyx_pybuffernd_grad.diminfo[2].shape;
- __Pyx_SET_CIMAG((__pyx_v_df[((2 * __pyx_v_j) + 1)]), (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_52, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_53, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_54, __pyx_pybuffernd_grad.diminfo[2].strides)).imag);
+ __pyx_t_54 = __pyx_v_isimplex;
+ __pyx_t_55 = __pyx_v_j;
+ if (__pyx_t_54 < 0) __pyx_t_54 += __pyx_pybuffernd_simplices.diminfo[0].shape;
+ if (__pyx_t_55 < 0) __pyx_t_55 += __pyx_pybuffernd_simplices.diminfo[1].shape;
+ __pyx_t_56 = (*__Pyx_BufPtrStrided2d(npy_int *, __pyx_pybuffernd_simplices.rcbuffer->pybuffer.buf, __pyx_t_54, __pyx_pybuffernd_simplices.diminfo[0].strides, __pyx_t_55, __pyx_pybuffernd_simplices.diminfo[1].strides));
+ __pyx_t_57 = __pyx_v_k;
+ __pyx_t_58 = 1;
+ if (__pyx_t_56 < 0) __pyx_t_56 += __pyx_pybuffernd_grad.diminfo[0].shape;
+ if (__pyx_t_57 < 0) __pyx_t_57 += __pyx_pybuffernd_grad.diminfo[1].shape;
+ if (__pyx_t_58 < 0) __pyx_t_58 += __pyx_pybuffernd_grad.diminfo[2].shape;
+ __pyx_t_59 = (*__Pyx_BufPtrStrided3d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_grad.rcbuffer->pybuffer.buf, __pyx_t_56, __pyx_pybuffernd_grad.diminfo[0].strides, __pyx_t_57, __pyx_pybuffernd_grad.diminfo[1].strides, __pyx_t_58, __pyx_pybuffernd_grad.diminfo[2].strides)).imag;
+ __Pyx_SET_CIMAG((__pyx_v_df[((2 * __pyx_v_j) + 1)]), __pyx_t_59);
}
__pyx_v_w = __pyx_f_5scipy_11interpolate_8interpnd__clough_tocher_2d_single_complex((&__pyx_v_info), __pyx_v_isimplex, __pyx_v_c, __pyx_v_f, __pyx_v_df);
+ __pyx_t_60 = __Pyx_CREAL(__pyx_v_w);
__pyx_t_26 = __pyx_v_i;
- __pyx_t_55 = __pyx_v_k;
+ __pyx_t_61 = __pyx_v_k;
if (__pyx_t_26 < 0) __pyx_t_26 += __pyx_pybuffernd_out.diminfo[0].shape;
- if (__pyx_t_55 < 0) __pyx_t_55 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_26, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_55, __pyx_pybuffernd_out.diminfo[1].strides)).real = __Pyx_CREAL(__pyx_v_w);
+ if (__pyx_t_61 < 0) __pyx_t_61 += __pyx_pybuffernd_out.diminfo[1].shape;
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_26, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_61, __pyx_pybuffernd_out.diminfo[1].strides)).real = __pyx_t_60;
- __pyx_t_56 = __pyx_v_i;
- __pyx_t_57 = __pyx_v_k;
- if (__pyx_t_56 < 0) __pyx_t_56 += __pyx_pybuffernd_out.diminfo[0].shape;
- if (__pyx_t_57 < 0) __pyx_t_57 += __pyx_pybuffernd_out.diminfo[1].shape;
- (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_56, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_57, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __Pyx_CIMAG(__pyx_v_w);
+ __pyx_t_60 = __Pyx_CIMAG(__pyx_v_w);
+ __pyx_t_62 = __pyx_v_i;
+ __pyx_t_63 = __pyx_v_k;
+ if (__pyx_t_62 < 0) __pyx_t_62 += __pyx_pybuffernd_out.diminfo[0].shape;
+ if (__pyx_t_63 < 0) __pyx_t_63 += __pyx_pybuffernd_out.diminfo[1].shape;
+ (*__Pyx_BufPtrStrided2d(__pyx_t_5numpy_complex_t *, __pyx_pybuffernd_out.rcbuffer->pybuffer.buf, __pyx_t_62, __pyx_pybuffernd_out.diminfo[0].strides, __pyx_t_63, __pyx_pybuffernd_out.diminfo[1].strides)).imag = __pyx_t_60;
}
__pyx_L6_continue:;
}
@@ -6670,8 +6691,8 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_grad.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__Pyx_ErrRestore(__pyx_type, __pyx_value, __pyx_tb);}
__Pyx_AddTraceback("scipy.interpolate.interpnd.CloughTocher2DInterpolator._evaluate_complex", __pyx_clineno, __pyx_lineno, __pyx_filename);
@@ -6681,15 +6702,15 @@ static PyObject *__pyx_pf_5scipy_11interpolate_8interpnd_26CloughTocher2DInterpo
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_grad.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_out.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_points.rcbuffer->pybuffer);
+ __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_simplices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_values.rcbuffer->pybuffer);
- __Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_vertices.rcbuffer->pybuffer);
__Pyx_SafeReleaseBuffer(&__pyx_pybuffernd_xi.rcbuffer->pybuffer);
__pyx_L2:;
__Pyx_XDECREF((PyObject *)__pyx_v_values);
__Pyx_XDECREF((PyObject *)__pyx_v_grad);
__Pyx_XDECREF((PyObject *)__pyx_v_out);
__Pyx_XDECREF((PyObject *)__pyx_v_points);
- __Pyx_XDECREF((PyObject *)__pyx_v_vertices);
+ __Pyx_XDECREF((PyObject *)__pyx_v_simplices);
__Pyx_XGIVEREF(__pyx_r);
__Pyx_RefNannyFinishContext();
return __pyx_r;
@@ -6867,8 +6888,10 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyArrayObject *__pyx_v_self, P
__pyx_v_f = NULL;
- __Pyx_INCREF(((PyObject *)__pyx_v_self->descr));
- __pyx_v_descr = __pyx_v_self->descr;
+ __pyx_t_4 = ((PyObject *)__pyx_v_self->descr);
+ __Pyx_INCREF(__pyx_t_4);
+ __pyx_v_descr = ((PyArray_Descr *)__pyx_t_4);
+ __pyx_t_4 = 0;
__pyx_v_hasfields = PyDataType_HASFIELDS(__pyx_v_descr);
@@ -6907,7 +6930,8 @@ static int __pyx_pf_5numpy_7ndarray___getbuffer__(PyArrayObject *__pyx_v_self, P
if (__pyx_t_1) {
- __pyx_v_t = __pyx_v_descr->type_num;
+ __pyx_t_5 = __pyx_v_descr->type_num;
+ __pyx_v_t = __pyx_t_5;
__pyx_t_1 = (__pyx_v_descr->byteorder == '>');
@@ -7999,13 +8023,13 @@ static __Pyx_StringTabEntry __pyx_string_tab[] = {
{&__pyx_n_s__rg, __pyx_k__rg, sizeof(__pyx_k__rg), 0, 0, 1, 1},
{&__pyx_n_s__self, __pyx_k__self, sizeof(__pyx_k__self), 0, 0, 1, 1},
{&__pyx_n_s__shape, __pyx_k__shape, sizeof(__pyx_k__shape), 0, 0, 1, 1},
+ {&__pyx_n_s__simplices, __pyx_k__simplices, sizeof(__pyx_k__simplices), 0, 0, 1, 1},
{&__pyx_n_s__start, __pyx_k__start, sizeof(__pyx_k__start), 0, 0, 1, 1},
{&__pyx_n_s__tol, __pyx_k__tol, sizeof(__pyx_k__tol), 0, 0, 1, 1},
{&__pyx_n_s__transpose, __pyx_k__transpose, sizeof(__pyx_k__transpose), 0, 0, 1, 1},
{&__pyx_n_s__tri, __pyx_k__tri, sizeof(__pyx_k__tri), 0, 0, 1, 1},
{&__pyx_n_s__values, __pyx_k__values, sizeof(__pyx_k__values), 0, 0, 1, 1},
{&__pyx_n_s__values_shape, __pyx_k__values_shape, sizeof(__pyx_k__values_shape), 0, 0, 1, 1},
- {&__pyx_n_s__vertices, __pyx_k__vertices, sizeof(__pyx_k__vertices), 0, 0, 1, 1},
{&__pyx_n_s__w, __pyx_k__w, sizeof(__pyx_k__w), 0, 0, 1, 1},
{&__pyx_n_s__warn, __pyx_k__warn, sizeof(__pyx_k__warn), 0, 0, 1, 1},
{&__pyx_n_s__warnings, __pyx_k__warnings, sizeof(__pyx_k__warnings), 0, 0, 1, 1},
@@ -8328,9 +8352,9 @@ static int __Pyx_InitCachedConstants(void) {
__Pyx_INCREF(((PyObject *)__pyx_n_s__points));
PyTuple_SET_ITEM(__pyx_k_tuple_56, 4, ((PyObject *)__pyx_n_s__points));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__points));
- __Pyx_INCREF(((PyObject *)__pyx_n_s__vertices));
- PyTuple_SET_ITEM(__pyx_k_tuple_56, 5, ((PyObject *)__pyx_n_s__vertices));
- __Pyx_GIVEREF(((PyObject *)__pyx_n_s__vertices));
+ __Pyx_INCREF(((PyObject *)__pyx_n_s__simplices));
+ PyTuple_SET_ITEM(__pyx_k_tuple_56, 5, ((PyObject *)__pyx_n_s__simplices));
+ __Pyx_GIVEREF(((PyObject *)__pyx_n_s__simplices));
__Pyx_INCREF(((PyObject *)__pyx_n_s__c));
PyTuple_SET_ITEM(__pyx_k_tuple_56, 6, ((PyObject *)__pyx_n_s__c));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__c));
@@ -8394,9 +8418,9 @@ static int __Pyx_InitCachedConstants(void) {
__Pyx_INCREF(((PyObject *)__pyx_n_s__points));
PyTuple_SET_ITEM(__pyx_k_tuple_58, 4, ((PyObject *)__pyx_n_s__points));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__points));
- __Pyx_INCREF(((PyObject *)__pyx_n_s__vertices));
- PyTuple_SET_ITEM(__pyx_k_tuple_58, 5, ((PyObject *)__pyx_n_s__vertices));
- __Pyx_GIVEREF(((PyObject *)__pyx_n_s__vertices));
+ __Pyx_INCREF(((PyObject *)__pyx_n_s__simplices));
+ PyTuple_SET_ITEM(__pyx_k_tuple_58, 5, ((PyObject *)__pyx_n_s__simplices));
+ __Pyx_GIVEREF(((PyObject *)__pyx_n_s__simplices));
__Pyx_INCREF(((PyObject *)__pyx_n_s__c));
PyTuple_SET_ITEM(__pyx_k_tuple_58, 6, ((PyObject *)__pyx_n_s__c));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__c));
@@ -8538,9 +8562,9 @@ static int __Pyx_InitCachedConstants(void) {
__Pyx_INCREF(((PyObject *)__pyx_n_s__points));
PyTuple_SET_ITEM(__pyx_k_tuple_66, 5, ((PyObject *)__pyx_n_s__points));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__points));
- __Pyx_INCREF(((PyObject *)__pyx_n_s__vertices));
- PyTuple_SET_ITEM(__pyx_k_tuple_66, 6, ((PyObject *)__pyx_n_s__vertices));
- __Pyx_GIVEREF(((PyObject *)__pyx_n_s__vertices));
+ __Pyx_INCREF(((PyObject *)__pyx_n_s__simplices));
+ PyTuple_SET_ITEM(__pyx_k_tuple_66, 6, ((PyObject *)__pyx_n_s__simplices));
+ __Pyx_GIVEREF(((PyObject *)__pyx_n_s__simplices));
__Pyx_INCREF(((PyObject *)__pyx_n_s__c));
PyTuple_SET_ITEM(__pyx_k_tuple_66, 7, ((PyObject *)__pyx_n_s__c));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__c));
@@ -8616,9 +8640,9 @@ static int __Pyx_InitCachedConstants(void) {
__Pyx_INCREF(((PyObject *)__pyx_n_s__points));
PyTuple_SET_ITEM(__pyx_k_tuple_68, 5, ((PyObject *)__pyx_n_s__points));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__points));
- __Pyx_INCREF(((PyObject *)__pyx_n_s__vertices));
- PyTuple_SET_ITEM(__pyx_k_tuple_68, 6, ((PyObject *)__pyx_n_s__vertices));
- __Pyx_GIVEREF(((PyObject *)__pyx_n_s__vertices));
+ __Pyx_INCREF(((PyObject *)__pyx_n_s__simplices));
+ PyTuple_SET_ITEM(__pyx_k_tuple_68, 6, ((PyObject *)__pyx_n_s__simplices));
+ __Pyx_GIVEREF(((PyObject *)__pyx_n_s__simplices));
__Pyx_INCREF(((PyObject *)__pyx_n_s__c));
PyTuple_SET_ITEM(__pyx_k_tuple_68, 7, ((PyObject *)__pyx_n_s__c));
__Pyx_GIVEREF(((PyObject *)__pyx_n_s__c));
@@ -8743,6 +8767,14 @@ PyMODINIT_FUNC PyInit_interpnd(void)
__pyx_m = PyModule_Create(&__pyx_moduledef);
#endif
if (unlikely(!__pyx_m)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ #if PY_MAJOR_VERSION >= 3
+ {
+ PyObject *modules = PyImport_GetModuleDict(); if (unlikely(!modules)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ if (!PyDict_GetItemString(modules, "scipy.interpolate.interpnd")) {
+ if (unlikely(PyDict_SetItemString(modules, "scipy.interpolate.interpnd", __pyx_m) < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
+ }
+ }
+ #endif
__pyx_b = PyImport_AddModule(__Pyx_NAMESTR(__Pyx_BUILTIN_MODULE_NAME)); if (unlikely(!__pyx_b)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
#if CYTHON_COMPILING_IN_PYPY
Py_INCREF(__pyx_b);
View
60 scipy/interpolate/interpnd.pyx.in
@@ -197,7 +197,7 @@ class LinearNDInterpolator(NDInterpolatorBase):
cdef np.ndarray[np.${DTYPE}_t, ndim=2] values = self.values
cdef np.ndarray[np.${DTYPE}_t, ndim=2] out
cdef np.ndarray[np.double_t, ndim=2] points = self.points
- cdef np.ndarray[np.npy_int, ndim=2] vertices = self.tri.vertices
+ cdef np.ndarray[np.npy_int, ndim=2] simplices = self.tri.simplices
cdef double c[NPY_MAXDIMS]
cdef ${CDTYPE} fill_value
cdef int i, j, k, m, ndim, isimplex, inside, start, nvalues
@@ -248,7 +248,7 @@ class LinearNDInterpolator(NDInterpolatorBase):
for j in xrange(ndim+1):
for k in xrange(nvalues):
- m = vertices[isimplex,j]
+ m = simplices[isimplex,j]
% if DTYPE == "double":
out[i,k] += c[j] * values[m,k]
% else:
@@ -559,20 +559,20 @@ cdef ${CDTYPE} _clough_tocher_2d_single_${DTYPE}(qhull.DelaunayInfo_t *d,
# XXX: optimize + refactor this!
- e12x = (+ d.points[0 + 2*d.vertices[3*isimplex + 1]]
- - d.points[0 + 2*d.vertices[3*isimplex + 0]])
- e12y = (+ d.points[1 + 2*d.vertices[3*isimplex + 1]]
- - d.points[1 + 2*d.vertices[3*isimplex + 0]])
+ e12x = (+ d.points[0 + 2*d.simplices[3*isimplex + 1]]
+ - d.points[0 + 2*d.simplices[3*isimplex + 0]])
+ e12y = (+ d.points[1 + 2*d.simplices[3*isimplex + 1]]
+ - d.points[1 + 2*d.simplices[3*isimplex + 0]])
- e23x = (+ d.points[0 + 2*d.vertices[3*isimplex + 2]]
- - d.points[0 + 2*d.vertices[3*isimplex + 1]])
- e23y = (+ d.points[1 + 2*d.vertices[3*isimplex + 2]]
- - d.points[1 + 2*d.vertices[3*isimplex + 1]])
+ e23x = (+ d.points[0 + 2*d.simplices[3*isimplex + 2]]
+ - d.points[0 + 2*d.simplices[3*isimplex + 1]])
+ e23y = (+ d.points[1 + 2*d.simplices[3*isimplex + 2]]
+ - d.points[1 + 2*d.simplices[3*isimplex + 1]])
- e31x = (+ d.points[0 + 2*d.vertices[3*isimplex + 0]]
- - d.points[0 + 2*d.vertices[3*isimplex + 2]])
- e31y = (+ d.points[1 + 2*d.vertices[3*isimplex + 0]]
- - d.points[1 + 2*d.vertices[3*isimplex + 2]])
+ e31x = (+ d.points[0 + 2*d.simplices[3*isimplex + 0]]
+ - d.points[0 + 2*d.simplices[3*isimplex + 2]])
+ e31y = (+ d.points[1 + 2*d.simplices[3*isimplex + 0]]
+ - d.points[1 + 2*d.simplices[3*isimplex + 2]])
e14x = (e12x - e31x)/3
e14y = (e12y - e31y)/3
@@ -672,13 +672,13 @@ cdef ${CDTYPE} _clough_tocher_2d_single_${DTYPE}(qhull.DelaunayInfo_t *d,
# Centroid of the neighbour, in our local barycentric coordinates
- y[0] = (+ d.points[0 + 2*d.vertices[3*itri + 0]]
- + d.points[0 + 2*d.vertices[3*itri + 1]]
- + d.points[0 + 2*d.vertices[3*itri + 2]]) / 3
+ y[0] = (+ d.points[0 + 2*d.simplices[3*itri + 0]]
+ + d.points[0 + 2*d.simplices[3*itri + 1]]
+ + d.points[0 + 2*d.simplices[3*itri + 2]]) / 3
- y[1] = (+ d.points[1 + 2*d.vertices[3*itri + 0]]
- + d.points[1 + 2*d.vertices[3*itri + 1]]
- + d.points[1 + 2*d.vertices[3*itri + 2]]) / 3
+ y[1] = (+ d.points[1 + 2*d.simplices[3*itri + 0]]
+ + d.points[1 + 2*d.simplices[3*itri + 1]]
+ + d.points[1 + 2*d.simplices[3*itri + 2]]) / 3
qhull._barycentric_coordinates(2, d.transform + isimplex*2*3, y, c)
@@ -809,7 +809,7 @@ class CloughTocher2DInterpolator(NDInterpolatorBase):
cdef np.ndarray[np.${DTYPE}_t, ndim=3] grad = self.grad
cdef np.ndarray[np.${DTYPE}_t, ndim=2] out
cdef np.ndarray[np.double_t, ndim=2] points = self.points
- cdef np.ndarray[np.npy_int, ndim=2] vertices = self.tri.vertices
+ cdef np.ndarray[np.npy_int, ndim=2] simplices = self.tri.simplices
cdef double c[NPY_MAXDIMS]
cdef ${CDTYPE} f[NPY_MAXDIMS+1]
cdef ${CDTYPE} df[2*NPY_MAXDIMS+2]
@@ -855,16 +855,16 @@ class CloughTocher2DInterpolator(NDInterpolatorBase):
for k in xrange(nvalues):
for j in xrange(ndim+1):
% if DTYPE == "double":
- f[j] = values[vertices[isimplex,j],k]
- df[2*j] = grad[vertices[isimplex,j],k,0]
- df[2*j+1] = grad[vertices[isimplex,j],k,1]
+ f[j] = values[simplices[isimplex,j],k]
+ df[2*j] = grad[simplices[isimplex,j],k,0]
+ df[2*j+1] = grad[simplices[isimplex,j],k,1]
% else:
- f[j].real = values[vertices[isimplex,j],k].real
- f[j].imag = values[vertices[isimplex,j],k].imag
- df[2*j].real = grad[vertices[isimplex,j],k,0].real
- df[2*j].imag = grad[vertices[isimplex,j],k,0].imag
- df[2*j+1].real = grad[vertices[isimplex,j],k,1].real
- df[2*j+1].imag = grad[vertices[isimplex,j],k,1].imag
+ f[j].real = values[simplices[isimplex,j],k].real
+ f[j].imag = values[simplices[isimplex,j],k].imag
+ df[2*j].real = grad[simplices[isimplex,j],k,0].real
+ df[2*j].imag = grad[simplices[isimplex,j],k,0].imag
+ df[2*j+1].real = grad[simplices[isimplex,j],k,1].real
+ df[2*j+1].imag = grad[simplices[isimplex,j],k,1].imag
% endif
w = _clough_tocher_2d_single_${DTYPE}(&info, isimplex, c,
View
3 scipy/spatial/SConscript
@@ -48,4 +48,5 @@ src = [join('qhull', 'src', s) for s in [
'random.c', 'rboxlib.c', 'stat.c', 'user.c', 'usermem.c',
'userprintf.c']]
-env.NumpyPythonExtension('qhull', source = ['qhull.c'] + src)
+env.NumpyPythonExtension('qhull', source = ['qhull.c'] + src,
+ CDEFINES={'qh_QHpointer': '1'})
View
55 scipy/spatial/__init__.py
@@ -12,19 +12,68 @@
cKDTree -- class for efficient nearest-neighbor queries (faster impl.)
distance -- module containing many different distance measures
-Delaunay triangulation:
+Delaunay triangulation, convex hulls and Voronoi diagrams:
.. autosummary::
:toctree: generated/
- Delaunay
- tsearch
+ Delaunay -- compute Delaunay triangulation of input points
+ ConvexHull -- compute a convex hull for input points
+ Voronoi -- compute a Voronoi diagram hull from input points
+
+Plotting helpers:
+
+.. autosummary::
+ :toctree: generated/
+
+ delaunay_plot_2d -- plot 2-D triangulation
+ convex_hull_plot_2d -- plot 2-D convex hull
+ voronoi_plot_2d -- plot 2-D voronoi diagram
+
+.. seealso:: :ref:`Tutorial <qhulltutorial>`
+
+
+Simplex representation
+======================
+
+The simplices (triangles, tetrahedra, ...) appearing in the Delaunay
+tesselation (N-dim simplices), convex hull facets, and Voronoi ridges
+(N-1 dim simplices) are represented in the following scheme::
+
+ tess = Delaunay(points)
+ hull = ConvexHull(points)
+ voro = Voronoi(points)
+
+ # coordinates of the j-th vertex of the i-th simplex
+ tess.points[tess.simplices[i, j], :] # tesselation element
+ hull.points[hull.simplices[i, j], :] # convex hull facet
+ voro.vertices[voro.ridge_vertices[i, j], :] # ridge between Voronoi cells
+
+For Delaunay triangulations and convex hulls, the neighborhood
+structure of the simplices satisfies the condition:
+
+ ``tess.neighbors[i,j]`` is the neighboring simplex of the i-th
+ simplex, opposite to the j-vertex. It is -1 in case of no
+ neighbor.
+
+Convex hull facets also define a hyperplane equation:
+
+ (hull.equations[i,:-1] * coord).sum() + hull.equations[i,-1] == 0
+
+Similar hyperplane equations for the Delaunay triangulation correspond
+to the convex hull facets on the corresponding N+1 dimensional
+paraboloid.
+
+The Delaunay triangulation objects offer a method for locating the
+simplex containing a given point, and barycentric coordinate
+computations.
"""
from kdtree import *
from ckdtree import *
from qhull import *
+from _plotutils import *
__all__ = filter(lambda s:not s.startswith('_'),dir())
__all__ += ['distance']
View
160 scipy/spatial/_plotutils.py
@@ -0,0 +1,160 @@
+import numpy as np
+from scipy.lib.decorator import decorator as _decorator
+
+__all__ = ['delaunay_plot_2d', 'convex_hull_plot_2d', 'voronoi_plot_2d']
+
+@_decorator
+def _held_figure(func, obj, ax=None, **kw):
+ import matplotlib.pyplot as plt
+
+ if ax is None:
+ fig = plt.figure()
+ ax = fig.gca()
+
+ was_held = ax.ishold()
+ try:
+ ax.hold(True)
+ return func(obj, ax=ax, **kw)
+ finally:
+ ax.hold(was_held)
+
+def _adjust_bounds(ax, points):
+ ptp_bound = points.ptp(axis=0)
+ ax.set_xlim(points[:,0].min() - 0.1*ptp_bound[0],
+ points[:,0].max() + 0.1*ptp_bound[0])
+ ax.set_ylim(points[:,1].min() - 0.1*ptp_bound[1],
+ points[:,1].max() + 0.1*ptp_bound[1])
+
+@_held_figure
+def delaunay_plot_2d(tri, ax=None):
+ """
+ Plot the given Delaunay triangulation in 2-D
+
+ Parameters
+ ----------
+ tri : scipy.spatial.Delaunay instance
+ Triangulation to plot
+ ax : matplotlib.axes.Axes instance, optional
+ Axes to plot on
+
+ Returns
+ -------
+ fig : matplotlib.figure.Figure instance
+ Figure for the plot
+
+ See Also
+ --------
+ Delaunay
+ matplotlib.pyplot.triplot
+
+ Notes
+ -----
+ Requires Matplotlib.
+
+ """
+ if tri.points.shape[1] != 2:
+ raise ValueError("Delaunay triangulation is not 2-D")
+
+ ax.plot(tri.points[:,0], tri.points[:,1], 'o')
+ ax.triplot(tri.points[:,0], tri.points[:,1], tri.simplices.copy())
+
+ _adjust_bounds(ax, tri.points)
+
+ return ax.figure
+
+@_held_figure
+def convex_hull_plot_2d(hull, ax=None):
+ """
+ Plot the given convex hull diagram in 2-D
+
+ Parameters
+ ----------
+ hull : scipy.spatial.ConvexHull instance
+ Convex hull to plot
+ ax : matplotlib.axes.Axes instance, optional
+ Axes to plot on
+
+ Returns
+ -------
+ fig : matplotlib.figure.Figure instance
+ Figure for the plot
+
+ See Also
+ --------
+ ConvexHull
+
+ Notes
+ -----
+ Requires Matplotlib.
+
+ """
+ if hull.points.shape[1] != 2:
+ raise ValueError("Convex hull is not 2-D")
+
+ ax.plot(hull.points[:,0], hull.points[:,1], 'o')
+ for simplex in hull.simplices:
+ ax.plot(hull.points[simplex,0], hull.points[simplex,1], 'k-')
+
+ _adjust_bounds(ax, hull.points)
+
+ return ax.figure
+
+@_held_figure
+def voronoi_plot_2d(vor, ax=None):
+ """
+ Plot the given Voronoi diagram in 2-D
+
+ Parameters
+ ----------
+ vor : scipy.spatial.Voronoi instance
+ Diagram to plot
+ ax : matplotlib.axes.Axes instance, optional
+ Axes to plot on
+
+ Returns
+ -------
+ fig : matplotlib.figure.Figure instance
+ Figure for the plot
+
+ See Also
+ --------
+ Voronoi
+
+ Notes
+ -----
+ Requires Matplotlib.
+
+ """
+ if vor.points.shape[1] != 2:
+ raise ValueError("Voronoi diagram is not 2-D")
+
+ ax.plot(vor.points[:,0], vor.points[:,1], '.')
+ ax.plot(vor.vertices[:,0], vor.vertices[:,1], 'o')
+
+ for simplex in vor.ridge_vertices:
+ simplex = np.asarray(simplex)
+ if np.all(simplex >= 0):
+ ax.plot(vor.vertices[simplex,0], vor.vertices[simplex,1], 'k-')
+
+ ptp_bound = vor.points.ptp(axis=0)
+
+ center = vor.points.mean(axis=0)
+ for pointidx, simplex in zip(vor.ridge_points, vor.ridge_vertices):
+ simplex = np.asarray(simplex)
+ if np.any(simplex < 0):
+ i = simplex[simplex >= 0][0] # finite end Voronoi vertex
+
+ t = vor.points[pointidx[1]] - vor.points[pointidx[0]] # tangent
+ t /= np.linalg.norm(t)
+ n = np.array([-t[1], t[0]]) # normal
+
+ midpoint = vor.points[pointidx].mean(axis=0)
+ direction = np.sign(np.dot(midpoint - center, n)) * n
+ far_point = vor.vertices[i] + direction * ptp_bound.max()
+
+ ax.plot([vor.vertices[i,0], far_point[0]],
+ [vor.vertices[i,1], far_point[1]], 'k--')
+
+ _adjust_bounds(ax, vor.points)
+
+ return ax.figure
View
3 scipy/spatial/bscript
@@ -2,4 +2,5 @@ from bento.commands import hooks
@hooks.pre_build
def pre_build(context):
- context.tweak_extension("qhull", use="FLAPACK CLIB")
+ context.tweak_extension("qhull", use="FLAPACK CLIB",
+ defines=['qh_QHpointer=1'])
View
29 scipy/spatial/generate_qhull.py
@@ -1,29 +0,0 @@
-#!/usr/bin/env python
-import tempfile
-import subprocess
-import os
-import sys
-import re
-import shutil
-
-tmp_dir = tempfile.mkdtemp()
-try:
- # Run Cython
- dst_fn = os.path.join(tmp_dir, 'qhull.c')
- ret = subprocess.call(['cython', '-o', dst_fn, 'qhull.pyx'])
- if ret != 0:
- sys.exit(ret)
-
- # Strip comments
- f = open(dst_fn, 'r')
- text = f.read()
- f.close()
-
- r = re.compile(r'/\*(.*?)\*/', re.S)
-
- text = r.sub('', text)
- f = open('qhull.c', 'w')
- f.write(text)
- f.close()
-finally:
- shutil.rmtree(tmp_dir)
View
28,919 scipy/spatial/qhull.c
21,015 additions, 7,904 deletions not shown because the diff is too large. Please use a local Git client to view these changes.
View
2 scipy/spatial/qhull.pxd
@@ -18,7 +18,7 @@ ctypedef struct DelaunayInfo_t:
int npoints
int nsimplex
double *points
- int *vertices
+ int *simplices
int *neighbors
double *equations
double *transform
View
1,328 scipy/spatial/qhull.pyx
@@ -16,7 +16,9 @@ cimport numpy as np
cimport cython
cimport qhull
-__all__ = ['Delaunay', 'tsearch']
+from numpy.compat import asbytes
+
+__all__ = ['Delaunay', 'ConvexHull', 'Voronoi', 'tsearch']
#------------------------------------------------------------------------------
# Qhull interface
@@ -31,6 +33,12 @@ cdef extern from "math.h":
double fabs(double x) nogil
double sqrt(double x) nogil
+cdef extern from "setjmp.h" nogil:
+ ctypedef struct jmp_buf:
+ pass
+ int setjmp(jmp_buf STATE) nogil
+ void longjmp(jmp_buf STATE, int VALUE) nogil
+
cdef extern from "qhull/src/qset.h":
ctypedef union setelemT:
void *p
@@ -40,7 +48,10 @@ cdef extern from "qhull/src/qset.h":
int maxsize
setelemT e[1]
-cdef extern from "qhull/src/qhull.h":
+ int qh_setsize(setT *set) nogil
+ void qh_setappend(setT **setp, void *elem) nogil
+
+cdef extern from "qhull/src/libqhull.h":
ctypedef double realT
ctypedef double coordT
ctypedef double pointT
@@ -56,16 +67,24 @@ cdef extern from "qhull/src/qhull.h":
unsigned id
setT *vertices
setT *neighbors
+ setT *ridges
+ setT *coplanarset
flagT simplicial
flagT flipped
flagT upperdelaunay
+ unsigned visitid
ctypedef struct vertexT:
vertexT *next
vertexT *previous
unsigned int id, visitid
pointT *point
- setT *neighbours
+ setT *neighbors
+
+ ctypedef struct ridgeT:
+ setT *vertices
+ facetT *top
+ facetT *bottom
ctypedef struct qhT:
boolT DELAUNAY
@@ -75,11 +94,17 @@ cdef extern from "qhull/src/qhull.h":
boolT NOerrexit
boolT PROJECTdelaunay
boolT ATinfinity
+ boolT UPPERdelaunay
+ boolT hasTriangulation
int normal_size
char *qhull_command
facetT *facet_list
facetT *facet_tail
+ vertexT *vertex_list
+ vertexT *vertex_tail
int num_facets
+ int num_points
+ int center_size
unsigned int facet_id
pointT *first_point
pointT *input_points
@@ -89,9 +114,10 @@ cdef extern from "qhull/src/qhull.h":
realT max_outside
realT MINoutside
realT DISTround
+ jmp_buf errexit
+ setT *other_points
-
- extern qhT qh_qh
+ extern qhT *qh_qh
extern int qh_PRINToff
extern int qh_ALL
@@ -109,15 +135,46 @@ cdef extern from "qhull/src/qhull.h":
void qh_checkpolygon() nogil
void qh_findgood_all() nogil
void qh_appendprint(int format) nogil
+ setT *qh_pointvertex() nogil
realT *qh_readpoints(int* num, int *dim, boolT* ismalloc) nogil
int qh_new_qhull(int dim, int numpoints, realT *points,
boolT ismalloc, char* qhull_cmd, void *outfile,
void *errfile) nogil
int qh_pointid(pointT *point) nogil
+ vertexT *qh_nearvertex(facetT *facet, pointT *point, double *dist) nogil
+ boolT qh_addpoint(pointT *furthest, facetT *facet, boolT checkdist) nogil
+ facetT *qh_findbestfacet(pointT *point, boolT bestoutside,
+ realT *bestdist, boolT *isoutside) nogil
+ void qh_setdelaunay(int dim, int count, pointT *points) nogil
+ void qh_restore_qhull(qhT **oldqh) nogil
+ qhT *qh_save_qhull() nogil
-# Qhull is not threadsafe: needs locking
-_qhull_lock = threading.Lock()
+cdef extern from "qhull/src/io.h":
+ ctypedef enum qh_RIDGE:
+ qh_RIDGEall
+ qh_RIDGEinner
+ qh_RIDGEouter
+
+ ctypedef void printvridgeT(void *fp, vertexT *vertex, vertexT *vertexA,
+ setT *centers, boolT unbounded)
+ int qh_eachvoronoi_all(void *fp, void* printvridge,
+ boolT isUpper, qh_RIDGE innerouter,
+ boolT inorder) nogil
+ void qh_order_vertexneighbors(vertexT *vertex) nogil
+ int qh_compare_facetvisit(void *p1, void *p2) nogil
+
+cdef extern from "qhull/src/geom.h":
+ pointT *qh_facetcenter(setT *vertices) nogil
+
+cdef extern from "qhull/src/poly.h":
+ void qh_check_maxout() nogil
+
+cdef extern from "qhull/src/mem.h":
+ void qh_memfree(void *object, int insize)
+
+from libc.string cimport memcpy
+from libc.stdlib cimport qsort
#------------------------------------------------------------------------------
# LAPACK interface
@@ -133,51 +190,272 @@ cdef extern from "qhull_blas.h":
#------------------------------------------------------------------------------
-# Delaunay triangulation using Qhull
+# Qhull wrapper
#------------------------------------------------------------------------------
-def _construct_delaunay(np.ndarray[np.double_t, ndim=2] points):
- """
- Perform Delaunay triangulation of the given set of points.
+# Qhull is not threadsafe: needs locking
+_qhull_lock = threading.Lock()
- """
+# Qhull has (swappable) global state: keep track which Qhull instance is active
+# and how many instances are alive
+cdef _Qhull _active_qhull = None
+cdef int _qhull_count = 0
+
+class QhullError(RuntimeError):
+ pass
+
+@cython.final
+cdef class _Qhull:
+ cdef qhT *_saved_qh
+ cdef list _point_arrays
+ cdef public bytes options
+ cdef public bytes mode_option
+ cdef public object furthest_site
+
+ cdef int numpoints, ndim, _is_delaunay
+ cdef np.ndarray _ridge_points
+
+ cdef list _ridge_vertices
+ cdef object _ridge_error
+ cdef int _nridges
+
+ cdef np.ndarray _ridge_equations
+
+ @cython.final
+ def __init__(self,
+ bytes mode_option,
+ np.ndarray[np.double_t, ndim=2] points,
+ bytes options=None,
+ bytes required_options=None,
+ furthest_site=False,
+ incremental=False):
+ global _active_qhull, _qhull_count
+ cdef int exitcode
+
+ points = np.ascontiguousarray(points, dtype=np.double)
+
+ self.numpoints = points.shape[0]
+ self.ndim = points.shape[1]
- # Run qhull with the options
- #
- # - d: perform delaunay triangulation
- # - Qbb: scale last coordinate for Delaunay
- # - Qz: reduces Delaunay precision errors for cospherical sites
- # - Qt: output only simplical facets (can produce degenerate 0-area ones)
- #
- cdef char *options = "qhull d Qz Qbb Qt"
- cdef int curlong, totlong
- cdef int dim
- cdef int numpoints
- cdef int exitcode
+ if self.numpoints <= 0:
+ raise ValueError("No points given")
+ if self.ndim < 2:
+ raise ValueError("Need at least 2-D data")
+
+ # Process options
+ option_set = set()
+ if options is not None:
+ option_set.update(options.split())
+ if furthest_site:
+ if b"Qz" in option_set:
+ option_set.remove(b"Qz")
+ option_set.add(b"Qu")
+ if required_options is not None:
+ required_option_set = set(required_options.split())
+ if b"QJ" in option_set and b"Qt" in required_option_set:
+ # safe to remove, QJ always produces simplical output
+ required_option_set.remove(b"Qt")
+ option_set.update(required_option_set)
+
+ if incremental:
+ incremental_bad_ops = set([b'Qbb', b'Qbk', b'QBk', b'QbB', b'Qz'])
+ bad_opts = []
+ for bad_opt in incremental_bad_ops:
+ if bad_opt in options:
+ bad_opts.append(bad_opt)
+ if bad_opts:
+ raise ValueError("Qhull options %r are incompatible "
+ "with incremental mode" % (bad_opts,))
+
+ if b"Qt" in option_set:
+ # Qhull wants this
+ option_set.add(b"Q11")
+
+ # We need to own the copy of the points in incremental mode
+ points = points.copy()
+
+ if mode_option in (b"d