37 点群基礎
パートVIIIでは現実世界を三次元データに変換する方法を扱った：レーザー三角測量（チャプター 33）、構造化照明、位相シフト、フォーカススタック——それぞれ独自の物理原理に基づくが、最終的には\((x,y,z)\)座標を持つ三次元点の集合という同じ出力を生み出す。チャプター 30 ではこの一連の撮像リンク全体を概観した。本章から始まるパートIXでは、取得済みの三次元データを処理する方法、具体的にはノイズ除去、位置合わせ（レジストレーション）、計測、フィッティングを主題とする。そしてこれらの処理対象は、圧倒的多数の場合において点群（point cloud）である。
点群は三次元世界における「画素」に相当するが、二次元画像よりもはるかに扱いが難しい。二次元画像は規則的な格子であり、各画素の近傍がどれであるか、いくつあるかは行と列のインデックスから直接決まる。一方、点群は順序付けられていない（unordered）点の集合であり、どの点同士が隣接しているかはデータ構造に記載されておらず、自力で計算する必要がある。点群は疎（sparse）であり、表面以外の空間には何の記録も存在しない。またノイズを含み（noisy）、各点は計測方向に沿って変動する。さらに外れ値（outlier）を含み、アーク、反射、マルチパスによって空中に浮かぶ「フライポイント（flyer）」が架空に生成されることがある。本章では合成シーンを用いてこれら4つの性質を示し（図 37.1、図 37.2）、これを基に最も基礎的な3つの事柄——三次元データの表現形式の種類、順序付けられていない点集合上で最近傍を高速に探索する方法、フライポイントを除去する方法——を解説する。
このシーン全体には20200個の点が含まれる：地面11000点、立方体上面2600点、4つの垂直側面3600点、球冠2800点、加えて200個の外れ値フライポイントである。すべての表面点には、各自の法線方向に沿って\(\sigma=0.10\) mmのガウス計測ノイズが重畳されており、高さ範囲は\(z\in[-0.43,\ 43.92]\) mmである（表面自体は最大で18 mmであり、それより高い部分はすべてフライポイントである）。このデータセットはcode/point_cloud_basics/のプロジェクトによって固定乱数シードを用いて決定的に生成されており、本章以降のすべての数値はこの実際の実行結果に基づいている。
37.1 三次元データの表現形式
同一の三次元物体でも、本質的に異なる3種類のデータ構造で記録することができ、工学的にはそれらを明確に区別する必要がある。
1つ目は点群であり、前述の順序付けられていない点集合\(\{\mathbf p_i=(x_i,y_i,z_i)\}\)である。各点には任意で色や法線を付加することができる。点群は最も汎用的な表現形式であり、任意のトポロジーの表面を表現でき、多視点の貼り合わせや自由曲面にも対応可能である。代償として、順序がなく不規則であるため、「近傍を探索する」操作には常に追加の空間インデックスが必要となる（次節の主題である）。
2つ目は距離画像（range image）であり、深度マップまたは2.5D画像とも呼ばれる。これは規則的な二次元格子であり、各格子セル\((u,v)\)に1つの高さ値\(z\)が格納される——本質的には「\(z\)を輝度値として扱う画像」である。最大の利点は規則性にある：近傍関係が格子インデックスから直接得られ、二次元画像のすべてのフィルタリング・モルフォロジー演算（チャプター 6）をそのまま利用でき、記憶領域も効率的である。
3つ目はメッシュ（mesh）であり、頂点と三角形パッチによって表面の接続関係を明示的に記録する。主にレンダリング、有限要素法、衝突検出に用いられ、産業計測では使用頻度が比較的低いため、本章では詳述しない。
「2.5D」という呼称は、距離画像の曖昧な位置を的確に表している。\(x,y\)平面と1つの\(z\)を持つため一見三次元的だが、各\((x,y)\)に対して1つの\(z\)しか持てないという制約がある。真の三次元表面では、垂直壁、オーバーハング面、空洞部では同一視線が複数の\(z\)を通過するが、これは距離画像では表現できない——距離画像は「ある固定視点から見て最も近い（または最も高い）層」しか記録できないのである。
点群と距離画像は相互に変換可能であるが、この変換は非可逆である。距離画像から点群への変換は単純であり、有効な各格子セルをそのインデックスと高さから三次元点に復元すればよい。逆に、点群から距離画像を生成するには、まず投影方向を選択し（ここでは上面視、つまり\(-z\)方向を見る）、\(XY\)平面を規則的な格子に分割してビン化し、各ビン内で代表的な1つの\(z\)だけを保持する。問題はまさにこの「1つだけ保持する」点にある。
これら20200個の点を\(200\times200\)の上面視格子にビン化し、各ビン内で最も高い\(z\)を採用した結果を@fig-pb-range に示す。格子化によるアーティファクトが明確に確認でき、地面はモザイク状に量子化され、球冠と立方体上面は残っているが、立方体の垂直側壁は完全に消滅しており、立方体の真下と球冠の真下にある遮蔽された地面部分も失われている。統計的には、元の20200点のうち、再投影後に生き残ったビンは12585個に過ぎず、7615点（37.7%）が失われた。これら失われた点はノイズではなく、実在する表面の実在するサンプル点である。それらが消失したのは、他の点と同一の\((x,y)\)を共有しつつ\(z\)が最も高くなかったために過ぎない——側壁点は上面点によって覆われ、遮蔽された地面点は上方の物体によって覆われたのである。
この37.7%という数値は無作為に選ばれたものではなく、「単一視点の本質」に関する教訓である。チャプター 33 でレーザー三角測量が直接出力するのは距離画像であり、その時点でレーザーの視点から照明可能な表面しか自然には見えていない——背面、側壁、空洞は最初から影の中にあるのである。距離画像の2.5D的な制限と、撮像方式の単一視点的な制限は、同一の事象に対する2つの表現に過ぎない。完全な三次元データを得るには、多視点でスキャンした後に位置合わせを行う（チャプター 39）か、最初から最後まで点群として扱い距離画像に縮退させないかのいずれかである。したがって実用的な工学的指針としては、点群を利用できる場合は早期に距離画像に変換しないことである。単一視点で十分であり、かつ二次元画像の成熟した演算子を利用したい場合に限り、その37.7%の情報と引き換えに距離画像の規則性を活用する価値が生まれる。
SciVisionでは、これら2つの形式はそれぞれSciPointCloud（点 + 任意の色 + 任意の法線、PLY/PCD/OBJの読み書きに対応）とSciRangeImage（ushort型の高さデータにresolutionX/Y/ZとoffsetZを付加し、16ビットPNGとして保存）に対応する。本章のシーンをSciPointCloudに読み込むとLength=20200となり、保存されたcloud_raw.plyには点ごとの高さに応じた色が付加されており、任意の点群ビューワーで直接開くことができる。
37.2 空間インデックス：KdTree
点群が順序付けられていないことによって生じる最初の、そして最も根本的な工学的問題は、「与えられた照会点に対し、その最近傍の点をいくつかどのようにして見つけるか」である。この最近傍（nearest neighbor）探索は、ほぼすべての三次元アルゴリズムの基盤となる——位置合わせでは点対を構築するために（チャプター 39）、フィルタリングでは近傍を定義するためにこれを利用し、法線推定、特徴記述、表面再構成も例外ではない。最も単純な手法は全探索である：各照会点に対し、\(N\)個すべての点との距離を計算してソートする。これは1回の探索で\(O(N)\)、\(M\)回行うと\(O(MN)\)の計算量となる。\(N\)と\(M\)がそれぞれ数万に達すると、数百億回の距離計算が必要となり、許容できないほど遅くなる。
kd木（kd-tree）はこの計算量を対数オーダーまで削減する。その思想は、座標軸に垂直な超平面で空間を再帰的に2分割することである：根ノードでは\(x\)座標の中央値を求め、値の小さいものを左部分木、大きいものを右部分木に振り分ける。次の層では\(y\)軸、その次の層では\(z\)軸に切り替え、このように軸をローテーションさせながら分割を続け、各葉ノードに含まれる点が少数になるまで続ける。構築後、1回の最近傍探索では、まず照会点が含まれる葉ノードまで下降し（この経路には\(O(\log N)\)回の比較しか必要ない）、候補となる近傍点の集合を得る。次に上向きにバックトラックし、各分割ノードで「照会点から分割超平面までの距離」が現在判明している第\(k\)近傍点までの距離よりも小さいかどうかを確認する——小さくない場合、超平面の反対側の部分木全体にはより近い点が存在し得ないため、一括して枝刈りして探索する必要がなくなる。この枝刈りこそが、平均探索計算量を\(O(\log N)\)まで削減する要因である。\(k\)近傍探索ではサイズ\(k\)の最大ヒープを用いて現在の最近傍\(k\)個を維持し、半径探索（radius search）では与えられた半径内にあるすべての点を収集し、バックトラック時の枝刈り判定基準を「第\(k\)近傍までの距離」から「半径」に置き換える。
kd木の対数計算量には、次元が高すぎないという前提条件がある。次元が数十から数百に上昇すると、「超平面までの距離が最近傍点までの距離よりもほぼ常に小さい」状態となり、枝刈りがほぼすべて無効になってkd木はほぼ全探索に縮退する——これが次元の呪い（curse of dimensionality）である。幸いにも三次元点群の次元は3しかなく、kd木はここでは非常に良く機能する。高次元特徴のマッチングでは、近似最近傍や局所鋭敏ハッシュを代わりに用いる。
本章のプロジェクトでは、暗黙的な配列型のkd木を実装し（std::nth_elementを用いて中央値で分割し、\(x/y/z\)軸をローテーションする）、全探索法との処理時間を比較した。20200点に対して8近傍探索を10000回実行したところ、kd木の処理時間は約60 ms、全探索法は約870 msであり、約14.5倍の高速化が達成された。さらに両者の結果はビット単位で一致している——kd木は近似ではなく精度を損なっておらず、不可能な領域を賢く回避しているだけで、正確な最近傍点を返している。この14.5倍の高速化は1回の探索では大きな差に見えないが、最近傍探索を繰り返し行う反復アルゴリズム（例えばICPでは毎回数万点の対応点を探索する）の内部では、数分と数十ミリ秒の差となって現れる——\(\log N\)が\(N\)に対して持つ工学的価値は、まさにこのような「数千万回呼び出される」内部ループにおいて実感されるのである。
37.3 ノイズと外れ値
三次元データに含まれる「汚れ」は、本質的に異なる2種類に分けられ、これらを混同すると誤ったツールを使用することになる。
1つ目は計測ノイズである。実在する各表面点には、計測方向（通常は表面法線）に沿って微小なランダムな変動が生じる。本章の合成データでは\(\sigma=0.10\) mmのガウスノイズとしている。このノイズの特徴は振幅が小さく、ゼロ平均で表面に密着していることであり、近傍内の点は互いに密集し密度は均一である。この種のノイズは平滑化フィルタ（パートIXの後続の章）で処理すべきであり、「点を削除する」方法で対処すべきではない。
2つ目は外れ値、すなわちフライポイントである。アーク、鏡面反射、マルチパス、センサーのアーティファクトによって、実在するいかなる表面上にも存在しない点が架空に生成され、空中に浮かんですべての正常な構造から離れた位置に存在する。本章の200個のフライポイントはまさにこのようにして、\(z\in[6,44]\) mmの空間領域に散布されている。その特徴は計測ノイズとは正反対で、表面から遠く、近傍が疎であることである——フライポイントの周囲にはほとんど他の点が存在せず、最近傍点との間に大きな隙間がある。
この「近傍が疎である」という特徴は、フライポイントを検出するためにそのまま利用できる。これが統計的外れ値除去（statistical outlier removal, SOR）である。アルゴリズムは非常に単純明快である：各点に対し、kd木を用いて\(k\)個の最近傍点を求め、それら\(k\)個の近傍点までの平均距離\(d_i\)を計算する。正常な表面点は密集しているため\(d_i\)は小さく、フライポイントは空中に孤立しているため\(d_i\)は大きくなる。全点の\(\{d_i\}\)を1つの分布と見なし、その平均\(\mu\)と標準偏差\(\sigma\)を求め、\(d_i>\mu+t\sigma\)となる点を外れ値と判定して除去する。これは実質的に、チャプター 2 で解説した\(3\sigma\)基準を「近傍平均距離」という量に適用したものであり、チャプター 26 の統計的異常判定と起源を同じくする——固定の距離閾値を事前に設定するのではなく、データ自身の分布によって閾値を定めるのである。
閾値を分布に基づいて定めるのは、なぜ固定のミリメートル値を設定するのではないのだろうか？「どれだけ離れていれば遠いとみなすか」は、点群の密度に完全に依存するからである。同じ0.5 mmの近傍平均距離でも、稠密なスキャンでは外れ値であるが、疎なスキャンでは正常である。\(\mu+t\sigma\)はスケールをデータ自身に委ねるため、別の点群に切り替えてもパラメータを再調整する必要がない——これがまさに、統計的閾値が硬い閾値に対して持つロバスト性である。
本章では\(k=8\)、\(t=2\)を用いる。このデータにおいて、近傍平均距離の\(\mu+2\sigma=3.23\) mmであり、これに基づいて190個の点が除去された。その190個はすべて真のフライポイントに一致し（200個のフライポイントに対する再現率は95%）、かつ誤って削除された点は0個であった——実在する表面点が誤って除去されることはなかった。図 37.4 は除去された190個の点を赤色の菱形でマークしており、それらは空間領域に整然と位置し、まさに浮遊するフライポイントであることがわかる。
では見逃された10個のフライポイントはどこへ行ったのだろうか？それらが除去されなかったのは、不運にも（あるいは幸運にも）実在する表面の近傍にランダムに落下したためである——フライポイントが地面や立方体のすぐ上に偶然漂着した場合、その8個の最近傍点に実在する表面点が混ざり込み、近傍平均距離が閾値以下に引き下げられるため、統計的な検出を免れるのである。これはアルゴリズムのバグではなく、統計的手法の固有の限界である：SORは「近傍密度の異常」を判定しており、正常な密度に偽装された外れ値を判別することはできない。この10個の見逃し点を正直に明記することは、100%の再現率を報告するよりも信頼できる——下流のアルゴリズム（ロバストフィッティングやICPなど）は、そもそも前処理ですべてを1回でクリーニングできると期待するのではなく、残留する少数の外れ値に対して耐性を持つべきである。
37.4 SciVision 実装
各データ形式の I/O は SciPointCloud で行うが、ここには正確に記載しなければならない落とし穴がある。SciPointCloud は点セットとして SciVector3dArray を受け入れるが、SciVector3dArray(float*, size) という「連続配列から一度に構築できる」という一見便利なコンストラクタは、本実行環境では動作しない。構築後に Length() を確認すると依然として 0 であり、点がまったく格納されていない。確実な方法は点ごとに Append することである：
SciVector3dArray pts;
for (int i = 0; i < n; ++i) {
SciVector3d v(xyz[3*i], xyz[3*i+1], xyz[3*i+2]);
pts.Append(v); // 点ごとに追加する；float* コンストラクタは本環境ではデータを格納しない
}
pc.SetPoints(pts); // この後 pc.Length() = n となる最近傍探索については、SDK の Sci3DKdTree は使用可能であり、その結果は自作の kd 木と一致する：
SCIMV::Sci3DKdTree kt;
long rc = kt.CreateKdTree(pc);
SciVector3d q(xyz[0], xyz[1], xyz[2]);
int num = 0; SciVector3dArray pos; SciIntArray ind; SciFloatArray dist;
kt.FindKNearestNeighbors(q, 8, &num, &pos, &ind, &dist);
// dist の最初の 4 要素 = 0.000 / 0.516 / 0.802 / 0.999 であり、自作 kd 木とビット単位で一致する最初の近傍の距離は 0.000 であり（照会点自身が見つかっている）、残りも自作実装と完全に一致する。これは本質的な検証材料となる：本章の計時と統計的外れ値除去（SOR）はいずれも自作の kd 木に基づいている（決定性と計時可能性のため）が、SDK の Sci3DKdTree はその正しさを裏付ける傍証として機能する。
データ変換を行う SciSv3DDataConvert について：順方向の ConvertPointCloudToRangeImage（点群から距離画像への変換）は使用可能であるが、逆方向の ConvertRangeImageToPointCloud は本環境では動作しない（生存点数が 0 になる）。そのため セクション 37.1 の 2.5D 損失実験は自作のビン分割で行い、決定性を制御可能としている。
SOR については、SDK の SciSv3DClean の SORFILTER 機能は本環境ではクラッシュする（実測によると終了コードはヒープレイアウトに応じて変動し、0xC0000005 のアクセス違反と 0xC0000409 のスタック保護エラーの両方が確認されている）。そのため本章では独立したサブプロセスのプローブでこの機能を検証し（クラッシュしてもメインの処理が巻き添えにならない）、実際の SOR 処理は自作実装を用いている：
// 自作 SOR：各点の k 近傍に対する平均距離を計算し、μ+2σ の閾値で削除
for (int i = 0; i < N; ++i) {
kd.knn(&xyz[3*i], K+1, &nb); // K+1：自身を含む
double sd = 0; int c = 0;
for (auto& p : nb) { // 自身をスキップし、K 個の近傍の距離を累積
if (p.second == i) continue;
sd += std::sqrt(p.first);
if (++c == K) break;
}
meanD[i] = sd / c; // 近傍平均距離 d_i
}
double mu = /* mean(meanD) */, sigma = /* std(meanD) */;
double thr = mu + 2.0 * sigma; // = 3.23 mm
for (int i = 0; i < N; ++i)
if (meanD[i] > thr) removed.push_back(i); // 190 点を削除この「SDK は I/O と交差検証に使用し、コアアルゴリズムは自作する」という役割分担は、パート IX の 3D モジュールの常識である。点群系のアルゴリズム（kd 木、SOR、ICP、PCA）はほとんどが自己完結的で自作が容易であり、自作版は決定性と計時可能性も保証できる。一方で本環境における SDK の 3D モジュールの成熟度は一様ではないため、動作する部分は利用し、クラッシュする部分は隔離する方針を採用している。完全な実行可能プロジェクトは code/point_cloud_basics/ に収録されている。
産業事例：レーザスキャン点群におけるフライポイントの嵐
溶接線のレーザスキャンはフライポイントの多発領域である。アーク光、スパッタ、母材の鏡面反射により、センサには虚偽のエコーが充満し、フライポイントが密に分布する。このような点群を用いて直接溶接余盛高さやアンダーカット深さを測定すると、フライポイントによって測定基準全体が歪められてしまう。正しい手順はまず SOR によるクリーニングを行い、その後に測定することであり、クリーニング後には溶接輪郭が直ちに安定する。
ただし SOR は諸刃の剣でもある。閾値を厳しく設定しすぎると、真の薄いエッジや鋭い特徴までフライポイントとともに誤削除してしまう。溶接線のバリの先端やアンダーカットの最深点は、それ自体が近傍の疎な「孤点」であり、統計的にはフライポイントと酷似している。過剰なクリーニングは測定対象である欠陥信号そのものを消去してしまうことになる。閾値は点群自身のノイズ分布に基づいて設定しなければならない（固定のミリメートル数ではなく \(\mu+t\sigma\) を用いる）。原則として：過剰なクリーニングで真の信号をノイズとして捨てるよりは、クリーニングを控えめにして、外れ値への耐性を下流の頑健なアルゴリズムに委ねるほうがよい。
37.5 まとめ
- 3 次元データには 3 つの形態があり、点群は最も汎用的、距離画像は最も規則的、メッシュは接続関係を管理する：点群は無順序点集合であり、任意のトポロジを表現できるが追加のインデックスが必要となる。距離画像は規則的な 2.5D グリッドであり、2 次元の演算子を直接再利用できるが、各 \((x,y)\) には 1 つの \(z\) しか格納されない。
- 点群から距離画像への変換は損失を伴い、その損失はシーンの「鉛直性」に比例する：本章の 20200 点を \(200\times200\) のグリッドに再投影したところ、残ったのはわずか 12585 点（37.7% の損失）であり、鉛直な側壁と遮蔽された領域がまとまって消失した。これは単一視点 2.5D の本質的な限界であり、レーザ三角測量の単一視点による取得方式と対応している。点群を使用できる場合は、早まって距離画像に縮退させてはならない。
- kd 木は無順序点群に対する最近傍探索の計算量を \(O(N)\) から \(O(\log N)\) に削減する：軸をローテーションしながら中央値で分割して木を構築し、バックトラッキング時に「超平面までの距離」で枝刈りを行う。本章の 10000 回の 8-NN 探索では、kd 木は約 60 ms、全探索は約 870 ms であり、14.5 倍の高速化とビット単位で正確な結果が得られた。これはすべての反復型 3D アルゴリズム（位置合わせ（レジストレーション）、フィルタリング、特徴抽出）の基盤である。
- 3 次元の「汚れ」には 2 種類があり、対策も異なる：測定ノイズは振幅が小さくゼロ平均で表面に近接しており、平滑化により処理できる。外れ値であるフライポイントは表面から離れており近傍が疎であるため、統計的外れ値除去（SOR）による点の削除で処理できる。両者を混同してはならない。
- SOR は「近傍平均距離」の \(\mu+t\sigma\) 閾値でフライポイントを検出し、閾値は分布に追従する：本章では \(k=8\)、\(\mu+2\sigma=3.23\) mm とすることで 190 点を削除し、真のフライポイントに対する適中率は 95%、誤削除は 0 であった。見逃された 10 点は表面近傍に位置し、正常な密度に擬装していた。これは統計的手法の本質的な限界であるため、下流の処理には残留外れ値に対する頑健性を持たせるべきである。
点群のデータ構造、空間インデックス、3 次元フィルタリングについてより体系的な解説は、Steger らの著作 (Steger, Ulrich, と Wiedemann 2018) を参照されたい。また、マシンビジョンの文脈における 3 次元幾何と距離画像処理の速習用参照としては (Szeliski 2022) がある。kd 木という空間インデックスの原初的な出典は、Bentley が多次元二分探索木を提案した古典的論文 (Bentley 1975) であり、これによって本章の最近傍探索の計算量が \(O(N)\) から対数オーダーに削減された。本章および以降の章で使用されるボクセルダウンサンプリング、統計的外れ値除去、法線推定などの点群演算子は、実用上は点群ライブラリ（PCL）の実装を参照することが多く、その全体像は Rusu と Cousins の論文 (Rusu と Cousins 2011) で解説されている。以降の チャプター 38 と チャプター 39 では、それぞれ 3 次元フィルタリング・ダウンサンプリングと点群の位置合わせ（レジストレーション）を詳しく展開する。