VTK 勉強会資料
大野
地球シミュレータセンター在籍時の勉強会資料 (2006 年ごろ)

Contents
1 基本と 2 次元データの可視化
1.1 勉強会の目標 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 VTK . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.3 データの可視化 序章 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.4 カラーコンター . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.5 等高線 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.6 同一の Window で複数の結果を表示 . . . . . . . . . . . . . . . . . . . . . . . . .
1.7 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3
3
3
4
4
9
12
17

2 基本と 2 次元データの可視化その 2
2.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 透視射影と正射影 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 凹凸つきスライス . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4 高さつき等高線 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.5 発展 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.6 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

19
19
19
26
33
37
41

3

VTK 形式のデータ
3.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 等間隔メッシュのデータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.3 Rectilinear Grid のデータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.4 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

44
44
44
48
48

4

3 次元スカラーデータの可視化 その 1
4.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 等値面 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3 データの切り出し : 断面切り . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.4 データの切り出し : Volume Of Interest (VOI) . . . . . . . . . . . . . . . . . . . .
4.5 データの間引き . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.6 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

49
49
49
52
55
62
68

5

3 次元スカラーデータの可視化 その 2
5.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2 ボリューム・レンダリング : レイ・キャスティング法 . . . . . . . . . . . . . . .
5.3 カラーコンター再び : テクスチャ・マッピングの利用 . . . . . . . . . . . . . . .
5.4 ボリューム・レンダリング : テクスチャ・マッピングの利用 . . . . . . . . . . .
5.5 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

69
69
69
74
80
86

1

6

3 次元ベクトルデータの可視化 その 1
88
6.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
6.2 HedgeHog(針鼠) によるベクトルデータの可視化 . . . . . . . . . . . . . . . . . . 88
6.3 矢印・コーンを使う . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
6.4 ベクトル版凹凸つきスライス . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96
6.5 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100

7

3 次元ベクトルデータの可視化 その 2
7.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7.2 流線 1 : 流線 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7.3 流線と矢印の組み合わせ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7.4 流線 2 : チューブ・リボン . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7.5 流線 3 : 流面 (Stream Surfaces) . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7.6 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

101
101
101
106
107
112
117

8

円柱座標・球座標
8.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8.2 円柱座標 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8.3 球座標 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8.4 VTK 形式のデータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
8.5 この章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

119
119
119
127
134
135

9

その他の事項
9.1 本章の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9.2 コマンドの組み込み . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9.3 時間発展のデータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9.4 ポリゴンデータ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9.5 アンチ・エイリアジング . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

136
136
136
141
148
152

2

Chapter 1
基本と 2 次元データの可視化
1.1

勉強会の目標多くのサンプルプログラムを通して VTK のプログラミング方法を習得し，自分のシミュレー
ション結果を自作ソフトで可視化できるようになる。
(「等値面を表示させてマウスでグリグリまわす」，このレベルのソフトなら作れるようになり
ます)

1.2

VTK

VTK (VISUALIZATION TOOLKIT) とは，フリーの可視化ライブラリである。これを利用すると，等値面生成などの可視化アルゴリズムや OpenGL などの CG プログラミングの知識がさほどなくとも，自分用の 3 次元可視化ソフトを容易に開発することができる。プログラミングは，C++, Tcl,
JAVA, Python でおこなえる（この勉強会では C++を使う）。また，この勉強会では，Version 4.2で作成したサンプルプログラムを使用する。

1.2.1

使用可能な可視化手法

カラーコンター，等高線，等値面，ボリューム・レンダリング，流線など，基本的なものはすべて使用可能。これ以外にも，チューブ状流線など面白い可視化手法が用意されている。

1.2.2

数値データ

型 Unsigned Char, Short, Int, Float, Double など
座標 等間隔，Rectilinear, StructuredGrid, UnstructuredGrid など

1.2.3

パイプライン

VTK を利用したプログラミングは，AVS のネットワークのように，VTK のオブジェクト同士を接続していくだけである。当然であるが，オブジェクト同士の接続においては，受け渡しす
るデータの型が一致していなければならない（接続可能かどうかは，オンライン・リファレン
スを調べるとわかる）。パイプラインの基本的な流れは，

3

Data → Filter → Mapper → Actor → Renderer → Window
Filter の部分で，等値面処理などを行う。

1.3

データの可視化 序章

本章では，ETOPO2 から日本近海を切り出した (VTK のデータ形式ではない) 生データ [zenkoku.dem: 900×900, short 型] をカラーコンターと等高線で可視化しつつ，VTK のプログラミングの雰囲
気をつかむ。

1.3.1

可視化の手順

等高線による可視化を行うプログラムを OpenGL を使って自作するとする。どのような作業をプログラムにさせたらよいだろうか？
次の手順は一つの例である (VTK のプログラムでも同様の手順をふむ)。

1. 可視化するデータを読み込む (サイズや座標の情報も取得する)
2. 等高線作成のアルゴリズムでポリゴンデータ (線分の集まり。カラーコンターなら三角形や四角形の集まり) を作る
3. 手順 2 で求めたポリゴンデータを OpenGL で描画する
注意すべき点は，手順 1 のカッコ内である。データがあっても，そのサイズ (200×150 など) や座標の情報 (等間隔なのか Rectilinear なのか? その刻み幅は? など) がなければ，可視化 (手順 2)は実行できない (AVS でもコントロールファイルに，それらの情報を記述しなければならない)。
VTK のプログラムでは，手順 1 でデータそのものと座標やサイズをあわせ持ったデータセッ
トを作らなければならない (あるいは，データを VTK 形式に変換しておかなければならない)。
データセットを作って初めて，手順 2 で等高線や等値面のポリゴンデータを生成するフィルタを使うことができる。

1.4

カラーコンター

等間隔メッシュを使って，日本近海の標高データを色で表現してみる。

1.4.1
1
2
3
4
5
6
7
8
9
10

サンプルプログラム 1/* ColorContour.cxx */#include <vtkImageData.h>#include <vtkShortArray.h>#include <vtkPointData.h>#include <vtkLookupTable.h>#include <vtkImageDataGeometryFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>
4

11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57

#include <unistd.h>#define MAP_X 900#define MAP_Y 900
short sdata[MAP_X * MAP_Y];
int main(int argc, char *argv[]){
int size[3] = { MAP_X, MAP_Y, 1 }; /* 2 次元データ */FILE *fpi;/* ハードディスクからデータを読み込む */if ((fpi = fopen("./zenkoku.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X * MAP_Y, fpi);
fclose(fpi);/* 読み込んだデータを sarray に代入 */vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X * MAP_Y, 1);/* 等間隔データを作成 */vtkImageData *imgData = vtkImageData::New();
imgData->SetScalarTypeToShort();imgData->SetDimensions(size);float sp = 0.1;
imgData->SetSpacing(sp, sp, 0);/* 等間隔メッシュの箱に，データを代入 */imgData->GetPointData()->SetScalars(sarray);/* カラーテーブルを作成 */vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.667, 0.0); /* 青 -> 赤 */lut->Build();/* imgData を基にポリゴンデータを作成*/
vtkImageDataGeometryFilter *igf
= vtkImageDataGeometryFilter::New();
igf->SetInput(imgData);
5

58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91

vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(igf->GetOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(-500.0, 3800.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren = vtkRenderer::New();
vren->AddActor(Actor);
vren->SetBackground(0.0, 0.0, 0.0); /* 背景を黒に設定 */vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer(vren);renWin->SetSize(600, 600); /*Window サイズを 600x600 に設定*//* 10 秒程度表示 */for(int i=0; i<10; i++){
renWin->Render();sleep(1);}sarray->Delete();imgData->Delete();lut->Delete();Mapper->Delete();Actor->Delete();vren->Delete();renWin->Delete();return 0;}

1.4.2

サンプルプログラム 1 の簡単な解説

このプログラムは，sdata に日本近海の標高データ (zenkoku.dem) を読み込み，それを VTK の等
間隔メッシュのデータにして，さらにそれを基にポリゴンを生成し，スカラーの値 (標高値) に応じて色づけして表示する (補足：2 次元データの可視化であるが，透視射影で表示している)。

1 - 12 : ヘッダファイルの読み込み
17 : 標高データを入れる short 型変数 “sdata” を宣言
21, 40 : Z 方向のサイズが “1” なので，このデータは 2 次元である
34 - 35 : VTK 用の short 型配列 (sarray) に sdata をセット
VTK のクラスは，new, delete 演算子は使わず，〜::New(), 〜->Delete() を使用する。35 行目のSetArray の 1 番目の引数はセットする配列データのポインタ，2 番目はサイズである。int, float,
6

double 型には，それぞれ vtkIntArray, vtkFloatArray, vtkDoubleArray が用意されている。
37 - 47 : VTK の等間隔メッシュデータ格納用のクラス (imgData) に sarray をセット
39 行でデータ型を，40 行でデータのサイズを，44 行でメッシュの幅を指定。原点の位置も指定できるが，このプログラムでは指定していない。
ようやく，(座標の情報を含んだ)VTK の等間隔メッシュのデータが完成。
49 - 52 : 色（青 → 赤）にセット
明度・彩度の指定や，自分でカラーテーブルを作ることもできる。
54 - 57 : imgData をポリゴン情報に変換
imgData(vtkImageData) をそのまま vtkPolyDataMapper に渡すことはできない。
59 - 63 : 色情報などと統合
このプログラムの可視化パイプラインを図で示すと，Figure 1.1 のようになる。

Figure 1.1: サンプルプログラム 1 のパイプライン構造
あるデータから別のデータに変換するオブジェクト (プロセス・オブジェクト) の出力は，a->GetOutput() と憶えて差し支えない (57 行，60 行を参考に)。であるから接続部分は，
入力側->SetInput(出力側->GetOutput());の形となる。
63 行で色づけするスカラーの範囲を指定しているが，-500.0(色づけ指定した範囲の最小値)より小さい値は青（最小値の色），3800.0（色づけ指定した範囲の最大値）より大きい値は，赤（最大値の色）になる。
70 : 背景の色を黒 (RGB = (0.0, 0.0, 0.0)) に指定
74 : Window の大きさを 600 × 600 に指定
カメラの位置などをまったく指定していないが，本サンプルプログラムの場合はデータ全体が表示されるように，VTK が自動的にカメラの位置・その他を決めている。もちろん指定することも可能である (たとえば下記リスト)。

1vren->GetActiveCamera()->SetPosition(30, 30, 100);7

Figure 1.2: カラーコンター ((左) カメラ位置指定なし, (右) カメラ位置指定)
2
3

vren->GetActiveCamera()->SetFocalPoint(30, 30, 0);vren->GetActiveCamera()->SetClippingRange(10, 200);

Figure 1.2 は，カメラ位置指定なしの場合と，(vtkCamera.h をインクルードして) 上記リス
トを 71 行目に挿入してカメラ位置を指定した場合に出力される画像。カメラ関係については，
次章でもう少し詳しく説明する。

1.4.3

画像保存

表示させた画像を保存しておきたい場合には，vtkRendererSource.h, vtkBMPWriter.h をインクルードして，下記リストを 81 行目，プログラムの最後に rs− >Delete(), bw− >Delete() を加
える。

1
2
3
4
5
6
7
8
9

vtkRendererSource *rs = vtkRendererSource::New();
rs->SetInput(vren);rs->WholeWindowOn();
rs->Modified();vtkBMPWriter *bw = vtkBMPWriter::New();
bw->SetInput(rs->GetOutput());
bw->SetFileName("Map.bmp");bw->Write();

このリストを挿入すると，Window に表示された画像が，BMP 形式でファイル名 “Map.bmp”
で保存される。vtkBMPWriter の類似のクラスとして，vtkJPEGWriter, vtkPNGWriter, vtkPostScriptWriter,
vtkTIFFWriter などがある。

8

1.5

等高線

次に，vtkContourFilter というフィルタを使って，等高線を描いてみる。このフィルタは，等間隔でも Rectiliner のデータでも使用可能で，2 次元データを渡すと等高線，3 次元データを渡すと等値面を生成する。今度は，Rectilinear のデータを使用してみる。

1.5.1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38

サンプルプログラム 2/* ContourLines.cxx */#include <vtkRectilinearGrid.h>#include <vtkFloatArray.h>#include <vtkShortArray.h>#include <vtkPointData.h>#include <vtkContourFilter.h>#include <vtkLookupTable.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <unistd.h>#define MAP_X 900#define MAP_Y 900#define NCONTLINES 10
short sdata[MAP_X * MAP_Y];
int main(int argc, char *argv[]){
int size[3] = { MAP_X, MAP_Y, 1 };FILE *fpi;
if ((fpi = fopen("./zenkoku.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X * MAP_Y, fpi);
fclose(fpi);
vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X * MAP_Y, 1);/* 座標データを格納するための配列 */vtkFloatArray *x_coord = vtkFloatArray::New();

9

39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

vtkFloatArray *y_coord = vtkFloatArray::New();
vtkFloatArray *z_coord = vtkFloatArray::New();
float z = 0.0;
float x[MAP_X];
float y[MAP_Y];
for (int i = 0; i < MAP_X; i++) {x[i] = (float) i *0.1;
y[i] = (float) i *0.1;}
x_coord->SetArray(x, MAP_X, 1);y_coord->SetArray(y, MAP_Y, 1);z_coord->SetArray(&z, 1, 1);
vtkRectilinearGrid *RectData= vtkRectilinearGrid::New();
RectData->SetXCoordinates(x_coord); /* 座標データをセット*/
RectData->SetYCoordinates(y_coord);
RectData->SetZCoordinates(z_coord);
RectData->SetDimensions(size);
RectData->GetPointData()->SetScalars(sarray);vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.667, 0.0);
lut->Build();/* 等高線を生成する */vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(RectData);/* 0.0 - 3800.0 の間で，NCONTLINES 本等高線を生成 */
contour->GenerateValues(NCONTLINES, 0.0, 3800.0);
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(contour->GetOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(-500.0, 3800.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren = vtkRenderer::New();
vren->AddActor(Actor);
vren->SetBackground(0.0, 0.0, 0.0);
10

86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109

vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer(vren);renWin->SetSize(600, 600);
for(int i=0; i<10; i++){
renWin->Render();sleep(1);}sarray->Delete();x_coord->Delete();y_coord->Delete();z_coord->Delete();RectData->Delete();lut->Delete();contour->Delete();Mapper->Delete();Actor->Delete();vren->Delete();renWin->Delete();return 0;}

1.5.2

サンプルプログラム 2 の簡単な説明

37 - 53 : Rectilinear の座標データ用の配列確保と，座標値の代入（結局等間隔にしています）
55 - 61 : 上記の座標データと sdata を Rectilinear のデータ (RectData) にセット。Rectilinear の 2
次元データ完成
68 - 72 : contour(vtkContourFilter) に RectData (vtkRectilinearGrid) を代入し，72 行で 0.0 〜 3800.0の間で等間隔に NCONTLINES 本 (10 本) 等高線を生成させる。等高線の値を直接指定することも可能であるこのサンプルプログラムの可視化パイプラインは，Figure 1.3 のようになる。
結果は Figure 1.4 。

1.5.3

パイプラインの接続可否

vtkContourFilter を例に考えてみる。
69 行 “contour->SetInput(RectData)”で，vtkContourFilter に RectilinearGrid のデータセット(のポインタ) を代入しているが，オンライン・リファレンスで SetInput の引数の定義を見てみると，SetInput(vtkRectilinearGrid *) ではなく，SetInput(vtkDataSet *) となっている。多重定
義されているわけでもない。これでもパイプラインが接続できる理由は，vtkRectilinearGrid がvtkDataSet の派生クラス (Figure 1.5) だからである。ちなみに vtkImageData も vtkDataSet の派
生クラスである。
このように，パイプラインの接続可否は，オンライン・リファレンスで該当クラスの継承関係を調べれるとわかる。§1.4.2 でも述べたが，例えば a と b を接続する場合，基本的に11

Figure 1.3: サンプルプログラム 2 のパイプライン構造
b->SetInput(a->GetOutput());という形をとるので，入力側の SetInput メソッドの引数 と出力側の GetOutput メソッドの返り値のデータ型に注目する。

1.6

同一の Window で複数の結果を表示

最後に「カラーコンター」，「等高線」，「その重ね合わせ」の 3 つの結果を同一の Window に表
示させる。

1.6.1
1
2
3
4
5
6
7
8
9
10
11
12
13
14

サンプルプログラム 3/* Contour3.cxx */#include <vtkRectilinearGrid.h>#include <vtkFloatArray.h>#include <vtkShortArray.h>#include <vtkPointData.h>#include <vtkContourFilter.h>#include <vtkLookupTable.h>#include <vtkRectilinearGridGeometryFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkCamera.h>#include <vtkActor.h>#include <vtkProperty.h>#include <vtkRenderer.h>12Figure 1.4: 等高線 ((左) カメラ位置指定なし, (右) カメラ位置指定)

Figure 1.5: クラスの継承
15
16
17
18
19
20
21
22
23
24
25
26
27

#include <unistd.h>#define MAP_X 900#define MAP_Y 900
short sdata[MAP_X*MAP_Y];
int main( int argc, char *argv[] ){int size[3] = {MAP_X, MAP_Y, 1};
FILE *fpi;

13

28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74

if ((fpi = fopen("./zenkoku.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X*MAP_Y, fpi);
fclose(fpi);
vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X*MAP_Y, 1);
vtkFloatArray *x_coord = vtkFloatArray::New();
vtkFloatArray *y_coord = vtkFloatArray::New();
vtkFloatArray *z_coord = vtkFloatArray::New();
vtkRectilinearGrid *RectData = vtkRectilinearGrid::New();
float z = 0.0;
float x[MAP_X];
float y[MAP_Y];
for(int i=0;i<MAP_X;i++){x[i] = (float)i*0.1;
y[i] = (float)i*0.1;}
x_coord->SetArray(x, MAP_X, 1);y_coord->SetArray(y, MAP_Y, 1);z_coord->SetArray(&z, 1, 1);
RectData->SetXCoordinates(x_coord);
RectData->SetYCoordinates(y_coord);
RectData->SetZCoordinates(z_coord);
RectData->SetDimensions(size);
RectData->GetPointData()->SetScalars(sarray);vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.667, 0.0);
lut->Build();
vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(RectData);
contour->SetValue(0, 0.0); /* スカラー値 0.0 の等高線を生成 */vtkRectilinearGridGeometryFilter *rgf
= vtkRectilinearGridGeometryFilter::New();
rgf->SetInput(RectData);
14

75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121

/* 等高線用の Mapper */vtkPolyDataMapper *Mapper1 = vtkPolyDataMapper::New();
Mapper1->SetInput(contour->GetOutput());/* スカラー値による色づけをしない */Mapper1->ScalarVisibilityOff();/* カラーコンター用の Mapper */vtkPolyDataMapper *Mapper2 = vtkPolyDataMapper::New();
Mapper2->SetInput(rgf->GetOutput());
Mapper2->SetLookupTable(lut);
Mapper2->SetColorModeToMapScalars();
Mapper2->SetScalarRange(-500.0, 3800.0);/* 等高線用の Actor */vtkActor *Actor1 = vtkActor::New();
Actor1->SetMapper(Mapper1);/* 等高線の色を白に設定 */Actor1->GetProperty()->SetColor(1.0, 1.0, 1.0);/* カラーコンター用の Actor */vtkActor *Actor2 = vtkActor::New();
Actor2->SetMapper(Mapper2);/* カラーコンターと等高線を同時に映す */vtkRenderer *vren= vtkRenderer::New();
vren->AddActor( Actor2 );vren->AddActor( Actor1 );vren->SetBackground( 0.0, 0.0, 0.0 );vren->GetActiveCamera()->SetPosition(30, 30, 100);
vren->GetActiveCamera()->SetFocalPoint(30, 30, 0);vren->GetActiveCamera()->SetClippingRange(10, 200);
vren->SetViewport(0.66, 0.0, 1.0, 1.0);/* 等高線を映す */vtkRenderer *vren1= vtkRenderer::New();
vren1->AddActor( Actor1 );vren1->SetBackground( 0.0, 0.0, 0.0 );vren1->SetViewport(0.33, 0.0, 0.66, 1.0);/* カラーコンターを映す */vtkRenderer *vren2= vtkRenderer::New();
vren2->AddActor( Actor2 );vren2->SetBackground( 0.0, 0.0, 0.0 );vren2->SetViewport(0.0, 0.0, 0.33, 1.0);
15

122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151

vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren1 );renWin->AddRenderer( vren2 );renWin->AddRenderer( vren );renWin->SetSize( 900, 300 );for(int i=0; i<10; i++){
renWin->Render();sleep(1);}sarray->Delete();x_coord->Delete();y_coord->Delete();z_coord->Delete();RectData->Delete();lut->Delete();contour->Delete();Mapper1->Delete();Mapper2->Delete();Actor1->Delete();Actor2->Delete();vren->Delete();vren1->Delete();vren2->Delete();renWin->Delete();return 0;}

1.6.2

サンプルプログラム 3 の簡単な説明

70 : 等高線の値を直接指定 (1 本目 (0 番) の等高線の値を 0.0 に指定。2 本目は 1 番, 3 本目は 2
番になる)。GenerateValues と違い，ピンポイントで指定できる72 - 74 : サンプルプログラム 1 の 54 - 57 行参照
81, 95 : スカラーによる色づけはせず，等高線を白 (RGB = (1.0, 1.0, 1.0)) で表示
103 - 104 : 同一の Renderer に Actor1(等高線) と Actor2(カラーコンター) をセットすることで，
カラーコンターと等高線を両方表示させる
124 - 126 : 同一の vtkRenderWindow に各 vtkRenderer をセットすることで，同一 Window に 3
つの結果を表示させる
109, 115, 121 : Window 内での位置を決定
WIndow の左下が (0.0, 0.0)，右上が (1.0, 1.0)。121 行では，Actor2 を映す範囲を Window の左側 1/3 使用するように SetViewport で設定している
127 : 3 つの結果を表示したいので，Window の大きさを 900 × 300 と横長に設定
Actor 以降のパイプラインは，Figure 1.6 のようになる。
16

Figure 1.6: サンプルプログラム 3 のパイプライン構造
結果と ViewPort の座標を Figure 1.7 に示す。データ同士の比較や，一つのデータをいろいろな角度から映した画像を並べるのに便利である。

Figure 1.7: 同一 Window 内での 3 パターンの可視化結果の表示（（左）カラーコンター，（中）等
高線，（右）組み合わせ）1.7

この章のまとめ• データの読み込み• 等間隔，Rectilinear• パイプラインの接続
17

本章で紹介したクラスは，Table 1.1 の通り。

Table 1.1: 本章て紹介したクラスおよびそのメソッド
配列

データセット

ポリゴン化
カラーテーブル等高線

MapperActor
Renderer
画像保存

Window

vtkShortArrayvtkFloatArray: SetArray(short *, vtkIdType, int) / vtkShortArrayvtkImageData, vtkRectilinearGrid: SetDimensions(int, int, int)
: SetSpacing(float, float, float): GetPointData()(->SetScalars(vtkDataArray*)): SetXCoordinates(vtkDataArray *) / vtkRectilinearGrid: SetYCoordinates(vtkDataArray *) / vtkRectilinearGrid: SetZCoordinates(vtkDataArray *) / vtkRectilinearGrid
vtkImageDataGeometryFilter
vtkRectilinearGridGeometryFilter
vtkLookupTable: SetHueRange(float, float): Build()
vtkContourFilter: GenerateValues(int N, float min, float max)
: SetValue(int, float)
vtkPolyDataMapper: SetColorModeToMapScalars(): SetLookupTable(vtkScalarsToColors *): SetScalarRange(float, float): ScalarVisibilityOff()vtkActor: GetProperty()(->SetColor(float, float, float))
vtkRenderer: SetViewPort(float, float, float, float)
vtkRendererSource
vtkBMPWriter, vtkPNGWriter, vtkJPEGWriter, vtkTIFFWriter
vtkRenderWindow

18Chapter 2
基本と 2 次元データの可視化その 2
2.1

本章の概要

前章に引き続き 2 次元の標高データを可視化するが，本章では 2 次元データを基に 3 次元オブジェクトを作成する。また，射影やカメラ，ポリゴンの法線ベクトルのことにも触れる。

2.2

透視射影と正射影

VTK は，透視射影と正射影のいずれもサポートしている。透視射影と正射影は，それぞれ以下のような特徴を持つ。
• 透視射影
物体の大きさが同じ場合，視点から近い物体ほど大きく，遠くにある物体ほど小さく画面に描く (遠近法で描いてくれる)• 正射影
近くとも遠くとも同じ大きさで画面に描く

2.2.1

透視射影

透視射影の場合は，Figure 2.1 で near と far の間の錐台の中が，画面に映し出される。この錘台の中に，映したいオブジェクトが入るようにカメラ関係のパラメータの設定をする。パラメー
タの設定は，vtkCamera と vtkRenderer で行える (サンプルプログラム 1 参照)。

• カメラの位置 (視点の位置) : SetPosition(double, double, double)• 焦点 : SetFocalPoint(double, double, double) [カメラの位置と焦点を結ぶ方向が，視線方向。
焦点より後ろが映らなくなるわけではない]• 上方向 : SetViewUp(double, double, double)• near, far : SetClippingRange(double near, double far) [ near と far はカメラからの距離]• ViewAngle : SetViewAngle(double)

19

Figure 2.1: 透視射影
上記は，vtkCamera のメソッドである。カメラの位置などを設定するには，まず vtkCamera を生成し，上記のメソッドでパラメータを設定，その後 vtkRenderer の SetActiveCamera(vtkCamera*) を介して，Renderer にカメラ情報を渡す (SetViewAngle を除いて多重定義されているので，
上記以外の設定方法もある)。カメラ位置などは，透視射影・正射影などかかわりなく，default
では，Figure 2.2 のように設定されている。

2.2.2

正射影

正射影の場合は，映し出される空間が錐台から直方体になる (Figure 2.3)。正射影では，映し出される画像は ViewAngle の影響を受けないが，映し出す範囲 (画面の縦の大きさ) を vtkCameraの下記メソッドで設定する必要がある。

• SetParallelScale(double scale)scale の値を変えることにより，映し出すオブジェクトに対する相対的な画面のサイズを変更できる。映し出すオブジェクトの高さが 6.0 の場合，画面いっぱいに映そうと思ったら，scale の値を 3.0(3 倍にする。画面の縦の長さは 2.0) にすればよい。正射影と透視射影の切り替えは，やはり vtkCamera の• ParallelProjectionOn()• ParallelProjectionOff()で行う。default では，透視射影である。

2.2.3

サンプルプログラム 1

このサンプルプログラムでは，vtkCamera と vtkRenderer でカメラ位置などを設定して，画面全体に可視化画像が表示されるようにしている。本章のサンプルデータ (MtFuji.dem) は，GTOPO30から富士山付近を切り出したサイズ 75 × 60 の Short 型バイナリデータである。

1/* ParallelPro.cxx */20

Figure 2.2: 座標
2
3
4
5
6
7
8
9
10
11
12
13

#include <vtkShortArray.h>#include <vtkImageData.h>#include <vtkPointData.h>#include <vtkContourFilter.h>#include <vtkLookupTable.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkCamera.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <unistd.h>

Figure 2.3: 正射影
21

14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60

#define MAP_X 75#define MAP_Y 60#define NCONTLINES 37
short sdata[MAP_X*MAP_Y];
int main( int argc, char *argv[] ){int size[3] = {MAP_X, MAP_Y, 1};
FILE *fpi;
if ((fpi =fopen("./MtFuji.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X*MAP_Y, fpi);
fclose(fpi);
vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X*MAP_Y, 1);
vtkImageData *imgData = vtkImageData::New();
imgData->SetSpacing(0.1,0.1,0);
imgData->SetDimensions(size);
imgData->GetPointData()->SetScalars(sarray);vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.667, 0.0);
lut->Build();
vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(imgData);
contour->GenerateValues(NCONTLINES, 0.0, 3776.0);
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(contour->GetOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(0.0, 3776.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren= vtkRenderer::New();
22

61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96

vren->AddActor( Actor );vren->SetBackground( 0.0, 0.0, 0.0 );/* カメラの設定 */vtkCamera *camera = vtkCamera::New();
camera->SetPosition(3.7, 2.95, 100); /* カメラ位置 */camera->SetFocalPoint(3.7, 2.95, 0); /* 焦点 */camera->SetViewUp(0, 1, 0);/* カメラの上方向 */camera->ParallelProjectionOn();/* 正射影 */camera->SetParallelScale(2.95);/* 画面サイズ*/
camera->SetClippingRange(50, 120);/* Near, Far *//* Renderer にカメラ情報を代入 */vren->SetActiveCamera(camera);
vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 740, 590 );for(int i=0;i<10;i++){
renWin->Render();sleep(1);}sarray->Delete();imgData->Delete();lut->Delete();contour->Delete();Mapper->Delete();Actor->Delete();vren->Delete();camera->Delete();renWin->Delete();return 0;}

2.2.4

サンプルプログラム 1 の簡単な説明

このプログラムは，カメラの X,Y 座標の位置をデータ (等高線) の中心 (XY 平面での) に据え，さらに画面の大きさもデータの大きさと同じにして，画面いっぱいに可視化像が映るようにしている。
カメラ関係の設定は，64 - 73 行目で行っている。

65 : camera(vtkCamera) の生成

23

Figure 2.4: カメラとオブジェクトの位置関係
66 : カメラの位置設定 1
67 : カメラの焦点設定
68 : カメラの上方向設定
69 : 射影方法を正射影に設定
70 : 画面 (表示するオブジェクト) の相対的なサイズを設定 (表示されるオブジェクトの高さは5.9 なので，引数に 2.95 を代入している)
71 : near, far を設定
74 : 66-71 行で設定した camera の内容を vtkRenderer に代入
カメラの位置関係は Figure 2.4，表示は Figure 2.5。

Figure 2.5: 富士山周辺の等高線：右下の山は箱根山

vtkImageData のメソッド SetOrigin(float, float, float) でデータの原点を (-3.7,-2.95,0.0) に，カメラと焦点の X,
Y 座標の値を 0.0 に設定しても，同じ画像が表示される。
124

2.2.5

カメラと焦点の位置の回転

FocalPoint の周りで，カメラの位置を簡単に回転させることができる。そのための vtkCameraのメソッドは，• Azimuth(double angle)• Elevation(double angle)
である。angle の単位は degree。これを使うと，カメラがオブジェクトの周りを周るアニメー

Figure 2.6: カメラ位置と焦点の回転方向
ションが簡単に作れる。サンプルプログラム 1 の 65 - 83 行目を下記のリストに差し替えると(vtkRenderSource.h と vtkBMPWriter.h のインクルードも忘れずに)，透視射影で等高線の周りを一週するアニメーションが Window に現れ，平行して各スナップショットの画像が (anim1 ディ
レクトリの下に)BMP 形式で保存される。この連番画像を使いムービーを作れば，VTK がなくともアニメーションを見ることができる。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16

vtkCamera *camera = vtkCamera::New();
camera->SetPosition(3.7, 2.95, 15);
camera->SetFocalPoint(3.7, 2.95, 0);
camera->SetViewUp(0, 1, 0);
camera->SetClippingRange(1, 20);
vren->SetActiveCamera(camera);
vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600);
vtkRendererSource *rs = vtkRendererSource::New();
rs->SetInput(vren);rs->WholeWindowOn();

25

17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35

vtkBMPWriter *bw = vtkBMPWriter::New();
bw->SetInput(rs->GetOutput());
char FileName[128];
for(int i=0;i<180;i++){sprintf(FileName,"anim1/FujiAnim%d.bmp",i);renWin->Render();rs->Modified();bw->SetFileName(FileName);bw->Write();/* Azimuth 方向に 2 度カメラを移動 */vren->GetActiveCamera()->Azimuth(2.0);}

カメラの設定は，vtkRenderer のメソッドを使って，直接変えることもできる (33 行)。24 - 33
行で，カメラの位置を 2 度ずつ回転させ，そのたびに表示される画像を BMP 形式で保存している。上記のリストを差し替えたバージョンのサンプルプログラムも RotateLines.cxx として保存されている。
33 行の Azimuth を Elevation に変えて実行すると，表示される画像がおかしくなる。これ
は，視線方向 (カメラと焦点を結ぶ方向) と ViewUp が垂直でなくなるからである。これを回避するには，34 行目に，

vren->GetActiveCamera()->OrthogonalizeViewUp();を加えて，常に垂直に保てばよい。
焦点の位置をカメラの位置を中心として回転させるメソッドもある。これを使うと，カメラの向き (視線方向) を変化させることができる。回転は，下記のメソッドで行う。

• Pitch(double angle)• Yaw(double angle)
Pitch でカメラの向きを上下方向に，Yaw で左右方向に変化させることができる。
カメラと焦点の位置は固定で，視線方向の周りでカメラを回転させるメソッド• Roll(double angle)
もある。

2.3

凹凸つきスライス

前セクションと同じデータ (MtFuji.dem:Short 型, 75×60) を 3 次元オブジェクトとして可視化しつつ，そのオブジェクトを Window の中でマウスを使って，回転・移動・ズームイン/アウトさせることができるプログラムを作る。同時に、カラーテーブルの作り方，ポリゴンの法線ベク
トルについても知識を増やす。

26

2.3.1

サンプルプログラム 2

前章のサンプルプログラム 1 は，標高 (スカラー値) に応じて色づけしただけであったが，本章は，ポリゴンの頂点を Z 方向に移動させることで，スライスに凹凸をつけている。さらに，Window 上でマウスのボタンを押しながらマウスカーソルを動かすと，オブジェクトが移動や回転をするようにしている。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

/* BumpySlice.cxx */#include <vtkShortArray.h>#include <vtkImageData.h>#include <vtkImageDataGeometryFilter.h>#include <vtkWarpScalar.h>#include <vtkPointData.h>#include <vtkLookupTable.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define MAP_X 75#define MAP_Y 60
short sdata[MAP_X*MAP_Y];
int main( int argc, char *argv[]){
int size[3] = {MAP_X, MAP_Y, 1};
FILE *fpi;
if ((fpi =fopen("./MtFuji.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X*MAP_Y, fpi);
fclose(fpi);
vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X*MAP_Y, 1);
vtkImageData *imgData = vtkImageData::New();
imgData->SetSpacing(0.1, 0.1, 0);
imgData->SetDimensions(size);
27

42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88

imgData->GetPointData()->SetScalars(sarray);vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();
vtkImageDataGeometryFilter *igf
= vtkImageDataGeometryFilter::New();
igf->SetInput(imgData);/* 頂点をスカラー値に応じて移動 */vtkWarpScalar *warp = vtkWarpScalar::New();
warp->SetInput(igf->GetOutput());
warp->SetScaleFactor(1.0/3776.0);
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(warp->GetPolyDataOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(0.0, 3776.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor( Actor );vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600 );/* Window をインタラクティブにする */vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();sarray->Delete();lut->Delete();igf->Delete();
28

89
90
91
92
93
94
95
96
97
98

warp->Delete();Mapper->Delete();Actor->Delete();vren->Delete();iwin->Delete();trackball->Delete();renWin->Delete();return 0;}(等間隔メッシュではなく，Rectilinear のデータの場合は，vtkImageDataGeometryFilter を
vtkRectilinearGridGeometryFilter に変えればよい。)

2.3.2

サンプルプログラム 2 の簡単な説明

このプログラムは，vtkImageDataGeometryFilter と vtkPolyDataMapper の間に，頂点を移動させるフィルタを入れて，スライスに凹凸をつける処理をしている。また，Window をインタラ
クティブにしている。

Figure 2.7: サンプルプログラム 2 のパイプライン構造
可視化パイプラインの構造を Figure 2.7 に示す。
このプログラムの第一のポイントは，52 - 55 行目である。

• vtkWarpScalar クラス
を使って，頂点をスカラー値に応じて Z 方向に移動させている。55 行目で，移動の度合を調整している。
第二のポイントは，74 - 84 行目である。renWin->Render() の代わりに，74 - 84 行目を書く
だけで，Window がインタラクティブになる ( 簡単な操作方法を，Table 2.1 にまとめた)。79 80 行目をみると，このプログラムは vtkInteractorStyleTrackballActor を使用しているのがわかる
が，vtkInteractorStyle は，これ以外にも数種類用意されている (オンライン・リファレンスでご確認を)。表示は，Figure 2.8(高さ方向は，かなり強調している)。

29

Table 2.1: Window の使い方
左ボタン + マウスカーソル上下左右
中ボタン + マウスカーソル上下左右
右ボタン + マウスカーソル上下
CTRL + 左ボタン + マウスカーソル上下左右
W, S
U
E, Q

オブジェクトの回転 (Elevation, Azimuth)
オブジェクトの移動
ズームイン・ズームアウト
オブジェクトの回転 (Roll)
W:ワイヤフレーム表示，S:面での表示
ユーザー定義コマンド (最終章参照)
終了

Figure 2.8: 富士山付近

2.3.3

ポリゴンの法線ベクトル

Figure 2.8 をみると，妙にポリゴンが目立っていることに気づくと思う。なぜこのようなことが起きているのか？ まず，ポリゴンへの “陰”のつけ方を紹介する。Figure 2.9 の拡散光と鏡面反射光に注目していただきたい。環境光だけでは，映しているものが 3 次元構造を持っていても，のっぺりとしか映らない。拡散光・鏡面反射光があってこそ，暗い部分や明るい部分ができて，3 次元的に見える。その拡散光と鏡面反射光の影響は，いずれもポリゴンの法線ベクトルの方向に依存している。サンプルプログラム 2 で，表面がガタガタして見える理由は，各ポリゴンの頂点に与えられている法線ベクトルが，同じ頂点を共有する隣のポリゴンの法線ベクトルと方向が異なっているので，陰のつけ方がポリゴンごとに不連続になるためである (Figure 2.10 の真ん中の図)。このガタガタを消すには，Figure 2.10 の右の図のような法線ベクトルの与え方をすればよい (このような法線ベクトルのセットの仕方をすれば，少ないポリゴンで綺麗な曲面を表現できる。)。これを実現してくれるのが，• vtkPolyDataNormals クラスである。vtkWarpScalar から直接 vtkPolyDataMapper にパイプラインをつながず，それらの間に，この Filter を入れればよい。サンプルプログラム 2 の 56 - 58 行目を下記リストで書き換える

30

Figure 2.9: 陰のつけ方(vtkPolyDataNormals.h のインクルードも忘れずに) とガタガタがかなりおさまる (Figure 2.11)。
1
2
3
4
5

vtkPolyDataNormals *PNormals = vtkPolyDataNormals::New();
PNormals->SetInput(warp->GetPolyDataOutput());
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(PNormals->GetOutput());

書き換えられているバージョンは，BumpySlice2.cxx として保存されている。(建物などを描いている場合は，スムーズにしないほうが良い)

2.3.4

光源

光源の位置 (照らす方向) を変えて，陰のつけ方を変えてみる。vtkLight と vtkRenderer を使う。

1
2
3
4
5
6
7
8
9
10
11

vtkLight *light1 = vtkLight::New();
light1->SetPosition(0.0, -1.0, 1.0);
light1->SetFocalPoint(0.0, 0.0, 0.0);
vtkLight *light2 = vtkLight::New();
light2->SetPosition(0.0, 0.0, -1.0);
light2->SetFocalPoint(0.0, 0.0, 0.0);
vtkRenderer *vren = vtkRenderer::New();
vren->AddActor(Actor);
vren->AddLight(light1); /* Light1 をセット*/31

Figure 2.10: ポリゴンの法線ベクトル : 左の図のような曲面を 2 枚のポリゴンで描く場合。真ん中の図のように，同じ頂点でもポリゴンにより法線ベクトルの方向が異なるとガタガタが目立つ。右の図のように法線ベクトルの方向が共通だと，スムーズになる。
vren->AddLight(light2); /* Light2 をセット*/12

BumpySlice.cxx や BumpySlice2.cxx の該当箇所を上記のリストで書き換えると，陰のつき方が変わる。
上記リストでは，光源の位置と焦点を与えて，光の方向のみ(Figure 2.12) を設定している。
光源の位置の後ろにあるオブジェクトにも，陰がつく。

2.3.5

カラーテーブル

カラーテーブルは，vtkLookupTable で簡単に指定できるが，自分で作りたい場合もある。下記
リストは，vtkLookupTable に自分で作ったカラーテーブルの渡し方の一例である。サンプルプ
ログラム 2 の 45 行目を下記のものと差し替えればよい。

1
2
3

lut->SetNumberOfColors(NCOLORS);
for(int i=0;i<NCOLORS;i++)lut->SetTableValue(i, FujiColor[i]);ただし，NCOLORS は使用する色の数，float FujiColor[NCOLORS][4] にはあらかじめ RGBA
が 0.0 - 1.0 の範囲で入っていなければならない。サンプルプログラム 2 の書き換えバージョンのBumpySlice3.cxx では，可視化されたスライスが，山っぽく見えるような色を set color 関数でFujiColor[NCOLORS][4] に入れるようにしてある。BumpySlice3.cxx を実行すると，Figure 2.13の左側の図が表示される。なお，右側の絵は，Actor のプロパティを変更する

Actor->GetProperty()->SetRepresentationToWireframe();の 1 行をプログラムに加えて，面でなくワイヤフレームで表示している 2 。

W キーでのワイヤフレーム表示と違い，アクターに入っているオブジェクトごとにワイヤフレームか否かを指定できる。
232

Figure 2.11: 富士山周辺：法線ベクトル修正

Figure 2.12: 光の方向

2.3.6

補足

Figure 2.13 のワイヤフレームを見るとわかるが，この凹凸つきスライスは，四角のポリゴンで構成されている。この場合，各ポリゴンの 4 つの頂点が同一平面上にあるとは限らない。ゆえに，vtkTriangleFilter (四角などのポリゴンを三角のポリゴンに変換してくれる) をパイプラインの中に入れて，四角のポリゴンをすべて三角のポリゴンに直しておくのが望ましいと思う。

2.4

高さつき等高線

VTK には，可視化した 3 次元オブジェクトを VRML2.0 形式や Inventor 形式で保存する便利な機能がある。高さつき等高線のプログラムを見ながら使い方をみてみる。

33

Figure 2.13: 富士山周辺：(左) 色づけ変更, (右) ワイヤフレーム版

2.4.1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

サンプルプログラム 3/* UnevenCLines.cxx */#include <vtkImageData.h>#include <vtkWarpScalar.h>#include <vtkShortArray.h>#include <vtkPointData.h>#include <vtkContourFilter.h>#include <vtkLookupTable.h>#include <vtkPolyData.h>#include <vtkVRMLExporter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define MAP_X 75#define MAP_Y 60#define NCONTLINES 37
short sdata[MAP_X*MAP_Y];
int main( int argc, char *argv[] ){FILE *fpi;
if ((fpi =fopen("./MtFuji.dem", "rb")) == NULL) {puts("cannot open");
34

31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77

exit(1);}
fread(sdata, sizeof(short), MAP_X*MAP_Y, fpi);
fclose(fpi);
vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X*MAP_Y, 1);
vtkImageData *imgData = vtkImageData::New();
imgData->SetDimensions(MAP_X, MAP_Y, 1);
imgData->SetSpacing(0.1, 0.1, 0);
imgData->GetPointData()->SetScalars(sarray);vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(imgData);
contour->GenerateValues(NCONTLINES, 0.0, 3776.0);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 頂点を移動して，高さをつける*/
vtkWarpScalar *warp = vtkWarpScalar::New();
warp->SetInput(contour->GetOutput());
warp->SetScaleFactor(1.0/3776.0);
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(warp->GetPolyDataOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(0.0, 3776.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor( Actor );vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);
35

78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106

vtkInteractorStyleTrackballCamera *trackball= vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();/* 高さつき等高線を VRML2.0 のファイルとして保存 */vtkVRMLExporter * vrml = vtkVRMLExporter::New();
vrml->SetFileName("ContourLines.wrl");
vrml->SetInput(renWin);
vrml->Write();sarray->Delete();imgData->Delete();lut->Delete();warp->Delete();contour->Delete();Mapper->Delete();Actor->Delete();vren->Delete();iwin->Delete();trackball->Delete();renWin->Delete();vrml->Delete();return 0;}

このプログラムを実行すると，Figure 2.14 のように高さのついた等高線が表示される。vtkImageData
を vtkRectlinearGrid に変えれば，特に何も直さずとも，Rectilinear のデータを扱える。

2.4.2

サンプルプログラム 3 の簡単な説明

等高線に高さをつけるには，vtkContourFilter で等高線を描いた後，vtkWarpScalar で頂点を移動させるだけでよい。パイプライン構造は，Figure 2.15。
VRML2.0 形式でオブジェクトを保存するのも簡単で，vtkVRMLExporter を使えばよい (86 90 行を加えるだけである。拡張子は wrl)。vtkBMPWriter は単に画像を保存するだけであるが，こちらは 3 次元構造も保存するので，学会その他の発表時に VRML ビューアでいろいろな角度から可視化オブジェクトを聴講者に見せることができる 3 。類似のクラスとして，vtkOBJExporter(Wavefront 形式)4 , vtkIVExporter (Inventor 形式), vtkRIBExporter (RenderMan 形式) がある。
Inventor 形式 (.iv) のファイルは，OpenGLPerformer を使うと，簡単に CAVE の VR 空間に表示すことがでる。
3
4

Cosmo Player や Cortona では見れなかったので，ホームページ掲載は少々問題ありかも。
色情報が保存されませんでした。

36

Figure 2.14: 高さつき等高線

2.5

発展

GTOPO30 を使って凹凸つきスライスを作るときに，• スカラー値と RGB を対応させつつカラーテーブルを作りたい• 凹凸と色付けを別々のスカラーデータで行いたい

2.5.1

サンプルプログラム 4

前述した 2 件について，次のサンプルプログラム 4 で方法を説明する。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

/* BumpySlice4.cxx */#include <vtkShortArray.h>#include <vtkFieldData.h>#include <vtkImageData.h>#include <vtkImageDataGeometryFilter.h>#include <vtkWarpScalar.h>#include <vtkPointData.h>#include <vtkPolyData.h>#include <vtkColorTransferFunction.h>#include <vtkTriangleFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>37

Figure 2.15: サンプルプログラム 3 のパイプライン構造
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47

#define MAP_X 75#define MAP_Y 60
short sdata[MAP_X*MAP_Y];
short sdata_m[MAP_X*MAP_Y];
int main( int argc, char *argv[]){
int i, j;
int size[3] = {MAP_X, MAP_Y, 1};
FILE *fpi;
if ((fpi =fopen("./MtFuji.dem", "rb")) == NULL) {puts("cannot open");exit(1);}
fread(sdata, sizeof(short), MAP_X*MAP_Y, fpi);
fclose(fpi);/* 色づけ用のデータを作る */for(j=0;j<MAP_Y;j++){for(i=0;i<MAP_X;i++){if(i < MAP_X/2) sdata_m[i+MAP_X*j] = - sdata[i+MAP_X*j];else sdata_m[i+MAP_X*j] = sdata[i+MAP_X*j];
}
}
38

48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94

vtkShortArray *sarray = vtkShortArray::New();
sarray->SetArray(sdata, MAP_X*MAP_Y, 1);
sarray->SetName("Picard"); /* Picard と名づける *//* 色づけ用のデータ */vtkShortArray *sarray2 = vtkShortArray::New();
sarray2->SetArray(sdata_m, MAP_X*MAP_Y, 1);
sarray2->SetName("Riker"); /* Riker と名づける *//* FieldData に 2 つの配列を代入 */vtkFieldData *fd = vtkFieldData::New();
fd->AddArray(sarray);fd->AddArray(sarray2);
vtkImageData *imgData = vtkImageData::New();
imgData->SetSpacing(0.1, 0.1, 0);
imgData->SetDimensions(size);
imgData->GetPointData()->ShallowCopy(fd); /* fd をセット *//* スカラーデータとして，Picard を使用 */imgData->GetPointData()->SetActiveScalars("Picard");/* スカラー値と RGB を対応させながらカラーテーブルを作る */vtkColorTransferFunction *ctf
= vtkColorTransferFunction::New();/*------------- スカラー値,
R,G,
B ---------*/
ctf->AddRGBPoint(-3800.0, 0.0, 0.0, 1.0);
ctf->AddRGBPoint( -500.0, 0.0, 1.0, 1.0);
ctf->AddRGBPoint(
0.0, 1.0, 1.0, 1.0);
ctf->AddRGBPoint( 500.0, 1.0, 1.0, 0.0);
ctf->AddRGBPoint( 1800.0, 1.0, 0.0, 0.0);
ctf->Build();
vtkImageDataGeometryFilter *igf
= vtkImageDataGeometryFilter::New();
igf->SetInput(imgData);
vtkTriangleFilter *tri = vtkTriangleFilter::New();
tri->SetInput(igf->GetOutput());
vtkWarpScalar *warp = vtkWarpScalar::New();
warp->SetInput(tri->GetOutput());
warp->SetScaleFactor(1.0/3776.0);
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(warp->GetPolyDataOutput());
39

95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139

Mapper->SetLookupTable(ctf);/* 色づけ用のデータは，FieldData の Riker を使う */Mapper->SetScalarModeToUsePointFieldData();
Mapper->ColorByArrayComponent("Riker", 0);vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor( Actor );vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();sarray->Delete();sarray2->Delete();imgData->Delete();fd->Delete();ctf->Delete();igf->Delete();tri->Delete();warp->Delete();Mapper->Delete();Actor->Delete();vren->Delete();iwin->Delete();trackball->Delete();renWin->Delete();return 0;}40

2.5.2

カラーテーブル作成

色付けに関しては，今まで使用してきた vtkLookupTable ではなく，vtkColorTransferFunction
を使えばよい。71 - 80 行目で，スカラー値と RGB 値を対応させつつカラーテーブルを作成している。たとえば，75 行目ではスカラー値 -3800.0 に RGB = (0.0, 0.0, 1.0) を割り当てている。このようにして作ったカラーテーブルを 95 行目で vtkLookupTable と同じ関数を使って，
vtkPolyDataMapper に渡している。
(作成したカラーテーブルは，負の値で青っぽく，正の値で赤っぽくなる。)

2.5.3

別のスカラーで色づけ

Figure 2.16: Field Data：複数のデータを持つことができる
凹凸をつけるスカラーデータとは別のスカラーデータで色づけするには，vtkFieldData というクラスを使えばよい。このクラスは，配列を複数まとめて持つことができる (FIgure 2.16)。
40 - 47 : 色づけ用のデータを元データから作成，sdata m に代入。このスカラーデータで色づけをする。実際には，ファイルから温度などのデータを読み込むことになる49 - 56 : 凹凸用データを sarray に，色づけ用データを sarray2 にセット。それぞれに “Picard”，“Riker”と SetName で名前をつけている
58 - 61 : sarray と sarray2 を fd(vtkFieldData) にセット
66 : 等間隔メッシュデータに，fd(sarray, sarray2 を保持) をセット
68 - 69 : SetActiveScalars で，ポリゴンデータを作るスカラーデータとして，“Picard”を使うように指示
97 - 99 : 色付けに fd にあるデータ・“Riker”を使うように指示
パイプライン構造は Figure 2.17，実行結果は Figure 2.18。X が小さいところと大きいところで，色 (青っぽい → 赤っぽい) が変わっているのが確認できる。これは，圧力の等値面を温度で色づけする，などに応用できる。

2.6

この章のまとめ• 透視射影・正射影・カメラ設定• インタラクティブな Window• 頂点移動・ポリゴンの法線ベクトル・光源の位置• VRML 2.0 形式での保存• カラーテーブル作成• ポリゴンデータをつくるのとは別のデータで，色づけする
本章で紹介したクラスは，Table 2.2 の通り。

41

Figure 2.17: サンプルプログラム 4 のパイプライン構造

Figure 2.18: 富士山周辺：地面の凹凸は “Picard”，色づけは “Riker”を使用

42

Table 2.2: 本章で紹介したクラスおよびそのメソッド
頂点移動
法線ベクトル計算
ポリゴンの三角化
照明設定

配列
配列関係
データセット

カラーテーブル

カラーテーブル

Mapper

カメラ設定

Window 関係
VRML 他

vtkWarpScalar
vtkPolyDataNormals
vtkTriangleFilter
vtkLight: SetPosition(float, float, float): SetFocalPoint(float, float, float)
vtkFloatArray: SetName(const char *)vtkFieldData: AddArray(vtkDataArray *)vtkImageData: GetPointData()(->ShallowCopy(vtkFieldData *)): GetPointData()(->SetActiveScalars(const char *))vtkLookupTable: SetNumberOfColors(int n): SetTableValue(vtkIdType, float [4])vtkColorTransferFunction: AddRGBPoint(float value, float r, float g, float b)
vtkPolyDataMapper: SetScalarModeToUsePointFieldData(): ColorByArrayComponent(char *, int)
vtkCamera: SetPosition(double, double, double): SetFocalPoint(double, double, double): SetViewUp(double,double, double): SetClippingRange(double, double): SetViewAngle(double) /degree で指定: ParallelProjectionOn(), 〜Off(): SetParallelScale(double): Azimuth(double) / degree で指定: Elevation(double) / degree で指定: Roll(double) / degree で指定: Pitch(double) / degree で指定: Yaw(double) / degree で指定: OrthogonalizeViewUp()
vtkRenderWindowInteractor, vtkInteractiveStyleTrackball
vtkVRMLExporter, vtkIVExporter, vtkOBJExporter, vtkRIBExporter

43

Chapter 3
VTK 形式のデータ
3.1

本章の概要

VTK のデータ形式の解説をする。この形式のデータを作ると，vtk〜Reader というクラスで読み込むことができるようになり，自ら C/C++ の関数を使いデータを読み込んむ部分を可視化プ
ログラムに組み込まなくて良い。また，データの格子サイズなどをプログラムの中に書き込まなくても良いので，プログラムの汎用性が高まる。
ここでは等間隔，Rectilinear Grid について，VTK 形式のデータの作り方を紹介する。
VTK 形式のデータは，アスキーやバイナリの生データの前に，ヘッダ (格子点のサイズなど) がついているだけである。

3.2

等間隔メッシュのデータ

等間隔メッシュのデータは，vtkImageData ではなく，その派生クラスの vtkStructuredPoints を使用する。ただ，このクラスは将来廃止される予定で，すでに独自のメソッドは，まったくないといってよい。
等間隔メッシュデータ (Float 型のバイナリデータ) のヘッダ部分は，下記のようになる。

1
2
3
4
5
6
7
8
9
10
11
12
13

# vtk DataFile Version 2.0
等間隔メッシュのデータ
BINARYDATASET STRUCTURED_POINTSDIMENSIONS 50 100 50
ORIGIN 0.0 0.0 0.0
SPACING 0.1 0.1 0.1
POINT_DATA 250000
SCALARS NameOfData_S float 1
LOOKUP_TABLE default
スカラーのバイナリデータ
VECTORS NameOfData_V float
ベクトルのバイナリデータ

2 行目はコメント，8 行目の数字は，格子点の総数である。このサンプルの場合，50×100×50=250000。
9 行目の NameOfData S はデータの名前 (前章の “Picard”, “Riker”にあたる)，“1”は要素数である。Float 型ではなく，Double 型，Int 型や Unsigned Char 型であれば，9 行目の float を double，
44

int や unsigned char と書けばよい。13 行目のベクトルデータは，X Y Z 成分が x1 y1 z1 x2 y2 z2 と順番に並んでいなければならない。データがバイナリ (Fortran の書式なしにあたる) ではなくアスキー (エディタで開ける形式) で与えられている場合は，3 行目を ASCII として，11, 13 行目の部分に数値を一つずつスペースで区切りながら，必要に応じて改行しつつ書いていく。
下記のプログラムは，スカラーのバイナリデータ (sample.dat) とベクトルのバイナリデータ(vector x,dat, vector y.dat, vector z.dat) を読み込み，VTK 形式のデータ (sample.vtk) を出力する
プログラムである。69, 112 行目を入れると，Big Endian のデータが Little Endian に変換され
る 1。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36

/* makevtkdata.cxx */#include <stdio.h>#include <stdlib.h>#include <vtkByteSwap.h>
typedef float DATA_TYPE;
char data_type[] = "float";char output[]= "./sample.vtk";
char comment[]= "SampleData";
int size[3]= { 50, 100, 50 };float spacing[3] = { 0.1, 0.1, 0.1 };float origin[3] = { 0.0, 0.0, 0.0 };char s_input[]= "./scalar.dat";char scalarname[] = "ScalarName";
int scalar_comp= 1;
DATA_TYPE *s_data;char v_input[3][128];
char v_input1[] = "./vector_x.dat";char v_input2[] = "./vector_y.dat";char v_input3[] = "./vector_z.dat";char vectorname[] = "VectorName";
DATA_TYPE *v_data;
DATA_TYPE *v_data_temp;
int main(int argc, char **argv){
FILE *fp, *fp2;
int total_size;
total_size = size[0] * size[1] * size[2];if ((fp = fopen(output, "wb")) == NULL) {puts("cannot open");exit(1);

vtkByteSwap には，逆の変換を行うメソッドも用意されている。また，Double 型などでは，別なメソッドを使用しなければならない。
1

45

37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83

}
fprintf(fp, "# vtk DataFile Version 2.0\n");fprintf(fp, "%s\n", comment);
fprintf(fp, "BINARY\n");fprintf(fp, "DATASET STRUCTURED_POINTS\n");fprintf(fp, "DIMENSIONS %d %d %d\n",size[0], size[1], size[2]);fprintf(fp, "ORIGIN %f %f %f\n",origin[0], origin[1], origin[2]);fprintf(fp, "SPACING %f %f %f\n",spacing[0], spacing[1], spacing[2]);fprintf(fp, "POINT_DATA %d\n", total_size);fprintf(fp, "SCALARS %s %s %d\n",scalarname, data_type, scalar_comp);
fprintf(fp, "LOOKUP_TABLE default\n");s_data = (DATA_TYPE *) malloc(sizeof(DATA_TYPE) * total_size);if ((fp2 = fopen(s_input, "rb")) == NULL) {puts("cannot open");free(s_data);
fclose(fp2);
exit(1);}
fread(s_data, sizeof(DATA_TYPE), total_size, fp2);
fclose(fp2);/* BE から LE に変更する場合に必要 */vtkByteSwap::Swap4BERange(s_data, total_size);fwrite(s_data, sizeof(DATA_TYPE), total_size, fp);
free(s_data);
fprintf(fp, "\n");fprintf(fp, "VECTORS %s %s\n", vectorname, data_type);
v_data = (DATA_TYPE *) malloc(
sizeof(DATA_TYPE) * total_size * 3);
v_data_temp = (DATA_TYPE *) malloc(
sizeof(DATA_TYPE) * total_size);
46

84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122

sprintf(v_input[0], "%s", v_input1);
sprintf(v_input[1], "%s", v_input2);
sprintf(v_input[2], "%s", v_input3);

for(int j=0;j<3;j++){if ((fp2 = fopen(v_input[j], "rb")) == NULL) {puts("cannot open");free(v_data);
free(v_data_temp);
fclose(fp);exit(1);}
fread(v_data_temp, sizeof(DATA_TYPE), total_size, fp2);
fclose(fp2);

for(int i=0;i<total_size;i++){v_data[3*i + j] = v_data_temp[i];}}/* BE から LE に変更する場合に必要 */vtkByteSwap::Swap4BERange(v_data, total_size*3);
fwrite(v_data, sizeof(DATA_TYPE), total_size*3, fp);
free(v_data);
free(v_data_temp);
fprintf(fp, "\n");fclose(fp);}

このようにして作成した VTK 形式のデータは，vtkStructuredPointsReader を利用して読み込むことができる。

47

3.3 Rectilinear Grid のデータ
RectilinearGrid のヘッダは，下記のようになる。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

# vtk DataFile Version 2.0
Rectilinear Grid のデータ
BINARYDATASET RECTILINEAR_GRIDDIMENSIONS 50 100 50
X_COORDINATES float 50
X 座標のバイナリデータ
Y_COORDINATES float 100
Y 座標のバイナリデータ
Z_COORDINATES float 50
Z 座標のバイナリデータ
POINT_DATA 250000
SCALARS NameOfData_S float 1
LOOKUP_TABLE default
スカラーのバイナリデータ
VECTORS NameOfData_V float
ベクトルのバイナリデータ

ORIGIN と SPACING が，X,Y, Z COORDINATES に置き換えられるだけで，あとは等間隔メッシュのデータと同じである。このようにして作ったデータは，vtkRectilinearGridReader を利用して読み込むことができる。

3.4

この章のまとめ• 等間隔，RectilinearGrid の VTK 形式のデータ• エンディアンの変換
本章で紹介したクラスは，Table 3.1 の通り。

Table 3.1: 本章で紹介したクラスおよびそのメソッド
リーダー
エンディアンの変換

vtkStructuredPointsReader
vtkRectilinearGridReader
vtkByteSwap: Swap4BERange(float *, int)

48

Chapter 4
3 次元スカラーデータの可視化 その 1
4.1

本章の概要

3 次元のスカラーデータを可視化する。同時に，データの切り出し・間引きの方法にも触れる。

4.2

等値面

3 次元スカラーデータを等値面を使って可視化する。前章までの知識で，なんの問題もなく可視化できる。

4.2.1

サンプルプログラム 1

地球磁気圏のプラズマの温度データ (plasma data.vtk : 等間隔メッシュ(間隔はすべての方向で0.1)，Float 型，サイズ 256×128×128) を読み込み，等値面を作成し，その等値面を 4 方向から映した画像を同じ Window 内に表示する。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

/* Isosurface.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkContourFilter.h>#include <vtkOutlineFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkCamera.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <unistd.h>
int main( int argc, char *argv[] ){/* カメラ位置 */
49

20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66

float position[4][3] ={ {-40.0, 6.35, 6.35},{12.75, -40.0, 6.35},{12.75, 6.35, 52.75},{20.0, -20.0, 23.75}};/* カメラの上方 */float viewup[4][3] ={{0.0, 0.0, 1.0}, {0.0, 0.0, 1.0},{0.0, 1.0, 0.0}, {0.0, 0.0, 1.0}};/* ビューポート */float viewport[4][4] = {{0.0, 0.0, 0.5, 0.5},{0.5, 0.0, 1.0, 0.5},{0.0, 0.5, 0.5, 1.0},{0.5, 0.5, 1.0, 1.0}};
int i;
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName("./plasma_data.vtk");
vtkImageData *imgData;
imgData = reader->GetOutput();/* 3 次元データから等値面を生成する */vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(imgData);
contour->SetValue(0, 12.0); /* 12.0 の等値面作成 */contour->ComputeNormalsOn(); /* 法線ベクトルを計算する */vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(contour->GetOutput());
Mapper->ScalarVisibilityOff();
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);Actor->GetProperty()->SetColor(0.5, 0.75, 0.6);/* 外枠 */
vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);
vtkRenderer *vren;
vtkRenderWindow *renWin = vtkRenderWindow::New();
50

67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100

for(i=0;i<4;i++){vren = vtkRenderer::New();
vren->AddActor(Actor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vren->GetActiveCamera()->SetPosition(position[i]);vren->GetActiveCamera()->SetFocalPoint(12.75, 6.35, 6.35);
vren->GetActiveCamera()->SetViewUp(viewup[i]);vren->GetActiveCamera()->OrthogonalizeViewUp();
vren->GetActiveCamera()->SetClippingRange(30.0, 200.0);
vren->SetViewport(viewport[i]);renWin->AddRenderer( vren );vren->Delete();}renWin->SetSize( 800, 600 );for(i=0;i<10;i++){
renWin->Render();sleep(1);}reader->Delete();contour->Delete();outline->Delete();Mapper->Delete();Actor->Delete();OLMapper->Delete();OLActor->Delete();renWin->Delete();return 0;}

4.2.2

サンプルプログラム 1 の簡単な説明

前章・前々章に使用した vtkContourFilter 1 に，3 次元のデータを代入しているだけであるし，カメラの設定は前章，ViewPort については前々章で触れているので，何の問題もないはずである。
新規箇所は，vtkOutlineFilter で，データの外枠を作るところくらいである。

19 - 31 : カメラ関係の設定値
44 : SetValue メソッドで，12.0 のレベルの等値面を作成する。GenerateValues と違い，ピンポ
イントで値を指定できる
1

3 次元等間隔メッシュデータに特化した vtkMarchingCubes という等値面生成専用クラスも存在する。

51

45 : 等値「面」であるので，法線も計算する
49 : ScalarVisibilityOff で，スカラー値による色づけはしないことを指示
53 : 等値面の色を指定 (上記と併用する)
56 - 63 : vtkOutlineFilter を使いデータの外枠を作る。作り方は簡単で，データを代入するだけで (57 行)，GetOutput() で外枠のポリゴンデータを取り出せる (60 行)
実行結果は，Figure 4.1。

Figure 4.1: 等値面：4 方向から

4.3

データの切り出し : 断面切りVTK を使うと，軸に垂直な面だけでなく，任意の断面で 3 次元データのカラーコンターや等高線を作ることができる。
使用するクラスは，vtkPlane と vtkCutter である。vtkPlane で平面を指定して，vtkCutter にターゲットとなるデータとともに vtkPlane を介して平面の情報を渡すと，断面のポリゴンデー
タを切り出してくれる。
さらに，カラーバーを表示させてみる。vtkScalarBarActor を使用すれば，容易に画面上にカラーバーを表示させることができる。また，このカラーバーは，オブジェクトを回転させても，静止したままである。

4.3.1

サンプルプログラム 2

斜め断面でカラーコンターを作る。また，vtkScalarBarActor を用いてカラーバーも同時に表示
させる。

1
2

/* CrossSection.cxx */#include <vtkImageData.h>
52

3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49

#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkPlane.h>#include <vtkCutter.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkScalarBarActor.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
int main( int argc, char *argv[] ){float range[2];char datafile[]="./plasma_data.vtk";
float plane_normal[3] = {-1.0, 0.5, 0.75};
float plane_origin[3] = {10.0, 6.35, 6.35};
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 任意の断面でスライスを作る *//* 面の情報 */vtkPlane *plane = vtkPlane::New();plane->SetOrigin(plane_origin); /* 面の原点 */plane->SetNormal(plane_normal); /* 面の法線ベクトル */vtkCutter *cutter = vtkCutter::New();
cutter->SetInput(imgData);
cutter->SetCutFunction(plane);/* 等高線の場合は，下記もパイプラインに入れる */
/*
53

50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96

vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(cutter->GetOutput());
contour->GenerateValues(15, range[0], range[1]);*/vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(cutter->GetOutput());
Mapper->SetLookupTable(lut);
Mapper->SetColorModeToMapScalars();
Mapper->SetScalarRange(range[0], range[1]);vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);/* カラーバーを表示する */vtkScalarBarActor *scalarbar = vtkScalarBarActor::New();scalarbar->SetTitle("Temperature");scalarbar->SetLookupTable(lut);
scalarbar->SetOrientationToVertical();scalarbar->SetWidth(0.075);
scalarbar->SetHeight(0.75);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(Actor);
vren->AddActor(OLActor);
vren->AddActor(scalarbar);vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =
54

97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117

vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();outline->Delete();Mapper->Delete();Actor->Delete();OLMapper->Delete();OLActor->Delete();scalarbar->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}

4.3.2

サンプルプログラム 2 の簡単な説明

40 - 42 : 平面の原点 (plane origin) と平面の法線ベクトル (plane normal) を与えて平面を決定する。plane origin と plane normal は，22 - 23 行で与えられている
44 - 46 : 上記の平面データとターゲットのデータを vtkCutter に渡す
48 - 53 : vtkCutter から直接 vtkPolyDataMapper につながず，vtkContourFilter を通してつなぐと，等高線を生成できる。コメントアウトした部分を参考に
74 - 80 : カラーバーの表示。vtkScalarBarActor という名のとおり，そのまま vtkRenderer につなげばよい (85 行)
76 : タイトル
77 : カラーテーブルを入力
78 : SetOrientationToHorizontal() で，バーを横にできる
79 - 80 : 引数の単位は，ViewPort のそれと (0.0 - 1.0) 同じである
データが Rectilinear でも，(リーダーを変えるだけで) 特に何も直さずとも動作する。
vtkPlane と vtkCutter のパイプライン接続は，Figure 4.2 のようになる。実行結果は，Figure 4.3。

4.4

データの切り出し : Volume Of Interest (VOI)

データが大きくなればなるほど，データ全体ではなく，可視化したい部分 (VOI : Volume Of
Interest) を切り出して可視化する。VTK には，データを切り出すためのクラスが用意されてい55

Figure 4.2: vtkPlane と vtkCutter の接続

Figure 4.3: 斜め断面とカラーバー
るので，自分で切り出す必要はない 2 。
データの切り出しは，vtkExtractVOI というクラスで行うことができる。このクラスをデータと可視化のためのクラスの間に挿入すればよい (Figure 4.4)。切り出しには，このクラスのメ
ソッド• SetVOI(i1 , i2 , j1 , j2 , k1 , k2 )で，格子点の番号 (0 以上の整数) で指定する (Figure 4.5)。たとえば，本章のサンプルデータで，X 方向の後半を半分切り出すなら• SetVOI(128, 255, 0, 127, 0,127)
y = 64 の平面を切り出すなら，• SetVOI(0, 255, 64, 64, 0,127)とすればよい。
2

計算機に乗せるにはあまりにも大きすぎる場合は，この限りではない。

56

Figure 4.4: vtkExtractVOI の接続位置

Figure 4.5: 格子点の指定
座標軸に垂直な面でのカラーコンターや等高線を描かせたい場合は，前セクションのように vtkPlane と vtkCutter を利用する方法よりも本セクションのように，平面を切り出して vtkContourFilter を使うほうが，処理時間は短くすむであろう。

4.4.1

サンプルプログラム 3

このプログラムは，サンプルデータを読み込み，そこから vtkExtractVOI で 2 次元および 3 次元のデータを切り出して，それぞれ等高線・等値面で可視化する。さらに，2 次元の面 (切り出す部分) を，オリジナルデータ内で移動させて，そのたびに可視化した画像を連番画像として保
存する。

1
2
3
4
5

/* ExtractData.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>
57

6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52

#include <vtkExtractVOI.h>#include <vtkLookupTable.h>#include <vtkContourFilter.h>#include <vtkOutlineFilter.h>#include <vtkAppendPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkCamera.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRendererSource.h>#include <vtkBMPWriter.h>#define NCONTLINES 15#define ISOSURFACE_VALUE 11.0
int main( int argc, char *argv[] ){float range[2];char datafile[] = "./plasma_data.vtk";
char FileName[128];/* 切り出す範囲 */int extent_xy[6] = {50, 200, 20, 107, 0, 0};
int extent_yz[6] = {0, 0, 20, 107, 20, 107};
int extent_zx[6] = {0, 255, 0, 0, 0, 127};
int extent_vol[6] = {75, 150, 15, 45, 15, 65};
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 元データから 2 次元データの切り出し */vtkExtractVOI *voi_xy = vtkExtractVOI::New();
voi_xy->SetInput(imgData);
vtkExtractVOI *voi_yz = vtkExtractVOI::New();
58

53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99

voi_yz->SetInput(imgData);
vtkExtractVOI *voi_zx = vtkExtractVOI::New();
voi_zx->SetInput(imgData);/* 切り出した 2 次元データの等高線生成 */vtkContourFilter *contour_xy = vtkContourFilter::New();
contour_xy->SetInput(voi_xy->GetOutput());
contour_xy->GenerateValues(NCONTLINES, range[0], range[1]);vtkContourFilter *contour_yz = vtkContourFilter::New();
contour_yz->SetInput(voi_yz->GetOutput());
contour_yz->GenerateValues(NCONTLINES, range[0], range[1]);vtkContourFilter *contour_zx = vtkContourFilter::New();
contour_zx->SetInput(voi_zx->GetOutput());
contour_zx->GenerateValues(NCONTLINES, range[0], range[1]);/* 3 次元データを切り出して，等値面を生成 */vtkExtractVOI *voi = vtkExtractVOI::New();
voi->SetInput(imgData);voi->SetVOI(extent_vol);vtkContourFilter *isosurf = vtkContourFilter::New();
isosurf->SetInput(voi->GetOutput());
isosurf->SetValue(0, ISOSURFACE_VALUE);
isosurf->ComputeNormalsOn();/* 等高線と等値面のポリゴンデータを一つにまとめる */vtkAppendPolyData *contours = vtkAppendPolyData::New();
contours->AddInput(contour_xy->GetOutput());
contours->AddInput(contour_yz->GetOutput());
contours->AddInput(contour_zx->GetOutput());
contours->AddInput(isosurf->GetOutput());
vtkPolyDataMapper *CMapper = vtkPolyDataMapper::New();
CMapper->SetInput(contours->GetOutput());
CMapper->SetLookupTable(lut);
CMapper->SetColorModeToMapScalars();
CMapper->SetScalarRange(range[0], range[1]);vtkActor *CActor = vtkActor::New();
CActor->SetMapper(CMapper);/* データの外枠 */vtkOutlineFilter *outline_xy = vtkOutlineFilter::New();
outline_xy->SetInput(voi_xy->GetOutput());
59

100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146

vtkOutlineFilter *outline_yz = vtkOutlineFilter::New();
outline_yz->SetInput(voi_yz->GetOutput());
vtkOutlineFilter *outline_zx = vtkOutlineFilter::New();
outline_zx->SetInput(voi_zx->GetOutput());
vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkOutlineFilter *outline_voi = vtkOutlineFilter::New();
outline_voi->SetInput(voi->GetOutput());/* 外枠のポリゴンデータを一つにまとめる */vtkAppendPolyData *outlines = vtkAppendPolyData::New();
outlines->AddInput(outline->GetOutput());
outlines->AddInput(outline_voi->GetOutput());
outlines->AddInput(outline_xy->GetOutput());
outlines->AddInput(outline_yz->GetOutput());
outlines->AddInput(outline_zx->GetOutput());
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outlines->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);vtkCamera *camera = vtkCamera::New();
camera->SetPosition(-35.0, -15.0, 20.0);
camera->SetFocalPoint( 10.0, 6.0, 3.0);
camera->SetViewUp( 0.0, 0.0, 1.0);
camera->SetClippingRange( 1.0, 100.0);
camera->OrthogonalizeViewUp();
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(CActor);
vren->AddActor(OLActor);
vren->SetActiveCamera(camera);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 750, 600 );vtkRendererSource *rs = vtkRendererSource::New();
rs->SetInput(vren);rs->WholeWindowOn();
60

147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193

vtkBMPWriter *bw = vtkBMPWriter::New();
bw->SetInput(rs->GetOutput());
for(int i=0;i<128;i++){/* 切り出す場所を変えて，等高線を作り直す */voi_xy->SetVOI(extent_xy[0], extent_xy[1], extent_xy[2],extent_xy[3], i, i);voi_yz->SetVOI(i*2, i*2, extent_yz[2], extent_yz[3],extent_yz[4], extent_yz[5]);voi_zx->SetVOI(extent_zx[0], extent_zx[1], i, i,
extent_zx[4], extent_zx[5]);renWin->Render();sprintf(FileName,"bmps/Anim%d.bmp",i);renWin->Render();rs->Modified();bw->SetFileName(FileName);bw->Write();}reader->Delete();lut->Delete();contour_xy->Delete();contour_yz->Delete();contour_zx->Delete();outline->Delete();outline_xy->Delete();outline_yz->Delete();outline_zx->Delete();CActor->Delete();OLActor->Delete();CMapper->Delete();OLMapper->Delete();
voi_xy->Delete();
voi_yz->Delete();
voi_zx->Delete();voi->Delete();isosurf->Delete();contours->Delete();outlines->Delete();vren->Delete();renWin->Delete();
61

194
195
196
197
198

rs->Delete();bw->Delete();return 0;}

4.4.2

サンプルプログラム 3 の簡単な説明

このプログラムは，vtkExtractVOI クラスを用いて，データから興味のある部分を切り出し，
それを vtkContourFilter で等高線や等値面で可視化している。なお，vtkExtractVOI は，等間隔メッシュのデータ専用のクラスで，Rectilinear のデータの場合は，vtkExtractRectilinearGrid,
Structured Grid のデータの場合は，vtkExtractGrid をそれぞれ使わなければならない。また，
vtkAppendPolyData と言うクラスを使うと，ポリゴンデータを一つにまとめることができて便利である。

29 - 33 : VOI の指定範囲を表す変数
42 : スカラーの範囲 (最大値・最小値) を調べる。range[2] に結果が入る
48 - 56 : 平面データ切り出し用に vtkExtractVOI を 3 つ生成
58 - 69 : vtkExtractVOI で切り出した平面データを vtkContourFilter に渡して，等高線を生成する
71 - 79 : vtkExtractVOI で 3 次元データを切り出し，等値面を生成
81 - 86 : vtkAppendPolyData で，等高線・等値面のポリゴンデータを一つにまとめる。こうす
ると，Mapper と Actor が一つで済む
97 - 111 : 切り出した範囲それぞれに vtkOutLineFilter で外枠を作る(162 : 画像は，bmps というディレクトリ下に保存される)
VOI の範囲の指定は，描画の直前 (152 - 158 行) で行っているのにもかかわらず，i=0 でもうま
くいく。これは，VTK が描画を指示したときに初めて計算を行うからである。i=1 以降では，すべてを計算しなおすのではなく，変化したところだけ再計算する。なお，41 行の「UpDate()」は，計算を強制的に実行させるためのメソッドである。
実行結果は，Figure 4.6。

4.5

データの間引き可視化したいデータが巨大な場合，あるいはポリゴン数が膨大になって，表示するのが困難な場合は，データあるいはポリゴンを間引くとよい。その方法を見てみる。

4.5.1

サンプルプログラム 4

このプログラムは，• 間引いたデータで等値面を作成• オリジナルデータで等値面を作成• オリジナルデータで等値面を作成してポリゴンデータを間引く

62

Figure 4.6: 2 次元 & 3 次元データの切り出しの 3 種類の等値面を同じ Window に表示する。データの格子点の間引きは vtkExtractVOI の
SetSampleRate というメソッドで，ポリゴンの間引きは vtkDecimatePro というクラスでできる。
これを例えば vtkContourFilter と vtkPolyDataMapper の間に挿入すれば (Figure 4.7)，ポリゴンを間引いた等値面が表示される。

Figure 4.7: vtkDecimatePro の接続位置
1
2
3
4
5
6
7
8
9
10

/* Decimate.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkExtractVOI.h>#include <vtkDecimatePro.h>#include <vtkContourFilter.h>#include <vtkOutlineFilter.h>#include <vtkAppendPolyData.h>
63

11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57

#include <vtkPolyDataMapper.h>#include <vtkPolyDataNormals.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define ISOSURF_VALUE 12.0
int main( int argc, char *argv[] ){int size[3];char datafile[] = "./plasma_data.vtk";
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetDimensions(size);/* 1/8 に間引いたデータで等値面を生成 */vtkExtractVOI *voi_all = vtkExtractVOI::New();
voi_all->SetInput(imgData);
voi_all->SetVOI(0, size[0]-1, 0, size[1]-1, 0, size[2]-1);
voi_all->SetSampleRate(2, 2, 2); // 各方向を 1/2 に間引く
vtkContourFilter *contour1 = vtkContourFilter::New();
contour1->SetInput(voi_all->GetOutput());
contour1->SetValue(0, ISOSURF_VALUE);
contour1->ComputeNormalsOn();
vtkPolyDataMapper *Mapper1 = vtkPolyDataMapper::New();
Mapper1->SetInput(contour1->GetOutput());
Mapper1->ScalarVisibilityOff();
vtkActor *Actor1 = vtkActor::New();
Actor1->SetMapper(Mapper1);
Actor1->GetProperty()->SetColor(0.75, 0.35, 0.35);/* 元データで等値面を生成 */vtkContourFilter *contour2 = vtkContourFilter::New();
contour2->SetInput(imgData);
64

58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104

contour2->SetValue(0, ISOSURF_VALUE);
vtkPolyDataNormals *surf2 = vtkPolyDataNormals::New();
surf2->SetInput(contour2->GetOutput());
vtkPolyDataMapper *Mapper2 = vtkPolyDataMapper::New();
Mapper2->SetInput(surf2->GetOutput());
Mapper2->ScalarVisibilityOff();
vtkActor *Actor2 = vtkActor::New();
Actor2->SetMapper(Mapper2);
Actor2->GetProperty()->SetColor(0.35, 0.75, 0.35);/* 元データで作った等値面のポリゴンを 7/8 間引く */vtkDecimatePro *decimate = vtkDecimatePro::New();
decimate->SetInput(contour2->GetOutput());
decimate->PreserveTopologyOn();
decimate->SetTargetReduction(7.0/8.0); // 7/8 間引く
vtkPolyDataNormals *surf3 = vtkPolyDataNormals::New();
surf3->SetInput(decimate->GetOutput());
vtkPolyDataMapper *Mapper3 = vtkPolyDataMapper::New();
Mapper3->SetInput(surf3->GetOutput());
Mapper3->ScalarVisibilityOff();
vtkActor *Actor3 = vtkActor::New();
Actor3->SetMapper(Mapper3);
Actor3->GetProperty()->SetColor(0.35, 0.35, 0.75);/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLAMapper = vtkPolyDataMapper::New();
OLAMapper->SetInput(outline->GetOutput());
vtkActor *OLAActor = vtkActor::New();
OLAActor->SetMapper(OLAMapper);/* 間引いたデータから作った等値面 */vtkRenderer *vren1= vtkRenderer::New();
vren1->AddActor(Actor1);
vren1->AddActor(OLAActor);
vren1->SetBackground( 0.0, 0.0, 0.0 );vren1->SetViewport( 0.0, 0.0, 0.333, 1.0);

65

105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151

/* 元データから作った等値面 */vtkRenderer *vren2= vtkRenderer::New();
vren2->AddActor(Actor2);
vren2->AddActor(OLAActor);
vren2->SetBackground( 0.0, 0.0, 0.0 );vren2->SetViewport( 0.333, 0.0, 0.666, 1.0 );/* 元データから作った等値面のポリゴンを間引いたもの */vtkRenderer *vren3= vtkRenderer::New();
vren3->AddActor(Actor3);
vren3->AddActor(OLAActor);
vren3->SetBackground( 0.0, 0.0, 0.0 );vren3->SetViewport( 0.666, 0.0, 1.0, 1.0);
vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren1 );renWin->AddRenderer( vren2 );renWin->AddRenderer( vren3 );renWin->SetSize( 1100, 500 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();surf1->Delete();surf2->Delete();surf3->Delete();contour1->Delete();contour2->Delete();decimate->Delete();Mapper1->Delete();Actor1->Delete();Mapper2->Delete();Actor2->Delete();Mapper3->Delete();Actor3->Delete();vren1->Delete();vren2->Delete();vren3->Delete();
66

152
153
154
155
156
157

iwin->Delete();trackball->Delete();renWin->Delete();return 0;}

Figure 4.8: 3 種類の等値面 ( (左):間引きデータによる等値面，(中):オリジナルデータによる等値
面，(右):中の結果のポリゴンを間引いたもの)

4.5.2

サンプルプログラム 4 の簡単な説明

36 - 40 : データ全体を取り出しているが，40 行目で，格子点を一つ置きに取るように指定している。だから，取り出したデータの格子点は，オリジナルの 1/8 である。SampleRate を 3 にすると，3 個に 1 個だけ格子点をとることになる 3 。このパラメータは，i, j, k 独立に与えることができる
71 - 75 : 55 - 58 行目で (オリジナルデータで) 作った等値面のポリゴンデータを間引きする。75
行目の SetTragetReduction で，ポリゴン数を 7/8 消去 (1/8 にする) するよう指定している。なお，
このクラスは，三角形のポリゴンしか受け付けないので，四角形ポリゴン等が入っている場合
は，vtkTriangleFilter などで，あらかじめ三角形ポリゴンに直しておかなければならない
119 - 134 : 3 つとも別々に回転させることができる
実行結果は，Figure 4.8。実行結果を見ると，元のデータを間引きして等値面を作ると多数の穴が開いているが (左)，オリジナルデータから作った等値面のポリゴンデータを間引きしたものには，そのようなものはない。ただし，計算時間は 3 つの中で最もかかる。
前セクションの内容と組み合わせると，データ全体を間引いて，切り出した部分はオリジナ
ルでそれぞれ等値面を作り，同じ Window に表示させることができる。
指定した範囲の境界の格子点を必ず拾うとは限らない。vtkExtractRectilinearGrid と vtkExtractGrid には，IncludeBoundaryOn() という，境界の格子点を必ず拾わせるメソッドがある。
367

4.5.3

補足 : vtkImageShrink3D

このセクションで紹介した vtkExtractVOI によるデータの間引きは，平均操作はせずに，単純に格子点を飛ばしながら拾っていく。しかし，等間隔メッシュのデータに限り，格子点に与えられているデータの平均を取りつつ，縮小するクラスがある。それが，vtkImageShrink3D である。サンプルプログラムは用意していないので，リファレンスを参照して欲しい。

4.6

この章のまとめ• データの外枠・カラーバー表示• 斜め断面• データの切り出し・間引き• ポリゴンの間引き
本章で紹介したクラスは，Table 4.1 の通り。

Table 4.1: 本章で紹介したクラスおよびそのメソッド
リーダー
データセット
外枠表示
断面表示

断面表示
カラーバー

ポリゴンをまとめる
データ切り出し・間引き

データ間引きポリゴンの間引き

vtkStructuredPointsReader: SetFileName(const char *)vtkImageData: GetScalarRange(float [2])vtkOutlineFilter
vtkPlane: SetOrigin(float o[3]): SetNormal(float n[3])vtkCutter: SetCutFunction(vtkImplicitFunction *)vtkScalarBarActor: SetTitle(const char *): SetOrientationTo〜(), Vertical, Horizontal: SetLookupTable(vtkScalarsToColors *): SetWidth(float): SetHeight(float)
vtkAppendPolyData
vtkExtractVOI (vtkExtractRectilinearGrid, vtkExtractGrid): SetVOI(int [6]): SetSampleRate(int, int, int)
vtkImageShrink3D
vtkDecimatePro: SetTargetReduction(float)

68Chapter 5
3 次元スカラーデータの可視化 その 2
5.1

本章の概要

ボリューム・レンダリングでスカラーデータを可視化する。また，テクスチャ・マッピングについても触れる。

5.2

ボリューム・レンダリング : レイ・キャスティング法

ボリューム・レンダリングとは，(スカラー) データ全体を半透明な形で可視化する手法である(等値面は，指定されたレベルのみ)。実行方法は，Figure 5.1 にあるとおりシンプルである。視点から視線を延ばして行き，データに入ったところで，サンプリング点における色をその点の不透明度 (伝達関数によって，スカラー値を色と不透明度に変換する) を重みとして加算して行き，投影面のピクセルの色を決める。この方法は，レイ・キャスティング法と呼ばれる。この方法はシンプルではあるが，計算量が多く，可視化像が出力されるまでに一般的には時間がかかる。
VTK には，このレイ・キャスティング法による方法と，後述するテクスチャ・マッピングを利用した方法が実装されている 1 。ただし，いずれの場合も，等間隔メッシュのデータで，Short型あるいは Unsigned Char 型のデータのみのサポートとなっている。

5.2.1

サンプルプログラム 1

このサンプル・プログラムは，レイ・キャスティング法によって，ボリューム・レンダリングを行うプログラムである。サンプルデータは，球状トカマクの圧力 (p128.vtk : Unsigned Char 型(0 - 255)，サイズ 128×128×128) である。

1
2
3
4
5
6
7
8

/* VolumeRendering_rc.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkPiecewiseFunction.h>#include <vtkColorTransferFunction.h>#include <vtkVolumeProperty.h>1

ボリュームプロという，ボリューム・レンダリング専用ハードウェアを利用するクラスも用意されている。

69

Figure 5.1: ボリューム・レンダリングの実行方法 (レイ・キャスティング法)
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36

#include <vtkVolumeRayCastCompositeFunction.h>#include <vtkVolumeRayCastMapper.h>#include <vtkVolume.h>#include <vtkCubeAxesActor2D.h>#include <vtkOutlineFilter.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkCamera.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
int main( int argc, char *argv[] ){char datafile[] = "./p128.vtk";
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();/* 色の伝達関数を作成 */vtkColorTransferFunction *tf4color= vtkColorTransferFunction::New();
tf4color->AddHSVPoint(0, 0.7,
70

1.0, 1.0);

37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83

tf4color->AddHSVPoint(63, 0.7*(1.0-63.0/255.0), 1.0, 1.0);
tf4color->AddHSVPoint(127,0.7*(1.0-127.0/255.0),1.0, 1.0);
tf4color->AddHSVPoint(190,0.7*(1.0-190.0/255.0),1.0, 1.0);
tf4color->AddHSVPoint(255,0.0,
1.0, 1.0);/* 不透明度の伝達関数を作成 */vtkPiecewiseFunction *tf4opacity= vtkPiecewiseFunction::New();
tf4opacity->AddPoint(0,0.0);
tf4opacity->AddPoint(1,0.01);
tf4opacity->AddPoint(128, 0.02);
tf4opacity->AddPoint(200, 0.02);
tf4opacity->AddPoint(255, 0.05);

//
//
//
//

/* 伝達関数を vtkVolumeProperty に代入 */vtkVolumeProperty *vp = vtkVolumeProperty::New();
vp->SetColor(tf4color);
vp->SetScalarOpacity(tf4opacity);
vp->SetInterpolationTypeToLinear();vp->ShadeOn(); /* 陰つけができる */vp->SetAmbient(0.7);
vp->SetDiffuse(0.3);
vp->SetSpecular(0.3);/* ボリュームレンダリング (Ray Casting) 用 Mapper */vtkVolumeRayCastMapper *vMapper= vtkVolumeRayCastMapper::New();
vtkVolumeRayCastCompositeFunction *cfunction= vtkVolumeRayCastCompositeFunction::New();
vMapper->SetVolumeRayCastFunction(cfunction);vMapper->SetSampleDistance(0.5);
vMapper->SetInput(imgData);

//
//
//
//

/* データの一部切り取りなどができる */vMapper->SetCroppingRegionPlanes(0.0, 90.0,
0.0, 127.0, 0.0, 127.0);
vMapper->SetCroppingRegionFlagsToSubVolume();vMapper->CroppingOn();/* ボリュームレンダリング用の Actor */vtkVolume *Volume = vtkVolume::New();
Volume->SetMapper(vMapper);Volume->SetProperty(vp);/* データの外枠 */
71

84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130

vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);/* 軸の生成 */vtkCubeAxesActor2D *Axes = vtkCubeAxesActor2D::New();
Axes->SetInput(imgData);
vtkCamera *camera = vtkCamera::New();
camera->SetPosition(-200.0, 63.5, 200.0);
camera->SetFocalPoint(63.5, 63.5, 63.5);
camera->SetViewUp(0.0, 1.0, 0.0);
camera->OrthogonalizeViewUp();
camera->SetClippingRange(30.0, 2000.0);
vtkRenderer *vren= vtkRenderer::New();
vren->AddVolume(Volume);
vren->AddActor (OLActor);
vren->AddActor (Axes);
vren->SetActiveCamera (camera);
vren->SetBackground( 0.0, 0.0, 0.0 );Axes->SetCamera(vren->GetActiveCamera()); /* これも忘れずに */vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 400, 300 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();tf4color->Delete();tf4opacity->Delete();
72

131
132
133
134
135
136
137
138
139
140
141
142
143
144
145

vp->Delete();cfunction->Delete();vMapper->Delete();Volume->Delete();outline->Delete();OLActor->Delete();Axes->Delete();iwin->Delete();trackball->Delete();vren->Delete();camera->Delete();renWin->Delete();return 0;}

Figure 5.2: 球状トカマクの圧力分布：右図は，一部を Cropping している

5.2.2

伝達関数

ボリューム・レンダリングを実行するには，伝達関数を定義しなければならない。伝達関数とは，スカラー値と色・不透明度の関係を表した関数である。伝達関数という名称は，本章で始めて出たが，カラーコンターや等高線などで，第 1 章から使っている。スカラー値が小さいところは青っぽく，大きいところは赤っぽくするよう vtkLookupTable で設定してきたが，まさにそれである。ただし，ボリューム・レンダリングには，色のほかに，スカラー値と不透明度 (アルファ値:RGBA の A) の関係も決めなければならない。不透明度の設定しだいで，得られる画像はかなり違ったものとなる。はじめは，消したいスカラー値の不透明度を 0.0 に近く，強調したい (重要な) スカラー値の不透明度を大きく (といっても 0.3 程度でよい) 設定してみるとよい。
VTK では，色の伝達関数を vtkColorTransferFunction で，不透明度を vtkPiecewiseFunctionで設定する。

73

5.2.3

サンプルプログラム 1 の簡単な説明

前章までのサンプルプログラムと大きく違うところは，• vtkPolyDataMapper → vtkVolumeRayCastMapper,• vtkActor → vtkVolume,となっているところである。ポリゴンデータを作って映す，という形ではないので，Mapper とActor も別物が用意されている (Figure 5.3)。
ボリューム・レンダリングに必要なパラメータとして，先述した伝達関数がある。伝達関数で，スカラー値と色・不透明度の関係を決める。色の伝達関数は 33 - 40 行で，不透明度の伝達関数は 42 - 49 行で作成している。そしてそれらの伝達関数を vtkVolumeProperty を介して
vtkVolume に渡している。また vtkCubeAxesActor2D を使って，軸を描いている。表示される
画像は，Figure 5.2

33 - 40 : vtkColorTransferFunction を使って，スカラー値と色の関係を決めている。前々章と違うのは，RGB ではなく HSV で色を指定しているところである。ちなみに 36 と 40 行だけでは，色の変化は青 → 紫 → 赤となる。理由は RGB の値が補間されるからである
42 - 49 : vtkPiecewiseFunction を使って，不透明度の伝達関数を作成している。AddPoint メソッ
ドでスカラー値と不透明度の関係を決めている (一つ目の引数がスカラー値，二つ目が不透明度)。消したいスカラー値の不透明度を 0.0 にして，強調したいスカラー値の不透明度を大きくする
53 - 55 : 伝達関数を vtkVolumeProperty にセットして，81 行で vtkVolume(Actor にあたる) にわ
たしている。
69 : レイの刻み幅を指定している。刻みが小さいほど絵は綺麗になるが，処理時間が大きくなる
93 - 95, 111 : これだけで，軸を出せる
VTK のボリューム・レンダリングは，陰つけやデータの切り出しも簡単にできる。例えば，
vtkVolumeProperty で ShadeOn() とするだけで陰がつくし (56 行)，73 - 76 行を入れるだけで，Figure 5.2 の右図のように一部を消すことができる (73 - 74 行での数値は，vtkExtractVOI のよ
うに格子番号ではなく，座標値を入れる)。
このボリューム・レンダリングは，計算負荷が高いので，PC で実行するときは，インタラ
クティブな Window で見るのではなく，カメラ設定で見たい角度からの画像を保存するという使い方の方が，良いかもしれない。

5.3

カラーコンター再び : テクスチャ・マッピングの利用

第 1 章で，カラーコンターの表示方法を紹介したが，データが等間隔の場合，別な方法もある。
おそらく，以下に紹介する方法の方が，(表示された後は) 高速である。

5.3.1

テクスチャ・マッピング

OpenGL などで建物を描くとき，単なる直方体ではなく，窓やドアを描きたい。それらを描く
とき，ポリゴンの頂点に色を指定して描いていくというのは，途方も無い作業であるし，ポリゴン数が膨大になり描画に多大な時間を要する。これを解決してくれるのが，テクスチャ・マッ
ピングである。これは，簡単に言うと，ポリゴンに画像を貼り付ける機能である。建物を描く
74

Figure 5.3: サンプルプログラム 1 のパイプライン構造
とき，建物の外形をポリゴンで作って，そのポリゴンにデジタルカメラで撮った写真 (画像) を貼り付ければ，簡単に窓やドアをつけれるし，質感も出せる (Figure 5.4)。この機能を使っても，
カラーコンターが可能なことが理解していただけると思う。
サンプルプログラム 2 を引数を “1”にして実行すると，スライス (y=127) のテクスチャを使ったカラーコンターが表示される。

Figure 5.4: 建物 (テクスチャマップされる前後)

5.3.2

サンプルプログラム 2

このプログラムは，テクスチャマッピングを使ったカラーコンターを表示するプログラムである。とりあえず，引数を “1”にして実行していただきたい (データは，前章とおなじ地球磁気圏

75

の温度データである)。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45

/* Colorcontour_again.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkImageDataGeometryFilter.h>#include <vtkPointData.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkImageActor.h>#include <vtkImageMapToColors.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define MAX_SLICES 120
int main( int argc, char *argv[] ){int size[3];float range[2];char datafile[] = "./plasma_data.vtk";
float color_Table[4];int i, nslices, dy;
float fdy;
if(argc != 2){
printf("Input the Number of slices (1-%d)\n", MAX_SLICES);
return 1;}
nslices = atoi(argv[1]);if(nslices < 1 || nslices > MAX_SLICES){
printf("Set 1 - %d\n", MAX_SLICES);
return 1;}
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;

76

46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92

imgData = reader->GetOutput();
imgData->Update();imgData->GetDimensions(size);
imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->SetNumberOfTableValues(256);
lut->SetRange(range[0], range[1]);lut->Build();/* lut の不透明度 (RGBA の A) のみ変更する */for(i=0;i<256;i++){lut->GetTableValue(i, color_Table);color_Table[3] = 1.0/nslices;lut->SetTableValue(i, color_Table);}/* スカラー値 -> 色 */vtkImageMapToColors *colordata;/* テクスチャを使用する場合の Actor */vtkImageActor *IActor;
if(nslices == 1) fdy = 0.0;
else fdy = (size[1]-1.0)/(nslices-1.0);
vtkRenderer *vren= vtkRenderer::New();
for(i=0;i<nslices;i++){colordata = vtkImageMapToColors::New();colordata->SetInput(imgData);
colordata->SetLookupTable(lut);
IActor = vtkImageActor::New();
IActor->SetInput(colordata->GetOutput());
dy = (int)(fdy* i);
IActor->SetDisplayExtent(0, size[0]-1,
size[1] -1 - dy, size[1]-1-dy, 0, size[2]-1);
vren->AddActor(IActor);
IActor->Delete();colordata->Delete();}/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
77

93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130

outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 800, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}

5.3.3

サンプルプログラム 2 の簡単な説明

引数を “1”にして実行したときの表示が，Figure 5.5 である。
このサンプルプログラムは，vtkImageMapToColors というクラスで，データを色情報に変換，vtkImageActor というテクスチャ・マッピングされたポリゴンを表示する Actor でカラーコ
ンター (画像) を表示している。表示するスライスは，83 - 84 行目で，指定している。

78

Figure 5.5: テクスチャ・マッピングを利用したカラーコンター

5.3.4

不透明度の設定

Figure 5.6: テクスチャ・マッピングを利用したボリューム・レンダリング
テクスチャ・マッピングでは，画像のピクセルに不透明度を設定することも可能である。さらに，半透明にしたスライスを前後に並べることで，混合された絵を表示することも可能である。これを利用すると，簡易的にボリューム・レンダリングができる。つまり，3 次元データをスライスしてそれぞれを画像化，不透明度を与えて前後に並べる (Figure 5.6)。
サンプルプログラム 2 で，これを実験できる。引数は，スライスの数である。また，不透明
度は，60 行目で (1.0/スライス数) になるように設定している。スライス数が上がるにつれ，ボリューム・レンダリングらしくなっていくのがわかる (Figure 5.7)。このプログラムは，Y(の負方向から見た) 方向からしか混合された絵を表示できないが，X, Z 方向用のテクスチャセットも作って，視点の位置に応じて切り替えるようにすれば，全方向に対応できる。これが，次に紹介する (2 次元) テクスチャ・マッピングを利用したボリューム・レンダリングである。

79

Figure 5.7: テクスチャ・マッピングを利用したボリューム・レンダリング実験; 左上からスライ
ス数が，2, 4, 8, 16, 32, 64 枚

5.4

ボリューム・レンダリング : テクスチャ・マッピングの利用

テクスチャ・マッピングを利用したボリュームレンダリングを行うためのクラスが，VTK には用意されている。一般的には，レイ・キャスティング法よりも高速である。

5.4.1

サンプルプログラム 3

テクスチャ・マッピング利用のボリュームレンダリングを使う場合は，レイ・キャスティングバージョンのプログラムの vtkVolumeRayCastMapper を vtkVolumeTextureMapper2D に変えるだけで良い。さらに，このプログラムで，Float 型 (Unsigned Char 型以外なら Float でなくても良い) のデータを読み込み，それを VTK のフィルタで Unsigned Char 型に直す方法も紹介する。
また，球体のポリゴンデータ生成やビルボード処理されたテキストの表示方法も紹介する。デー
タは，引き続き地球磁気圏の温度 (ただし間隔は 1.0) である。

1
2
3
4
5
6
7
8
9
10
11
12
13

/* VolumeRendering_2dt.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkSphereSource.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkImageShiftScale.h>#include <vtkPiecewiseFunction.h>#include <vtkColorTransferFunction.h>#include <vtkVolumeProperty.h>#include <vtkVolumeTextureMapper2D.h>
80

14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60

#include <vtkVolume.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkCamera.h>#include <vtkScalarBarActor.h>#include <vtkVectorText.h>#include <vtkFollower.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
int main( int argc, char *argv[] ){int i;
float range[2];char datafile[] = "./tempe200.vtk";
float shift, scale;float opacity_table[256];
float color_table[3*256], temp_color[4];float earth_position[3] = {100.5, 63.5, 63.5};
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);/* Float 型のデータを Unsigned Char 型に変換する */shift = -range[0];scale = 255.0/(range[1]-range[0]);vtkImageShiftScale *f2uc = vtkImageShiftScale::New();
f2uc->SetShift(shift); // 最小値を 0 にする
f2uc->SetScale(scale); // 値域を 0 - 255 にする
f2uc->SetOutputScalarTypeToUnsignedChar();f2uc->SetInput(imgData);/* vtkLookupTable を利用して色の伝達関数を作る */vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->SetNumberOfTableValues(256);
81

61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107

lut->SetRange(range);
lut->Build();/* vtkLookupTable から色を取り出す */for(i=0;i<256;i++){lut->GetTableValue(i, temp_color);
color_table[i*3]= temp_color[0];color_table[i*3+1] = temp_color[1];color_table[i*3+2] = temp_color[2];}vtkColorTransferFunction *tf4color= vtkColorTransferFunction::New();
tf4color->SetColorSpaceToRGB();
tf4color->BuildFunctionFromTable(0, 255, 256, color_table);/* 不透明度の伝達関数作成 */for(i=0;i<256;i++){if(i > 187) opacity_table[i] = i/255.0 * 0.025;
else opacity_table[i] = 0.0;}
vtkPiecewiseFunction *tf4opacity= vtkPiecewiseFunction::New();
tf4opacity->BuildFunctionFromTable(
0, 255, 256, opacity_table, 1);/* vtkVolumeProperty に伝達関数を代入 */vtkVolumeProperty *vp = vtkVolumeProperty::New();
vp->SetColor(tf4color);
vp->SetScalarOpacity(tf4opacity);
vp->SetInterpolationTypeToLinear();/* ボリュームレンダリング (2D Texture Technique) 用 Mapper*/vtkVolumeTextureMapper2D *vMapper= vtkVolumeTextureMapper2D::New();
vMapper->SetInput(f2uc->GetOutput());
vMapper->SetMaximumNumberOfPlanes(256);vMapper->SetTargetTextureSize(256, 128);
vtkVolume *Volume = vtkVolume::New();
Volume->SetMapper(vMapper);Volume->SetProperty(vp);/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
82

108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154

vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);/* 青い球体を地球の位置に表示する */vtkSphereSource *earth = vtkSphereSource::New();
earth->SetRadius(3.3);
vtkPolyDataMapper *EMapper = vtkPolyDataMapper::New();
EMapper->SetInput(earth->GetOutput());
vtkActor *EActor = vtkActor::New();
EActor->SetMapper(EMapper);
EActor->SetPosition(earth_position);
EActor->GetProperty()->SetColor(0.0, 0.0, 1.0);/* ビルボード処理されたテキストを表示する *//* EARTH という文字のポリゴンデータ生成 */vtkVectorText *EText = vtkVectorText::New();
EText->SetText("
EARTH");
vtkPolyDataMapper *TMapper = vtkPolyDataMapper::New();
TMapper->SetInput(EText->GetOutput());/* ビルボード処理用 Actor */vtkFollower *TActor = vtkFollower::New();
TActor->SetMapper(TMapper);TActor->SetScale(5.0, 5.0, 5.0);
TActor->AddPosition(earth_position);/* カラーバー */vtkScalarBarActor *scalarbar = vtkScalarBarActor::New();scalarbar->SetTitle("Temperature");scalarbar->SetLookupTable(lut);
scalarbar->SetOrientationToVertical();scalarbar->SetWidth(0.075);
scalarbar->SetHeight(0.75);
vtkCamera *camera = vtkCamera::New();
camera->SetPosition(-200.0, -200.0, 237.5);
camera->SetFocalPoint(127.5, 63.5, 63.5);
camera->SetViewUp(0.0, 0.0, 1.0);
camera->OrthogonalizeViewUp();
camera->SetClippingRange(30.0, 2000.0);
83

155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201

vtkRenderer *vren= vtkRenderer::New();
vren->AddVolume(Volume);
vren->AddActor( EActor);
vren->AddActor( OLActor);
vren->AddActor( TActor);
vren->AddActor( scalarbar);vren->SetActiveCamera (camera);
vren->SetBackground( 0.0, 0.0, 0.0 );TActor->SetCamera(vren->GetActiveCamera()); /* これも忘れずに */vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 800, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();f2uc->Delete();lut->Delete();tf4color->Delete();tf4opacity->Delete();vp->Delete();vMapper->Delete();Volume->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();earth->Delete();EMapper->Delete();EActor->Delete();EText->Delete();TMapper->Delete();TActor->Delete();scalarbar->Delete();iwin->Delete();trackball->Delete();
84

202
203
204
205
206
207

vren->Delete();camera->Delete();renWin->Delete();return 0;}

5.4.2

サンプルプログラム 3 の簡単な説明

ボリューム・レンダリングに関しては，Mapper がレイ・キャスティングバージョンと違うだけである。(新規事項ではないが，色の伝達関数作成を工夫した (57 - 75 行))

95 - 96 : vtkVolumeRayCastMapper ではなく，vtkVolumeTextureMapper2D を使う
97 : Unsigned Char 型に直したデータをセット
98 : スライス数の最大値設定
99 : テクスチャの画像解像度設定。Default では，512×512。ターゲットのデータサイズは，256×128×128 であるので，256×128 に設定している。なお，このテクスチャのサイズは，2 のべき乗 (おそらく 64 以上の) でなければならない。
以上が，レイ・キャスティングバージョンからの主な変更点であるが，このプログラムには以下の新しい事項が含まれている。

• Float 型から Unsigned Char 型への変換• 球体のポリゴン生成• ビルボード処理されたテキスト
以下，順に説明する。
Float 型から Unsigned Char 型への変換 (47 - 55 行)
vtkImageShiftScale を使うと良い (クラス名からわかるが，等間隔データ用である)。
48, 52 : 最小値が 0.0 になるように，データのスカラー値を “Shift” (range[0] には，データの最小値が代入されている)。データの値域は，[0.0, Max - Min] となっている
49, 53 : 値が，[0.0, 255.0] になるように Scale を変更
54 : 出力データを Unsigned Char 型に指定
球体のポリゴン生成 (115 - 117 行)
vtkSphereSource で球体のポリゴン (正確には球面のポリゴン) が生成できるので，それを vtkPlyDataMapper, vtkActor とつなげれば，簡単に球体を映すことができる。さらに，VTK には球体だ
けでなく，立方体，コーン，矢印などのポリゴンを生成するクラスも用意されている (Figure 9.3
参照。左上から vtkConeSource, vtkArrowSource, vtkCubeSource, vtkSphereSource, vtkCylinderSource, vtkDiskSource を使ってポリゴンデータを生成した)。本プログラムでは，これを利用して，地球の位置 (正確でないかもしれません) に青い球体を表示するようにしている。次章のベクトル場の可視化では，矢印を使う予定である
ビルボード処理されたテキスト (127 - 139, 165 行)
常に画面に平行な面にテキストを表示させる方法がある。vtkVectorText, vtkPolyDataMapper,

85

Figure 5.8: ポリゴンの生成
vtkFollower(vtkActor の代わり) を利用すればよい。vtkVectorText は，テキストのポリゴンデー
タを生成するだけなので，vtkFollower ではなく vtkActor でも表示できるが，それだと，カメラの視線方向が変わっても，テキストは回転しない。ちなみに，vtkFollower は，vtkActor の派生クラス。本プログラムでは，青い球体の横に “EARTH”と表示させている

このプログラムを実行すると，Figure 5.9 のような画像が表示される。

Figure 5.9: 地球磁気圏の温度分布

5.5

この章のまとめ• レイ・キャスティング法によるボリュームレンダリング86• テクスチャマッピングを利用したカラーコンター• テクスチャマッピングを利用したボリューム・レンダリング• Float 型から Unsigned Char 型への変換• 球体 (面) の表示，軸の表示，ビルボード処理されたテキスト表示
本章で紹介したクラスは，Table 5.1 の通り。

Table 5.1: 本章で紹介したクラスおよびそのメソッド
データセット
データ変換

カラーテーブル関係

カラーテーブル関係
球体その他の表示
軸
テクスチャ用 Actor
テキスト表示ボリュームのプロパティ(VR 用 Mapper)(VR 用 Actor)

vtkImageData: GetDimensions(int [3])vtkImageShiftScale: SetShift(float): SetScale(float): SetOutputScalarTypeTo〜()
vtkImageMapToColors
vtkPiecewiseFunction: AddPoint(float, float): BuildFunctionFromTable(float min, float max, int size, float *, 1)
vtkColorTransferFunction: BuildFunctionFromTable(float min, float max, int size, float *)vtkSphereSource, vtkConeSource, vtkArrowSource
vtkCubeSource, vtkCylinderSource, vtkDiskSource
vtkCubeAxesActor
vtkImageActor
vtkVectorText, vtkFollower
vtkVolumeProperty: SetColor(vtkColorTransferFunction *): SetOpacity(vtkPiecewiseFunction *): SetInterpolationTypeToLinear(): ShadeOn()
vtkVolumeRayCastMapper
vtkVolumeTextureMapper2D———————: SetSampleDistance(float) / RayCast: SetMaximumNumberOfPlanes(int ) / Texture: SetTargetTextureSize(int, int) / Texture: SetCroppingRegionPlanes(float, float, float, float, float, float): SetCroppingRegionFlagsToSubVolume(): CroppingOn()
vtkVolume: SetProperty(vtkVolumeProperty *)87

Chapter 6
3 次元ベクトルデータの可視化 その 1
6.1

本章の概要

ベクトルデータを矢印やスライスを使って可視化する。

6.2

HedgeHog(針鼠) によるベクトルデータの可視化

ここでは，VTK 用のベクトルデータの生成などを中心に話を進める。本章のサンプルデータ
は，私が自作したものである。Float 型・サイズ 300×300×100 で，x, y, z 成分が別々のファイル ( vector x.dat, vector y.dat, vector z.dat) になっている。また，ベクトルの絶対値のファイル(vector abs.dat) も用意した。
このサンプルプログラムは，データ内からランダムに点を選び出して，そこでのベクトルを基に線分や矢印を描くことでベクトル場を可視化する。ベクトル場の可視化で良く用いられる手法である。

6.2.1
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

サンプルプログラム 1/* viz_vect.cxx */#include <vtkImageData.h>#include <vtkFloatArray.h>#include <vtkPointData.h>#include <vtkMaskPoints.h>#include <vtkHedgeHog.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define SIZE_X 300
88

18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64

#define SIZE_Y 300#define SIZE_Z 100
float vector[SIZE_X*SIZE_Y*SIZE_Z][3];float vector_abs[SIZE_X*SIZE_Y*SIZE_Z];
char datafile_x[]= "./vector_x.dat";char datafile_y[]= "./vector_y.dat";char datafile_z[]= "./vector_z.dat";char datafile_abs[] = "./vector_abs.dat";
float *vector_x;
float *vector_y;
float *vector_z;
void load_data();
int main( int argc, char *argv[] ){int size[3] = {SIZE_X, SIZE_Y, SIZE_Z};
float range[2];load_data(); /* ベクトル・スカラーデータの読み込み */vtkFloatArray *varray = vtkFloatArray::New();
varray->SetNumberOfComponents(3); /* 要素数は 3 *//* 格子点の総数 */varray->SetNumberOfTuples(SIZE_X* SIZE_Y* SIZE_Z);
varray->SetArray(vector[0], SIZE_X* SIZE_Y* SIZE_Z* 3, 1);
vtkFloatArray *sarray = vtkFloatArray::New();
sarray->SetArray(vector_abs, SIZE_X* SIZE_Y* SIZE_Z, 1);
vtkImageData *imgData = vtkImageData::New();
imgData->SetDimensions(size);
imgData->SetSpacing(1.0, 1.0, 1.0);
imgData->SetScalarTypeToFloat();
imgData->GetPointData()->SetScalars(sarray);/* ベクトルデータの代入 */imgData->GetPointData()->SetVectors(varray);
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();
89

65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111

/* ランダムに格子点を取り出す */vtkMaskPoints *mask = vtkMaskPoints::New();
mask->SetInput(imgData);
mask->SetOnRatio(1000);
mask->RandomModeOn();/* 上記で取り出した格子点のところに，
ベクトルデータを基に線分のポリゴンデータを生成 */vtkHedgeHog *hedgehog = vtkHedgeHog::New();
hedgehog->SetInput(mask->GetOutput());
hedgehog->SetScaleFactor(5.0/range[1]);vtkPolyDataMapper * vecMapper = vtkPolyDataMapper::New();
vecMapper->SetInput(hedgehog->GetOutput());
vecMapper->SetScalarRange(range[0], range[1]);vecMapper->SetLookupTable(lut);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(vecMapper);/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(Actor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();

90

112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158

iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();varray->Delete();sarray->Delete();imgData->Delete();lut->Delete();mask->Delete();hedgehog->Delete();vecMapper->Delete();Actor->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}
void load_data(){
int i,j,k;
FILE *fpi;
vector_x = (float *)malloc(sizeof(float)*SIZE_X*SIZE_Y*SIZE_Z);
vector_y = (float *)malloc(sizeof(float)*SIZE_X*SIZE_Y*SIZE_Z);
vector_z = (float *)malloc(sizeof(float)*SIZE_X*SIZE_Y*SIZE_Z);
if ((fpi = fopen(datafile_x, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_x, sizeof(float), SIZE_X*SIZE_Y*SIZE_Z, fpi);
fclose(fpi);
if ((fpi = fopen(datafile_y, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_y, sizeof(float), SIZE_X*SIZE_Y*SIZE_Z, fpi);
91

159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200

fclose(fpi);
if ((fpi = fopen(datafile_z, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_z, sizeof(float), SIZE_X*SIZE_Y*SIZE_Z, fpi);
fclose(fpi);
if ((fpi = fopen(datafile_abs, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_abs, sizeof(float), SIZE_X*SIZE_Y*SIZE_Z, fpi);
fclose(fpi);for(k=0;k<SIZE_Z;k++){for(j=0;j<SIZE_Y;j++){for(i=0;i<SIZE_X;i++){vector[i + j* SIZE_X + k*SIZE_X*SIZE_Y][0]= vector_x[i+j*SIZE_X + k*SIZE_X*SIZE_Y];
vector[i + j* SIZE_X + k*SIZE_X*SIZE_Y][1]= vector_y[i+j*SIZE_X + k*SIZE_X*SIZE_Y];
vector[i + j* SIZE_X + k*SIZE_X*SIZE_Y][2]= vector_z[i+j*SIZE_X + k*SIZE_X*SIZE_Y];}}}
free(vector_x);free(vector_y);free(vector_z);
}

6.2.2

サンプルプログラム 1 の簡単な説明 1: VTK 用のベクトルデータ

vtkImageData や vtkRectilinearGrid にベクトル場をセットしなければ，可視化のためのフィルタを使うことはできない。VTK 形式のデータ作成の項で，VTK 形式のベクトルデータの作り方
92

はすでに述べたので，ここでは，マニュアルでベクトルデータをセットする方法をサンプルプ
ログラムを通して説明する。
VTK のベクトルデータは，X, Y, Z 成分のデータがメモリ上に X Y Z の順番で並んでいなければならない。具体的にサンプルプログラム 1 の例で言うと，float vector[SIZE X*SIZE Y*SIZE Z][3]の中に，X, Y, Z の各成分が入ったものを用意して (load data 関数)，それを vtkFloatArray を介して，imgData(等間隔メッシュ) に代入してベクトルデータを作っている。
サンプルプログラム 1 の load data 関数 (135 - 200 行) は，別々のファイルから，それぞれ x,
y, z 成分のデータを読み込み，float vector[SIZE X*SIZE Y*SIZE Z][3] にまとめている。その後，42 - 47 行で，varray (vtkFloatArray) にこのデータを渡している。スカラーのときは，
47 : varray->SetArray(vector[0], SIZE X*SIZE Y*SIZE Z*3, 1);とするだけであったが，今回はベクトルなので，
43 : varray->SetNumberOfComponent(3); /* 成分の数を指定 */46 : varray->SetNumberOfTuples(SIZE X*SIZE Y*SIZE Z); /* 格子点の総数を指定 */としている (Figure 6.1)。こうして用意した VTK の配列を，58 行目で等間隔メッシュのデータの箱に

imgData->GetPointData()->SetVectors(varray);として代入している (このサンプルプログラム 1 では，ベクトルデータと同時に，スカラーデー
タ (ベクトルの絶対値のデータ) も同じ箱 (imgData) に代入している)。
この imgData を使うことで，ようやく，VTK でのベクトル場の可視化が可能になる。

Figure 6.1: Component と Tuple

6.2.3

サンプルプログラム 1 の簡単な説明 2: 可視化部分

ベクトル場に，線分や矢印をちりばめるには，一部の点をデータ内から取り出さなければならない (すべての点で線分を作ると，何がなんだかわからなくなるので)。それには，• vtkMaskPoints クラス 1を使えばよい。vtkMaskPoints は，ランダムあるいは周期的に，データの格子点を取り出してくれる。また，• vtkHedgeHog クラス
このクラスのほかに，vtkThresholdPoints というスカラー値の閾値を与えて，点を取り出すクラスもある。例
えば，高温なところの点のみを取り出すなどに使える。vtkMaskPoints と組み合わせても良い。
193

は，入力された点に，その点のベクトルデータを考慮して，線分のポリゴンデータを生成する。
このクラスに，vtkMaskPoints で取り出した点を入力すればよい。
このサンプルプログラムの可視化処理部分は，vtkMaskPoints を使い全格子点の 1/1000 の点をランダムに取り出し，その点で線分のポリゴンデータを vtkHedgeHog を使い生成して表示する。色づけは，スカラーデータを使っている ( このプログラムでは，ベクトルの絶対値をスカ
ラーデータとして持っているので，ベクトルの絶対値が大きいところでは赤っぽくなり，小さいところでは青っぽくなる)。

66 - 70 : imgData から，1/1000 の割合で格子点を取り出す。
69 : SetOnRatio で，1/1000 の割合で格子点を取り出すよう指定
70 : RandomModeOn() で，ランダムに点を選び出す。default では，周期的に選び出す
74 - 76 : 取り出した点でのベクトルの値を基に，線分のポリゴンデータを作成する
その他は，前章までの知識で，理解していただけると思う。パイプラインを Figure 6.2 に示す。

Figure 6.2: サンプルプログラム 1 のパイプライン構造

Figure 6.3: 線分・矢印によるベクトル場表示
表示は，Figure 6.3(左)。

94

6.3

矢印・コーンを使う

HedgeHog による可視化では，ベクトルの向きまではわからない。そこで，線分の代わりに，矢印やコーンを使って，方向までわかるようにしてみる。それには，• vtkGlyph3D クラスという，SetInput で入力された点に，SetSource で入力されたポリゴンデータを方向も考慮に入れて配置していくクラスを vtkHedgeHog と差し替えればよい。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

/* 矢印のポリゴンデータ生成 */vtkArrowSource *arrow = vtkArrowSource::New();
arrow->SetTipResolution(4.0);
arrow->SetShaftResolution(4.0);/* mask で取り出した点に矢印のポリゴンデータを配置する */vtkGlyph3D *glyphs = vtkGlyph3D::New();
glyphs ->SetInput (mask ->GetOutput());
glyphs ->SetSource (arrow ->GetOutput());
glyphs ->SetScaleFactor (2.0);
glyphs ->SetScaleModeToScaleByVector();glyphs ->SetColorModeToColorByScalar();vtkPolyDataMapper * vecMapper = vtkPolyDataMapper::New();
vecMapper->SetInput(glyphs->GetOutput());
vecMapper->SetScalarRange(range[0], range[1]);vecMapper->SetLookupTable(lut);

サンプルプログラム 1 の 72 - 81 行を上記のリストと入れ替えると，線分が矢印に変わる (入れ替えたプログラムは，viz vect2.cxx として保存されている)。

1 - 4 : 矢印のポリゴンデータ生成
8 : 取り出した点を入力
9 : 矢印のポリゴンデータを入力
11 : 矢印を表示するとき，大きさはベクトルの大きさと連動するようにする。SetScaleModeToScalingOff() で，ベクトルの大きさにかかわらず，同じ大きさの矢印を使うことができる
パイプラインを Figure 6.4 に示す。表示は，Figure 6.3(右)。矢印ではなく，コーンを使いたい
ときは，vtkConeSource でコーンのポリゴンデータを生成して，9 行目の SetSource でそれを入力すればよい。ここまでの知識を組み合わせると，Figure 6.5 のように，等値面上の点でのベクトルを矢印で表現することもできる。気圧の等値面を描き，それを温度で色付けして，その面上での風速を矢印で表示，なども可能である。
この 3 次元の矢印やコーンを利用すると，計算機によっては，描画が遅くなる場合がある。このような場合は，2 次元の矢印を使うと良い。2 次元の矢印のポリゴンデータは，vtkGlyphSource2Dで生成することができる。vtkGlyphSource2D を利用すると，Figure 6.6 のようなポリゴンデー
タが得られる。
これ以降本章のサンプルプログラムでは，リーダーを使いデータを読み込む方法をとる。

95

Figure 6.4: 矢印を使う場合のパイプライン構造

Figure 6.5: ベクトルを矢印で表現：等値面上の点を選んである
補足
言うまでもないが，vtkHedgeHog や vtkGlyph3D に入力する点には，ベクトルデータが付随していなければならない。

6.4

ベクトル版凹凸つきスライス

スカラーデータの可視化で (第 2 章)，vtkWarpScalar というクラスを使い，頂点をスカラーの値に応じて移動させることで，スライスに凹凸をつけて富士山を可視化した。この vtkWarpScalarのベクトル版で， vtkWarpVector というクラスがある。これは，ベクトルの方向と大きさで頂点を移動させるクラスである。

6.4.1
1

サンプルプログラム 2/*Vec_BumpySlice.cxx */96

Figure 6.6: 2 次元グリフ
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36

#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkImageDataGeometryFilter.h>#include <vtkWarpVector.h>#include <vtkPolyDataNormals.h>#include <vtkTriangleFilter.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
int main( int argc, char *argv[] ){float range[2];char datafile[]= "./vector.vtk";
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);

97

37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83

lut->Build();
vtkImageDataGeometryFilter *igf
= vtkImageDataGeometryFilter::New();
igf->SetInput(imgData);
igf->SetExtent(75, 225, 75, 225, 80, 80);
vtkTriangleFilter *tri = vtkTriangleFilter::New();
tri->SetInput(igf->GetOutput());/* ベクトルデータを基に頂点を移動 */vtkWarpVector *v_carpet = vtkWarpVector::New();
v_carpet->SetInput(tri->GetOutput());
v_carpet->SetScaleFactor(10.0);
vtkPolyDataNormals *normals = vtkPolyDataNormals::New();
normals->SetInput(v_carpet->GetPolyDataOutput());
vtkPolyDataMapper * vecMapper = vtkPolyDataMapper::New();
vecMapper->SetInput(normals->GetOutput());
vecMapper->SetScalarRange(range[0], range[1]);vecMapper->SetLookupTable(lut);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(vecMapper);vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);vtkPolyDataMapper *gridMapper = vtkPolyDataMapper::New();
gridMapper->SetInput(igf->GetOutput());
gridMapper->ScalarVisibilityOff();
vtkActor *GActor = vtkActor::New();
GActor->SetMapper(gridMapper);
GActor->GetProperty()->SetColor(1.0,1.0,1.0);
GActor->GetProperty()->SetAmbient(0.5);
GActor->GetProperty()->SetDiffuse(0.5);
GActor->GetProperty()->SetRepresentationToWireframe();vtkRenderer *vren= vtkRenderer::New();
98

84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123

vren->AddActor(Actor);
vren->AddActor(GActor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();igf->Delete();tri->Delete();v_carpet->Delete();normals->Delete();vecMapper->Delete();Actor->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();gridMapper->Delete();GActor->Delete();vren->Delete();renWin->Delete();iwin->Delete();trackball->Delete();return 0;}

6.4.2

サンプルプログラム 2 の簡単な説明

vtkWarpVector (47 - 50) は，vtkWarpScalar から使い方を類推できると思う。実行結果は，Figure 6.7。下のほうに見える白い網は，移動する前の頂点の位置を表している。2 つを見比べると，渦巻いているのがわかる。

99

Figure 6.7: 凹凸つきスライスでベクトル場の一部を表示

6.5

この章のまとめ• 線分・矢印によるベクトル場の可視化• 凹凸つきスライスによるベクトル場の可視化
本章で紹介したクラスは，Table 6.1 の通り。

Table 6.1: 本章で紹介したクラスおよびそのメソッド
配列

点の抽出

点の抽出
ベクトル場の線分作成
2 次元の矢印
ポリゴンデータ配置
頂点移動

vtkFloatArray: SetNumberOfComponents(int): SetNumberOfTuples(const vtkIdType)vtkMaskPoints: SetOnRatio(int): RandomModeOn()
vtkThresholdPoints
vtkHedgeHog
vtkGlyphSource2D
vtkGlyph3D
vtkWarpVector

100Chapter 7
3 次元ベクトルデータの可視化 その 2
7.1

本章の概要

ベクトルデータを流線や流線と他のフィルタを組み合わせて可視化する。

7.2

流線 1 : 流線

当然であるが，VTK には，流線追跡機能がある。単純に曲線 (実際は細かい線分の集まり) を描
かせることもできるし，他のフィルタと組み合わせて，面白い可視化をすることもできる。
流線計算には，vtkStreamLine (破線版は，vtkDashedStreamLine) を使う。このクラスに，デー
タと出発点，Integrator を代入して刻み幅などのパラメータを設定すれば，流線のポリゴンを出力してくれる。

7.2.1

サンプルプログラム 1

下記が，複数の流線を描くプログラムである。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

/* StreamLines.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkPoints.h>#include <vtkPointSet.h>#include <vtkRungeKutta4.h>#include <vtkStreamLine.h>#include <vtkDashedStreamLine.h>#include <vtkSphereSource.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>
101

19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65

#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define NSP 10 /* 流線の出発点の数 */int main( int argc, char *argv[] ){int i;
float range[2];/* 流線の出発点 */float startpoints[NSP][3]={{100.0, 100.0, 15.0},{110.0, 110.0, 15.0},{120.0, 120.0, 15.0},{130.0, 130.0, 15.0},{140.0, 140.0, 15.0},{150.0, 150.0, 15.0},{160.0, 160.0, 15.0},{170.0, 170.0, 15.0},{180.0, 180.0, 15.0},{190.0, 190.0, 15.0}};
char datafile[]= "./vector.vtk";

vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 出発点のセッティング */vtkPoints *points = vtkPoints::New();
points->Allocate(NSP);
for (i = 0; i < NSP; i++) points->InsertPoint(i, startpoints[i]);vtkPolyData *pset = vtkPolyData::New();
pset->SetPoints(points);
102

66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112

/* 流線の計算 */vtkStreamLine *slines = vtkStreamLine::New();
vtkRungeKutta4 *integ = vtkRungeKutta4::New();slines->SetInput(imgData);
slines->SetSource(pset);slines->SetIntegrator(integ);slines->SetMaximumPropagationTime(300);
slines->SetIntegrationStepLength(0.1);
slines->SetStepLength(0.5);
vtkPolyDataMapper * slineMapper = vtkPolyDataMapper::New();slineMapper->SetInput(slines->GetOutput());
slineMapper->SetScalarRange(range[0], range[1]);slineMapper->SetLookupTable(lut);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(slineMapper);/* 出発点に赤い球体を配置 */vtkSphereSource *sphere = vtkSphereSource::New();
sphere->SetRadius(1.0);
vtkPolyDataMapper *SMapper = vtkPolyDataMapper::New();
SMapper->SetInput(sphere->GetOutput());
vtkActor *SActor[NSP];
vtkRenderer *vren= vtkRenderer::New();
for(i=0;i<NSP;i++){SActor = vtkActor::New();
SActor->SetMapper(SMapper);SActor->GetProperty()->SetColor(1.0, 0.0, 0.0);
SActor->SetPosition(startpoints[i]);vren->AddActor(SActor);
SActor->Delete();}/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
103

113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153

OLActor->SetMapper(OLMapper);vren->AddActor(Actor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();points->Delete();pset->Delete();slines->Delete();integ->Delete();slineMapper->Delete();Actor->Delete();sphere->Delete();SMapper->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}

7.2.2

サンプルプログラム 1 の簡単な説明

新規な部分は，58 - 76 行である。58 - 65 では，startpoints[NSP][3] に代入しておいた流線の出発点を vtkPolyData の点データとして pset に保持させている (もっとスマートな方法があるかも104しれません)。67 - 76 行で，4 次のルンゲ・クッタ法を用いて，積分計算している。vtkStreamLine は，基本的に，速度場 を可視化するためのクラスです。なお，vtkStreameLine の代わりに
vtkDashedStreameLine を 使うと，破線が得られる。

23 : NSP = 出発点の数
31 - 41 : float startpoints[NSP][3] に出発点の座標値をセット
58 - 62 : startpoints[NSP][3] を points(vtkPoints) に代入
64 - 65 : points を pset(vtkPolyData) に代入 (もう少しスマートな方法があるかもしれないです)71 : ベクトルのデータをセット
72 : 出発点 (pset) をセット
73 : 4 次のルンゲ・クッタを使用する (2 次も用意されている)
74 : 積分を続ける長さを時間で指定。これが大きいほど，流線は長くなるが，計算時間もかかる (引数が同じ値でも，入力されているベクトル場により，流線の長さは変わる)
75 : 積分の細かさを指定する。小さいほど正確になる。セルの大きさにたいする割合 (0.1 なら
セルの 1/10 刻みで流線計算する)
76 : ポリゴン (線分の集まり) の頂点を出力する時間間隔。小さいほど多くの点を出力して綺麗なカーブを描くが，描画が遅くなる。また，速度が小さいところほど，多くの点を出力する
86 - 100 : vtkSphereSource を使い，出発点に赤い球体を表示させる

Figure 7.1: 流線による可視化
結果は，Figure 7.1。パイプラインを Figure 7.2 に示す。
「積分の細かさ」「点を出力する時間間隔」では，意味がよく分からないと思うので，IntegrationStepLength と StepLength を大きくした場合，流線にどのような変化がおきるか，具体的な結果を Figure 7.3 に示す。色つきが，サンプルプログラム 1 と同じパラメータで，白い流線がそれぞれの StepLength を大きくしたものである。IntegrationSteLength を大きくした Figure 7.3(左)をみると，元の線とずれてしまっている (ポリゴンの頂点の数は，減っていない)。理由は，積分が荒くなったためである。StepLength を大きくした (右) をみると，ポリゴンの点自体は，元の線上に存在するのがわかる。
速度場以外を可視化する場合は特に，「時間」や「刻み幅」の設定を気をつけなければならない。
積分する方向は，SetIntegrationDirectionTo〜(Forward(), Backward, BothDirections()) で，指
定できる。

105

Figure 7.2: 流線のパイプライン構造

Figure 7.3: StepLength：(左)IntegratonStepLength を大きくした場合，(右)StepLength を大きくした場合
補足
前章の vtkHedgeHog や vtkGlypho3D の場合とは違い，流線計算のためのベクトルデータは別に入力するので，vtkStreamLine に入力する点にはベクトルデータが付随している必要はない。

7.3

流線と矢印の組み合わせ

前章のサンプルプログラム 2 とこのサンプルプログラム 1 を組み合わせると，Figure 7.4(左) にあるように，流線上に矢印を配置することもできるし，(右) のように流線上の点に別のベクト
ル場の矢印を配置することもできる。
前者は，流線のポリゴンデータを vtkMaskPoints に代入した後に，その出力を vtkGlyph3Dに入力すればよい。
後者は，少々工夫が必要である。 vtkGlyph3D に入力する点には，ベクトルデータが付随
していなければならないが，流線計算の出力をそのまま代入したのでは，流線と同じベクトル
データが付随しているので，別のベクトル場の矢印を表示することはできない。この問題を解決するのが，vtkProbeFilter である。このフィルタは，点とデータセットを入力すると，入力点にデータセットから (格子点上になければ補間して) スカラーやベクトルのデータを付随させた

106

Figure 7.4: 流線上に矢印を配置
上で出力してくれる。具体的に，プログラムの一部を示す。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20

vtkStreamLine *sline = vtkStreamLine::New();
vtkRungeKutta4 *integ = vtkRungeKutta4::New();sline->SetInput(imgData);
sline->SetSource(pset);sline->SetIntegrator(integ);sline->SetMaximumPropagationTime(300);
sline->SetIntegrationStepLength(0.1);
sline->SetStepLength(0.5);
vtkMaskPoints *mask = vtkMaskPoints::New();
mask->SetInput(sline->GetOutput());
mask->SetOnRatio(50);
vtkProbeFilter *probe = vtkProbeFilter::New();
probe->SetInput(mask->GetOutput());
probe->SetSource(imgData2);
vtkGlyph3D *glyphs = vtkGlyph3D::New();
glyphs ->SetInput (probe ->GetOutput());
glyphs ->SetSource (arrow ->GetOutput());

1-8 : imgData のベクトルデータで流線を計算
10 - 12 : 流線のポリゴンデータの頂点を一部を抜き取る
14 - 16 : 15 行で入力された点に imgData2 のベクトルやスカラーデータを付随させて，出力する (19 行)7.4

流線 2 : チューブ・リボン

VTK には，線 (分) をチューブで覆う vtkTubeFilter というクラスやリボン上に変換する vtkRibbonFilter というクラスがある。これを使うと，単なる線での表示よりも，奥行き感などがでる。
107

いずれのクラスも，線 (分) のポリゴンデータを入力すると，チューブやリボンのポリゴンデー
タを出力する。これらのクラスと流線を組み合わせると，単なる曲線での表示よりも構造を把握しやすくなる。

7.4.1

サンプルプログラム 2

サンプルプログラム 1 と違い，vtkPointSource を使い出発点をランダムに決めている。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39

/* StreamLines2.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkPointSet.h>#include <vtkRungeKutta4.h>#include <vtkStreamLine.h>#include <vtkSphereSource.h>#include <vtkPointSource.h>#include <vtkTubeFilter.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define NSP 10/* 出発点の数 *//* 球体の中からランダムに NSP 点選び出発点とする */#define Bradius 5.0 /* 球体の半径 */float Bposition[3] = {160, 149.0, 10.0}; /* 球体の中心 */int main( int argc, char *argv[] ){float range[2];char datafile[]= "./vector.vtk";
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();

108

40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86

imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 球体の中から NSP 点ランダムに点を選ぶ */vtkPointSource *pset = vtkPointSource::New();
pset->SetCenter(Bposition);pset->SetRadius(Bradius);
pset->SetNumberOfPoints(NSP);/* 流線計算 */vtkStreamLine *slines = vtkStreamLine::New();
vtkRungeKutta4 *integ = vtkRungeKutta4::New();slines->SetInput(imgData);
slines->SetSource(pset->GetOutput());
slines->SetIntegrator(integ);slines->SetMaximumPropagationTime(300);
slines->SetIntegrationStepLength(0.1);
slines->SetIntegrationDirectionToForward();slines->SetStepLength(0.5);/* 流線をチューブで囲む */vtkTubeFilter *tubes = vtkTubeFilter::New();
tubes->SetInput(slines->GetOutput());
tubes->SetRadius(0.5);
tubes->SetRadiusFactor(5.0);
tubes->SetVaryRadiusToVaryRadiusByScalar();tubes->SetNumberOfSides(8);vtkPolyDataMapper * tubeMapper = vtkPolyDataMapper::New();
tubeMapper->SetInput(tubes->GetOutput());
tubeMapper->SetScalarRange(range[0], range[1]);tubeMapper->SetLookupTable(lut);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(tubeMapper);/* 球体をワイヤフレームで表示 */vtkSphereSource *sphere = vtkSphereSource::New();
sphere->SetRadius(Bradius);
sphere->SetCenter(Bposition);
vtkPolyDataMapper *SMapper = vtkPolyDataMapper::New();
109

87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133

SMapper->SetInput(sphere->GetOutput());
vtkActor *SActor = vtkActor::New();
SActor->SetMapper(SMapper);SActor->GetProperty()->SetColor(1.0, 1.0, 1.0);
SActor->GetProperty()->SetRepresentationToWireframe();/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(Actor);
vren->AddActor(SActor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();pset->Delete();slines->Delete();integ->Delete();tubes->Delete();tubeMapper->Delete();Actor->Delete();outline->Delete();
110

134
135
136
137
138
139
140
141
142

OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}

7.4.2

サンプルプログラム 2 の簡単な説明

新規な部分は，47 - 51 行と 65 - 71 行である。47 - 51 では，vtkPointSource を使い出発点をランダムに選び，65 - 71 行では，vtkStreamLine の出力結果 (線分) を vtkTubeFilter を使いチュー
ブで囲んでいる。

47 - 51 : 「BPosition を中心とした半径 BRadius の球の中から，ランダムに NSP 個の点を選び出せ」，という意味。メソッドの名前がわかりやすいので，これだけで理解いただけると思う
67 : 流線のポリゴンを vtkTubeFikter に代入
68 : チューブの半径の最小値
69 : チューブの半径の最大値を最小値の何倍にするか
70 : スカラーの大きさに応じて，半径を変化させる。....ByVector() も存在する。また，SetVaryRadiusToVaryRadiusOff() で半径一定モードになる
71 : 引数が N だと，チューブが N 角形になる。この数字を大きくすると，綺麗になるが，ポリゴン数が増えるので，描画は遅くなる
81 - 92 : 点を取り出した球体をワイヤフレームで表示
表示は，Figure 7.5(左)。vtkTubeFilter の代わりに vtkRibbonFilter を使うと，Figure 7.5(右) のよ
うな可視化結果となる。この vtkTubeFilter は，いうまでもないが，入力のポリゴンデータは，

Figure 7.5: (左) チューブ，(右) リボン
流線の結果でなくてもよい。等高線をチューブで囲むこともできる。

111

7.5

流線 3 : 流面 (Stream Surfaces)

流線を沢山引くのも良いが，流線と流線の間に平面を入れると，少ない流線でもベクトル場の構造がわかりやすくなる。前セクションの “リボン”と違うのは，前者は単独の流線を面状にするだけであるが，こちらは，流線間に面を張るところである。
流線間に面を張るので，ある意味流線間にある流線の構造も見えている。ベクトル場全体の把握には，こちらの方が適しているかもしれない。

7.5.1

サンプルプログラム 3

指定した線分上の点をいくつか拾い上げ，そこから流線を引き，vtkRuledSurfaceFilter で流線と流線の間に面を張る。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34

/* StreamLines3.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkPoints.h>#include <vtkPointSet.h>#include <vtkRungeKutta4.h>#include <vtkStreamLine.h>#include <vtkLineSource.h>#include <vtkRuledSurfaceFilter.h>#include <vtkTubeFilter.h>#include <vtkSphereSource.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>#define NSP 56 /* 解像度 */float pos1[3] = {100.0, 140.0, 5.0};
float pos2[3] = {200.0, 200.0, 15.0};
int main( int argc, char *argv[] ){float range[2];char datafile[]= "./vector.vtk";

112

35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81

vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update();imgData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 解像度 NSP の線分を生成 */vtkLineSource *pole = vtkLineSource::New();
pole->SetPoint1(pos1);
pole->SetPoint2(pos2);
pole->SetResolution(NSP);/* 流線計算。上記線分の頂点が出発点となる */vtkStreamLine *slines = vtkStreamLine::New();
vtkRungeKutta4 *integ = vtkRungeKutta4::New();slines->SetInput(imgData);
slines->SetSource(pole->GetOutput());
slines->SetIntegrator(integ);slines->SetMaximumPropagationTime(300);
slines->SetIntegrationStepLength(0.1);
slines->SetIntegrationDirectionToForward();slines->SetStepLength(0.5);/* 流線の間に面を生成する */vtkRuledSurfaceFilter *ssurfaces = vtkRuledSurfaceFilter::New();
ssurfaces->SetInput(slines->GetOutput());
ssurfaces->SetOnRatio(2); /* 一つおき */ssurfaces->SetOffset(0);ssurfaces->SetDistanceFactor(10.0);
ssurfaces->SetRuledModeToResample();ssurfaces->SetResolution(100, 1);
vtkPolyDataMapper * ssMapper = vtkPolyDataMapper::New();
ssMapper->SetInput(ssurfaces->GetOutput());
ssMapper->SetScalarRange(range[0], range[1]);ssMapper->SetLookupTable(lut);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(ssMapper);
113

82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128

/* 流線そのものも表示する */vtkPolyDataMapper * slinesMapper = vtkPolyDataMapper::New();slinesMapper->SetInput(slines->GetOutput());
slinesMapper->SetScalarRange(range[0], range[1]);slinesMapper->SetLookupTable(lut);
vtkActor *Actor2 = vtkActor::New();
Actor2->SetMapper(slinesMapper);/* 出発点を発生させた線分をチューブで囲み表示 */vtkTubeFilter *tubes = vtkTubeFilter::New();
tubes->SetInput(pole->GetOutput());
tubes->SetRadius(1.0);
tubes->SetNumberOfSides(8);vtkPolyDataMapper *PMapper = vtkPolyDataMapper::New();
PMapper->SetInput(tubes->GetOutput());
vtkActor *PActor = vtkActor::New();
PActor->SetMapper(PMapper);/* 線分の始点と終点に球体を表示 */vtkSphereSource *sphere = vtkSphereSource::New();
sphere -> SetRadius(2.0);
vtkPolyDataMapper *sMapper = vtkPolyDataMapper::New();sMapper->SetInput(sphere->GetOutput());
vtkActor *SActor1 = vtkActor::New();
SActor1->SetMapper(sMapper);SActor1->SetPosition(pos1);
vtkActor *SActor2 = vtkActor::New();
SActor2->SetMapper(sMapper);SActor2->SetPosition(pos2);/* データの外枠 */vtkOutlineFilter *outline = vtkOutlineFilter::New();
outline->SetInput(imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);

114

129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175

vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(Actor);
vren->AddActor(Actor2);
vren->AddActor(PActor);
vren->AddActor(OLActor);
vren->AddActor(SActor1);
vren->AddActor(SActor2);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();reader->Delete();lut->Delete();pole->Delete();slines->Delete();integ->Delete();ssurfaces->Delete();ssMapper->Delete();Actor->Delete();slinesMapper->Delete();Actor2->Delete();tubes->Delete();PMapper->Delete();PActor->Delete();sphere->Delete();sMapper->Delete();SActor1->Delete();SActor2->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();
115

176
177
178
179

renWin->Delete();return 0;}

7.5.2

サンプルプログラム 3 の簡単な説明

Figure 7.6: vtkLineSource で線分のポリゴンデータを生成
新規な部分は，48 - 52 行と 66 - 73 行である。48 - 52 行では，vtkLineSource を使い，線分のポリゴンを生成している。流線計算の出発点は，線分の頂点がセットされる (Figure 7.6)。66 73 行では，vtkRuledSurfaceFilter を使い，隣り合う流線 (実際には入力された順番が「隣り合った」流線) に面を張っている。
50 - 51 : SetPoint1, 2 で，線分の始点と終点を設定
52 : SetResolution で線分の解像度を指定。引数を N とすると，上記の線分は，(N+1) 個の頂点を持つ (頂点は線分上に等間隔に配置される)
59 : vtkLineSource の出力を vtkStreamLine に代入。上記の頂点が，流線計算の出発点になる
69 : vtkRuledSurfaceFilter のメソッド SetOnRatio で，面のはり方を設定する。2 で一つおき。1
でまったく飛ばさずに面を張っていく。例えば，流線が 6 本あったとすると，引数を 1 とすると 1-2, 2-3,... の流線の間に面を張る。2 とすると，1-2, 3-4, 5-6 番目の流線の間，3 とすると 1-2,
5-6 番目の流線の間に面を張る (Figure 7.7)
71 : 引数を d とすると，流線上の点と点の距離が，出発点同士の距離の d 倍以上のところでは，面を張らなくなる
72 - 73 : 流線上の点を等間隔に 100(73 行の一つ目の引数) 個リサンプルして，面を作る。73 行の二つ目の引数は，線と線の間にいくつポリゴンを作るか。72 行で，SetRuledModeToPointwalk()とすると，流線上の点をそのまま用いて面を張る (ポリゴンが密なところ疎なところができる)
83 - 90 : 流線も表示させる
92 - 102 : 流線の出発点を取り出した線分を，チューブで表示
104 - 117 : 上記線分の始点終点に，球体を表示
表示は，Figure 7.8。

116

Figure 7.7: 面を張る位置

Figure 7.8: Stream Surfaces

7.6

この章のまとめ• 流線によるベクトル場の可視化• 点データの操作 (流線と矢印の組み合わせ)• チューブ・リボン• 流面 (Stream Surfaces)
本章で紹介したクラスは，Table 7.1 の通り。

117

Table 7.1: 本章で紹介したクラスおよびそのメソッド
点データの操作
点データ
ポリゴンデータ
流線計算

vtkProbeFilter
vtkPoints
vtkPolyData
vtkStreamLine, vtkDashedStreamLine: SetIntegrator(vtkInitialValueProblemSolver *): SetMaximumPropagationTime(float): SetIntegrationStepLength(float): SetStepLength(float): SetIntegrationDirectionTo〜(), Forward, Backward, BothDirections
流線計算 (Integrator) vtkRungeKutta4
点データ生成vtkPointSource: SetCenter(float [3]): SetRadius(float): SetNumberOfPoints(vtkIdType)
線データ生成vtkLineSource: SetPoint1,2(float [3]): SetResolution(int)チューブ
vtkTubeFilter: SetRadius(float): SetRadiusFactor(float): SetVaryRadiusToVaryRadius〜(), ByScalar,ByVector, Off: SetNumberOfSides(int)リボン
vtkRibbonFilter
面
vtkRuledSurfaceFilter: SetOnRatio(int): SetOffSet(int): SetDistanceFactor(float): SetRuledModeTo〜(), Resample, Pointwalk: SetResolution(int, int)

118

Chapter 8
円柱座標・球座標
8.1

本章の概要

StructuredGrid で，円柱座標や球座標のデータを可視化する。

8.2

円柱座標

vtkStructuredGrid を使うと，円柱座標や球座標のデータも取り扱うことができる。vtkStructuredGrid で円柱や球座標のデータさえ作ってしまえば，それ以降の可視化パイプラインは，前章までの知識で組むことができる。

8.2.1

サンプルプログラム 1

ひとまず，サンプルプログラムを見ながら，円柱座標データの作り方を解説する。本章のサン
プルデータは，前章 (ベクトル場の可視化) と同様のもので，スカラー場とベクトル場いずれも
存在する。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

/* Cylindrical Coordinates *//* cylind1.cxx */#include <vtkStructuredGrid.h>#include <vtkFloatArray.h>#include <vtkPointData.h>#include <vtkMaskPoints.h>#include <vtkHedgeHog.h>#include <vtkContourFilter.h>#include <vtkLookupTable.h>#include <vtkStructuredGridOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
119

19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65

#define SIZE_T 300#define SIZE_Z 300#define SIZE_R 100
float vector[SIZE_T*SIZE_Z*SIZE_R][3];float vector_abs[SIZE_T*SIZE_Z*SIZE_R];
char datafile_t[]= "./vector_x.dat";char datafile_z[]= "./vector_y.dat";char datafile_r[]= "./vector_z.dat";char datafile_abs[] = "./vector_abs.dat";
float *vector_t;
float *vector_z;
float *vector_r;
float points[SIZE_T*SIZE_Z*SIZE_R][3];const double PI = 3.14159;
const double dt = 180.0/299.0; /* Theta の刻み幅 */void load_data();/* データ読み込み */void cylind_coords(); /* 格子点のデータ作成 */int main( int argc, char *argv[] ){int size[3] = {SIZE_T, SIZE_Z, SIZE_R};
float range[2];cylind_coords();
load_data();vtkFloatArray *varray = vtkFloatArray::New();
varray->SetNumberOfComponents(3);varray->SetNumberOfTuples(SIZE_T* SIZE_Z* SIZE_R);
varray->SetArray(vector[0], SIZE_T* SIZE_Z* SIZE_R*3, 1);
vtkFloatArray *sarray = vtkFloatArray::New();
sarray->SetArray(vector_abs, SIZE_T* SIZE_Z* SIZE_R, 1);/* 格子点の座標データ */vtkFloatArray *carray = vtkFloatArray::New();
carray->SetNumberOfComponents(3);carray->SetNumberOfTuples(SIZE_T* SIZE_Z* SIZE_R);
carray->SetArray(points[0], SIZE_T* SIZE_Z* SIZE_R*3, 1);
vtkPoints *coordinates = vtkPoints::New();
120

66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112

coordinates->SetDataTypeToFloat();coordinates->SetData(carray);
vtkStructuredGrid *CylindData = vtkStructuredGrid::New();
CylindData->SetDimensions(size);/* 格子点の座標データをセット */CylindData->SetPoints(coordinates);
CylindData->GetPointData()->SetVectors(varray);
CylindData->GetPointData()->SetScalars(sarray);
CylindData->Update();CylindData->GetScalarRange(range);
vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
lut->Build();/* 以下 HedgeHog と 等値面のパイプライン *//* HedgeHog */vtkMaskPoints *mask = vtkMaskPoints::New();
mask->SetInput(CylindData);
mask->SetOnRatio(1000);
mask->RandomModeOn();
vtkHedgeHog *hedgehog = vtkHedgeHog::New();
hedgehog->SetInput(mask->GetOutput());
hedgehog->SetScaleFactor(5.0/range[1]);vtkPolyDataMapper * vecMapper = vtkPolyDataMapper::New();
vecMapper->SetInput(hedgehog->GetOutput());
vecMapper->SetScalarRange(range[0], range[1]);vecMapper->SetLookupTable(lut);
vtkActor *vecActor = vtkActor::New();
vecActor->SetMapper(vecMapper);/* 等値面 */vtkContourFilter *contour = vtkContourFilter::New();
contour->SetInput(CylindData);
contour->SetValue(0, 3.0);
contour->ComputeNormalsOn();
vtkPolyDataMapper * isoMapper = vtkPolyDataMapper::New();
isoMapper->SetInput(contour->GetOutput());
isoMapper->ScalarVisibilityOff();
vtkActor *isoActor = vtkActor::New();
121

113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159

isoActor->SetMapper(isoMapper);/* データの外枠 */vtkStructuredGridOutlineFilter *outline= vtkStructuredGridOutlineFilter::New();
outline->SetInput(CylindData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);
vtkRenderer *vren= vtkRenderer::New();
vren->AddActor(vecActor);
vren->AddActor(isoActor);
vren->AddActor(OLActor);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 800, 600 );vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();varray->Delete();sarray->Delete();coordinates->Delete();CylindData->Delete();lut->Delete();mask->Delete();hedgehog->Delete();contour->Delete();vecMapper->Delete();vecActor->Delete();isoMapper->Delete();isoActor->Delete();outline->Delete();
122

160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206

OLMapper->Delete();OLActor->Delete();iwin->Delete();trackball->Delete();vren->Delete();renWin->Delete();return 0;}
void cylind_coords(){
int i,j,k;double theta, r, z;/* 格子点の位置 */for(k=0; k<SIZE_R; k++){for(j=0; j<SIZE_Z; j++){for(i=0; i<SIZE_T; i++){theta = PI/180.0*i*dt;
r = 50.0 + 0.5*k; /* 半径方向は 50, 99.5*/
z = j;
points[i + SIZE_T*j + k*SIZE_T*SIZE_Z][0]= r * cos(theta);
points[i + SIZE_T*j + k*SIZE_T*SIZE_Z][1]= r * sin(theta);
points[i + SIZE_T*j + k*SIZE_T*SIZE_Z][2]
= z;
}
}
}
}
void load_data(){
int i, j, k;double theta;
float temp_t, temp_z, temp_r;
FILE *fpi;
vector_t = (float *)malloc(sizeof(float)*SIZE_T*SIZE_Z*SIZE_R);
123

207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253

vector_z = (float *)malloc(sizeof(float)*SIZE_T*SIZE_Z*SIZE_R);
vector_r = (float *)malloc(sizeof(float)*SIZE_T*SIZE_Z*SIZE_R);

if ((fpi = fopen(datafile_t, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_t, sizeof(float), SIZE_T*SIZE_Z*SIZE_R, fpi);
fclose(fpi);

if ((fpi = fopen(datafile_z, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_z, sizeof(float), SIZE_T*SIZE_Z*SIZE_R, fpi);
fclose(fpi);

if ((fpi = fopen(datafile_r, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_r, sizeof(float), SIZE_T*SIZE_Z*SIZE_R, fpi);
fclose(fpi);
if ((fpi = fopen(datafile_abs, "rb")) == NULL) {puts("cannot open");exit(1);}
fread(vector_abs, sizeof(float), SIZE_T*SIZE_Z*SIZE_R, fpi);
fclose(fpi);/********************************//* x 成分 -> theta 成分*//* y 成分 -> z 成分*//* z 成分 -> radial 成分*//********************************/
124

254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282

for(k=0;k<SIZE_R;k++){for(j=0;j<SIZE_Z;j++){for(i=0;i<SIZE_T;i++){temp_t = vector_t[i+j*SIZE_T + k*SIZE_T*SIZE_Z];
temp_z = vector_z[i+j*SIZE_T + k*SIZE_T*SIZE_Z];
temp_r = vector_r[i+j*SIZE_T + k*SIZE_T*SIZE_Z];
theta = PI/180.0*i*dt;
vector[i + j* SIZE_T + k*SIZE_T*SIZE_Z][0]= -temp_t*sin(theta) + temp_r * cos(theta);
vector[i + j* SIZE_T + k*SIZE_T*SIZE_Z][1]= temp_t*cos(theta) + temp_r * sin(theta);
vector[i + j* SIZE_T + k*SIZE_T*SIZE_Z][2]= temp_z;}}}free(vector_t);free(vector_z);
free(vector_r);
}

8.2.2

サンプルプログラム 1 の簡単な説明 : 格子点・データの作り方

格子点の作り方は，ベクトルデータの作り方に似ている。このサンプルプログラム 1 では，cylind coords() 関数 (170 - 196) で，全格子点の座標を float points[SIZE T*SIZT Z*SIZE R][3]に代入している。その後，points を carray(vtkFloatArray:60 - 63 行) と coordinates(vtkPoints) を介して，CylindData (vtkStructuredGrid) に代入している。この points に入っているデータは，座標値なので成分も 3 つある。それ故ベクトルデータのときと同様に，

carray->SetNumberOfComponents(3);carray->SetNumberOfTuples(SIZE T* SIZE Z* SIZE R);としている。これは，残念ながら直接 CylindData に代入できなく，一旦 coordinates(vtkPoints:65- 67 行) に代入してからCylidData->SetPoints(coordinates);

125

とする。
データ (スカラーやベクトル) の CylindData への代入方法は，等間隔メッシュや Rectilinearの場合と同様である。このサンプルプログラムでは，前章と同じデータの X 成分 → Theta 成
分 (0 度から 0.602 度刻み)，Y 成分 →z 成分，Z 成分 →R 成分 (50.0 から 0.5 刻み) となるように，
load data 関数 (198 - 282 行) で変換している。
格子点の位置やベクトルデータは，X Y Z で与えることに注意。

20 - 22 : Theta 方向 (SZIE T), Z 方向 (SIZE Z) の解像度はいずれも 300, R 方向 (SIZE R) は 100
以下 cylind coords 関数
170 - 196 :Theta は [0, π, Z は [0,300)，R は [50, 100) になるようにしている
以下 load data 関数
198 - 282 : データを vector, vector abs に読み込む。そのとき，sin, cos を使って，X 成分 → Theta
成分，Y 成分 → Z 成分，Z 成分 → R 成分となるように変換している (VTK とは関係ない)
59 - 63 : carray (vtkFloatArray) に，格子点の位置データ (float points[...][3]) を代入
65 - 67 : carray を coordinates (vtkPoints) に代入
72 - 73 : 円柱座標を作るために用意した箱 (CylindData) に，coordinates をセット (晴れて円柱座
標の箱になる)74 - 75 : スカラー及びベクトルデータを CylindData に代入
Figure8.1 にパイプラインを示す。

Figure 8.1: 円柱座標の作り方

8.2.3

サンプルプログラムの簡単な説明:可視化処理部

このプログラムは，ベクトル場を HedgeHog で，スカラー場を等値面で可視化する。可視化パ
イプラインは，前章までと同様である。唯一違うのが，外枠を作るクラスである (115 - 118 行)。
vtkOutlineFilter ではなく，vtkStructuredGridOutlineFilter を使っている。vtkOutlineFilter を使用
してもエラーや警告はでないが，表示される外枠は，円柱形の枠ではなく，データ全体が入る
直方体の枠になる。表示される画像は，Figure 8.2(左)。ここではソースリストを紹介しないが，もう一つ円柱座標のサンプルプログラム (cylind2.cxx) を用意した。カラーコンターと流線で可視化している。実行結果のみ，Figure 8.2(右) に示す。
参考のため，等間隔メッシュを使い，サンプルプログラム 1 と同じ可視化パイプラインをつないだ結果を，Figure 8.3 に示す。

126

Figure 8.2: (左)HedgeHog と 等値面による可視化，(右) カラーコンターと流線による可視化

Figure 8.3: HedgeHog と 等値面による可視化 (等間隔メッシュ)

8.3

球座標

円柱座標のときと同様に，球座標の格子点のデータさえ作ってしまえば，等値面でも流線でも使えることは理解いただけていると思う。ここでは，ボリューム・レンダリングに挑戦する。
(ここで紹介する方法は，StructuredGrid の他に，Rectilinear でも使えるはずである)

8.3.1

サンプルプログラム 2

座標は，X 方向 →Theta 方向 (解像度 300, [π/4, 3π/4])，Y 方向 →Phi 方向 (解像度 300, [0, π/2])，
Z 方向 →R 方向 (解像度 100, [50, 75)) となるようにしている (spheri coords 関数 202 - 228 行)。
座標やデータの設定方法は，サンプルプログラム 1 を参考にしていただければ，理解していただけると思う。なお，半径 50.0 の球体も同時に表示するようにしている。

1
2
3/* Spherical Coordinates *//* spheri3.cxx */#include <vtkStructuredGrid.h>

127

4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50

#include <vtkStructuredGridGeometryFilter.h>#include <vtkFloatArray.h>#include <vtkPointData.h>#include <vtkSphereSource.h>#include <vtkDataSetMapper.h>#include <vtkPlane.h>#include <vtkCutter.h>#include <vtkPolyDataNormals.h>#include <vtkLookupTable.h>#include <vtkStructuredGridOutlineFilter.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkCamera.h>#include <vtkRendererSource.h>#include <vtkBMPWriter.h>#include <unistd.h>#define SIZE_T 300#define SIZE_P 300#define SIZE_R 100
float vector_abs[SIZE_T*SIZE_P*SIZE_R];
char datafile_abs[] = "./vector_abs.dat";
float points[SIZE_T*SIZE_P*SIZE_R][3];const double PI = 3.14159;
const double dt = 90.0/299.0; /* Theta, Phi の刻み幅 */void load_data();void spheri_coords();
int main( int argc, char *argv[] ){int i;
int size[3] = {SIZE_T, SIZE_P, SIZE_R};
float range[2];float theta = PI/180.0*SIZE_T/2.0*dt + PI/4.0;
float r = 50.0 + 0.25*SIZE_R/2.0;
float phi = PI/180.0*dt * SIZE_P/2.0;
float Color[4];float center[3], normal[3];

128

51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97

spheri_coords();
load_data();vtkFloatArray *sarray = vtkFloatArray::New();
sarray->SetArray(vector_abs, SIZE_T* SIZE_P* SIZE_R, 1);/* 格子点の座標データ */vtkFloatArray *carray = vtkFloatArray::New();
carray->SetNumberOfComponents(3);carray->SetNumberOfTuples(SIZE_T* SIZE_P* SIZE_R);
carray->SetArray(points[0], SIZE_T* SIZE_P* SIZE_R*3, 1);
vtkPoints *coordinates = vtkPoints::New();
coordinates->SetDataTypeToFloat();coordinates->SetData(carray);
vtkStructuredGrid *SphericalData = vtkStructuredGrid::New();
SphericalData->SetDimensions(size);/* 格子点の座標データをセット */SphericalData->SetPoints(coordinates);
SphericalData->GetPointData()->SetScalars(sarray);SphericalData->Update();SphericalData->GetScalarRange(range);/* カラーテーブル作成 */vtkLookupTable *lut = vtkLookupTable::New();
lut->SetNumberOfColors(256);
lut->SetHueRange(0.677, 0.0);
lut->Build();for(i=0;i<256;i++){lut->GetTableValue(i, Color);
Color[3] = i/255.0;
if(i > 100 && i < 150) Color[3] = 0.0;
lut->SetTableValue(i, Color);}/* スライス */center[0] = r * sin(theta) * cos(phi);
center[1] = r * sin(theta) * sin(phi);
center[2] = r * cos(theta);
vtkPlane *plane = vtkPlane::New();plane->SetOrigin(center);
129

98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144

vtkCutter *cutter = vtkCutter::New();
cutter->SetInput(SphericalData);
cutter->SetCutFunction(plane);cutter->SetSortByToSortByCell();cutter->GenerateCutScalarsOff();
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(cutter->GetOutput());
Mapper->SetLookupTable(lut);
Mapper->SetScalarRange(range);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);/* 半径 50 の球体を表示する */vtkSphereSource *ball = vtkSphereSource::New();ball->SetRadius(50.0);
ball->SetPhiResolution(100);ball->SetThetaResolution(100);
vtkPolyDataNormals *norm = vtkPolyDataNormals::New();
norm->SetInput(ball->GetOutput());
vtkPolyDataMapper *ballMapper = vtkPolyDataMapper::New();
ballMapper->SetInput(norm->GetOutput());
vtkActor *ballActor = vtkActor::New();
ballActor->SetMapper(ballMapper);ballActor->GetProperty()->SetColor(1.0, 0.75, 0.75);
vtkCamera *camera = vtkCamera::New();
camera->SetPosition(150, 300, 0.0);
camera->SetFocalPoint(center[0], center[1], center[2]);camera->SetClippingRange(10, 800);
camera->SetViewUp(0,0,1);/* 投影面の法線ベクトルを求め，平面の法線ベクトルとする */camera->ComputeViewPlaneNormal();camera->GetViewPlaneNormal(normal);plane->SetNormal(normal);cutter->GenerateValues(64, -15., 15.);/* データの外枠 */vtkStructuredGridOutlineFilter *outline= vtkStructuredGridOutlineFilter::New();
130

145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191

outline->SetInput(SphericalData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New();
OLMapper->SetInput(outline->GetOutput());
vtkActor *OLActor = vtkActor::New();
OLActor->SetMapper(OLMapper);
vtkRenderer *vren = vtkRenderer::New();
vren->AddActor(ballActor);
vren->AddActor(OLActor);
vren->AddActor(Actor);
vren->SetActiveCamera(camera);
vren->SetBackground( 0.0, 0.0, 0.0 );vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer( vren );renWin->SetSize( 500, 375 );char FileName[]="volren.bmp";
vtkRendererSource *rs = vtkRendererSource::New();
rs->SetInput(vren);rs->WholeWindowOn();
vtkBMPWriter *bw = vtkBMPWriter::New();
bw->SetInput(rs->GetOutput());
renWin->Render();rs->Modified();bw->SetFileName(FileName);bw->Write();sarray->Delete();coordinates->Delete();SphericalData->Delete();lut->Delete();plane->Delete();cutter->Delete();Mapper->Delete();Actor->Delete();outline->Delete();OLMapper->Delete();OLActor->Delete();ball->Delete();ballMapper->Delete();
131

192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238

ballActor->Delete();camera->Delete();vren->Delete();renWin->Delete();rs->Delete();bw->Delete();return 0;}
void spheri_coords(){
int i,j,k;double theta, phi, r;/* 格子点の位置 */for(k=0; k<SIZE_R; k++){for(j=0; j<SIZE_P; j++){for(i=0; i<SIZE_T; i++){theta = PI/180.0*i*dt + PI/4.0;
r = 50.0 + 0.25*k;
phi = PI/180.0*dt * j;points[i + SIZE_T*j + k*SIZE_T*SIZE_P][0]= r * sin(theta) * cos(phi);points[i + SIZE_T*j + k*SIZE_T*SIZE_P][1]= r * sin(theta) * sin(phi);points[i + SIZE_T*j + k*SIZE_T*SIZE_P][2]= r * cos(theta);}}}}
void load_data(){
FILE *fpi;
if ((fpi = fopen(datafile_abs, "rb")) == NULL) {puts("cannot open");exit(1);
}

132

239
240
241
242
243

fread(vector_abs, sizeof(float), SIZE_T*SIZE_P*SIZE_R, fpi);
fclose(fpi);
}

8.3.2

サンプルプログラム 2 の簡単な説明

VTK のボリューム・レンダリングは，等間隔メッシュに限られているが，工夫次第では，様々な座標系で実行できる。テクスチャ・マッピングを利用したボリューム・レンダリングを思い出していただきたい。不透明度を含んだ画像を貼り付けたポリゴンを前後に並べてブレンドする方法である。テクスチャ機能が使えなくても，第 1 章で解説したカラーコンターの方法で代用できる。しかもデータは等間隔メッシュである必要がない。つまり，頂点に色だけでなく不
透明度も割り当てたカラーコンターのスライスを前後に並べるのである。スライスの引き出しには，vtkPlane と vtkCutter が使える。これを使えば，3 次元テクスチャを利用したボリュー
ム・レンダリングのように，投影面に平行なスライスを使ってボリューム・レンダリングできる (Figure 8.4)。

Figure 8.4: 投影面と平行な面でスライス
以下，スライス関係の説明

91 - 103 : データの中心を通る平面 (plane，原点のみ設定) と cutter を生成している
135 - 137 : 投影面の法線ベクトルを計算して，float normal[3] に代入 (規格化済み)
139 : 上で求めた法線ベクトルを plane の法線ベクトルに設定
140 : -15.0 〜 15.0 の間で，64 枚スライスを生成。意味するところは，第 3 章に登場したときは，
F (x, y, z) = 0.0 の平面でデータを切ったが，本章では，F (x, y, z) = An ( An は，-15.0〜15.0 の間から 64 個等間隔に値を割り当てられる) の (64 枚の) 平面で切ったスライスを表示させる
色関係は，76 - 89 行で設定している。カラーテーブルを一旦作って (80 行) から，不透明度を変えている (82 - 89 行)。これで，球座標のボリューム・レンダリングが実行できる。表示画像は，
Figure 8.5。
ただ，この方法は，処理時間が非常に長いという欠点がある。Figure 8.5 を出すのに，venus
で 20 分強費やしている。スライス数を増やすと，さらに処理時間がかかると思われる。

133

Figure 8.5: 球座標データのボリューム・レンダリング

8.3.3

その他の球座標サンプル

spheri1.cxx, spheri2.cxx という 2 つのサンプルプログラムを用意した。実行結果は，それぞれFigure 8.6 の (左) と (右) である。

Figure 8.6: (左) 矢印 と 等値面による可視化，(右) 等高線と流線 (チューブ) による可視化

8.4

VTK 形式のデータ

Structured Grid のデータを VTK 形式にするときのヘッダを紹介する。Rectilinear Grid の座標を
格子点の座標に変えるだけである。
1
2
3
4
5

# vtk DataFile Version 2.0
Structured Grid のデータ
BINARYDATASET STRUCTURED_GRIDDIMENSIONS 50 100 50
134

6
7
8
9
10
11
12
13

POINTS 250000 float
格子点のバイナリデータ
POINT_DATA 250000
SCALARS NameOfData_S float 1
LOOKUP_TABLE default
スカラーのバイナリデータ
VECTORS NameOfData_V float
ベクトルのバイナリデータ
このようにして作成したデータは，vtkStructuredGridReader で読み込むことができる。

8.5

この章のまとめ• vtkStructuredGrid を使っての円柱・球座標データの作成と可視化• ポリゴンスライスによる球座標データのボリュームレンダリング• vtkStructuredGrid の VTK 形式のデータ
本章で紹介したクラスは，Table 8.1 の通り。

Table 8.1: 本章で紹介したクラスおよびそのメソッド
データセット
リーダー
カメラ関係

断面
外枠

vtkStructuredGrid: SetPoint(vtkDataArray *)vtkStructuredGridReader
vtkCamera: ComputeViewPlaneNormal(): GetViewPlaneNormal(float [3])vtkCutter: GenerateValues(int, float, float)
vtkStructuredGridOutlineFilter

135Chapter 9
その他の事項
9.1

本章の概要

可視化手法以外で，知っていると便利なことなどなどを紹介する。

9.2

コマンドの組み込み

vtkCommand クラスを使うと，インタラクティブ Window に，自分で定義したコマンドを組み込むことができます (キーボードの「U」ボタンがユーザー定義コマンド用に用意されています)。例えば，「U」ボタンを押すと，表示されている画像を BMP ファイルとして保存するなどできるようになります。また，等値面を表示するプログラムで，インタラクティブに等値面レベルの変更ができるようになると，非常に便利ではないだろうか。vtkCommand は，このようなことを可能にする。

9.2.1

サンプルプログラム 1

これは，第 3 章のサンプルプログラム “ExtractData.cxx”を改造したものである。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

/* ExtractData_k.cxx */#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkExtractVOI.h>#include <vtkLookupTable.h>#include <vtkContourFilter.h>#include <vtkOutlineFilter.h>#include <vtkAppendPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkCommand.h>#include <vtkRendererSource.h>
136

18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64

#include <vtkBMPWriter.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>
vtkExtractVOI *voi;
vtkContourFilter *contour;
vtkRenderer *vren;
vtkRenderWindowInteractor *iwin;
int extent_vol[6];float range[2];float iso_value;class OurCommand:public vtkCommand{
public:static OurCommand *New (){return new OurCommand;}
virtual void Execute (vtkObject * caller, unsigned long, void *){/* 「U」が押されると実行される */static int num = 0;char FileName[128];
char c;vtkRendererSource *rs = vtkRendererSource::New ();vtkBMPWriter *bw = vtkBMPWriter::New ();printf ("Command\n");printf ("[s]napshot\n");printf ("e[x]tract data\n");printf ("[i]sosurface value :: ");scanf ("%c", &c);
fflush (stdin);switch (c){case ’s’:/* 表示されている画像を BMP ファイルとして保存 */rs->SetInput (vren);rs->WholeWindowOn ();rs->Modified ();bw->SetInput (rs->GetOutput ());
137

65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111

sprintf (FileName, "shot%d.bmp", num);bw->SetFileName (FileName);bw->Write ();num++;
break;
case ’x’:/* データの一部を切り出す */puts ("Input i1 i2 j1 j2 k1 k2");
scanf ("%d %d %d %d %d %d", &extent_vol[0], &extent_vol[1],&extent_vol[2], &extent_vol[3], &extent_vol[4],&extent_vol[5]);printf ("[%d - %d] [%d - %d] [%d - %d]\n", extent_vol[0],extent_vol[1], extent_vol[2], extent_vol[3],extent_vol[4], extent_vol[5]);voi->SetVOI (extent_vol);fflush (stdin);break;
case ’i’:/* 等値面・等高線のレベルを変える */printf ("Input Isosurface value (%f - %f) :: ",range[0], range[1]);scanf ("%f", &iso_value);
printf ("value = %f\n", iso_value);
contour->SetValue (0, iso_value);
fflush (stdin);break;
default:
puts ("No Command");break;}
iwin->Render();rs->Delete ();bw->Delete ();
}
};
138

112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158

int main (int argc, char *argv[]){
char datafile[] = "./plasma_data.vtk";
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->Update ();imgData->GetScalarRange (range);
vtkLookupTable *lut = vtkLookupTable::New ();lut->SetHueRange (0.7, 0.0);
lut->Build ();voi = vtkExtractVOI::New ();voi->SetInput (imgData);voi->SetVOI (0, 255, 0, 127, 0, 127);
contour = vtkContourFilter::New ();contour->SetInput (voi->GetOutput ());
contour->SetValue (0, 10);
contour->ComputeNormalsOn ();vtkPolyDataMapper *CMapper = vtkPolyDataMapper::New ();CMapper->SetInput (contour->GetOutput ());
CMapper->SetLookupTable (lut);
CMapper->SetColorModeToMapScalars ();CMapper->SetScalarRange (range[0], range[1]);vtkActor *CActor = vtkActor::New ();CActor->SetMapper (CMapper);vtkOutlineFilter *outline = vtkOutlineFilter::New ();outline->SetInput (imgData);
vtkOutlineFilter *outline_voi = vtkOutlineFilter::New ();outline_voi->SetInput (voi->GetOutput ());
vtkAppendPolyData *outlines = vtkAppendPolyData::New ();outlines->AddInput (outline->GetOutput ());
outlines->AddInput (outline_voi->GetOutput ());
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New ();
139

159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204

OLMapper->SetInput (outlines->GetOutput ());
vtkActor *OLActor = vtkActor::New ();OLActor->SetMapper (OLMapper);vren = vtkRenderer::New ();vren->AddActor (CActor);
vren->AddActor (OLActor);
vren->SetBackground (0.0, 0.0, 0.0);
vtkRenderWindow *renWin = vtkRenderWindow::New ();renWin->AddRenderer (vren);renWin->SetSize (500, 375);
renWin->LineSmoothingOn();
iwin = vtkRenderWindowInteractor::New ();iwin->SetRenderWindow (renWin);/* 定義したコマンドのクラスを生成して，iwin に登録 */OurCommand *oc = OurCommand::New ();iwin->AddObserver (vtkCommand::UserEvent, oc);
oc->Delete ();vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New ();iwin->SetInteractorStyle (trackball);iwin->Initialize ();iwin->Start ();reader->Delete ();lut->Delete ();outline->Delete ();CActor->Delete ();OLActor->Delete ();CMapper->Delete ();OLMapper->Delete ();voi->Delete ();contour->Delete ();outlines->Delete ();vren->Delete ();renWin->Delete ();trackball->Delete ();iwin->Delete ();return 0;

140

205

}

9.2.2

コマンド組み込みの方法

コマンドを組み込むには，vtkCommand を継承したクラスを自作することが必要である。といっても，難しいことはなく，下記の 10 行のところに，「U」が押されたら実行させたい処理を書くだけである。もちろん，クラス名は OurCommand でなくても構わない。

1
2
3
4
5
6
7
8
9
10
11
12

class OurCommand : public vtkCommand{
public:static OurCommand *New (){return new OurCommand;}
virtual void Execute (vtkObject * caller, unsigned long, void *){/* ここに，実行させたい処理を C/C++で書く */}};このようにして作ったクラスをサンプルプログラムの 177 - 180 行のようにして，vtkRednerWindowInteractor に登録すると，「U」ボタンを押されたときに，上記 10 行のところに書いた処理が実行される。
サンプルプログラム 1 では，「U」ボタンが押されると，ターミナルで 3 種類のコマンド (“s”,“x”, “i”) の入力待ち (scanf を使ってます) になり，“s”でスナップショット，“x”でデータ切り出
し (さらに切り出す範囲を入力する)，“i”で等値面のレベル変更 (レベルをさらに入力する) ができるようにしてある (Figure 9.1)。
これを使えば，かなり高機能なインタラクティブ可視化ソフトを自作できる (はず)。(bf 105
行の iwin->Render() は，再描画の指示を Window に送る。これをいれないとマウスで Window
をなぞるまで等値面の再計算などは行われない)

9.3

時間発展のデータ

いままでは，一つのデータだけを使い可視化してきた (アニメーションを作るときも)。何も工夫しなくても，時間発展のデータを扱うことができる。

9.3.1

サンプルプログラム 2:Tiny MovieMaker

このサンプルプログラムは，第 4 章のテクスチャ・マッピングを利用したボリューム・レンダリ
ングのプログラム (VolumeRendering 2d.cxx) を改造したものである。実行すると，カメラアン
グル (Azimuth 方向に 3 度ずつ) を変えつつデータも変えていく (時間発展させる)。同時に，デー
タ 1 つにつき，3 枚の画像を保存するようにしてある。サンプルデータは，同じく，地球磁気圏の温度で，本章では 10 個用意した。

1/* tinymm.cxx */141

Figure 9.1: コマンド待ちの様子
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28

#include <vtkImageData.h>#include <vtkStructuredPoints.h>#include <vtkStructuredPointsReader.h>#include <vtkPointData.h>#include <vtkSphereSource.h>#include <vtkLookupTable.h>#include <vtkOutlineFilter.h>#include <vtkImageShiftScale.h>#include <vtkPiecewiseFunction.h>#include <vtkColorTransferFunction.h>#include <vtkVolumeProperty.h>#include <vtkVolumeTextureMapper2D.h>#include <vtkVolume.h>#include <vtkPolyData.h>#include <vtkPolyDataMapper.h>#include <vtkRenderWindow.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkProperty.h>#include <vtkCamera.h>#include <vtkScalarBarActor.h>#include <vtkVectorText.h>#include <vtkFollower.h>#include <vtkRendererSource.h>#include <vtkBMPWriter.h>#define CAMERA_STEP 30
142

29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75

#define DATA_NUM

10float range[2] = { 4.2, 13.3 };int
main (int argc, char *argv[]){
int i, file_num = 0, pic_num = 0;char datafile[128];
float shift, scale;float opacity_table[256];
float color_table[3 * 256], temp_color[4];float earth_position[3] = { 100.5, 63.5, 63.5 };char FileName[128];
sprintf (datafile, "./tempe%d.vtk", 200);
vtkStructuredPointsReader *reader= vtkStructuredPointsReader::New();
reader->SetFileName(datafile);
vtkImageData *imgData;
imgData = reader->GetOutput();
imgData->SetSpacing (1.0, 1.0, 1.0);
shift = -range[0];scale = 255.0 / (range[1] - range[0]);vtkImageShiftScale *f2uc = vtkImageShiftScale::New ();f2uc->SetShift (shift);f2uc->SetScale (scale);f2uc->SetOutputScalarTypeToUnsignedChar ();f2uc->SetInput (imgData);
vtkLookupTable *lut = vtkLookupTable::New ();lut->SetHueRange (0.7, 0.0);
lut->SetNumberOfTableValues (256);
lut->SetRange (range);
lut->Build ();for (i = 0; i < 256; i++){lut->GetTableValue (i, temp_color);
color_table[i * 3] = temp_color[0];color_table[i * 3 + 1] = temp_color[1];color_table[i * 3 + 2] = temp_color[2];
}
143

76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122

vtkColorTransferFunction *tf4color= vtkColorTransferFunction::New ();tf4color->SetColorSpaceToRGB ();tf4color->BuildFunctionFromTable (0, 255, 256, color_table);for (i = 0; i < 256; i++){
if (i > 187)
opacity_table[i] = i / 255.0 * 0.025;
else
opacity_table[i] = 0.0;}
vtkPiecewiseFunction *tf4opacity= vtkPiecewiseFunction::New ();tf4opacity->BuildFunctionFromTable (0, 255, 256, opacity_table, 1);
vtkVolumeProperty *vp = vtkVolumeProperty::New ();vp->SetColor (tf4color);
vp->SetScalarOpacity (tf4opacity);
vp->SetInterpolationTypeToLinear ();vtkVolumeTextureMapper2D *vMapper= vtkVolumeTextureMapper2D::New ();vMapper->SetInput (f2uc->GetOutput ());
vMapper->SetMaximumNumberOfPlanes (256);vMapper->SetTargetTextureSize (256, 128);
vtkVolume *Volume = vtkVolume::New ();Volume->SetMapper (vMapper);Volume->SetProperty (vp);
vtkOutlineFilter *outline = vtkOutlineFilter::New ();outline->SetInput (imgData);
vtkPolyDataMapper *OLMapper = vtkPolyDataMapper::New ();OLMapper->SetInput (outline->GetOutput ());
vtkActor *OLActor = vtkActor::New ();OLActor->SetMapper (OLMapper);
vtkSphereSource *earth = vtkSphereSource::New ();earth->SetRadius (3.3);

144

123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169

vtkPolyDataMapper *EMapper = vtkPolyDataMapper::New ();EMapper->SetInput (earth->GetOutput ());
vtkActor *EActor = vtkActor::New ();EActor->SetMapper (EMapper);
EActor->SetPosition (earth_position);
EActor->GetProperty ()->SetColor (0.0, 0.0, 1.0);
vtkVectorText *EText = vtkVectorText::New ();EText->SetText ("EARTH");
vtkPolyDataMapper *TMapper = vtkPolyDataMapper::New ();TMapper->SetInput (EText->GetOutput ());
vtkFollower *TActor = vtkFollower::New ();TActor->SetMapper (TMapper);TActor->SetScale (5.0, 5.0, 5.0);
TActor->AddPosition (earth_position);
vtkScalarBarActor *scalarbar = vtkScalarBarActor::New ();scalarbar->SetTitle ("Temperature");scalarbar->SetLookupTable (lut);
scalarbar->SetOrientationToVertical ();scalarbar->SetWidth (0.075);
scalarbar->SetHeight (0.75);
vtkCamera *camera = vtkCamera::New ();camera->SetPosition (-200.0, -200.0, 237.5);
camera->SetFocalPoint (127.5, 63.5, 63.5);
camera->SetViewUp (0.0, 0.0, 1.0);
camera->OrthogonalizeViewUp ();camera->SetClippingRange (30.0, 2000.0);
vtkRenderer *vren = vtkRenderer::New ();vren->AddVolume (Volume);
vren->AddActor (EActor);
vren->AddActor (OLActor);
vren->AddActor (TActor);
vren->AddActor (scalarbar);vren->SetActiveCamera (camera);
vren->SetBackground (0.0, 0.0, 0.0);
TActor->SetCamera (vren->GetActiveCamera ());
vtkRenderWindow *renWin = vtkRenderWindow::New ();renWin->AddRenderer (vren);renWin->SetSize (500, 375);
145

170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216

vtkRendererSource *rs = vtkRendererSource::New ();vtkBMPWriter *bw = vtkBMPWriter::New ();for (i = 0; i < CAMERA_STEP; i++){
if (i % (CAMERA_STEP / DATA_NUM) == 0 && i != 0){/* 次のデータを読み込む */file_num++;
sprintf (datafile, "./tempe%d.vtk", 200 + file_num * 5);reader->SetFileName(datafile);
reader->Update();printf ("Data Read %d\n", 200 + file_num * 5);}
renWin->Render ();rs->SetInput (vren);rs->WholeWindowOn ();rs->Modified ();sprintf (FileName, "tm%d.bmp", pic_num);bw->SetInput (rs->GetOutput ());
bw->SetFileName (FileName);bw->Write ();pic_num++;
camera->Azimuth (3.0); /* カメラ位置を 3 度移動 */}reader->Delete ();f2uc->Delete ();lut->Delete ();tf4color->Delete ();tf4opacity->Delete ();vp->Delete ();vMapper->Delete ();Volume->Delete ();outline->Delete ();OLMapper->Delete ();
146

217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232

OLActor->Delete ();earth->Delete ();EMapper->Delete ();EActor->Delete ();EText->Delete ();TMapper->Delete ();TActor->Delete ();scalarbar->Delete ();vren->Delete ();camera->Delete ();renWin->Delete ();rs->Delete ();bw->Delete ();return 0;}

9.3.2

サンプルプログラム 2 の簡単な解説

177 - 189 行で，データを更新しているだけで，なにも新しいことはない (その他，(時系列すべてを通した) スカラー値の最大値・最小値も予め与えてある)。これだけのプログラムで，MovieMaker(3)の真似事をさせることができた。表示は，Figure 9.2 (185 行を忘れずに。手動でデータを作る

Figure 9.2: 時間発展するデータの可視化
ときは，例えば farray( vtkFloatArray ) にデータを入れているのであれば，データを読み込んだ
あとに，farray->Modified() を実行して，VTK にデータが変化したことを知らせなければなら
ない)

147

9.4

ポリゴンデータ

可視化像と同時に線などを表示したい場合がある ( 例えば，気圧の等高線と同時に，低気圧の(これまでの) 進路などを表示させる)。矢印や円柱など，もともと VTK に用意されているものについては，それを使えばよいが，用意されていない場合，vtkPolyData を利用し自作すること
もできる。

9.4.1

ポリゴンデータの作成

点・連結した線分・ポリゴン・連結した三角形のデータ (Figure 9.3) を格納できる。これらは，OpenGL の GL POINTS, GL LINE STRIP, GL POLYGON, GL TRIANGLE STRIP に相当する。

Figure 9.3: ポリゴンの種類
次のサンプルプログラムは，ポリゴンデータ (vtkPolyData) をマニュアル作成する例である。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

/* pyramid.cxx */#include <vtkCellArray.h>#include <vtkFloatArray.h>#include <vtkPolyData.h>#include <vtkLookupTable.h>#include <vtkPoints.h>#include <vtkPointData.h>#include <vtkPolyDataMapper.h>#include <vtkActor.h>#include <vtkRenderer.h>#include <vtkRenderWindow.h>#include <vtkRenderWindowInteractor.h>#include <vtkInteractorStyleTrackballCamera.h>const int n_points = 9;
const int n_polys = 5;
const int n_lines = 6;
const int dp = 8; /* dummy point id */int main( int argc, char *argv[] ){int i;
float vpoints[n_points][3]={
148

25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71

{-2.5, -2.5, -2.5}, { 2.5, -2.5, -2.5},{ 2.5, 2.5, -2.5}, {-2.5, 2.5, -2.5},{-2.5, -2.5, 2.5}, { 2.5, -2.5, 2.5},{ 2.5, 2.5, 2.5}, {-2.5, 2.5, 2.5},{ 0.0, 0.0, 0.0}};
int poly_nverts[n_polys] = {3, 3, 3, 3, 4};
vtkIdType poly_pts[n_polys][4] = {{8, 0, 1, dp}, {8, 1, 2, dp},
{8, 2, 3, dp}, {8, 3, 0, dp},
{0, 3, 2, 1}};
int line_nverts[n_lines] = {5, 5, 2, 2, 2, 2};
vtkIdType line_pts[n_lines][5]={{0, 1, 2, 3, 0}, {4, 5, 6, 7, 4},{0, 4, dp, dp, dp}, {1, 5, dp, dp, dp},
{2, 6, dp, dp, dp}, {3, 7, dp, dp, dp}};

float scalars[n_points] = {5.0, 0.0, 5.0, 0.0, 7.5, 7.5, 7.5, 7.5, 10.0};

vtkPoints *points = vtkPoints::New();for (i=0;i<n_points;i++)points->InsertPoint(i, vpoints[i]);vtkCellArray *polys = vtkCellArray::New();for (i=0;i<n_polys;i++)polys->InsertNextCell(poly_nverts[i], poly_pts[i]);vtkCellArray *lines = vtkCellArray::New();for (i=0; i<n_lines;i++)lines->InsertNextCell(line_nverts[i], line_pts[i]);vtkFloatArray *farray = vtkFloatArray::New();
farray->SetArray(scalars, n_points, 1);
vtkPolyData *pyramid = vtkPolyData::New();
pyramid->SetPoints(points);
pyramid->SetPolys(polys);pyramid->SetLines(lines);
pyramid->GetPointData()->SetScalars(farray);vtkLookupTable *lut = vtkLookupTable::New();
lut->SetHueRange(0.7, 0.0);
149

72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115

lut->Build();
vtkPolyDataMapper *Mapper = vtkPolyDataMapper::New();
Mapper->SetInput(pyramid);
Mapper->SetLookupTable(lut);
Mapper->SetScalarRange(0.0, 10.0);
vtkActor *Actor = vtkActor::New();
Actor->SetMapper(Mapper);
vtkRenderer *ren = vtkRenderer::New();
ren->AddActor(Actor);
ren->SetBackground(0,0,0);
vtkRenderWindow *renWin = vtkRenderWindow::New();
renWin->AddRenderer(ren);renWin->SetSize(500, 500);
vtkInteractorStyleTrackballCamera *trackball =vtkInteractorStyleTrackballCamera::New();
vtkRenderWindowInteractor *iwin= vtkRenderWindowInteractor::New();
iwin->SetRenderWindow(renWin);iwin->SetInteractorStyle(trackball);iwin->Initialize();iwin->Start();points->Delete();polys->Delete();lines->Delete();farray->Delete();pyramid->Delete();lut->Delete();Mapper->Delete();Actor->Delete();ren->Delete();renWin->Delete();trackball->Delete();iwin->Delete();return 0;}
このサンプルプログラムの実行結果は，Figure 9.4 である
ポリゴンデータを作成するには，以下 (の作業) が必要である。

150

Figure 9.4: マニュアル作成したデータ ((左)ScalarVisibilityOff())
1. 頂点のデータ
2. ポリゴンを構成する頂点の結び方の指定 (4 種類ある:Figure 9.3 参照)
3. 頂点に付随するスカラー，ベクトル，法線ベクトル (なくとも良い)
以下，順に見ていく。
1 : 本サンプルプログラムでは，頂点のデータの設定は，24 - 29, 49 - 51, 65 行でおこなっている。vtkPoints を介して，vtkPolyData に代入している。これは，流線の出発点の設定のところでもおこなっている。こうして格納した頂点のデータは，Figure 9.5 のように，代入した順番に番号がふられる。
2 : 上で格納した点をどのように結びポリゴンを作るのか，を指定する。例えば，0 番と 4 番を結び線を引く，0 番・1 番・8 番の頂点を結び 3 角形のポリゴンを作る，など。サンプルプログラムでの関係箇所は，31- 41, 53 - 59, 66 - 67 行。ポリゴンのデータ (頂点の結び方) を作るには，vtkCellArrayを使用する。vtkCellArray のデータ構造は，[n0 , pid00 , pid01 , pid02 ..., n1 , pid10 , pid11 , pid12 , ..] のような整数の羅列である (nx , pidy ... は，すべて整数)。n はセルの配列の数 (ポリゴンの頂点の数)，「0 番と 4 番を結び線を引き」たい場合は [2, 0, 4]，「0 番・
pid は対応する頂点の id。上の例で，1 番・8 番の頂点を結び 3 角形のポリゴンを作」りたい場合は，[3, 0, 1, 8] という整数が入った
vtkCellArray を作ればよい。53 - 59 行で，この作業を行っている。InsertNextCell の一つ目の引数が n であり，二つ目の引数は id が入った配列のポインタである。
このようにして作成した vtkCellArray を 66 - 67 行で，vtkPolyData に代入している。代入す
るときに，ポリゴンの種類によって，SetVerts, SetLines, SetPolys, SetStrips を使い分ける。
サンプルプログラムでは，8・0・1，8・1・2，8・2・3，8・3・0 を結び三角形を，0・1・2・3(Figure 9.5) を結び四角形のポリゴンをを作っている。また，このポリゴンを囲むように線を引いている。
3 : 頂点に付随するスカラーやベクトル，法線ベクトルもセットできる。サンプルプログラム
では，スカラーをセットしている。44 - 46, 61 - 62, 68 行でこの作業を行っている。方法は，
vtkImageData や vtkRectilinearGrid などの場合と，まったく同様である 1 。
1

法線ベクトルの場合は，ベクトルの大きさを 1.0 にしなければならない。

151

Figure 9.5: 頂点と頂点の番号

このようにして作成したポリゴンデータは，可視化のためのフィルタからの出力のポリゴンデー
タと同様に vtkPolyDataMapper を通して表示できる。

9.5

アンチ・エイリアジング

サンプルプログラム 1 の 172 行目のように，...SmoothingOn() とするとアンチ・エイリアジン
グがかかり，画像が少し綺麗になる。点・線・ポリゴンについて，それぞれ下記のメソッドが
vtkRenderWindow に用意されている。

• PointSmoothingOn();• LineSmoothingOn();• PolygonSmoothingOn();

あとがき
VTK の可視化ライブラリを使えば，可視化アルゴリズムも OpenGL もほとんど知らなくとも，
かなり良くできた 3 次元可視化ツールを開発することができます (実際，本章では MovieMakerの真似事をやりました)。これまでの内容を基に，自分だけの可視化ツールを開発してください。
健闘を祈ります。なお，本勉強会で使用したサンプルプログラムの作成に当たり，Vtisualization
Toolkit (3rd Edition)， VTK User’s Guide や附属 CD に収録されているサンプルプログラム (主に Tcl) を参考にしました。
長い間，おつき合いありがとうございました。

152

Acknowledgments
サンプルデータとして使用した球状トカマクの圧力データは，核融合科学研究所 理論・データ
解析センター 水口直紀博士から提供していただいたシミュレーションデータを加工したものです。また，地球磁気圏の温度データは，名古屋大学 太陽環境研究所 荻野龍樹教授から提供していただいたシミュレーションデータを加工したものです。

参考文献(1) Will Schroeder, Ken Martin and Bill Lorensen: Visualization Toolkit(3rd Edition), Kitware (2002)(2) VTK User’s Guide, Kitware (2003)(3) H. UEHARA, S. KAWAHARA, N. OHNO, M. FURUICHI, F. ARAKI and A. KAGEYAMA:”MovieMaker: A Parallel Movie-Making Software for Large Scale Simulations”, J. Plasma Physics
Vol.72, part 6, pp.841-844 (2006)

153