多視点画像で三次元形状を認識する

本記事では三次元形状認識を題材とした、外観検査アルゴリズムコンテスト2020「X線CTによる工業製品の内部検査～ボクセルデータからの形状と材質の分類～」[^1]への提出作品[^2]の解説を行います[^注1]。PythonでのVTK（Visualization Toolkit）による簡単なボリュームレンダリングの実装例も紹介します。

## 本記事の要約
- 三次元形状の認識を高速に行う手法としてMVCNN（Multi-view CNN）[^3]を試しました。形状が既知の三次元形状の識別は、多視点画像を用いた手法で比較的高速に精度よく処理可能なことが分かりました。
- 未知形状への対応のため、クラス識別タスクではなく、距離学習による照合タスクとしました。
- 画像処理による事前の部品分離精度、および学習データの品質が正解率に寄与していると考えられました。
- 1データセットから250個程度の部品を抽出し、12種類の部品個数の計数にかかった時間はCPU実行で2分半ほどでした。
- VTKを利用して、ボリュームレンダリングによる可視化結果から多視点画像を作成しました。ボリュームレンダリングの方法もご紹介します。

## 外観検査アルゴリズムコンテスト2020課題
![課題](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/4e7dc029-b8b6-3fe0-484c-325260a78fe4.png "外観検査アルゴリズムコンテスト2020課題")
2020年の課題は「X線CTによる工業製品の内部検査～ボクセルデータからの形状と材質の分類～」でした。 ビーカーのような容器の中に、図の左上のような部品が多数入っており、それをX線CT撮影した断層画像群が5セット、提供されました。さらに、容器の中に入っている部品のうち、12種類については、部品単体でX線CT撮影した断層画像群も提供されました。 課題は、12種類の部品が、容器の中にそれぞれ幾つ入っているかを計量するものでした。容器の中には、12種類の部品以外の部品も入っており、12種類の中には剛体ではない（形状が変わる）ものも含まれていました。

## 事前検討

### データの観察
まず、事前検討として、提供されたデータを眺めて、以下の考察を得ました。
- 部品入り容器をX線CT撮影した各提供データセットにおいて、画像の輝度値は各々で正規化されていて、同じ素材や空気でも、輝度値が異なっていました。すなわち、輝度値から材質を推定するのは困難と考えられました。
![輝度値の違い](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/149be431-b077-c619-7dde-50e9df193fa4.png "データセットごとの輝度値の違い")

- ひときわ明るく映っている金属部品があり、その周囲に金属由来のアーチファクトが発生している影響で周囲の部品の境界が分かりづらくなっていました。（下図中白丸内）
![アーチファクト](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/02be11ce-5b56-b77b-c545-d7852f50e599.png "アーチファクト例")

- 1つの閾値で輝度値の低い部品と背景を分けるのは難しそうでした。以下ヒストグラムを見ると、一番大きいピークの部分は背景（空気部分）だと考えられますが、低輝度値の前景（部品部分）のピークの1つが隣接しています。
![ヒストグラム](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/19ce926e-73ca-c311-e48d-bf9162e215c5.png "輝度値ヒストグラム")

- 1つのデータセットあたり、概ねタテ1000px、ヨコ1000px程度の画像が1500枚程であり、全量をメモリ上に読み込むとメモリが貧弱なPCでは処理が難しそうでした。評価用PCスペックとして公開されていたものがメモリ16GBだったのですが、浮動小数点型を扱う三次元ガウシアンフィルタのような処理も適用したため、メモリ16GBのPCではプログラム製作途中にbad_alloc例外で頻繁に落ちてしまいました。余談ではありますが、省メモリな処理とするために単精度浮動小数点型を使用、配列を使いまわす等実装上の苦労もありました。

### 三次元形状マッチングの手法
本課題における三次元形状マッチングの手法の候補として、以下のようなものを事前に検討しました。前述のとおり、データセットごとに異なる輝度値の正規化が行われていたため、本作品では材質を元にした部品の識別は考慮せず、形状のみを考慮した識別を行うことを考えました。
1. 断層画像のマッチング
2. ボリュームデータのマッチング
3. 深度画像（Depth map）でのマッチング
4. 多視点ボリュームレンダリング画像でのマッチング
5. 多視点深度画像（Depth map）でのマッチング
6. ボリュームデータ表面の点群同士のマッチング


人間が容器の中身を取り出して部品の個数を数える速さを（勝手に）目標にしていましたので、出来るだけ高速な手法が望ましいです。
人間が部品の種類別に個数を数える場合、外側の形状を複数の方向から見て部品の種類を判断します。同様の発想で形状のマッチングを行うとすれば、複数視点からの深度情報付の画像を利用することが良さそうであると考えました。加えて、画像であればボクセルデータほどの大きな情報量を取り扱わなくとも、マッチングが可能であり、計算機リソース面、計算時間面でも利点があります。

![部品の識別](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/beadee9c-2570-4234-09a8-e1cbc62159bd.png "部品の識別")

当初5、または4と5の併用として検討していましたが、手元の開発環境の中でもVTKで深度画像が作れる環境と作れない環境があり、評価環境での動作に不安があったため最終的には4の**多視点ボリュームレンダリング画像でのマッチング**としました。

### 多視点画像
本作品ではボクセルデータをボリュームレンダリングして作成した多視点画像を利用しました。VTKを利用した簡単な多視点画像の作り方について本記事の後ろの方でご紹介（記事内リンク：[ボリュームレンダリングと多視点画像作成](#ボリュームレンダリングと多視点画像作成)）しています。

## 手法概要
下図に、本作品の処理の流れを示しています。まず、X線CTによるボクセルデータを部品ごとに切り離したボクセルデータにした後、部品単位で形状特徴量を計算し、特徴量を比較することで参照部品と同じ部品かどうかを判断します。

![手法概要](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/94f2eb32-5ed4-4f4c-a461-bcc1fe443db6.png "手法概要")

### 1. 部品分離
二値化画像処理と三次元ラベリング処理で部品ごとの分離を行いました。 通常の二値化処理を行っただけでは多くの部品がくっついた状態になってしまいましたので、できるだけ部品同士を切り離すために三次元ガウシアンフィルタでぼかして閾値処理を行い、部品間の接合部分が小さいものは切り離す処理を行いました。この処理で切り離すことが難しかったもの、あるいは金属由来のアーチファクトの影響で部品形状が一部欠けてしまったものが、最終的に誤認識の原因の１つになってしまいました。

### 2. 三次元形状特徴量計算
三次元特徴量として、多視点画像を利用した三次元形状認識手法である、MVCNNを利用して形状特徴量の計算を行いました。本課題の提供データ中の部品の中には大小さまざまなものが含まれていましたが、多視点画像はすべて同じサイズで作成しています。このため、本処理で抽出する特徴量は、大きさの情報が失われ、形状のみに依存するものになっています。

#### MVCNN
MVCNNのモデル構造を下図に示します。部品ごとのボクセルデータから作成した、16視点の多視点画像を入力として、最終層の出力を128次元の特徴量としました。

![MVCNN](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/c9927e9b-9905-c3e4-6a1e-82ef536c998b.png "MVCNNモデル構造")
MVCNNは元々、多視点画像による三次元形状認識（クラス識別）手法です。本課題において提供された12個の参照部品に関して試しに学習したところ、部品分離処理が成功したものについてのクラス識別ではかなり高精度での識別が可能でした。 ただし、本課題で提供されたデータの中に含まれていない未知の部品があった場合には、既知のどれかのクラスとして識別するような誤識別が発生しうるため、本作品では距離学習を行い、形状特徴量の比較によって部品種別を判定することにしました。

出力層を128次元とし、損失関数にTripletLoss[^4]を採用しました。 提供データから切り出した部品は全部で30種類強ありましたが、筆者の目で見て形状の近いものをまとめて21クラスに分類しました。学習は128次元の特徴量空間上で、同じクラス同士のユークリッド距離は小さく、異なるクラス同士のユークリッド距離は大きくなるよう学習を行いました。 データ拡張として、多視点画像の視点違いと部品の色味（濃淡）違いを採用しました。

### 3. 特徴量比較
ある部品について、参照部品と形状が似通っているものほど二つの形状特徴量のユークリッド距離が小さくなるため、形状特徴量間のユークリッド距離に（経験的な）閾値を設けて形状を判断しました。 加えて、形状特徴量は大きさの情報を失っていますので、似通った形状で大きさが違うものは判別がつかないことがあるために、部品ごとに切り出した時点での体積についても加味して部品の判定を行いました。

## 結果と考察
150～250個程度の部品を含む、提供された5つのデータセットについて、CPU（Corei7-6700K）で1分～2分半程度の時間で処理できました。 判定を誤ったものとしては、部品の分離に失敗しているもの、学習データ作成時に筆者に見分けがつかなかったもの等がありました。
画像処理による部品の分離にはまだ改善の余地があったと思います。精度向上のためには深層学習手法の適用も視野に入れた方が良いと思いました。また、学習データの品質は大事ですので、適切に分類しましょう、ということではありますが、機械学習モデルの出力をそのまま利用するのではなく、後処理で体積を加味した判定を行っていることで、アノテーション誤り等によるモデル（学習データ）の不出来をカバーできていた面もありますので、前処理、後処理はとても重要です。

![誤判定の例](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/fc12ddd0-1eac-b1a3-0346-5de7f95d84dc.png "誤判定の例")

## ボリュームレンダリングと多視点画像作成
ここまでにご説明した作品では、ボクセルデータから多視点画像を生成しています。 これはVTKを用いて実現しています。PythonやC++、最近ではJavaScriptなどからVTKを使って、ボリュームレンダリングによる可視化や多視点画像作成が手軽に出来ます。本記事では、Pythonでの実装例をご紹介します。[^注2]手元ではPython 3.7.4、vtk 9.0.1、Windows10で動作していました。（2020年当時のバージョンです）
参考までに、インストールコマンドを下記に示します。サンプルスクリプトで使用しているので、OpenCVとNumpyも一緒にインストールするものです。

~~~ console
> pip install opencv-python numpy vtk
~~~

### ボリュームレンダリングによる可視化
VTKWithNumpyのサンプル[^5]を元にしています。
200px*200pxの画像を200枚読み込んで、VTKでボリュームレンダリングする例です。render関数はちょっと長いですが、スクリプト中のコメントをご参照ください。

~~~ python
import vtk
import numpy as np
import cv2
from vtk.util.numpy_support import vtk_to_numpy

def render(volume_data):
    """
    volume_dataは三次元のnumpy.ndarray(type:numpy.uint8)
    """
    # Numpyの配列をVTKの画像形式へ変換
    dataImporter = vtk.vtkImageImport()
    data_string = volume_data.tobytes()
    dataImporter.CopyImportVoidPointer(data_string, len(data_string))

    dataImporter.SetNumberOfScalarComponents(1)
    # 各次元のサイズを設定
    dataImporter.SetDataExtent(0, volume_data.shape[2]-1, 0, volume_data.shape[1]-1, 0, volume_data.shape[0]-1)
    dataImporter.SetWholeExtent(0, volume_data.shape[2]-1, 0, volume_data.shape[1]-1, 0, volume_data.shape[0]-1)

    # 可視化時の透け具合を設定
    # 輝度値50以下は透明
    alphaChannelFunc = vtk.vtkPiecewiseFunction()
    alphaChannelFunc.AddPoint(50, 0.0)
    alphaChannelFunc.AddPoint(100, 0.05)
    alphaChannelFunc.AddPoint(200, 0.2)

    # 可視化の色を設定
    # 輝度値80以下が赤、120が緑、170以上は青
    colorFunc = vtk.vtkColorTransferFunction()
    colorFunc.AddRGBPoint(80, 1.0, 0.0, 0.0)
    colorFunc.AddRGBPoint(120, 0.0, 1.0, 0.0)
    colorFunc.AddRGBPoint(170, 0.0, 0.0, 1.0)

    volumeProperty = vtk.vtkVolumeProperty()
    volumeProperty.SetColor(colorFunc)
    volumeProperty.SetScalarOpacity(alphaChannelFunc)

    volumeMapper = vtk.vtkFixedPointVolumeRayCastMapper()
    # 多視点画像生成で深度画像を生成したい場合は以下
    # volumeMapper = vtk.vtkGPUVolumeRayCastMapper()

    volumeMapper.SetInputConnection(dataImporter.GetOutputPort())

    # ボリュームデータにMapperとPropertyを反映
    volume = vtk.vtkVolume()
    volume.SetMapper(volumeMapper)
    volume.SetProperty(volumeProperty)

    # 可視化のウィンドウを設定
    renderer = vtk.vtkRenderer()
    renderWin = vtk.vtkRenderWindow()
    renderWin.AddRenderer(renderer)
    renderInteractor = vtk.vtkRenderWindowInteractor()
    renderInteractor.SetRenderWindow(renderWin)

    renderer.AddVolume(volume)
    renderer.SetBackground(0,0,0)
    renderWin.SetSize(400, 400)

    # 表示開始
    renderInteractor.Initialize()
    renderWin.Render()
    renderInteractor.Start()

if __name__ == '__main__':
    # 画像のサイズ、枚数、ファイル名を設定
    img_width = 200
    img_height = 200
    num_img = 200
    img_name_prefix = 'img'
    volume_data = np.zeros((num_img, img_height, img_width), np.uint8)
    # OpenCV等で画像を読み込み
    for i in range(num_img):
        path = img_name_prefix + str(i) + '.png'
        img = cv2.imread(path)
        volume_data[i,:,:] = img[:,:,0]
    # 可視化！
    render(volume_data)
~~~

これだけです。下図のようなウィンドウが現れて、ぐるぐる回して見られます！

![VTK可視化ウィンドウ1](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/bac7ae69-f46a-3d8f-2bcc-a0ac324e944a.png "VTK可視化ウィンドウ1") ![VTK可視化ウィンドウ2](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/1d023ebe-cde4-f383-a6d8-e241d260120f.png "VTK可視化ウィンドウ2")

### 多視点画像として保存
多視点画像は、上記の可視化ウィンドウの中に見えている可視化結果を、カメラの視点を変えて保存するだけで作成可能です。 下のスクリプトは1視点分の画像保存の例です。
筆者が本作品で採用をあきらめた、深度画像生成についてもコメントアウト部分を活かせば動きます。

~~~ python
######################################################
# 上のスクリプトのrender関数内の”#　表示開始”コメント以降の3行を以下で置き換える
######################################################
    volumeMapper.RenderToImageOn()

    # カメラの角度を設定
    el_angle = 10
    azi_angle = 0

    renderer.ResetCamera()
    renderer.UseShadowsOn()
    camera = renderer.GetActiveCamera()
    camera.Elevation(el_angle)
    camera.Azimuth(azi_angle)

    renderWin.Render()

    img = vtk.vtkImageData()
    volumeMapper.GetColorImage(img)
    # 深度画像も生成できる（はず）
    # volumeMapper.GetDepthImage(img)
    rows, cols, _ = img.GetDimensions()
    sc = img.GetPointData().GetScalars()
    # VTKImageをnumpy.arrayに変換
    out_img = vtk_to_numpy(sc)
    out_img = out_img.reshape(rows, cols, -1)
    out_img = out_img[:,:,[2,1,0]] # RGB->BGR

    cv2.imwrite("out.png", out_img)
~~~
上記スクリプト中の「カメラの角度を設定」のコメント部分の角度の値を変更することで、別視点の画像となります。下の図はぐるっと10視点分の画像で作成しました。
今回はボクセルデータなのでVTKを利用した多視点画像作成ですが、点群データやメッシュデータならPCL[^6]かOpen3D[^7]を利用するほうが便利かもしれません。

![10視点アニメ](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2990701/c9bc30a9-fb7b-cfd6-3ffc-c11f45f58c14.gif "10視点アニメ")

## おわりに
お読みいただきありがとうございました。何かのお役に立てていただければ幸いです。
以上、情報通信研究部　水谷麻紀子でした。

[^1]:http://alcon.itlab.org/
[^2]:MV-CNNを用いた三次元類似形状検索～ボクセルデータを多視点画像により形状比較することで高速・省メモリに形状比較～, 水谷麻紀子, 外観検査アルゴリズムコンテスト2020
[^3]:Su, Hang, et al. "Multi-view convolutional neural networks for 3d shape recognition." Proceedings of the IEEE international conference on computer vision. 2015.
[^4]:Schroff, Florian, Dmitry Kalenichenko, and James Philbin. "Facenet: A unified embedding for face recognition and clustering." Proceedings of the IEEE conference on computer vision and pattern recognition. 2015.
[^5]:https://kitware.github.io/vtk-examples/site/Python/Utilities/VTKWithNumpy/
[^6]:https://pointclouds.org/
[^7]:http://www.open3d.org/

[^注1]:本作品は、2020年の外観検査アルゴリズムコンテストで最優秀賞をいただきました。
[^注2]:なおコンテスト作品では、MVCNNの学習プログラムはPythonで実装したものの、提出作品は部品分離、多視点画像作成から推論処理含めてすべてC++で実装しました。C++でもほぼ同じような実装で多視点画像が作成できます。