我有一个二维数据和二维四边形网格，描述了细分为补丁的域。数据在每个网格节点处定义。数据中的不连续性存在于补丁边界处，即数据在同一位置处被多重定义。
如何使用 Python 绘制这些数据，并在节点之间进行线性插值并正确表示沿每个面片面的不连续值？
下面是三个示例元素或补丁，每个元素或补丁都有六个节点值。
节点位置和值数据可能存储在[Kx3x2]数组，其中 K 是元素数量。例如，
x = np.array( [
[ [0.0, 1.0], [0.0, 1.0], [0.0, 1.0] ], #element 0
[ [1.0, 2.0], [1.0, 2.0], [1.0, 2.0] ], #element 1
[ [2.0, 3.0], [2.0, 3.0], [2.0, 3.0] ], #element 2
] )
y = np.array( [
[ [0.0, 0.0], [0.5, 0.5], [1.0, 1.0] ], #element 0
[ [0.0, 1.0], [0.5, 1.5], [1.0, 2.0] ], #element 1
[ [1.0, 1.0], [1.5, 1.5], [2.0, 2.0] ], #element 2
] )
z = np.array( [
[ [0.0, 0.5], [0.0, 0.8], [0.0, 1.0] ], #element 0
[ [0.3, 1.0], [0.6, 1.2], [0.8, 1.3] ], #element 1
[ [1.2, 1.5], [1.3, 1.4], [1.5, 1.7] ], #element 2
] )
我考虑过pyplot.imshow()。这无法同时考虑整个域并且仍然表示多值不连续节点。打电话可能有用imshow()分别为每个补丁。但是，我如何在同一轴上绘制每个补丁图像？imshow()对于非矩形补丁也是有问题的，这是我的一般情况。
我考虑过pyplot.pcolormesh()，但它似乎只适用于以细胞为中心的数据。
一种选择是通过对所有元素进行三角测量，然后使用 matplotlib 进行绘图tripcolor() http://matplotlib.org/api/axes_api.html?highlight=tri#matplotlib.axes.Axes.tripcolor我现在发现的功能。两个有用的演示是here http://matplotlib.org/examples/pylab_examples/tripcolor_demo.html and here http://matplotlib.org/examples/pylab_examples/triplot_demo.html?highlight=triangulation.
Auto-triangulation of my global domain can be problematic, but Delaunay triangulation of a single quadrilateral works very well:
I create a global triangulation by appending the triangulation of each element. This means shared nodes are actually duplicated in the position arrays and value arrays. This allows for the discontinuous data at element faces.
Drawing with linear interpolation and discontinuities as desired can be achieved with the tripcolor() function, supplying the node locations, and values for each node.
I was a little concerned how contour plotting might work, since element faces are no longer logically connected. tricontour() still works as expected. (shown here with triangulation overlaid)
用以下代码重现：
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.tri as tri
x = np.array( [
[ [0.0, 1.0], [0.0, 1.0], [0.0, 1.0] ], #element 0
[ [1.0, 2.0], [1.0, 2.0], [1.0, 2.0] ], #element 1
[ [2.0, 3.0], [2.0, 3.0], [2.0, 3.0] ], #element 2
] )
y = np.array( [
[ [0.0, 0.0], [0.5, 0.5], [1.0, 1.0] ], #element 0
[ [0.0, 1.0], [0.5, 1.5], [1.0, 2.0] ], #element 1
[ [1.0, 1.0], [1.5, 1.5], [2.0, 2.0] ], #element 2
] )
z = np.array( [
[ [0.0, 0.5], [0.0, 0.8], [0.0, 1.0] ], #element 0
[ [0.3, 1.0], [0.6, 1.2], [0.8, 1.3] ], #element 1
[ [1.2, 1.5], [1.3, 1.4], [1.5, 1.7] ], #element 2
] )
global_num_pts = z.size
global_x = np.zeros( global_num_pts )
global_y = np.zeros( global_num_pts )
global_z = np.zeros( global_num_pts )
global_triang_list = list()
offset = 0;
num_triangles = 0;
#process triangulation element-by-element
for k in range(z.shape[0]):
points_x = x[k,...].flatten()
points_y = y[k,...].flatten()
z_element = z[k,...].flatten()
num_points_this_element = points_x.size
#auto-generate Delauny triangulation for the element, which should be flawless due to quadrilateral element shape
triang = tri.Triangulation(points_x, points_y)
global_triang_list.append( triang.triangles + offset ) #offseting triangle indices by start index of this element
#store results for this element in global triangulation arrays
global_x[offset:(offset+num_points_this_element)] = points_x
global_y[offset:(offset+num_points_this_element)] = points_y
global_z[offset:(offset+num_points_this_element)] = z_element
num_triangles += triang.triangles.shape[0]
offset += num_points_this_element
#go back and turn all of the triangle indices into one global triangle array
offset = 0
global_triang = np.zeros( (num_triangles, 3) )
for t in global_triang_list:
global_triang[ offset:(offset+t.shape[0] )] = t
offset += t.shape[0]
plt.figure()
plt.gca().set_aspect('equal')
plt.tripcolor(global_x, global_y, global_triang, global_z, shading='gouraud' )
#plt.tricontour(global_x, global_y, global_triang, global_z )
#plt.triplot(global_x, global_y, global_triang, 'go-') #plot just the triangle mesh
plt.xlim((-0.25, 3.25))
plt.ylim((-0.25, 2.25))
plt.show()