プログラミング
Poisson Disk Sampling
Poisson Disk Sampling (stripeacross.com)
要約
この記事では、コンピュータグラフィックスやシミュレーションで一般的に発生する「ランダムだが近すぎないように配置する」問題に対する、Robert Bridsonによる効率的なアルゴリズム「Poisson Disk Sampling」を解説しています。アルゴリズムの基本原理、グリッド分割、アクティブリストを用いたサンプリング手法に加え、親子の関係を利用した2次元での改善策や、距離分布を操作する高次元でも適用可能な改善策についても詳細に説明しています。
全文翻訳
Poisson Disk Sampling
2024年、9人の数学者チームが、幾何学的ラングランズ予想に関する1000ページ近い壮大な証明を発表しました。これは純粋数学における頂点的な成果であり、私は彼らが証明した命題のステートメントさえ理解できないだろうと諦めています。それとは全く対照的に、2007年にRobert Bridsonは1ページの論文を発表し、それは1000件近い引用があり、完全に理解するのに10分もかかりません。それは、コンピュータグラフィックスやシミュレーションで一般的に発生する問題、つまりランダムに配置するが、近すぎないように配置する、という問題に対する簡単な解決策を提示しています。
例えば、手続き的に森を生成しようとして、木の配置方法が必要だとします。単純なランダムサンプリングの問題点は明白です。木が重なり合ってしまうでしょう!私たちが求めているのは、任意の2本の木の間に最小距離を設定できる能力です。このルールに従う木の分布をポアソンディスク分布と呼びます。ランダムなダーツを投げ、他の点から最小距離内に落ちた点をすべて却下する、という単純な却下サンプリングアプローチを試すことができますが、より巧妙なデータ構造がない場合、各サンプルに対して衝突をチェックするのに線形時間かかり、却下率は急速に1に近づきます。Bridsonのアルゴリズムは、これを効率的に行う方法を提供します。
Bridsonのアルゴリズム
目的とする点間の最小距離をrとし、d次元空間で作業するとします。Bridsonのアルゴリズムは以下の通りです。
空間を、辺の長さが r/√d のグリッドに分割します。これにより、各グリッドセルには最大1つの点しか含まれないことが保証されます。
アクティブリストに、空間から一様に選ばれたランダムな点を1つ追加して初期化します。
アクティブリストが空でない間:
アクティブリストからランダムな要素pを均一に選択します。
pを中心とする半径rから2rの間のドーナツ(アニュラス)から、最大k回、均一にサンプリングします。有効なポアソンディスクサンプルが見つかった場合、グリッドを使用して効率的に衝突検出を行い、アクティブリストに追加し、新しいpを選択します。k回の試行で有効な点が見つからない場合は、pをアクティブリストから削除します。Bridsonはk=30を推奨しています。
ドーナツを均一にサンプリングする最も簡単な方法は、ランダムな単位ベクトル v∈ℝd と、区間 [1/2^d, 1) から一様に選ばれた数xを生成し、最終的なサンプルを 2rx^(1/d)v とすることです。数年前に、これがなぜ機能するのかを説明するビデオを作成しました。2次元では、単位ベクトルを選択することは、角度 θ∈[0, 2π) を選択することと同等です。高次元では、各成分が正規分布からサンプリングされたベクトルを正規化することができます。
改善策
Bridsonのアルゴリズムには、同じ数の点を生成するために必要な反復回数を劇的に減らすことができる、2つの簡単な改善策があります。1つ目は2次元で機能しますが、2つ目は高次元でも機能します。
まず2次元の改善策から始めましょう。アルゴリズムが点pを配置し、次に新しい点qを得るためにそのドーナツをサンプリングする場合を考えます。pをqの親と呼びます。これらの点間の関係には貴重な情報が格納されています。必然的にqを中心とするドーナツをサンプリングする際に、pに近すぎる点を除外する必要がある角度の範囲全体が存在します。この範囲は、以下の図の破線で表されます。
pq |p - q| = 1.50r
視覚的な直感は把握しやすいですが、それを数式に落とし込むのは面倒な三角法の演習です。詳細については省略し、証明なしで主張しますが、破線で形成される円錐は角度αを中心に持ち、その幅は2βです。
α = atan2(py - qy, px - qx)
β = min(arccos(|p - q|^2 + 3r^2) / (4r * |p - q|), arccos(|p - q|) / (2r))
この数式で興味深いのは、βの計算式に現れる最小値です。これは、ドーナツの内円または外円のいずれかが円錐を境界付ける可能性があるという事実を考慮しており、点間の距離によって異なります。円錐とドーナツの交差が完全に円内にあることを保証するために、最小値を選択する必要があります。点間の距離が√3 * r を超えると、円錐の境界点が外円から内円にジャンプすることに注意してください。
これを実装するには、各点の親を保存するだけで十分です。その後、円錐の角度を計算し、新しいサンプルの角度θを円錐の外側の範囲で生成できます。以下のグラフは、親の最適化ありとなしでBridsonのアルゴリズムによって生成された点の数を示しています。各データポイントは、辺の長さℓ=100、r=1のグリッドで100回の試行の平均です。
この改善策は高次元にも一般化できる可能性がありますが、ドーナツの接触ベクトルの保存にさらに多くのスペースが必要になり、高次元ではドーナツと球の交差の体積が比例して重要でなくなるため、効果は低下する可能性があります。同様に、親に加えて点の子供を保存して、ドーナツのさらに多くのセクションを排除することを妨げるものはありませんが、これもより多くのストレージが必要になり、θの選択がはるかに遅くなります。
次のサンプルの角度の選択方法を変更するのではなく、2番目の改善策は次のサンプルの距離の選択方法を変更します。ドーナツ内の各点から中心までの距離の分布を考えます。その累積分布関数(CDF)は、区間 [r, 2r] で x^d に比例します。この説明については、私のビデオを参照してください。しかし、指数をd以外の定数cに変更した場合はどうなるでしょうか?そうすれば、点を中心に近づけたり遠ざけたりすることができます。c≠0の場合、正確なCDFは以下のようになります。
Fc(x) = {
0 if x ≤ r,
(x^c - r^c) / ((2^c - 1)r^c) if r < x ≤ 2r,
1 if x > 2r.
}
以下の図は、cが変化したときのCDFとドーナツ内の500個のランダムサンプルを示しています。c=2は一様分布を与えることを覚えておいてください。
0
r
2r
0.0
0.5
1.0
distance
CDFc = 2.00
スライダーでc=0を設定できることに気づいたかもしれませんが、F0(x)はゼロ除算のため未定義です。これを修正するために、F0(x) = lim(c→0) Fc(x) と定義します。
F0(x) = {
0 if x ≤ r,
log2(x) - log2(r) if r < x ≤ 2r,
1 if x > 2r.
}
任意のcの値を使用して半径をサンプリングするには、逆変換サンプリングを適用できます。c=0の場合、半径は r * 2^x となるべきです。ここでxは区間 [0, 1) からの一様乱数です。それ以外の場合は、yが区間 [1, 1/2^c) からの一様乱数である場合、r * 2 * y^(1/c) を使用します。区間の境界は、cが正か負かによって入れ替わります。
以下のヒートマップは、cがBridsonのアルゴリズムによって生成される点の数に与える影響を示しています。辺の長さℓ=100、r=1のグリッドで100回の試行という同じ実験設定を使用しています。
これは、cを非常に負の値に設定するか、cが負の無限大に近づく極限を取り、各点が親から正確にrの距離を持つように強制すべきであることを示唆しているようです。しかし、それは...