事实上，plot_trisurf似乎应该是完美的这项任务！此外，你可以使用tri.UniformTriRefiner来得到一个带有更小三角形的Triangulation：

import numpy
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib import tri, cm

x = numpy.array([0, 0, 1, 1])
y = numpy.array([0.5, 0.75, 1, 0.5])
z = numpy.array([0, 0.5, 1, 0])

triang = tri.Triangulation(x, y)
refiner = tri.UniformTriRefiner(triang)
new, new_z = refiner.refine_field(z, subdiv=4)

norm = plt.Normalize(vmax=abs(y).max(), vmin=-abs(y).max())
kwargs = dict(triangles=new.triangles, cmap=cm.jet, norm=norm, linewidth=0.2)

fig = plt.figure()
ax = Axes3D(fig)
pt = ax.plot_trisurf(new.x, new.y, new_z, **kwargs)
plt.show()

产生以下图像：

三角网格细化最近才添加到matplotlib中，因此您需要版本1。3使用它。如果你被版本1卡住了。2如果注解掉import matplotlib.tri.triinterpolate行和所有的refine_field方法，你也应该能够直接使用source from Github。然后，您需要使用refine_triangulation方法并使用griddata来插值新的相应Z值。

**编辑：**上面的代码使用三次插值来确定新三角形的Z值，但对于线性插值，您可以替换/添加以下行：

interpolator = tri.LinearTriInterpolator(triang, z)
new, new_z = refiner.refine_field(z, interpolator, subdiv=4)

或者，使用scipy.interpolate.griddata进行插值：

from scipy.interpolate import griddata

new = refiner.refine_triangulation(subdiv = 4)
new_z = griddata((x,y),z, (new.x, new.y), method='linear')

这个问题归结为如何在matplotlib中进行曲面的插值着色，即。也就是说，相当于Matlab的shading('interp')特性。简短的回答是：这是不可能的。它本身并不受支持，所以最好的办法就是手工完成，这也是目前为止所提出的解决方案的目标。
几年前我也走上了这条路，当时我也对Matlab的shading('interp')感到沮丧：它的工作原理是简单地在每个四边形上插入4个角颜色，这意味着颜色梯度的方向在相邻四边形上可以不同。我想要的是，每个颜色带将正好在z轴上的两个明确定义的值之间，相邻单元格之间没有视觉中断。
做三角测量绝对是个好主意。但是，我的方法不是简单地细化网格并希望达到相邻三角形的颜色在视觉上无法区分的点（而不是达到首先出现伪影的点），而是计算三角测量上的轮廓带，然后将它们绘制成3D。
当我第一次实现它的时候，matplotlib不支持三角测量轮廓。现在它通过_tri.TriContourGenerator实现。如果这也提供了提取的多边形顶点的z值，我们就完成了。不幸的是，在Python级别上无法访问它们，因此我们需要尝试通过比较create_filled_contours()和create_contours()的输出来重建它们，这在以下代码中完成：

import numpy as np
from mpl_toolkits.mplot3d import Axes3D
from mpl_toolkits.mplot3d.art3d import Poly3DCollection
import matplotlib.pyplot as plt
from matplotlib import _tri, tri, cm

def contour_bands_3d(x, y, z, nbands=20):
    # obtain the contouring engine on a triangulation
    TRI = tri.Triangulation(x, y)
    C = _tri.TriContourGenerator(TRI.get_cpp_triangulation(), z)

    # define the band breaks
    brks = np.linspace(z.min(), z.max(), nbands+1)

    # the contour lines
    lines = [C.create_contour(b) for b in brks]

    # the contour bands
    bands = [C.create_filled_contour(brks[i], brks[i+1]) for i in xrange(nbands)]

    # compare the x, y vertices of each band with the x, y vertices of the upper
    # contour line; if matching, z = z1, otherwise z = z0 (see text for caveats)
    eps = 1e-6
    verts = []
    for i in xrange(nbands):
        b = bands[i][0]
        l = lines[i+1][0]
        z0, z1 = brks[i:i+2]
        zi = np.array([z1 if (np.abs(bb - l) < eps).all(1).any() else z0 for bb in b])
        verts.append(np.c_[b, zi[:,None]])
    return brks, verts

x = np.array([0, 0, 1, 1])
y = np.array([0.5, 0.75, 1, 0.5])
z = np.array([0, 0.5, 1,0])

fig = plt.figure()
ax = Axes3D(fig)
verts = [zip(x, y, z)]
PC = Poly3DCollection(verts)
ax.add_collection3d(PC)

# calculate the 3d contour bands
brks, verts = contour_bands_3d(x, -y, z)

cmap = cm.get_cmap('jet')
norm = plt.Normalize(vmax=abs(y).max(), vmin=-abs(y).max())

PC = Poly3DCollection(verts, cmap=cmap, norm=norm, edgecolors='none')
PC.set_array(brks[:-1])
ax.add_collection(PC)
ax.set_ylim((-1, 1))
plt.show()

这就是结果：

请注意，z值的重建并不完全正确，因为我们还需要检查x，y顶点是否实际上是原始数据集的一部分，在这种情况下，必须采用其原始z值。然而，修改轮廓绘制算法的C++代码以跟踪z值将容易得多。这将是一个很小的变化，而试图在Python中覆盖所有情况无异于一场噩梦。
至于效率，我们正试图在Python级别上完成图形卡的工作，所以这将是可怕的。但这对所有mplot3d都是一样的。如果需要一个性能实现，我推荐VTK中的BandedContourFilter()。它的工作速度非常快，也可以在Python中使用。