tl;dr: 我认为您不应该为此尝试使用 matplotlib,这会很困难并且效果不佳。我建议使用专用的 vtk 库,或者是裸 vtk、更高级别的 mayavi.mlab 或我最近获得的最喜欢的 pyvista。我会详细说明这一切。
数据
首先,这是您输入数据的一个小型、独立版本(因为您在问题中链接的数据既太大,而且迟早可能会成为断开的链接),采用旧版 VTK 文件格式。我已将您的数据减少到三个不同大小的rectangular cuboids 以接近您的数字。
# vtk DataFile Version 3.1
MCVE VTK file
ASCII
DATASET UNSTRUCTURED_GRID
POINTS 16 float
0. 0. 0.
0. 0. 3.
0. 2. 0.
0. 2. 3.
4. 0. 0.
4. 0. 3.
4. 2. 0.
4. 2. 3.
5. 0. 0.
5. 0. 3.
5. 2. 0.
5. 2. 3.
13. 0. 0.
13. 0. 3.
13. 2. 0.
13. 2. 3.
CELLS 3 27
8 0 1 3 2 4 5 7 6
8 4 5 7 6 8 9 11 10
8 8 9 11 10 12 13 15 14
CELL_TYPES 3
12 12 12
CELL_DATA 3
SCALARS elem_val float
LOOKUP_TABLE default
1
2
3
让我们讨论一下这个文件代表什么。标头指定它是一个非结构化网格。这意味着它可以包含以任意方式排列的点。基本上是一袋积分。你可以找到一些关于文件格式here的解释。
第一个块POINTS包含16 float行,每行对应一个点在3d中的坐标,共16个点。
第二个块CELLS 定义了 3 行,每行对应一个单元格(较小的单位,在本例中为体积),该单元格根据点的从 0 开始的索引定义。第一个数字(8)表示给定单元格中的顶点数,后面的数字是相应顶点的点索引。上述示例数据文件中的所有三个单元都包含 8 个顶点,因为我们要绘制的每个长方体都有 8 个顶点。 CELLS 行的第二个数字是这个块中数字的总数,即3 * (8+1),即 27。
第三个块CELL_TYPES 定义了每个3 单元格的类型。在这种情况下,它们都是12 类型,对应于“六面体”。借用from Fig. 2 and 3 of the already linked examples的资料图:
这些列出了主要的细胞类型及其各自的索引。
最后一个块SCALARS 包含每个单元格的标量(数字),稍后将根据它着色。标量 1 到 3 将映射到颜色图上,为您提供图中所示的红到蓝过渡。
为什么不用 matplotlib?
我不熟悉meshio,但我怀疑它可以让您访问VTK 文件中的上述块。您显示的 mesh.cells 属性表明它识别出每个单元格都是“六面体”,并列出了每个单元格及其各自的 8 个顶点索引。 mesh.points 属性可能是一个形状为 (n,3) 的数组,在这种情况下,mesh.points[cell_inds, :] 为您提供由其 8 长度数组 cell_inds 定义的给定单元格的 (8,3) 形状的坐标。
你会如何用 matplotlib 可视化这个?首先,您的实际数据非常庞大,它包含 84480 个单元格,即使从远处看,它们看起来很像我上面的示例数据。所以你必须
- 想出一种方法将所有这些单元坐标转换为要使用 matplotlib 绘制的表面,这并不容易,
- 然后意识到 80k 表面将导致 matplotlib 中巨大的内存和 CPU 开销,最后
- 请注意,matplotlib 具有 2d 渲染器,因此复杂(阅读:不相交、互锁)表面的 3d 可视化often goes awry。
考虑到所有这些因素,我绝对不会尝试为此使用 matplotlib。
然后呢?
使用 ParaView 在后台使用的东西:VTK!您仍然可以通过低级 vtk 模块或高级(er)级 mayavi.mlab 模块以编程方式使用机器。还有 mayavi 相关的 tvtk 模块,它是一种中间立场(出于这些目的,它仍然是低级别的 VTK,但具有对 python 更友好的 API),但我将把它留作练习给读者。
1。 vtk
使用 vtk 读取和绘制非结构化网格有点复杂(就像使用裸 vtk 一样,因为您必须自己组装管道),但使用 this archived wiki page 加上 these corrections that have changed since 可以管理(这实际上使用了遗留阅读器,并且代码可以针对更现代的输入格式进行相当多的现代化改造,请参阅this newer example):
from vtk import (vtkUnstructuredGridReader, vtkDataSetMapper, vtkActor,
vtkRenderer, vtkRenderWindow, vtkRenderWindowInteractor)
file_name = "mesh_mcve.vtk" # minimal example vtk file
# Read the source file.
reader = vtkUnstructuredGridReader()
reader.SetFileName(file_name)
reader.Update() # Needed because of GetScalarRange
output = reader.GetOutput()
output_port = reader.GetOutputPort()
scalar_range = output.GetScalarRange()
# Create the mapper that corresponds the objects of the vtk file
# into graphics elements
mapper = vtkDataSetMapper()
mapper.SetInputConnection(output_port)
mapper.SetScalarRange(scalar_range)
# Create the Actor
actor = vtkActor()
actor.SetMapper(mapper)
# Create the Renderer
renderer = vtkRenderer()
renderer.AddActor(actor)
renderer.SetBackground(1, 1, 1) # Set background to white
# Create the RendererWindow
renderer_window = vtkRenderWindow()
renderer_window.AddRenderer(renderer)
# Create the RendererWindowInteractor and display the vtk_file
interactor = vtkRenderWindowInteractor()
interactor.SetRenderWindow(renderer_window)
interactor.Initialize()
interactor.Start()
请注意,我只对原始 wiki 版本进行了最低限度的更改。可以在here 找到使用旧版阅读器的更复杂示例。这是视口旋转后的输出:
实际颜色取决于默认颜色图和标量的缩放比例。 vtk 模块的上述默认值似乎默认使用jet 颜色映射,并且它对标量进行规范化,以便将值映射到完整的颜色范围。
2。 mayavi.mlab
就个人而言,我发现vtk 使用起来非常痛苦。它涉及大量搜索,而且往往是在库中定义的子模块和类的迷宫中挖掘。这就是为什么我总是尝试通过mayavi.mlab 的更高级别功能来使用vtk。当您不使用 VTK 文件时(例如尝试可视化在 numpy 数组中定义的数据时),此模块特别有用,但在这种情况下它恰好为我们节省了大量工作,因为好吧,同时提供额外的功能。这是使用mlab 的相同可视化:
from mayavi import mlab
from mayavi.modules.surface import Surface
file_name = "mesh_mcve.vtk" # minimal example vtk file
# create a new figure, grab the engine that's created with it
fig = mlab.figure()
engine = mlab.get_engine()
# open the vtk file, let mayavi figure it all out
vtk_file_reader = engine.open(file_name)
# plot surface corresponding to the data
surface = Surface()
engine.add_filter(surface, vtk_file_reader)
# block until figure is closed
mlab.show()
工作量少很多!我们将整个 VTK 解析怪物推送到 mayavi,以及映射器、演员和渲染器的混乱......
它的外观如下:
以上是最小化、最省力的可视化,但您当然可以从这里开始更改任何您喜欢的内容以使其符合您的需求。你可以改变背景,改变颜色图,以奇怪的方式操作数据,你可以命名它。请注意,这里的颜色与vtk 的情况相比是相反的,因为默认颜色图或标量到颜色图(查找表)的映射是不同的。你越是偏离 mlab 的高级 API,它就越脏(因为你越来越接近引擎盖下的裸 VTK),但你通常仍然可以使用@987654385 节省大量工作和混淆代码@。
最后,mayavi 的图形窗口支持各种宝石:管道和场景的交互式修改、坐标轴等注释、切换正交投影,甚至能够在自动生成中记录您以交互方式更改的任何内容蟒蛇脚本。我肯定会建议尝试使用 mayavi 实现您想要做的事情。如果你知道你会用 ParaView 做什么,那么通过使用它的交互式会话记录功能将它移植到 mayavi 是相当容易的。
3。 pyvista
最近有人向我指出pyvista,这是一个建立在vtk 之上的令人愉快的多功能且强大的库。尽管它的 API 需要一些时间来适应,但有很多示例和详尽的 API 参考 in the documentation。开始使用该库有一些学习曲线,但您可能会发现使用它的高级界面和“按我的意思做”机制会更有效率。我特别欣赏它的开源特性,包括公共问题跟踪器和响应迅速的核心开发人员。
那么我们如何使用pyvista 读取和绘制网格?如下:
import pyvista as pv
# read the data
grid = pv.read('mesh_mcve.vtk')
# plot the data with an automatically created Plotter
grid.plot(show_scalar_bar=False, show_axes=False)
这会产生
正如你所看到的颜色非常不同:这是因为pyvista 使用了matplotlib 的感知统一viridis colormap,这对于数据可视化非常有用!如果您坚持使用更奇怪的颜色,您可以将cmap='jet' 传递给grid.plot 的调用。这里有很多要说的默认照明和阴影有何不同,但我建议仔细阅读有关过滤和绘制数据集的所有选项和方法的文档。