使用开源工具进行 3D 数据可视化:VTK 使用教程
已发表: 2022-03-11在 Toptal 博客上的最新文章中,熟练的数据科学家 Charles Cook 写了关于使用开源工具进行科学计算的文章。 他的教程重点介绍了开源工具以及它们在轻松处理数据和获取结果方面可以发挥的作用。
但是一旦我们解决了所有这些复杂的微分方程,另一个问题就会出现。 我们如何理解和解释来自这些模拟的大量数据? 我们如何可视化潜在的千兆字节数据,例如大型模拟中具有数百万网格点的数据?
在我为硕士论文研究类似问题的过程中,我接触到了可视化工具包(VTK)——一个专门用于数据可视化的强大图形库。
在本教程中,我将快速介绍 VTK 及其管道架构,并继续讨论使用叶轮泵中模拟流体数据的真实 3D 可视化示例。 最后我会列出图书馆的长处,以及我遇到的弱点。
数据可视化和 VTK 管道
开源库 VTK 包含具有许多复杂可视化算法的可靠处理和渲染管道。 然而,它的功能并不止于此,随着时间的推移,图像和网格处理算法也被添加进来。 在我目前与一家牙科研究公司合作的项目中,我正在使用 VTK 在基于 Qt 的类似 CAD 的应用程序中执行基于网格的处理任务。 VTK 案例研究显示了广泛的适用应用。
VTK 的架构围绕着一个强大的管道概念。 这个概念的基本轮廓如下所示:
- 资源处于管道的最开始,并创造“无中生有”。 例如,一个
vtkConeSource
创建一个 3D 锥体,一个vtkSTLReader
读取*.stl
3D 几何文件。 - 过滤器将源或其他过滤器的输出转换为新的东西。 例如,
vtkCutter
使用隐式函数(例如平面)在算法中切割前一个对象的输出。 VTK 附带的所有处理算法都以过滤器的形式实现,并且可以自由链接在一起。 - 映射器将数据转换为图形基元。 例如,它们可用于指定用于着色科学数据的查找表。 它们是指定显示内容的抽象方式。
- Actor表示场景中的一个对象(几何图形加上显示属性)。 此处指定颜色、不透明度、阴影或方向等内容。
- 渲染器和窗口最终以独立于平台的方式将渲染描述到屏幕上。
典型的 VTK 渲染管道从一个或多个源开始,使用各种过滤器将它们处理成多个输出对象,然后使用映射器和演员分别渲染。 这个概念背后的力量是更新机制。 如果过滤器或源的设置发生更改,所有相关的过滤器、映射器、演员和渲染窗口都会自动更新。 另一方面,如果管道下游的对象需要信息以执行其任务,则它可以轻松获取信息。
另外,不需要直接处理像OpenGL这样的渲染系统。 VTK 以独立于平台和(部分)渲染系统的方式封装所有低级任务; 开发人员在更高的水平上工作。
带有转子泵数据集的代码示例
让我们看一个数据可视化示例,该示例使用 2011 年 IEEE 可视化竞赛中旋转叶轮泵中的流体流动数据集。数据本身是计算流体动力学模拟的结果,与 Charles Cook 的文章中描述的非常相似。
特色泵的压缩模拟数据大小超过 30 GB。 它包含多个部分和多个时间步长,因此尺寸很大。 在本指南中,我们将使用其中一个时间步长的转子部分,其压缩大小约为 150 MB。
我选择使用 VTK 的语言是 C++,但还有其他几种语言的映射,如 Tcl/Tk、Java 和 Python。 如果目标只是单个数据集的可视化,则根本不需要编写代码,而是可以使用 Paraview,这是 VTK 大部分功能的图形前端。
数据集以及为什么需要 64 位
我从上面提供的 30 GB 数据集中提取了转子数据集,方法是在 Paraview 中打开一个时间步并将转子部分提取到一个单独的文件中。 它是一个非结构化的网格文件,即由点和 3D 单元组成的 3D 体积,如六面体、四面体等。 每个 3D 点都有相关的值。 有时单元格也有关联的值,但在这种情况下没有。 该培训将专注于点的压力和速度,并尝试在其 3D 环境中可视化这些。
压缩文件大小约为 150 MB,使用 VTK 加载时内存大小约为 280 MB。 然而,通过在 VTK 中处理它,数据集在 VTK 管道中被多次缓存,我们很快就达到了 32 位程序的 2 GB 内存限制。 使用 VTK 时有一些方法可以节省内存,但为了简单起见,我们将在 64 位中编译和运行该示例。
致谢:数据集由德国克劳斯塔尔大学应用力学研究所提供(Dipl. Wirtsch.-Ing. Andreas Lucius)。
目标
我们将使用 VTK 作为工具实现的是下图所示的可视化。 作为 3D 上下文,数据集的轮廓使用部分透明的线框渲染来显示。 然后使用数据集的左侧部分使用表面的简单颜色编码来显示压力。 (我们将跳过此示例中更复杂的体积渲染)。 为了可视化速度场,数据集的右侧部分填充了流线,这些流线按照速度的大小进行了颜色编码。 这种可视化选择在技术上并不理想,但我想让 VTK 代码尽可能简单。 此外,这个例子是可视化挑战的一部分,即流动中的大量湍流是有原因的。
一步步
我将逐步讨论 VTK 代码,展示渲染输出在每个阶段的外观。 完整的源代码可以在培训结束时下载。
让我们首先从 VTK 中包含我们需要的所有内容并打开 main 函数。
#include <vtkActor.h> #include <vtkArrayCalculator.h> #include <vtkCamera.h> #include <vtkClipDataSet.h> #include <vtkCutter.h> #include <vtkDataSetMapper.h> #include <vtkInteractorStyleTrackballCamera.h> #include <vtkLookupTable.h> #include <vtkNew.h> #include <vtkPlane.h> #include <vtkPointData.h> #include <vtkPointSource.h> #include <vtkPolyDataMapper.h> #include <vtkProperty.h> #include <vtkRenderer.h> #include <vtkRenderWindow.h> #include <vtkRenderWindowInteractor.h> #include <vtkRibbonFilter.h> #include <vtkStreamTracer.h> #include <vtkSmartPointer.h> #include <vtkUnstructuredGrid.h> #include <vtkXMLUnstructuredGridReader.h> int main(int argc, char** argv) {
接下来,我们设置渲染器和渲染窗口以显示我们的结果。 我们设置背景颜色和渲染窗口大小。
// Setup the renderer vtkNew<vtkRenderer> renderer; renderer->SetBackground(0.9, 0.9, 0.9); // Setup the render window vtkNew<vtkRenderWindow> renWin; renWin->AddRenderer(renderer.Get()); renWin->SetSize(500, 500);
使用此代码,我们已经可以显示一个静态渲染窗口。 相反,我们选择添加一个vtkRenderWindowInteractor
以交互地旋转、缩放和平移场景。
// Setup the render window interactor vtkNew<vtkRenderWindowInteractor> interact; vtkNew<vtkInteractorStyleTrackballCamera> style; interact->SetRenderWindow(renWin.Get()); interact->SetInteractorStyle(style.Get());
现在我们有一个正在运行的示例,显示一个灰色的空渲染窗口。
接下来,我们使用 VTK 附带的众多读取器之一加载数据集。
// Read the file vtkSmartPointer<vtkXMLUnstructuredGridReader> pumpReader = vtkSmartPointer<vtkXMLUnstructuredGridReader>::New(); pumpReader->SetFileName("rotor.vtu");
简短介绍 VTK 内存管理:VTK 使用一个方便的自动内存管理概念,围绕引用计数。 然而,与大多数其他实现不同的是,引用计数保存在 VTK 对象本身中,而不是智能指针类。 这样做的好处是可以增加引用计数,即使 VTK 对象是作为原始指针传递的。 有两种主要的方法来创建托管的 VTK 对象。 vtkNew<T>
和vtkSmartPointer<T>::New()
,主要区别在于vtkSmartPointer<T>
可以隐式转换为原始指针T*
,并且可以从函数返回。 对于vtkNew<T>
的实例,我们必须调用.Get()
来获取原始指针,并且我们只能通过将其包装到vtkSmartPointer
中来返回它。 在我们的示例中,我们从不从函数返回,并且所有对象都一直存在,因此我们将使用简短的vtkNew
,只有上述例外用于演示目的。
此时,尚未从文件中读取任何内容。 我们或更下游的过滤器必须调用Update()
才能真正发生文件读取。 让 VTK 类自己处理更新通常是最好的方法。 但是,有时我们希望直接访问过滤器的结果,例如获取此数据集中的压力范围。 然后我们需要手动调用Update()
。 (我们不会通过多次调用Update()
来降低性能,因为结果会被缓存。)
// Get the pressure range pumpReader->Update(); double pressureRange[2]; pumpReader->GetOutput()->GetPointData()->GetArray("Pressure")->GetRange(pressureRange);
接下来,我们需要使用vtkClipDataSet
提取数据集的左半部分。 为了实现这一点,我们首先定义一个定义拆分的vtkPlane
。 然后,我们将第一次看到 VTK 管道是如何连接在一起的: successor->SetInputConnection(predecessor->GetOutputPort())
。 每当我们从clipperLeft
请求更新时,此连接现在将确保所有前面的过滤器也是最新的。
// Clip the left part from the input vtkNew<vtkPlane> planeLeft; planeLeft->SetOrigin(0.0, 0.0, 0.0); planeLeft->SetNormal(-1.0, 0.0, 0.0); vtkNew<vtkClipDataSet> clipperLeft; clipperLeft->SetInputConnection(pumpReader->GetOutputPort()); clipperLeft->SetClipFunction(planeLeft.Get());
最后,我们创建我们的第一个演员和映射器来显示左半边的线框渲染。 请注意,映射器连接到其过滤器的方式与过滤器彼此连接的方式完全相同。 大多数时候,渲染器本身正在触发所有actor、映射器和底层过滤器链的更新!
唯一不是不言自明的行可能是leftWireMapper->ScalarVisibilityOff();
- 它禁止通过设置为当前活动数组的压力值对线框进行着色。
// Create the wireframe representation for the left part vtkNew<vtkDataSetMapper> leftWireMapper; leftWireMapper->SetInputConnection(clipperLeft->GetOutputPort()); leftWireMapper->ScalarVisibilityOff(); vtkNew<vtkActor> leftWireActor; leftWireActor->SetMapper(leftWireMapper.Get()); leftWireActor->GetProperty()->SetRepresentationToWireframe(); leftWireActor->GetProperty()->SetColor(0.8, 0.8, 0.8); leftWireActor->GetProperty()->SetLineWidth(0.5); leftWireActor->GetProperty()->SetOpacity(0.8); renderer->AddActor(leftWireActor.Get());
此时,渲染窗口终于显示了一些东西,即左侧部分的线框。

右侧部分的线框渲染以类似的方式创建,通过将(新创建的) vtkClipDataSet
的平面法线切换到相反方向并稍微改变(新创建的)映射器和演员的颜色和不透明度。 请注意,这里我们的 VTK 管道从同一个输入数据集分成两个方向(右和左)。
// Clip the right part from the input vtkNew<vtkPlane> planeRight; planeRight->SetOrigin(0.0, 0.0, 0.0); planeRight->SetNormal(1.0, 0.0, 0.0); vtkNew<vtkClipDataSet> clipperRight; clipperRight->SetInputConnection(pumpReader->GetOutputPort()); clipperRight->SetClipFunction(planeRight.Get()); // Create the wireframe representation for the right part vtkNew<vtkDataSetMapper> rightWireMapper; rightWireMapper->SetInputConnection(clipperRight->GetOutputPort()); rightWireMapper->ScalarVisibilityOff(); vtkNew<vtkActor> rightWireActor; rightWireActor->SetMapper(rightWireMapper.Get()); rightWireActor->GetProperty()->SetRepresentationToWireframe(); rightWireActor->GetProperty()->SetColor(0.2, 0.2, 0.2); rightWireActor->GetProperty()->SetLineWidth(0.5); rightWireActor->GetProperty()->SetOpacity(0.1); renderer->AddActor(rightWireActor.Get());
输出窗口现在按预期显示两个线框部分。
现在我们准备可视化一些有用的数据! 为了将压力可视化添加到左侧部分,我们不需要做太多事情。 我们创建了一个新的映射器并将其连接到clipperLeft
,但这次我们通过压力数组着色。 也是在这里,我们最终利用了我们上面推导出的pressureRange
范围。
// Create the pressure representation for the left part vtkNew<vtkDataSetMapper> pressureColorMapper; pressureColorMapper->SetInputConnection(clipperLeft->GetOutputPort()); pressureColorMapper->SelectColorArray("Pressure"); pressureColorMapper->SetScalarRange(pressureRange); vtkNew<vtkActor> pressureColorActor; pressureColorActor->SetMapper(pressureColorMapper.Get()); pressureColorActor->GetProperty()->SetOpacity(0.5); renderer->AddActor(pressureColorActor.Get());
输出现在看起来如下图所示。 中间的压力非常低,将物料吸入泵中。 然后,这种材料被运送到外面,迅速获得压力。 (当然应该有一个带有实际值的彩色地图图例,但我将其省略以使示例更短。)
现在更棘手的部分开始了。 我们想在右侧绘制速度流线。 流线是通过在源点的矢量场内积分生成的。 向量场已经是“速度”向量数组形式的数据集的一部分。 所以我们只需要生成源点。 vtkPointSource
生成一个随机点球体。 我们将生成 1500 个源点,因为它们中的大多数无论如何都不会位于数据集中,并且会被流跟踪器忽略。
// Create the source points for the streamlines vtkNew<vtkPointSource> pointSource; pointSource->SetCenter(0.0, 0.0, 0.015); pointSource->SetRadius(0.2); pointSource->SetDistributionToUniform(); pointSource->SetNumberOfPoints(1500);
接下来我们创建流跟踪器并设置其输入连接。 “等等,多个连接?”,你可能会说。 是的 - 这是我们遇到的第一个具有多个输入的 VTK 过滤器。 正常输入连接用于向量场,源连接用于种子点。 由于“Velocity”是clipperRight
中的“活动”向量数组,我们不需要在这里明确指定它。 最后,我们指定应该从种子点向两个方向进行积分,并将积分方法设置为 Runge-Kutta-4.5。
vtkNew<vtkStreamTracer> tracer; tracer->SetInputConnection(clipperRight->GetOutputPort()); tracer->SetSourceConnection(pointSource->GetOutputPort()); tracer->SetIntegrationDirectionToBoth(); tracer->SetIntegratorTypeToRungeKutta45();
我们的下一个问题是通过速度大小为流线着色。 由于向量的大小没有数组,我们将简单地将大小计算到一个新的标量数组中。 正如您所猜到的,这个任务也有一个 VTK 过滤器: vtkArrayCalculator
。 它接受一个数据集并将其输出不变,但只添加一个从一个或多个现有数组计算而来的数组。 我们将这个数组计算器配置为获取“Velocity”向量的大小并将其输出为“MagVelocity”。 最后,我们再次手动调用Update()
,以得出新数组的范围。
// Compute the velocity magnitudes and create the ribbons vtkNew<vtkArrayCalculator> magCalc; magCalc->SetInputConnection(tracer->GetOutputPort()); magCalc->AddVectorArrayName("Velocity"); magCalc->SetResultArrayName("MagVelocity"); magCalc->SetFunction("mag(Velocity)"); magCalc->Update(); double magVelocityRange[2]; magCalc->GetOutput()->GetPointData()->GetArray("MagVelocity")->GetRange(magVelocityRange);
vtkStreamTracer
直接输出折线, vtkArrayCalculator
将它们原封不动地传递。 因此,我们可以直接使用新的映射器和参与者显示magCalc
的输出。
相反,在本次培训中,我们选择通过显示功能区来使输出更好一些。 vtkRibbonFilter
生成 2D 单元格以显示其输入的所有折线的色带。
// Create and render the ribbons vtkNew<vtkRibbonFilter> ribbonFilter; ribbonFilter->SetInputConnection(magCalc->GetOutputPort()); ribbonFilter->SetWidth(0.0005); vtkNew<vtkPolyDataMapper> streamlineMapper; streamlineMapper->SetInputConnection(ribbonFilter->GetOutputPort()); streamlineMapper->SelectColorArray("MagVelocity"); streamlineMapper->SetScalarRange(magVelocityRange); vtkNew<vtkActor> streamlineActor; streamlineActor->SetMapper(streamlineMapper.Get()); renderer->AddActor(streamlineActor.Get());
现在仍然缺少的,实际上也是生成中间渲染所需要的,是实际渲染场景和初始化交互器的最后五行。
// Render and show interactive window renWin->Render(); interact->Initialize(); interact->Start(); return 0; }
最后,我们得到了完成的可视化,我将在这里再次展示:
可以在此处找到上述可视化的完整源代码。
黄金三镖客
我将列出我个人对 VTK 框架的优缺点列表来结束这篇文章。
优点:积极开发:VTK 正在由几个贡献者积极开发,主要来自研究社区。 这意味着一些尖端算法可用,许多 3D 格式可以导入和导出,错误被积极修复,问题通常在讨论板上有现成的解决方案。
缺点:可靠性:然而,将来自不同贡献者的许多算法与 VTK 的开放式管道设计相结合,可能会导致不寻常的过滤器组合出现问题。 为了弄清楚为什么我的复杂过滤器链没有产生预期的结果,我不得不多次进入 VTK 源代码。 我强烈建议以允许调试的方式设置 VTK。
优点:软件架构:VTK 的管道设计和通用架构似乎经过深思熟虑,并且很高兴与之合作。 几行代码可以产生惊人的结果。 内置的数据结构易于理解和使用。
缺点:微架构:一些微架构设计决策超出了我的理解。 常量正确性几乎不存在,数组作为输入和输出传递,没有明显区别。 我通过放弃一些性能并使用我自己的
vtkMath
包装器为我自己的算法减轻了这一点,它利用了自定义 3D 类型,如typedef std::array<double, 3> Pnt3d;
.Pro : Micro Documentation : 所有类和过滤器的 Doxygen 文档内容广泛且可用,wiki 上的示例和测试用例对理解过滤器的使用方式也有很大帮助。
缺点:宏文档:网络上有几个很好的 VTK 教程和介绍。 但是据我所知,没有大的参考文档来解释具体的事情是如何完成的。 如果你想做一些新的事情,期望在一段时间内搜索如何去做。 此外,很难找到任务的特定过滤器。 但是,一旦找到它,Doxygen 文档通常就足够了。 探索 VTK 框架的一个好方法是下载并试用 Paraview。
优点:隐式并行化支持:如果您的源可以分成几个可以独立处理的部分,那么并行化就像在处理单个部分的每个线程中创建一个单独的过滤器链一样简单。 大多数大型可视化问题通常属于这一类。
缺点:没有显式并行化支持:如果您没有幸免于大的、可分割的问题,但您想利用多个内核,那么您只能靠自己了。 您必须弄清楚哪些类是线程安全的,甚至可以通过反复试验或阅读源代码来重新进入。 我曾经追踪到一个 VTK 过滤器的并行化问题,该过滤器使用静态全局变量来调用某个 C 库。
Pro : Buildsystem CMake : 多平台元构建系统 CMake 也是由 Kitware(VTK 的制造商)开发的,并在 Kitware 之外的许多项目中使用。 它与 VTK 很好地集成在一起,使得为多个平台设置构建系统变得不那么痛苦。
Pro :平台独立性、许可和长寿: VTK 是开箱即用的平台独立性,并在非常宽松的 BSD 风格许可下获得许可。 此外,对于需要专业支持的重要项目,我们还提供专业支持。 Kitware 得到了许多研究机构和其他公司的支持,并将持续一段时间。
遗言
总的来说,VTK 是解决我喜欢的各种问题的最佳数据可视化工具。 如果您遇到需要可视化、网格处理、图像处理或类似任务的项目,请尝试使用输入示例启动 Paraview,并评估 VTK 是否适合您。