オープンソースツールを使用した3Dデータの視覚化：VTKを使用したチュートリアル

公開: 2022-03-11

Toptalのブログに関する最近の記事で、熟練したデータサイエンティストのCharles Cookは、オープンソースツールを使用した科学計算について書いています。 彼のチュートリアルでは、オープンソースツールと、データを簡単に処理して結果を取得する際にツールが果たすことができる役割について重要なポイントを示しています。

しかし、これらの複素微分方程式をすべて解くとすぐに、別の問題が発生します。 これらのシミュレーションから得られる膨大な量のデータをどのように理解し、解釈するのでしょうか。 大規模なシミュレーション内で数百万のグリッドポイントを持つデータなど、潜在的なギガバイトのデータをどのように視覚化するのでしょうか。

3Dデータ視覚化ツールに関心のあるデータサイエンティスト向けのデータ視覚化トレーニング。

修士論文の同様の問題に取り組んでいるときに、データの視覚化に特化した強力なグラフィックライブラリであるVisualization Toolkit（VTK）に接触しました。

このチュートリアルでは、VTKとそのパイプラインアーキテクチャを簡単に紹介し、インペラーポンプでシミュレートされた流体からのデータを使用した実際の3D視覚化の例について説明します。 最後に、ライブラリの長所と、遭遇した短所をリストします。

データの視覚化とVTKパイプライン

オープンソースライブラリVTKには、多くの高度な視覚化アルゴリズムを備えた堅実な処理およびレンダリングパイプラインが含まれています。 ただし、時間の経過とともに画像およびメッシュ処理アルゴリズムも追加されているため、この機能はそれだけではありません。 歯科研究会社との現在のプロジェクトでは、QtベースのCADのようなアプリケーション内のメッシュベースの処理タスクにVTKを利用しています。 VTKのケーススタディは、幅広い適切なアプリケーションを示しています。

VTKのアーキテクチャは、強力なパイプラインの概念を中心に展開しています。 この概念の基本的な概要を次に示します。

これは、VTKデータ視覚化パイプラインがどのように見えるかです。

  • ソースはパイプラインの最初にあり、「何もないところから何か」を作成します。 たとえば、 vtkConeSourceは3Dコーンを作成し、vtkSTLReaderは*.stl vtkSTLReaderジオメトリファイルを読み取ります。
  • フィルタは、ソースまたは他のフィルタの出力を新しいものに変換します。 たとえば、 vtkCutterは、平面などの陰関数を使用して、アルゴリズム内の前のオブジェクトの出力をカットします。 VTKに付属するすべての処理アルゴリズムはフィルターとして実装されており、自由にチェーンできます。
  • マッパーはデータをグラフィックスプリミティブに変換します。 たとえば、科学データを着色するためのルックアップテーブルを指定するために使用できます。 これらは、何を表示するかを指定する抽象的な方法です。
  • アクターは、シーン内のオブジェクト（ジオメトリと表示プロパティ）を表します。 ここでは、色、不透明度、陰影、向きなどを指定します。
  • レンダラーとWindowsは、最終的にプラットフォームに依存しない方法で画面へのレンダリングを記述します。

典型的なVTKレンダリングパイプラインは、1つ以上のソースから始まり、さまざまなフィルターを使用してそれらをいくつかの出力オブジェクトに処理し、マッパーとアクターを使用して別々にレンダリングします。 この概念の背後にある力は、更新メカニズムです。 フィルタまたはソースの設定が変更されると、依存するすべてのフィルタ、マッパー、アクター、およびレンダリングウィンドウが自動的に更新されます。 一方、パイプラインのさらに下流にあるオブジェクトがそのタスクを実行するために情報を必要とする場合、それは簡単に取得できます。

さらに、OpenGLのようなレンダリングシステムを直接扱う必要はありません。 VTKは、すべての低レベルのタスクをプラットフォームおよび（部分的に）レンダリングシステムに依存しない方法でカプセル化します。 開発者ははるかに高いレベルで作業します。

ローターポンプデータセットを使用したコード例

IEEE Visualization Contest 2011の回転インペラーポンプ内の流体の流れのデータセットを使用したデータ視覚化の例を見てみましょう。データ自体は、CharlesCookの記事で説明されているものとよく似た計算流体力学シミュレーションの結果です。

注目のポンプのzip形式のシミュレーションデータのサイズは30GBを超えています。 複数のパーツと複数のタイムステップが含まれているため、サイズが大きくなります。 このガイドでは、これらのタイムステップの1つで、圧縮サイズが約150MBのローター部分を試してみます。

私がVTKを使用するために選択した言語はC++ですが、Tcl / Tk、Java、Pythonなどの他のいくつかの言語のマッピングがあります。 ターゲットが単一のデータセットの視覚化だけである場合、コードを記述する必要はまったくなく、代わりに、VTKのほとんどの機能のグラフィカルフロントエンドであるParaviewを利用できます。

データセットと64ビットが必要な理由

Paraviewで1つのタイムステップを開き、ローターパーツを別のファイルに抽出することにより、上記の30GBデータセットからローターデータセットを抽出しました。 これは非構造格子ファイルです。つまり、六面体、四面体などのポイントと3Dセルで構成される3Dボリュームです。 各3Dポイントには、関連付けられた値があります。 セルに値が関連付けられている場合もありますが、この場合はそうではありません。 このトレーニングでは、ポイントでの圧力と速度に焦点を当て、これらを3Dコンテキストで視覚化しようとします。

VTKをロードした場合、圧縮ファイルサイズは約150 MB、メモリ内サイズは約280MBです。 ただし、VTKで処理することにより、データセットはVTKパイプライン内で複数回キャッシュされ、32ビットプログラムの2GBのメモリ制限にすぐに到達します。 VTKを使用するときにメモリを節約する方法はいくつかありますが、簡単にするために、64ビットでサンプルをコンパイルして実行します。

謝辞：データセットは、ドイツのクラウスタール大学応用力学研究所（Dipl。Wirtsch.-Ing。Andreas Lucius）の厚意により提供されています。

ターゲット

VTKをツールとして使用して達成するのは、下の画像に示す視覚化です。 3Dコンテキストとして、データセットのアウトラインは、部分的に透明なワイヤーフレームレンダリングを使用して表示されます。 次に、データセットの左側を使用して、表面の単純な色分けを使用して圧力を表示します。 （この例では、より複雑なボリュームレンダリングをスキップします）。 速度フィールドを視覚化するために、データセットの右側の部分は、速度の大きさによって色分けされた流線で埋められています。 この視覚化の選択は技術的には理想的ではありませんが、VTKコードをできるだけ単純に保ちたいと思いました。 さらに、この例が視覚化の課題の一部である理由があります。つまり、流れに多くの乱流があります。

これは、VTKチュートリアルの例から得られた3Dデータの視覚化です。

ステップバイステップ

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に付属している多くのリーダーの1つを使用してデータセットをロードします。

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

VTKメモリ管理への短いエクスカーション：VTKは、参照カウントを中心に展開する便利な自動メモリ管理の概念を使用しています。 ただし、他のほとんどの実装とは異なり、参照カウントは、スマートポインタークラスではなく、VTKオブジェクト自体に保持されます。 これには、VTKオブジェクトがrawポインターとして渡された場合でも、参照カウントを増やすことができるという利点があります。 マネージドVTKオブジェクトを作成する主な方法は2つあります。 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());

最後に、最初のアクターとマッパーを作成して、左半分のワイヤーフレームレンダリングを表示します。 マッパーは、フィルターが相互に接続されているのとまったく同じ方法でフィルターに接続されていることに注意してください。 ほとんどの場合、レンダラー自体がすべてのアクター、マッパー、および基盤となるフィルターチェーンの更新をトリガーしています。

自明ではない唯一の行は、おそらく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());

この時点で、レンダリングウィンドウは最終的に何か、つまり左側のワイヤーフレームを表示しています。

これは、VTKツールからの3Dデータ視覚化の結果の例でもあります。

右側のワイヤフレームレンダリングも同様の方法で作成されます。これは、（新しく作成された） vtkClipDataSetの平面法線を反対方向に切り替え、（新しく作成された）マッパーとアクターの色と不透明度をわずかに変更することです。 ここで、VTKパイプラインが同じ入力データセットから2つの方向（右と左）に分割されていることに注意してください。

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

予想どおり、出力ウィンドウに両方のワイヤーフレームパーツが表示されます。

VTKの例に従って、データ視覚化出力ウィンドウに両方のワイヤーフレームパーツが表示されるようになりました。

これで、いくつかの有用なデータを視覚化する準備が整いました。 左側の部分に圧力の視覚化を追加するために、多くのことを行う必要はありません。 新しいマッパーを作成し、それを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フィルターです。 通常の入力接続はベクトル場に使用され、ソース接続はシードポイントに使用されます。 「Velocities」はclipperRightの「アクティブな」ベクトル配列であるため、ここで明示的に指定する必要はありません。 最後に、シードポイントから両方向に積分を実行するように指定し、積分方法をルンゲクッタ4.5に設定します。

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

次の問題は、速度の大きさで流線を着色することです。 ベクトルの大きさの配列がないため、大きさを計算して新しいスカラー配列にします。 ご想像のとおり、このタスクにはvtkArrayCalculatorというVTKフィルターもあります。 データセットを取得して変更せずに出力しますが、既存の1つ以上の配列から計算された配列を1つだけ追加します。 この配列計算機は、「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());

現在も欠落しており、中間レンダリングを生成するためにも実際に必要なのは、シーンを実際にレンダリングしてインタラクターを初期化するための最後の5行です。

 // Render and show interactive window renWin->Render(); interact->Initialize(); interact->Start(); return 0; }

最後に、完成した視覚化に到達します。これをもう一度ここで紹介します。

VTKトレーニングの演習により、この完全な視覚化の例が得られます。

上記の視覚化の完全なソースコードは、ここにあります。

良い、悪い、そして醜い

この記事を、VTKフレームワークの個人的な長所と短所のリストで締めくくります。

  • プロ：活発な開発：VTKは、主に研究コミュニティ内から、いくつかの貢献者から活発に開発されています。 これは、いくつかの最先端のアルゴリズムが利用可能であり、多くの3D形式をインポートおよびエクスポートでき、バグが積極的に修正され、問題は通常、ディスカッション掲示板で既成の解決策を持っていることを意味します。

  • 短所：信頼性：さまざまな貢献者からの多くのアルゴリズムをVTKのオープンパイプライン設計と組み合わせると、異常なフィルターの組み合わせで問題が発生する可能性があります。 複雑なフィルターチェーンが目的の結果を生成しない理由を理解するために、VTKソースコードを数回調べなければなりませんでした。 デバッグが可能な方法でVTKを設定することを強くお勧めします。

  • プロ：ソフトウェアアーキテクチャ：VTKのパイプライン設計と一般的なアーキテクチャはよく考えられているようで、一緒に作業するのは楽しいことです。 数行のコード行で驚くべき結果が得られます。 組み込みのデータ構造は、理解と使用が簡単です。

  • 短所：マイクロアーキテクチャ：いくつかのマイクロアーキテクチャ設計の決定は、私の理解を逃れます。 Const-correctnessはほとんど存在せず、配列は明確な区別なしに入力と出力として渡されます。 パフォーマンスを放棄し、 typedef std::array<double, 3> Pnt3d;のようなカスタム3Dタイプを利用するvtkMathの独自のラッパーを使用することで、独自のアルゴリズムでこれを軽減しました。 。

  • Pro ：マイクロドキュメンテーション：すべてのクラスとフィルターのDoxygenドキュメンテーションは広範で使用可能です。また、wikiの例とテストケースは、フィルターの使用方法を理解するのに非常に役立ちます。

  • 短所：マクロドキュメント：Web上にVTKの優れたチュートリアルと紹介がいくつかあります。 しかし、私が知る限り、特定のことがどのように行われるかを説明する大きなリファレンスドキュメントはありません。 あなたが何か新しいことをしたいのなら、しばらくの間それをする方法を探すことを期待してください。 さらに、タスクの特定のフィルターを見つけるのは困難です。 ただし、それを見つけたら、通常はDoxygenのドキュメントで十分です。 VTKフレームワークを探索する良い方法は、Paraviewをダウンロードして実験することです。

  • Pro ：暗黙的な並列化のサポート：ソースを個別に処理できる複数の部分に分割できる場合、並列化は、単一の部分を処理する各スレッド内に個別のフィルターチェーンを作成するのと同じくらい簡単です。 ほとんどの大規模な視覚化の問題は通常、このカテゴリに分類されます。

  • 短所：明示的な並列化のサポートなし：大きくて分割可能な問題に恵まれていないが、複数のコアを利用したい場合は、自分で行ってください。 どのクラスがスレッドセーフであるか、または試行錯誤によって、またはソースを読み取ることによって再入可能であるかを把握する必要があります。 私はかつて、Cライブラリを呼び出すために静的グローバル変数を使用するVTKフィルターへの並列化の問題を追跡しました。

  • Pro ： Buildsystem CMake ：マルチプラットフォームのmeta-build-system CMakeもKitware（VTKのメーカー）によって開発され、Kitware以外の多くのプロジェクトで使用されています。 VTKと非常にうまく統合され、複数のプラットフォーム用のビルドシステムのセットアップがはるかに簡単になります。

  • Pro ：プラットフォームの独立性、ライセンス、および寿命：VTKは、すぐに使用できるプラットフォームに依存せず、非常に寛容なBSDスタイルのライセンスの下でライセンスされます。 さらに、それを必要とする重要なプロジェクトに対して専門的なサポートを利用できます。 Kitwareは多くの研究機関や他の企業に支えられており、しばらくの間存在するでしょう。

最後の言葉

全体として、VTKは、私が好きな種類の問題に最適なデータ視覚化ツールです。 視覚化、メッシュ処理、画像処理、または同様のタスクを必要とするプロジェクトに遭遇した場合は、入力例を使用してParaviewを起動し、VTKがツールになり得るかどうかを評価してみてください。