sampling

A module for generating points with a user specified distribution.

rbf.pde.sampling.poisson_discs(rfunc, domain, seeds=10, ntests=50, rmax_factor=1.5)

Generates Poisson disc points within the domain defined by vert and smp. Poisson disc points are tightly packed but are no closer than a user specified value. This algorithm is based on [1], and it has been modified to allow for spatially variable spacing. This works for two and three spatial dimension.

Parameters:

  • rfunc (callable) – A function that takes a (n, d) array of points as input and returns the desired minimum nearest neighbor distance for those points.

  • domain ((p, d) float array and (q, d) int array) – The vertices of the domain and the connectivity of the vertices

  • seeds (int) – The number of initial points, which are generated from a Halton sequence.

  • ntests (int, optional) – The maximum number of attempts at finding a new neighbor for a point before giving up. Increasing this generally results in tighter packing.

  • rmax_factor (float, optional) –

Returns:

(n, d) float array

References

[1] Bridson, R., Fast Poisson Disk Sampling in Arbitrary Dimensions.

rbf.pde.sampling.rejection_sampling(size, rho, domain, start=0, max_sample_size=1000000)

Returns points within the boundaries defined by vert and smp and with density rho. The nodes are generated by rejection sampling.

Parameters:

  • size (int) – Number of points to return

  • rho (callable) – A function that takes an (n, d) array of points and returns the density at those points. This should be normalized so that the density is between 0 and 1.

  • domain ((p, d) float array and (q, d) int array) – The vertices of the domain and the connectivity of the vertices

  • start (int, optional) – The starting index for the Halton sequence, which is used to propose new points. Setting this value is akin to setting the seed for a random number generator.

  • max_sample_size (int, optional) – max number of nodes allowed in a sample for the rejection algorithm. This prevents excessive RAM usage

Returns:

(size, d) float array

Examples

'''
This script demonstrates how to generate poisson discs with variable
radii.
'''
import time
import numpy as np
from scipy.spatial import cKDTree

import matplotlib.pyplot as plt
from matplotlib.patches import Circle

from rbf.pde.sampling import poisson_discs

def radius(x):
    '''disc radius as a function of position'''
    r = np.linalg.norm(x - 0.5, axis=1)
    return 0.01 + 0.1*r

# define the domain
vert = np.array([[0.0, 0.0], [1.0, 0.0], [1.0, 1.0], [0.0, 1.0]])
smp = np.array([[0, 1], [1, 2], [2, 3], [3, 0]])

# generate the discs
start = time.time()
centers = poisson_discs(radius, (vert, smp))
runtime = time.time() - start
print('generated %s samples in %s s' % (len(centers), runtime))

# plot the domain and discs
fig, ax = plt.subplots()
for s in smp:
    ax.plot(vert[s, 0], vert[s, 1], 'k-')

for c, r in zip(centers, radius(centers)):
    ax.plot(c[0], c[1], 'C0.')
    ax.add_artist(Circle(c, r, color='C0', alpha=0.1))

ax.set_aspect('equal')
ax.set_title('Poisson discs with variable radii')
plt.tight_layout()
plt.savefig('../figures/sampling.a.1.png')

# plot the nearest neighbor distance divided by the disc radii. This
# verifies that all nodes are spaced sufficiently far away
fig, ax = plt.subplots()
dist = cKDTree(centers).query(centers, 2)[0][:, 1]
ax.hist(dist/radius(centers), 20)
ax.set_xlabel('(nearest neighbor distance) / (disc radius)')
ax.set_ylabel('count')
ax.grid(ls=':', color='k')
plt.tight_layout()
plt.savefig('../figures/sampling.a.2.png')
plt.show()
_images/sampling.a.1.png _images/sampling.a.2.png