使用开源工具进行3D数据可视化：使用VTK的教程

本文概述

• 数据可视化和VTK管道
• Rotor Pump数据集的代码示例
• 优点和缺点
• 总结

熟练的数据科学家查尔斯·库克(Charles Cook)在最近在srcmini博客上发表的文章中, 谈到了使用开源工具进行科学计算。他的教程重点介绍了开源工具及其在轻松处理数据和获取结果中可以发挥的作用。

但是, 一旦我们解决了所有这些复杂的微分方程, 就会出现另一个问题。我们如何理解和解释这些模拟产生的大量数据？我们如何可视化潜在的千兆字节数据, 例如大型仿真中具有数百万个网格点的数据？

在为硕士论文解决类似问题的过程中, 我接触了Visualization Toolkit或VTK-一个功能强大的图形库, 专门用于数据可视化。

在本教程中, 我将快速介绍VTK及其管线体系结构, 并继续讨论使用叶轮泵中模拟流体数据得出的真实3D可视化示例。最后, 我将列出库的优点以及遇到的缺点。

数据可视化和VTK管道

开源库VTK包含具有许多复杂可视化算法的可靠处理和渲染管道。但是, 它的功能不止于此, 因为随着时间的流逝, 还添加了图像和网格处理算法。在我目前与一家牙科研究公司合作的项目中, 我正在利用VTK在基于Qt的类似CAD的应用程序中执行基于网格的处理任务。 VTK案例研究显示了广泛的适合应用。

VTK的架构围绕着强大的管道概念展开。此概念的基本轮廓如下所示：

• 来源位于管道的最开始, 创建”一无所有的东西”。例如, vtkConeSource创建一个3D圆锥体, 而vtkSTLReader读取* .stl 3D几何文件。
• 筛选器将源或其他筛选器的输出转换为新的内容。例如, vtkCutter使用隐式函数(例如平面)在算法中剪切先前对象的输出。 VTK随附的所有处理算法均实现为过滤器, 并且可以自由链接在一起。
• 映射器将数据转换为图形基元。例如, 它们可以用于指定用于为科学数据着色的查找表。它们是指定显示内容的抽象方法。
• 角色代表场景中的一个对象(几何形状和显示属性)。此处指定了颜色, 不透明度, 阴影或方向之类的内容。
• 渲染器和Windows最终以独立于平台的方式在屏幕上描述了渲染。

典型的VTK渲染管道从一个或多个源开始, 使用各种过滤器将它们处理成几个输出对象, 然后使用映射器和角色分别进行渲染。这个概念背后的力量是更新机制。如果更改了滤镜或源的设置, 则会自动更新所有相关的滤镜, 映射器, actor和渲染窗口。另一方面, 如果位于管道下游的对象需要信息以执行其任务, 则可以轻松获取它。

另外, 不需要直接处理像OpenGL这样的渲染系统。 VTK以平台和(部分)与系统无关的方式封装了所有低级任务。开发人员的工作水平更高。

Rotor Pump数据集的代码示例

让我们看一个数据可视化示例, 该示例使用来自IEEE Visualization Contest 2011的旋转叶轮泵中的流体流数据集。数据本身是计算流体动力学仿真的结果, 非常类似于Charles Cook的文章中所述。

特色泵的压缩模拟数据大小超过30 GB。它包含多个部分和多个时间步长, 因此尺寸较大。在本指南中, 我们将介绍这些时间步之一的转子部分, 压缩后的大小约为150 MB。

我选择使用VTK的语言是C ++, 但是还有其他几种语言的映射, 例如Tcl / Tk, Java和Python。如果目标仅是单个数据集的可视化, 则完全不需要编写代码, 而可以利用Paraview(VTK的大多数功能的图形前端)。

数据集以及为什么需要64位

通过在Paraview中打开一个时间步并将转子零件提取到一个单独的文件中, 我从上面提供的30 GB数据集中提取了转子数据集。它是非结构化的网格文件, 即由点和3D单元(如六面体, 四面体等)组成的3D体积。每个3D点都有关联的值。有时单元格也具有关联的值, 但在这种情况下则没有。该培训将集中于这些点的压力和速度, 并尝试在其3D上下文中将其可视化。

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

致谢：数据集由德国克劳斯塔尔大学应用​​力学研究所(Dipl。Wirtsch.-Ing。Andreas Lucius)提供。

目标

使用VTK作为工具, 我们将实现的是下图所示的可视化效果。作为3D上下文, 使用部分透明的线框渲染显示数据集的轮廓。然后, 数据集的左侧部分通过简单的表面颜色编码用于显示压力。 (在此示例中, 我们将跳过更复杂的体积渲染)。为了可视化速度场, 数据集的右侧填充了流线, 流线通过其速度大小进行颜色编码。这种可视化选择在技术上并不理想, 但我想保持VTK代码尽可能简单。另外, 该示例有理由成为可视化挑战的一部分, 即, 流动中存在许多湍流。

一步步

我将逐步讨论VTK代码, 以显示渲染输出在每个阶段的外观。完整的源代码可以在培训结束时下载。

让我们首先包括VTK所需的一切, 然后打开主要功能。

#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和mapper来显示左半部分的线框渲染。请注意, 映射器的连接方式与过滤器彼此完全相同。大多数时候, 渲染器本身会触发所有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的平面法线切换到相反的方向, 并略微更改(新创建的)映射器和actor的颜色和不透明度, 以类似的方式创建右侧部分的线框渲染。注意, 这里我们的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);

接下来, 我们创建streamtracer并设置其输入连接。你可能会说：”等等, 多个连接？”是的-这是我们遇到的多个输入的第一个VTK过滤器。普通输入连接用于矢量场, 源连接用于种子点。由于”速度”是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不变地传递折线。因此, 我们可以直接使用新的映射器和actor直接显示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的管道设计和通用体系结构似乎经过深思熟虑, 并且很高兴一起使用。几行代码会产生惊人的结果。内置的数据结构易于理解和使用。

• 缺点：微体系结构：一些微体系结构设计决策使我无法理解。 const正确性几乎不存在, 数组作为输入和输出传递, 没有明显的区别。我通过放弃一些性能并将自己的包装器用于vtkMath来减轻我自己算法的负担, 该包装器使用了诸如typedef std :: array <double, 3> Pnt3d;之类的自定义3D类型。

• 优点：微型文档：所有类和过滤器的Doxygen文档都是广泛且可用的, Wiki上的示例和测试用例也有助于理解如何使用过滤器。

• 缺点：宏文档：Web上有一些很好的VTK教程和介绍。但是据我所知, 还没有大量的参考文档来说明特定操作的完成方式。如果你想做一些新的事情, 请期待一段时间后再做些什么。此外, 很难找到任务的特定过滤器。但是, 一旦找到它, Doxygen文档通常就足够了。探索VTK框架的一个好方法是下载并试用Paraview。

• 优点：隐式并行支持：如果你的源代码可以分为几个部分, 可以独立处理, 那么并行化就像在处理单个部分的每个线程中创建单独的过滤器链一样简单。大多数大型可视化问题通常都属于此类。

• 缺点：没有明确的并行化支持：如果你没有遇到大型的, 可分割的问题, 但是你想利用多个内核, 则只能靠自己了。你必须通过反复试验或阅读源代码找出哪些类是线程安全的, 甚至可以重入。我曾经追踪到VTK过滤器的并行化问题, 该过滤器使用静态全局变量来调用某些C库。

• 优点：Buildsystem CMake：多平台元构建系统CMake也由Kitware(VTK的制造商)开发, 并在Kitware之外的许多项目中使用。它与VTK很好地集成在一起, 使为多个平台设置构建系统的工作变得轻松得多。

• 优点：平台独立性, 许可和寿命：VTK是开箱即用的平台独立性, 并根据非常宽松的BSD风格许可进行许可。此外, 专业支持可用于需要它的那些重要项目。 Kitware得到了许多研究机构和其他公司的支持, 并将持续一段时间。

总结

总体而言, VTK是解决我喜欢的各种问题的最佳数据可视化工具。如果你遇到需要可视化, 网格处理, 图像处理或类似任务的项目, 请尝试使用输入示例启动Paraview并评估VTK是否适合你。