scipy.interpolate.griddata在numpy.linspace()生成的x/y对中的问题

Question

我正在玩一个叫做 plot_surface 的函数，这个函数是用来在 matplotlib 中绘制三维表面的。

我有一个简单的脚本，用来绘制一个表面，x、y 和 z 坐标都是基于随机生成的布朗噪声序列：

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import scipy.interpolate
np.random.seed(1000)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

num = 50
x_y_z = [[],[],[]]
for i in range(len(x_y_z)):
    x_y_z[i] = np.random.randn(num).cumsum()

xs, ys, zs = x_y_z
#xs = np.linspace(1, 25, num)
#ys = np.linspace(1, 25, num)
X, Y = np.meshgrid(xs, ys)
print 'X:', X
print 'Y:', Y
Z = scipy.interpolate.griddata((xs, ys), zs, (X, Y))

ax.plot_surface(X,Y,Z)

plt.show()

这段代码生成了下面的图形（这只是实验，并没有特别的意义，我只是想试试 plot_surface 函数）：

在上面的脚本中，我有两行被注释掉了：

#xs = np.linspace(1, 25, num)
#ys = np.linspace(1, 25, num)

这两行代码在 scipy.interpolate.griddata 这部分出错了，我不太明白为什么。基本上，我只是想把布朗噪声的 x 和 y 坐标换成一个简单的线性序列。出现的错误信息是：

QH6154 qhull precision error: initial facet 1 is coplanar with the interior point
ERRONEOUS FACET:
- f1
    - flags: bottom simplicial upperDelaunay flipped
    - normal:    0.7071  -0.7071        0
    - offset:         -0
    - vertices: p50(v2) p49(v1) p0(v0)
    - neighboring facets: f2 f3 f4

While executing:  | qhull d Qz Qbb Qt
Options selected for Qhull 2010.1 2010/01/14:
  run-id 1327050034  delaunay  Qz-infinity-point  Qbbound-last  Qtriangulate
  _pre-merge  _zero-centrum  Pgood  _max-width 24  Error-roundoff 3.5e-14
  _one-merge 2.4e-13  _near-inside 1.2e-12  Visible-distance 6.9e-14
  U-coplanar-distance 6.9e-14  Width-outside 1.4e-13  _wide-facet 4.2e-13

precision problems (corrected unless 'Q0' or an error)
      2 flipped facets

The input to qhull appears to be less than 3 dimensional, or a
computation has overflowed.

Qhull could not construct a clearly convex simplex from points:
- p32(v3):    17    17   9.7
- p50(v2):    13    13    24
- p49(v1):    25    25    22
- p0(v0):     1     1     0

The center point is coplanar with a facet, or a vertex is coplanar
with a neighboring facet.  The maximum round off error for
computing distances is 3.5e-14.  The center point, facets and distances
to the center point are as follows:

center point    13.92    13.92    13.87

facet p50 p49 p0 distance=    0
facet p32 p49 p0 distance=    0
facet p32 p50 p0 distance=    0
facet p32 p50 p49 distance=    0

These points either have a maximum or minimum x-coordinate, or
they maximize the determinant for k coordinates.  Trial points
are first selected from points that maximize a coordinate.

The min and max coordinates for each dimension are:
  0:         1        25  difference=   24
  1:         1        25  difference=   24
  2:         0        24  difference=   24

If the input should be full dimensional, you have several options that
may determine an initial simplex:
  - use 'QJ'  to joggle the input and make it full dimensional
  - use 'QbB' to scale the points to the unit cube
  - use 'QR0' to randomly rotate the input for different maximum points
  - use 'Qs'  to search all points for the initial simplex
  - use 'En'  to specify a maximum roundoff error less than 3.5e-14.
  - trace execution with 'T3' to see the determinant for each point.

If the input is lower dimensional:
  - use 'QJ' to joggle the input and make it full dimensional
  - use 'Qbk:0Bk:0' to delete coordinate k from the input.  You should
    pick the coordinate with the least range.  The hull will have the
    correct topology.
  - determine the flat containing the points, rotate the points
    into a coordinate plane, and delete the other coordinates.
  - add one or more points to make the input full dimensional.
Traceback (most recent call last):
  File "time_dist_unreachable_contour.py", line 21, in <module>
    Z = scipy.interpolate.griddata((xs, ys), zs, (X, Y))
  File "/usr/lib/python2.7/dist-packages/scipy/interpolate/ndgriddata.py", line 182, in griddata
    ip = LinearNDInterpolator(points, values, fill_value=fill_value)
  File "interpnd.pyx", line 192, in interpnd.LinearNDInterpolator.__init__ (scipy/interpolate/interpnd.c:2524)
  File "qhull.pyx", line 917, in scipy.spatial.qhull.Delaunay.__init__ (scipy/spatial/qhull.c:4030)
  File "qhull.pyx", line 170, in scipy.spatial.qhull._construct_delaunay (scipy/spatial/qhull.c:1269)
RuntimeError: Qhull error

有没有人能解释一下，为什么 scipy.interpolate.griddata 能处理噪声序列，但处理线性序列时却出问题？

Answer 1

从一条线上给定的数据，无法进行插值（或者有意义的外推）到二维空间。

如果所有的点都在一条线上，那么在这两个点之间的二维坐标集合的大小是零。