点群データが熱い

これからはセンサーの時代です！（？）
ということでPythonでOpen3Dで遊んでみます。

本当にこの Open3D というライブラリすごいです。
こういうのを作り公開されてくださる方々に心から感謝します……。

Open3Dを使って手持ちの点群で遊んでみる

ライブラリの準備

ライブラリ（Open3D）をインストールする。

> pip install open3d

データの準備

今回用意したのははシンプルな通称 .xyz 形式のファイルです。
X Y Z が並んでいるだけです。

0.365193 0.609847 0.426132
0.320558 0.615929 0.409627
0.171359 0.759521 0.636285
    :

最初の一歩で何も考えずに動かしてみる

最小サンプルを動かして確認してみます！

import open3d as o3d

pcd = o3d.io.read_point_cloud("vase.xyz")
print(pcd)

o3d.visualization.draw_geometries([pcd])

これだけで点群が3D表示されて楽しい！
マウスで画面をぐりぐりできて、Ctrl押しながらで平行移動、Shiftおしながらで回転。Escで終了。

クラスタリングとかしてみる

少し欲が出てきたからクラスタリング（点の集まりを見つける）してみる。
上手くできない場合は eps と min_points の値をいじると良いです。

import open3d as o3d
import numpy as np
import random

pcd = o3d.io.read_point_cloud("test.xyz")
print(pcd)

with o3d.utility.VerbosityContextManager(o3d.utility.VerbosityLevel.Debug) as cm:
    labels = np.array( pcd.cluster_dbscan(eps=0.7, min_points=5, print_progress=True))
print(labels)

max_label = labels.max()
print(f"point cloud has {max_label + 1} clusters")

colors_array = []
for i in range(max_label+1):
    colors_array.append((random.random(),random.random(),random.random()))

colors = []
for i in range(len(labels)):
    if labels[i]==-1:
        colors.append((0.9,0.9,0.9))
    else:
        colors.append(colors_array[labels[i]])
colors = np.array(colors)

pcd.colors = o3d.utility.Vector3dVector(colors)

o3d.visualization.draw_geometries([pcd])

花瓶の中の何かがクラスタリングされて浮き出て色がかわっています！

面を割り出してみる

点群だけだと面白くないから面を見つけたい。
そのために法線（面に対して垂直にのびる線）を見つける。

コードを実行したら黒い線がちゃんと面の法線方向に向いてますか？

import open3d as o3d
import numpy as np
import random

pcd = o3d.io.read_point_cloud("vase.xyz")
print(pcd)

size = np.abs((pcd.get_max_bound() - pcd.get_min_bound())).max() / 10
pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius=size, max_nn=50))

o3d.visualization.draw_geometries([pcd],point_show_normal=True)

わあ、、、ちょっと気持ち悪い、、、
でもちゃんと法線を見つけてくれています。

法線が見つかったら面を作る。

import open3d as o3d
import numpy as np
import random

pcd = o3d.io.read_point_cloud("vase.xyz")
print(pcd)

size = np.abs((pcd.get_max_bound() - pcd.get_min_bound())).max() / 10
pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius=size, max_nn=50))

oboxes = pcd.detect_planar_patches(
    normal_variance_threshold_deg=60,
    coplanarity_deg=75,
    outlier_ratio=0.75,
    min_plane_edge_length=0,
    min_num_points=0,
    search_param=o3d.geometry.KDTreeSearchParamKNN(knn=100))

print("Detected {} patches".format(len(oboxes)))

geometries = []
for obox in oboxes:
    mesh = o3d.geometry.TriangleMesh.create_from_oriented_bounding_box(obox, scale=[1, 1, 0.0001])
    mesh.paint_uniform_color(obox.color)
    geometries.append(mesh)
    geometries.append(obox)
geometries.append(pcd)

o3d.visualization.draw_geometries(geometries)

そりゃあ局面の物体から面を作ったらこうなりますよね、、、

メッシュを作ってみる

面ではなくメッシュ（いわゆるポリゴン）を割り出してみます。

import open3d as o3d
import numpy as np

pcd = o3d.io.read_point_cloud("vase.xyz")
print(pcd)

size = np.abs((pcd.get_max_bound() - pcd.get_min_bound())).max() / 10
pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius=size, max_nn=50))

pcd.orient_normals_consistent_tangent_plane(10)

distances = pcd.compute_nearest_neighbor_distance()
avg_dist = np.mean(distances)
radius = 2*avg_dist   
radii = [radius, radius * 2]
recMeshBPA = o3d.geometry.TriangleMesh.create_from_point_cloud_ball_pivoting(pcd, o3d.utility.DoubleVector(radii))
    
o3d.visualization.draw_geometries([recMeshBPA])

すごい……。

ひと通り遊んでみました！

簡単です。
スゴイです。

手元では手持ちのセンサーでとったデータを使っていたのでノイズもあったり等こんなチュートリアルみたいに簡単にはいきませんが、そういうの込みでこれからが熱い界隈だと感じました！！