在他最近在Toptal博客上的一篇文章中， 熟练的数据科学家 查尔斯·库克写道 使用开源工具进行科学计算. 他的教程对开源工具及其在轻松处理数据和获取结果方面可以发挥的作用提出了重要的观点.

但是一旦我们解完这些复杂的微分方程另一个问题就来了. 我们如何理解和解释这些模拟产生的大量数据? 我们如何可视化潜在的千兆字节数据, 如数据与数百万个网格点内的大型模拟?

在我的工作中遇到了类似的问题 Master’s Thesis, 我接触到了可视化工具包, 或VTK -一个强大的图形库，专门用于数据可视化.

在本教程中，我将快速介绍VTK及其管道架构, 并继续讨论现实生活中的3D可视化示例，使用来自叶轮泵中的模拟流体的数据. 最后，我将列出这个库的优点，以及我遇到的缺点.

数据可视化和VTK管道

开源库 VTK 包含一个坚实的处理和渲染管道与许多复杂的可视化算法. It’s capabilities, however, don’t stop there, 随着时间的推移，图像和网格处理算法也被添加. 在我目前与一家牙科研究公司的项目中, 我利用VTK基于网格的处理任务在一个 Qt基于cad的应用程序. The VTK case studies 显示广泛的适用范围.

VTK的架构围绕着一个强大的管道概念. 这个概念的基本轮廓如下:

• Sources 都在管道的最开始，创造“无中生有”的东西. For example, a vtkConeSource 创建一个三维锥体 vtkSTLReader reads *.stl 3D geometry files.
• Filters 将源或其他过滤器的输出转换为新的内容. For example a vtkCutter 使用隐式函数e切割算法中前一个对象的输出.g., a plane. VTK附带的所有处理算法都作为滤波器实现，可以自由地链接在一起.
• Mappers 将数据转换为图形原语. 例如，它们可用于指定用于着色的查找表 scientific data. 它们是指定显示内容的抽象方式.
• Actors 表示场景中的一个对象(几何图形加上显示属性). 颜色、不透明度、阴影或方向等都在这里指定.
• Renderers & Windows 最后以一种与平台无关的方式描述在屏幕上的渲染.

典型的VTK渲染管道从一个或多个源开始, 使用各种过滤器将它们处理成几个输出对象, 然后分别使用映射器和演员渲染. 这个概念背后的力量是更新机制. 如果更改了过滤器或源的设置, 所有相关过滤器, mappers, 演员和渲染窗口自动更新. If, on the other hand, 管道下的对象需要信息来执行其任务, 它很容易得到它.

此外，不需要直接处理像OpenGL这样的渲染系统. VTK encapsulates all the low level task in a platform- and (partially) rendering system-independent way; the developer works on a much higher level.

代码示例与转子泵数据集

让我们看一个数据可视化的例子，使用来自2011年IEEE可视化竞赛的旋转叶轮泵流体流动数据集. 数据本身是计算流体动力学模拟的结果, 就像查尔斯·库克的文章中描述的那样.

特色泵的压缩仿真数据大小超过30gb. 它包含多个部分和多个时间步长，因此尺寸较大. In this guide, 我们将使用其中一个时间步长的转子部分, 它的压缩大小约为150 MB.

我使用VTK的语言选择是c++, 但是也有一些其他语言的映射，比如Tcl/Tk, Java, and Python. 如果目标只是单个数据集的可视化, 人们根本不需要编写代码，而是可以利用 Paraview，一个图形前端的大部分VTK的功能.

数据集和为什么64位是必要的

我从上面提供的30 GB数据集中提取了转子数据集, 通过在Paraview中打开一个时间步骤并将转子部分提取到单独的文件中. 它是一个非结构化网格文件，i.e.由点和三维单元组成的三维体，如六面体、四面体等. 每个3D点都有相关的值. 有时单元格也有相关的值，但在本例中没有. 这个训练将集中在压力和速度在点，并试图可视化这些在他们的3D环境.

压缩文件大小约为150mb，内存大小约为280mb，加载VTK时. However, 通过在VTK中处理, 数据集在VTK管道内缓存多次，我们很快达到32位程序的2 GB内存限制. 使用VTK时，有很多方法可以节省内存, 但为了简单起见，我们只编译并运行64位的例子.

Acknowledgements:数据集由应用力学研究所提供, Clausthal大学, Germany (Dipl. Wirtsch.-Ing. Andreas Lucius).

The Target

使用VTK作为工具，我们将实现如下图所示的可视化效果. 作为3D上下文，数据集的轮廓使用部分透明的线框渲染来显示. 然后，数据集的左侧部分使用简单的表面颜色编码来显示压力. (在这个例子中，我们将跳过更复杂的体渲染). 为了使速度场可视化，数据集的右侧填充了 streamlines，它们以速度的大小用颜色编码. 这种可视化选择在技术上并不理想, 但我想保持VTK代码尽可能简单. 此外，这个示例作为可视化挑战的一部分还有一个原因.e.气流中有很多湍流.

Step by Step

我将一步一步地讨论VTK代码，展示每个阶段的渲染输出. 完整的源代码可以在培训结束时下载.

让我们首先从VTK中包含我们需要的所有内容并打开main函数.

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

Int main(Int argc, char** argv)
{

接下来，我们设置渲染器和渲染窗口，以便显示我们的结果. 我们设置了背景颜色和渲染窗口大小.

  //设置渲染器
  vtkNew renderer;
  renderer->SetBackground(0.9, 0.9, 0.9);

  //设置渲染窗口
  vtkNew renWin;
  renWin->AddRenderer(renderer.Get());
  renWin->SetSize(500, 500);

使用这段代码，我们已经可以显示一个静态渲染窗口. 相反，我们选择加a vtkRenderWindowInteractor 为了交互式旋转，缩放和平移场景.

  //设置渲染窗口交互器
  vtkNew interact;
  vtkNew style;
  interact->SetRenderWindow(renWin.Get());
  interact->SetInteractorStyle(style.Get());

现在我们有一个运行的例子，显示一个灰色的、空的渲染窗口.

接下来，我们使用VTK附带的众多阅读器之一加载数据集.

  // Read the file
  vtkSmartPointer pumpReader = vtkSmartPointer::New();
  pumpReader->SetFileName("rotor.vtu");

对VTK内存管理的简短介绍: VTK使用一个方便的自动内存管理概念，围绕引用计数. 但是与大多数其他实现不同, 引用计数保存在VTK对象本身中, 而不是智能指针类. 这样做的好处是可以增加引用计数, 即使VTK对象作为原始指针传递. 创建托管VTK对象有两种主要方法. vtkNew and vtkSmartPointer::New()，主要区别是a vtkSmartPointer 隐式可转换为原始指针吗 T*，并且可以从函数返回. For instances of vtkNew 我们会打电话的 .Get() 来获得一个原始指针，我们只能通过将它包装到 vtkSmartPointer. 在我们的例子中, 我们永远不会从函数中返回，所有的对象在整个过程中都是存在的, 因此我们使用缩写 vtkNew，但上述例外情况仅作示范用途.

此时，还没有从文件中读取任何内容. 我们或者链下的过滤器必须调用 Update() 让文件读取实际发生. 让VTK类自己处理更新通常是最好的方法. However, 有时我们希望直接访问过滤器的结果, 例如，获取该数据集中的压力范围. 那我们就得打电话 Update() manually. (我们不会因为调用而损失性能 Update() 多次，因为结果被缓存.)

  //获取压力范围
  pumpReader->Update();
  双pressureRange [2];
  pumpReader->GetOutput()->GetPointData()->GetArray("Pressure")->GetRange(pressureRange);

接下来，我们需要提取数据集的左半部分，使用 vtkClipDataSet. 要实现这一点，我们首先定义a vtkPlane 这就定义了分割. 然后，我们将首次看到VTK管道是如何连接在一起的: successor->SetInputConnection(predecessor->GetOutputPort()). 当我们请求更新时 clipperLeft 这个连接现在将确保所有前面的过滤器也是最新的.

  //从输入中截取左边的部分
  vtkNew planeLeft;
  planeLeft->SetOrigin(0.0, 0.0, 0.0);
  planeLeft->SetNormal(-1.0, 0.0, 0.0);

  vtkNew clipperLeft;
  clipperLeft->SetInputConnection(pumpReader->GetOutputPort());
  clipperLeft->SetClipFunction(planeLeft.Get());

最后，我们创建了第一个演员和映射器来显示左半部分的线框渲染. Notice, 映射器以与过滤器之间完全相同的方式连接到它的过滤器. Most of the time, 渲染器本身会触发所有角色的更新, 映射器和底层过滤器链!

唯一不能解释的行是“可能” leftWireMapper->ScalarVisibilityOff(); -禁止按压力值对线框着色, 哪些被设置为当前活动的数组.

  //创建左边部分的线框图
  vtkNew leftWireMapper;
  leftWireMapper->SetInputConnection(clipperLeft->GetOutputPort());
  leftWireMapper->ScalarVisibilityOff();

  vtkNew 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());

此时，渲染窗口终于显示了一些东西，i.e.，左边部分的线框图.

以类似的方式创建右侧部分的线框渲染, 通过切换(新创建的)的平面法线 vtkClipDataSet 到相反的方向，稍微改变(新创建的)映射器和角色的颜色和不透明度. Notice, 在这里，我们的VTK管道从相同的输入数据集分成两个方向(右和左).

  //从输入中截取右边的部分
  vtkNew planeRight;
  planeRight->SetOrigin(0.0, 0.0, 0.0);
  planeRight->SetNormal(1.0, 0.0, 0.0);

  vtkNew clipperRight;
  clipperRight->SetInputConnection(pumpReader->GetOutputPort());
  clipperRight->SetClipFunction(planeRight.Get());

  //为右边的部分创建线框图
  vtkNew rightWireMapper;
  rightWireMapper->SetInputConnection(clipperRight->GetOutputPort());
  rightWireMapper->ScalarVisibilityOff();

  vtkNew 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 我们在上面推导过.

  //创建左边部分的压力表示
  vtkNew pressureColorMapper;
  pressureColorMapper->SetInputConnection(clipperLeft->GetOutputPort());
  pressureColorMapper->SelectColorArray("Pressure");
  pressureColorMapper->SetScalarRange(pressureRange);

  vtkNew pressureColorActor;
  pressureColorActor->SetMapper(pressureColorMapper.Get());
  pressureColorActor->GetProperty()->SetOpacity(0.5);
  renderer->AddActor(pressureColorActor.Get());

输出现在看起来如下图所示. 中间的压力很低，将物料吸入泵内. 然后，这种材料被运送到外面，迅速获得压力. (当然应该有一个带有实际值的彩色图例, 但为了让例子简短，我省略了它.)

现在棘手的部分开始了. 我们想在右边画出速度流线. 流线是由源点在矢量场内的积分产生的. 矢量场已经以“velocity”矢量数组的形式成为数据集的一部分. 所以我们只需要生成源点. vtkPointSource 生成一个随机点的球体. 我们将生成1500个源点, 因为它们中的大多数无论如何都不会位于数据集中，并且会被流跟踪器忽略.

  //为流线创建源点
  vtkNew pointSource;
  pointSource->SetCenter(0.0, 0.0, 0.015);
  pointSource->SetRadius(0.2);
  pointSource->SetDistributionToUniform();
  pointSource->SetNumberOfPoints(1500);

接下来，我们创建流跟踪程序并设置其输入连接. “Wait, multiple connections?”, you might say. 是的，这是我们遇到的第一个具有多个输入的VTK滤波器. 普通输入连接用于向量场, 源连接用于种子点. 由于“velocity”是中的“活动”向量数组 clipperRight，我们不需要在这里明确地指定它. 最后，我们规定了积分应该从种子点向两个方向进行, 将积分方法设为 Runge-Kutta-4.5.

  vtkNew tracer;
  tracer->SetInputConnection(clipperRight->GetOutputPort());
  tracer->SetSourceConnection(pointSource->GetOutputPort());
  tracer->SetIntegrationDirectionToBoth();
  tracer->SetIntegratorTypeToRungeKutta45();

下一个问题是用速度大小给流线上色. 因为没有数组表示向量的大小, 我们将简单地将大小计算到一个新的标量数组中. 正如你所猜到的，这个任务也有一个VTK过滤器: vtkArrayCalculator. 它接受一个数据集并输出它, 但是只添加一个数组，该数组是从一个或多个现有数组中计算出来的. 我们将这个数组计算器配置为获取Velocity矢量的大小并将其输出为MagVelocity. Finally, we call Update() 再次手动，以导出新数组的范围.

  //计算速度大小并创建带状
  vtkNew magCalc;
  magCalc->SetInputConnection(tracer->GetOutputPort());
  magCalc->AddVectorArrayName("Velocity");
  magCalc->SetResultArrayName("MagVelocity");
  magCalc->SetFunction("mag(Velocity)");

  magCalc->Update();
  双magVelocityRange [2];
  magCalc->GetOutput()->GetPointData()->GetArray("MagVelocity")->GetRange(magVelocityRange);

vtkStreamTracer 直接输出折线和 vtkArrayCalculator 将它们原封不动地传递下去. 因此我们可以直接显示的输出 magCalc 直接使用新的映射器和actor.

相反，在这次训练中，我们选择通过显示缎带来使输出更好一些. vtkRibbonFilter 生成2D单元格以显示其输入的所有折线的带.

  //创建并渲染彩带
  vtkNew ribbonFilter;
  ribbonFilter->SetInputConnection(magCalc->GetOutputPort());
  ribbonFilter->SetWidth(0.0005);

  vtkNew streamlineMapper;
  streamlineMapper->SetInputConnection(ribbonFilter->GetOutputPort());
  streamlineMapper->SelectColorArray("MagVelocity");
  streamlineMapper->SetScalarRange(magVelocityRange);

  vtkNew streamlineActor;
  streamlineActor->SetMapper(streamlineMapper.Get());
  renderer->AddActor(streamlineActor.Get());

现在还缺少什么, 实际上也需要生成中间效果图, 最后五行是用来渲染场景和初始化交互器的吗.

  //渲染和显示交互式窗口
  renWin->Render();
  interact->Initialize();
  interact->Start();
  return 0;
}

最后，我们到达了完成的可视化，我将在这里再次展示:

可以找到上述可视化的完整源代码 here.

好，坏，丑

我将以我个人对VTK框架的优缺点来结束这篇文章.

• Pro: 积极开发: VTK正在几个贡献者的积极发展中, 主要来自研究领域. 这意味着一些尖端的算法是可用的, 许多3d格式可以导入和导出, bug被积极修复, 问题通常在讨论区有现成的解决方案.

• Con: Reliability将来自不同贡献者的许多算法与VTK的开放管道设计相结合, 不寻常的过滤器组合会导致问题吗. 我不得不进入VTK源代码几次，以弄清楚为什么我复杂的过滤器链没有产生预期的结果. 我强烈建议以允许调试的方式设置VTK.

• Pro: 软件架构: VTK的管道设计和总体架构似乎经过深思熟虑，很高兴与之合作. 几行代码就可以产生惊人的结果. 内置的数据结构易于理解和使用.

• Con: 微体系结构一些微建筑设计决策超出了我的理解范围. const正确性几乎不存在, 数组作为输入和输出传递，没有明显的区别. 为了减轻这个问题，我放弃了一些性能，并使用了我自己的包装器 vtkMath 它利用自定义3D类型，如 typedef std::array Pnt3d;.

• Pro: 微的文档氧的所有类和过滤器的文档是广泛的和可用的, wiki上的示例和测试用例对理解如何使用过滤器也有很大帮助.

• Con: 宏的文档:网上有一些很好的VTK教程和介绍. 然而，据我所知, 没有大的参考文档来解释具体的事情是如何完成的. 如果你想尝试新事物，就要花些时间去寻找方法. 此外，很难为任务找到特定的过滤器. 然而，一旦您找到了它，氧气文档通常就足够了. 探索VTK框架的一个好方法是下载并试用Paraview.

• Pro: 隐式并行化支持:如果您的源可以分成几个部分，可以独立处理, 并行化就像在处理单个部件的每个线程中创建单独的过滤器链一样简单. 大多数大型可视化问题通常属于这一类.

• Con: 没有显式并行化支持如果你没有大的福气, 可分的问题, 但你想利用多核, 你得靠自己了. 您必须弄清楚哪些类是线程安全的, 甚至可以通过反复试验或阅读源代码来重新进入. 我曾经追踪到一个并行化问题的VTK过滤器，它使用一个静态全局变量来调用一些C库.

• Pro: Buildsystem CMake:多平台元构建系统 CMake 也是由Kitware (VTK的制造商)开发的，并用于Kitware之外的许多项目. 它与VTK集成得非常好，并使为多个平台建立构建系统变得不那么痛苦.

• Pro: 平台独立性、许可和寿命: VTK是独立于平台的开箱即用的，并根据a 非常宽松的bsd风格许可证. 此外，对于那些需要它的重要项目，可以提供专业支持. Kitware得到了许多研究机构和其他公司的支持，并将在一段时间内出现.

Last Word

总的来说，对于我喜欢的问题，VTK是最好的数据可视化工具. 如果你遇到一个需要可视化的项目, mesh processing, 图像处理或类似的任务, 尝试用输入示例启动Paraview，并评估VTK是否适合您.