快速泊松碟采樣算法實作(Fast Poisson Disc Sampling)
前(fei)言(hua)
最近在看一些随機地圖生成算法,涉及到生成Voronoi圖,這需要提前在一個平面内随機生成一堆的點,這些點還要滿足随機而且盡量均勻分布在平面上。一般文章都提到采用Lloyd Relaxation算法,不過這個算法比較複雜,消耗也比較大,後來看到這個快速泊松碟采樣算法,也是用于生成一堆均勻分布的點的,而且算法複雜度在 O ( N ) O(N) O(N)。
算法的說明論文可以參考這裡。個人覺得描述得很不(fan)清(ren)晰(lei),後來我看了别人的視訊說明才真正了解這個算法,其實算法本身很簡單,也很直覺,下面會講解下。完整代碼我放在個人Github上。
算法說明
下面以2維平面為例:
假設我們需要在一個寬高為 ( w i d t h , h e i g h t ) (width,height) (width,height)的平面内平均生成一堆的點,且這些點之間的距離不能小于 r r r。
我們可以先從平面内随機選一個點,然後在這個點附近随機找一些點,并判斷這些點是否合法,合法的話則在這些點附近繼續随機尋找,直到找不到合法點為止。
算法簡單說起來就是這樣,下面詳細說下細節。
- 為了保證能盡量填滿整個平面, 随機找點時,采用與中心點距離為 [ r , 2 r ) [r,2r) [r,2r)的圓環内找,這個距離能保證找的點距離自身大于 r r r,且不會離得太遠,能填滿整個平面。
- 怎麼判斷一個點的附近已經找不到合法點了?算法定義了一個常量 k k k,對于每個點,我們嘗試在它附近随機找 k k k次,如果都找不到,那麼就認為這個點附近已經沒有合法點。
- 怎麼快速判定随機找出的點是否合法?這個是算法的關鍵,可以采用一些空間劃分方法來做(遊戲場景也會經常用到),首先将平面劃分成 m m m行 n n n列的格子,每個格子都儲存了格子内部的點。這樣當我需要判斷一個點是否合法時,我隻要和附近的格子内的點做判斷即可。
-
那怎麼确定每個格子的大小?我們盡量讓每個格子内部最多隻能有1個點,這樣資料結構就會簡單很多。怎麼做到呢?我們假設每個格子都是正方形,那正方形内部距離最遠的點就是對角線的2個點,是以我們隻要保證正方形的對角線長度大于等于 r r r,則正方形内部任意2個點之間的距離肯定小于 r r r,進而保證每個正方形内部肯定最多隻能有1個點。假設正方形邊長為 a a a,對角線長度為 r r r,那麼有: a 2 + a 2 = r 2 a^2+a^2=r^2 a2+a2=r2,那麼 a = r 2 a=\frac{r}{\sqrt{2}} a=2
r。
代碼
public static class Algorithm
{
public static List<Vector2> Sample2D(float width, float height, float r, int k = 30)
{
return Sample2D((int)DateTime.Now.Ticks, width, height, r, k);
}
public static List<Vector2> Sample2D(int seed, float width, float height, float r, int k = 30)
{
// STEP 0
// 次元,平面就是2維
var n = 2;
// 計算出合理的cell大小
// cell是一個正方形,為了保證每個cell内部不可能出現多個點,那麼cell内的任意點最遠距離不能大于r
// 因為cell内最長的距離是對角線,假設對角線長度是r,那邊長就是下面的cell_size
var cell_size = r / Math.Sqrt(n);
// 計算出有多少行列的cell
var cols = (int)Math.Ceiling(width / cell_size);
var rows = (int)Math.Ceiling(height / cell_size);
// cells記錄了所有合法的點
var cells = new List<Vector2>();
// grids記錄了每個cell内的點在cells裡的索引,-1表示沒有點
var grids = new int[rows, cols];
for (var i = 0; i < rows; ++i) {
for (var j = 0; j < cols; ++j) {
grids[i, j] = -1;
}
}
// STEP 1
var random = new Random(seed);
// 随機選一個起始點
var x0 = new Vector2(random.Range(width), random.Range(height));
var col = (int)Math.Floor(x0.x / cell_size);
var row = (int)Math.Floor(x0.y / cell_size);
var x0_idx = cells.Count;
cells.Add(x0);
grids[row, col] = x0_idx;
var active_list = new List<int>();
active_list.Add(x0_idx);
// STEP 2
while (active_list.Count > 0) {
// 随機選一個待處理的點xi
var xi_idx = active_list[random.Range(active_list.Count)]; // 區間是[0,1),不用擔心溢出。
var xi = cells[xi_idx];
var found = false;
// 以xi為中點,随機找與xi距離在[r,2r)的點xk,并判斷該點的合法性
// 重複k次,如果都找不到,則把xi從active_list中去掉,認為xi附近已經沒有合法點了
for (var i = 0; i < k; ++i) {
var dir = random.insideUnitCircle();
var xk = xi + (dir.normalized * r + dir * r); // [r,2r)
if (xk.x < 0 || xk.x >= width || xk.y < 0 || xk.y >= height) {
continue;
}
col = (int)Math.Floor(xk.x / cell_size);
row = (int)Math.Floor(xk.y / cell_size);
if (grids[row, col] != -1) {
continue;
}
// 要判斷xk的合法性,就是要判斷附近有沒有點與xk的距離小于r
// 由于cell的邊長小于r,是以隻測試xk所在的cell的九宮格是不夠的(考慮xk正好處于cell的邊緣的情況)
// 正确做法是以xk為中心,做一個邊長為2r的正方形,測試這個正方形覆寫到的所有cell
var ok = true;
var min_r = (int)Math.Floor((xk.y - r) / cell_size);
var max_r = (int)Math.Floor((xk.y + r) / cell_size);
var min_c = (int)Math.Floor((xk.x - r) / cell_size);
var max_c = (int)Math.Floor((xk.x + r) / cell_size);
for (var or = min_r; or <= max_r; ++or) {
if (or < 0 || or >= rows) {
continue;
}
for (var oc = min_c; oc <= max_c; ++oc) {
if (oc < 0 || oc >= cols) {
continue;
}
var xj_idx = grids[or, oc];
if (xj_idx != -1) {
var xj = cells[xj_idx];
var dist = (xj - xk).magnitude;
if (dist < r) {
ok = false;
goto end_of_distance_check;
}
}
}
}
end_of_distance_check:
if (ok) {
var xk_idx = cells.Count;
cells.Add(xk);
grids[row, col] = xk_idx;
active_list.Add(xk_idx);
found = true;
break;
}
}
if (!found) {
active_list.Remove(xi_idx);
}
}
return cells;
}
}
// 測試代碼
class Program
{
static void Main(string[] args)
{
var width = 1024;
var height = 1024;
var r = 50f;
var points = Algorithm.Sample2D(width, height, r);
var image = new Bitmap(width, height);
using (var graphics = Graphics.FromImage(image)) {
graphics.FillRectangle(Brushes.Black, 0f, 0f, width, height);
var dot_r = 3f;
var pen = new Pen(Color.DarkRed, 2f);
foreach (var p in points) {
graphics.FillEllipse(Brushes.Yellow, p.x - dot_r, p.y - dot_r, 2f * dot_r, 2f * dot_r);
graphics.DrawEllipse(pen, p.x - r / 2f, p.y - r / 2f, r, r);
}
}
image.Save("out.png");
}
}
輸出的圖檔:可以看到黃點分布随機且比較平均,且任意2個黃點之間的距離都小于 r r r(紅色圓的半徑是 r / 2 r/2 r/2)