点群の法線をPCAで求めたら37%が裏返っていた。スタンフォードバニーでMST伝播による向き揃えを試す
約14分で読めます

点群の法線をPCAで求めたら37%が裏返っていた。スタンフォードバニーでMST伝播による向き揃えを試す


CG(コンピュータグラフィックス)の世界には「hello world」に当たる点群データがある。スタンフォードバニー、通称うさぎだ。1994年にスタンフォード大学でウサギの置物をレーザースキャンして作られたこのデータは、以来30年間、点群処理・メッシュ復元・法線推定といったアルゴリズムのベンチマークとして使われ続けている。

今回はこのバニーの生の点群(XYZ座標だけ、メッシュのつながりは使わない)から、表面の向き(法線)を推定する実験をした。やり方自体は教科書的——各点の近傍をPCA(主成分分析)にかけて、分散が一番小さい方向を法線とする——なのだが、これには昔からよく知られた落とし穴がある。符号の不定性だ。実際にやってみると、本当に隣り合う点同士で法線があっさり裏返った。

先に結果を見てもらうのが早い。

法線のxyz成分をそのままRGBの色に変換して(法線マップと同じ考え方)点群に塗ってみた図だ。表面はなめらかに湾曲しているはずなのに、色がまだら模様に荒れている。これが「符号の不定性」の見た目だ。

そもそも法線とは何か、何ができたら「成功」なのか

具体的な手順に入る前に、この記事の土台になっている言葉を整理しておく。法線とは、表面のある場所に垂直に立てた矢印のことだ。地面に旗を立てたときのポールの向きを想像すると近い。メッシュ(三角形の集まり)なら各三角形の向きから法線は自明に計算できるが、今回扱うのはXYZ座標が8171個並んでいるだけの生の点群で、「どの点とどの点が同じ面を作っているか」という情報がそもそも存在しない。だから、各点のまわりの点の散らばり方から「ここの表面はだいたいこの向きだろう」と推定してやる必要がある。

なぜそんな矢印が欲しいのかというと、法線は点群処理のほとんどの下流タスクの前提だからだ。陰影をつけてリアルに表示する(ライティング計算は法線がないと始まらない)、点群から面のつながったメッシュを復元する(定番のPoisson再構成は「向きの揃った法線」が入力に必須)、形状の表と裏を判定する——どれも法線が要る。そしてこのとき重要なのは、矢印が表面に垂直なだけでなく、みんなが一貫して同じ側(普通は外側)を向いていることだ。半分が外向き・半分が内向きに混ざった法線は、下流ではほぼ使いものにならない。

というわけで今回の「成功」の定義はシンプルで、隣り合う点の法線が、同じ側を向いて揃っていること。表面がなめらかなら隣接2点の法線はほぼ平行になるはずなので、「隣接ペアのうち、法線同士の内積が負(=90度以上ズレて互いに逆半球を向いている)になっている割合」を反転率と呼び、この数字が低いほど成功、とする。

図の見方も先に補足しておくと、法線は(x, y, z)の3成分を持つ長さ1のベクトルなので、これをそのまま(R, G, B)に割り当てて色として塗れる。この変換では、符号が反転した法線はちょうど補色(RGBをすべて反転した色)になる。つまり冒頭の図で隣同士の点の色が激しく入り乱れて見えるのは、「隣なのに法線が逆向き」という状態をそのまま色にしたものだ。

データの準備: 8171点、そのまま使える大きさだった

点群は raw.githubusercontent.com/mikedh/trimesh が配布しているASCII PLY形式のバニーを使った(1本目の候補URLでそのまま取得できた)。ヘッダをパースすると

  • 頂点数: 8171点
  • x範囲: -0.0946 〜 0.0609
  • y範囲: 0.0334 〜 0.1866(高さ方向)
  • z範囲: -0.0619 〜 0.0587

単位はメートルで、実物のバニー像とだいたい同じ高さ(15cmほど)のスキャンだった。8000点程度ならk近傍探索も現実的な時間で終わるので、間引きせずそのまま全点を使うことにした。

素朴にPCAをやると、37.1%の隣接ペアが裏返る

各点について近傍15点をとり、共分散行列の最小固有値に対応する固有ベクトルを法線とした。

nn = NearestNeighbors(n_neighbors=k + 1).fit(pts)
_, idx = nn.kneighbors(pts)
for i in range(len(pts)):
    nbrs = pts[idx[i]]
    mean = nbrs.mean(axis=0)
    cov = (nbrs - mean).T @ (nbrs - mean) / len(nbrs)
    eigvals, eigvecs = np.linalg.eigh(cov)
    normals[i] = eigvecs[:, 0]  # 最小固有値の固有ベクトル

eigh が返す固有ベクトルの符号は数学的には「+でも−でもよい」——アルゴリズムの内部事情でどちらが出るかが決まる。つまり、隣り合う2点の法線が同じ面の同じ向きを指しているとは限らない。実際に隣接点ペアで法線の内積を測ってみると、

  • 45,529 / 122,565ペア(37.1%)が、負の内積(=互いに逆半球)を向いていた

分母の122,565という半端な数字は、8171点それぞれについて「その点と、その近傍15点」のペアを数えたもの(8171×15)だ。面がなめらかに湾曲しているなら理屈上は0%に近いはずが、4割近くが裏返っている。「符号の不定性」という言葉は知っていても、こうして具体的な割合で見ると想像以上だった。

スタンフォードバニーの側面図。左は素朴なPCA法線をRGB色にした散布図で、隣接点なのに色がまだらに入り乱れている箇所が2つ赤い矢印で示されている。右はMST伝播で向きを揃えた後の散布図で、同じ面がなめらかに同じ色調でグラデーションしており、37.1%から0.6%へ反転率が下がったことが注記されている

直す: 最小全域木で「向きを伝播」させる

法線推定の古典的な対処法(Hoppe et al., 1992)を実装した。考え方はシンプルで、

  1. 全点の法線について、近傍点とのペアごとに「向きの食い違いコスト」1 - |dot(n_i, n_j)| を辺の重みとしたグラフを作る(絶対値を取っているので、符号が違っても法線の”向いている軸”が近ければコストは小さい)
  2. そのグラフの**最小全域木(MST)**を求める(scipy.sparse.csgraph.minimum_spanning_tree)。MSTは「符号を無視した向きの近さ」が高い辺を優先的につなぐので、木をたどれば向きの近い点同士が隣り合う
  3. 適当な1点(今回は一番高い場所=耳の先端)を根にしてBFSで木を辿り、親の法線と内積が負なら子の法線の符号を反転する
graph = csr_matrix((weights, (rows, cols)), shape=(n, n))
graph = graph.minimum(graph.T)
mst = minimum_spanning_tree(graph)
order, predecessors = breadth_first_order(mst + mst.T, i_start=root, directed=False)
for node in order:
    parent = predecessors[node]
    if np.dot(oriented[parent], oriented[node]) < 0:
        oriented[node] *= -1

結果、BFSは根から8171点全点に到達した(=点群のk近傍グラフは1つの連結成分にまとまっていた)。向きを伝播させた後の反転率を同じ指標で測り直すと、

  • 694 / 122,565ペア(0.6%)まで低下

回転させて見比べると一目瞭然だった。最初のインタラクティブ図がまだら模様だったのに対して、こちらは背中も耳もお腹も、なめらかな色のグラデーションになっている。

色だけだと少し抽象的なので、法線を実際の矢印として描いた図も作った。胴体の真ん中あたり(z座標が中央値±2mmの範囲)で点群を厚さ4mmに輪切りにすると、ウサギの横顔の輪郭線が242点分だけ現れる。その輪郭上の点に、推定した法線の矢印を立てたのが下の図だ(矢印は見やすさのため94本に間引いてある)。

スタンフォードバニーをz座標の中央値付近で厚さ4mmに輪切りにした断面(242点)に、推定した法線を矢印として描いた上下2パネルの図。背景には全点群の側面シルエットが薄いグレーで敷かれている。上のパネルは素朴なPCA法線(赤い矢印)で、輪郭線上の矢印が外向きと内向きにでたらめに混ざり、頭のてっぺんや胸、足元では隣り合う点同士なのに矢印が正反対を向いている。下のパネルはMST伝播後(緑の矢印)で、ほぼすべての矢印が輪郭の外側を向いて整列しており、脚の付け根のくぼみなどにだけ向きの乱れがわずかに残っている

上の素朴なPCA版では、すぐ隣の点なのに矢印が正反対を向いている場所がそこら中にある。MST伝播後は、ほぼすべての矢印がきれいに外側を向いて揃った。「37.1%が裏返る」「0.6%まで直る」という数字の中身は、要するにこの矢印の反転のことだ。

正直な結果: 0.6%は「ゼロ」ではなかった

ここで終わらせず、残った0.6%がどこにあるかを確認した。反転ペアに関わっていた点は268点(全体の3.3%)あり、その高さ(y座標)を見ると 0.034〜0.187とほぼ全域に散らばっていた。1箇所の”継ぎ目”に集中しているわけではなく、耳の付け根や脚の裏など、局所的に法線がねじれやすい(曲率が急に変わる)場所にぽつぽつと小さな矛盾が残っていた形だ。

さらに、重心から各点への方向を「外向きの目安」として比べると、オリエンテーション後の法線のうち91.0%が大まかに外向きだった。バニーは凸形状ではない(耳の内側や股の間は凹んでいる)ので100%にはならないのが当然だが、逆に言うと残り9%は「MSTの伝播としては局所的に一貫しているが、直感的な外向きとはズレている」領域がそれなりにあるということでもある。MST伝播は”隣同士のつじつま”は合わせてくれるが、“全体としてどちらが外か”までは保証しない——アルゴリズムの性質として腑に落ちる結果だった。

kを振ってみる: 小さすぎても大きすぎてもダメ

近傍点数kは、法線推定の唯一のハイパーパラメータだ。k=3から130まで11段階で振って、法線場がどう変わるか見てみた。

スタンフォードバニーの側面図を3パネル並べたもの。k=5ではノイズで色が細かくざらつき耳の輪郭も乱れている。k=20では耳の稜線までなめらかに保たれている。k=130では耳が本体の色に溶け込むように塗りつぶされ輪郭が消えている、という3つの結論が矢印付きで示されている

kを3から130まで変化させながら法線カラーの点群が変わっていくアニメーション。小さいkでは色がざらつき、大きいkでは耳の色が本体に溶け込んでいく

定量的にも確認した。2つの指標を使う。

  • ガタつき指標: 各点の法線と、常に固定した近傍6点の法線との内積のズレ(kとは独立な評価用グラフで測定)。値が小さいほど法線場がなめらか
  • 耳の曲率推定値: 耳の先端付近(上位5%の高さの点)での「最小固有値/固有値の和」の平均。大きいほど「そこに本当のエッジ(耳のヘリ)がある」と正しく検出できていることを意味する

横軸をkの対数スケールにした折れ線グラフ。青線(法線のガタつき)はkが大きくなるほど単調に減少する。赤線(耳の稜線の曲率推定値)はk=3から増加してk=60付近でピークを迎えたあと、k=90・130で減少に転じる。k=60あたりが山であることが赤字で注記されている

2本の線は違う形をしていた。ガタつきはkを大きくするほど単調に下がり続ける(近傍を増やせば増やすほど、平均化されて滑らかになるのは当然)。一方で耳の曲率推定値はk=60付近で山を作り、そこから先は下がっていく。k=3の時点(0.0059)から徐々に上がってk=45〜60で0.14付近に達したあと、k=90で0.1382、k=130で0.1295まで下がった。近傍が大きくなりすぎると耳という薄い構造の”裏側”の点まで巻き込んでしまい、本当は鋭いはずのヘリが平らに均されてしまうということだ。

「ノイズを消したいなら大きく、細部を残したいなら小さく」という単純なトレードオフだと思っていたが、実際にはノイズ低減は単調に効き続ける一方で、細部の検出には最適な範囲があってそこを超えると悪化するという、非対称な形の関係だった。片方だけ見ていたら「大きいほど良い」と誤解していたと思う。

手を動かして意外だったこと

一番の収穫は、教科書の1行で済まされがちな「符号は不定」という注意書きの重みを、37.1%という具体的な数字で体感できたことだった。アルゴリズムの実装としてはeighを呼ぶだけの1行なのに、その裏で起きている符号の揺れが、可視化するとここまで露骨に破綻して見えるとは思っていなかった。

もう1つは、MST伝播が「万能の修理」ではなかったこと。0.6%まで下げられたのは十分実用的な数字だが、ゼロにはならなかったし、その残りは1箇所の欠陥ではなく全体にまばらに散った小さなねじれだった。木構造で伝播させる以上、根本的に「局所的には正しくても大域的には矛盾する」ケースを完全には防げない——これも実際に手を動かして初めて実感できた限界だった。うまくいった数字だけでなく、うまくいかなかった0.6%まで含めて報告する方が、この実験の実態に近いと思う。

まとめ

  • スタンフォードバニーの点群(8171点、ASCII PLY)にk近傍PCAで法線を推定した
  • 素朴に固有ベクトルの符号をそのまま使うと、隣接点ペアの37.1%が逆半球を向く(符号の不定性が実際にどれくらい起きるかを定量化できた)
  • 近傍グラフのMSTを使った向き伝播(Hoppe法)で0.6%まで低減。ただし残った矛盾は1箇所に集中せず全体にまばらに散っており、伝播だけでは大域的な一貫性まで保証されないことが分かった
  • 重心基準の「外向きらしさ」は91.0%——凹凸のある形状ではMSTの局所整合と直感的な外向きが完全には一致しない
  • 近傍点数kを3〜130で振ると、ガタつきは単調に減少する一方、耳の輪郭を捉える曲率推定値はk=60付近をピークに減少に転じた。ノイズ低減と細部保持は同じ方向のトレードオフではなく、非対称な形をしていた