Linear-Time Poisson-Disk Patterns

There is a long-standing interest in Poisson-disk patterns in the graphics community, primarily for their use in sampling [Yel83,Coo86,Mit87]. There have been many algorithms for generating such patterns. Direct implementation of "dart-throwing" [Mit87,MF92] produces true Poisson-disk patterns, but is slow to converge. Approximations from relaxation [Llo83] or tiling [ODJ04,HDK01] can produce patterns similar to Poisson-disk patterns more efficiently. Recently, exact methods taking log-linear time (O(N log N ) where N is the total number of points) have been described [DH06,Jon06], as well as a method with empirical O(N ) speed, but lacking a rigorous proof of this performance [WCE07].

We present an algorithm with provable O(N ) performance. The algorithm maintains two data structures: a grid of regions in which points might still be inserted, and a bucket (i.e., an unordered set) of regions where a point will be generated (a subset of the grid). At each step of the algorithm, a region is taken from the bucket, a new point is inserted in that region, and nearby regions are updated and possibly added to the bucket. The bucket is only empty when no more points can be added (i.e., the pattern is maximal). The work for each iteration is O(1), for a total cost of O(N ).

Python source code is included in the ancillary data with this paper.

The Poisson Disk distribution can be defined as the limit of a uniform two-dimensional Poisson process with a minimum-distance rejection criterion. Successive points are independently drawn from the uniform distribution on [0, 1] 2 . If a new point is at least distance R from all points already accepted, it is also accepted. Otherwise, it is rejected. We call this the naïve algorithm. The choice of R controls the minimum distance between points (for N points in the unit square, πR 2 N/4 ≈ 0.548 as R → 0 [DWJ91]).

Efficient algorithms for Poisson-disk patterns rely on generating new points in regions where they are guaranteed (or highly probable) to be accepted [DH06,Jon06,WCE07]. In order to guarantee equivalence of results with the naïve algorithm, these methods have used O(log N ) area-weighted binary search to find where to insert a new point [DH06,Jon06], or weighted spatial indexing [WCE07] with theoretical O(log N ) but empirical O(1) cost.

Our algorithm can be seen as an optimization of the naïve algorithm using a spatial data structure. We Figure 1: A point p is shown with its neighbor grid squares in grey. If p’s arrival time is earlier than any of its neighbors, it will be accepted and added to the output. The free regions of p’s neighbors will then be updated to the dark gray areas. This may result in points such as q being invalidated. In q’s case, a replacement point q will be generated in its new free region, A q . The time of arrival of q will be t q + t + , where t + is drawn from an exponential distribution parameterized by the size of the updated free region. store a grid with spacing ≤ R/ √ 2 such that no more than one point can land in any grid square in the final pattern. We model (implicitly) a uniform 2D Poisson process on [0, 1] 2 , with rate λ = 1, by storing at each grid square the location and time of the earliest point landing in that square. These points and their arrival times will be updated as the algorithm progresses.

Each grid square has three associated pieces of data: the free region within that square where new points might be generated, a random point within that free region, and a time-of-arrival for that point under a the Poisson process. Initially, the free region for each grid square is the entire square, each point is chosen uniformly within its square, and the times of arrival are drawn from A 0 e -A0t , where A 0 is the area of a grid square.

We also define a neighbor relationship from points to grid squares, where the neighbors of a point are any grid squares within R of the point (see figure 1).

The first insight of our paper is that any point p that has time-of-arrival t p lower than any of its neighbors can be added to the output immediately, as this indicates that p arrives before any other point that could prevent it from being accepted. On acceptance, the free regions of p’s neighbors are updated (see figure 1).

It is possible that accepting p will invalidate a point q from another grid square with ||p -q|| 2 < R and t q > t p , in which case q is removed from the grid and a new point q in the updated free region is created with a new and later t q (see figure 1).

The second key insight of our algorithm is that the new t q should be t q plus a random variable drawn from the exponential distribution parameterized by the area of the updated free region A q , i.e., t q = t q + t + , where t + is drawn from A q e -A q t+ . The logic is as follows. The points p and q represent the first arrivals in their respective grid squares, at times t p and t q (in the naïve algorithm). When q is invalidated by p, the time until anot

This content is AI-processed based on open access ArXiv data.