Finding random x,y coordinates in a figure where you know the boundary

Question

The goal is to find coordinates in a figure with an unknown shape. What IS known is a list of coordinates of the boundary of that figure, for example:

boundary = [(0,0),(1,0),(2,0),(3,0),(3,1),(3,2),(3,3),(2,3),(2,2),(1,2),(1,3),(0,3),(0,2),(0,1]

which would look something like this:

Square with a gab

This is a very basic example and i'd like to do it with very larg lists of very different kinds of figures.

The question is how to get a random coordinate that lies within the figure WITHOUT hardcoding the anything about the shape of the figure, because this will be unknown at the beginning? Is there a way to know for certain or is making an estimate the best option? How would I implement an estimate like that?

Answer 1

Here is tentative answer. You sample numbers in two steps.

Before, do preparation work - split your figure into simple elementary objects. In your case you split it into rectangles, often people triangulate and split it into triangles.

So you have number N of simple objects, each with area of A_i and total area A = Sum(A_i).

First sampling step - select which rectangle you pick point from. In some pseudocode

r = randomU01(); // random value in [0...1) range
for(i in N) {
    r = r - A_i/A;
    if (r <= 0) {
        k = i;
        break;
    }
}

So you picked up one rectangle with index k, and then just sample point uniformly in that rectangle

x = A_k.dim.x * randomU01();
y = A_k.dim.y * randomU01();

return (x + A_k.lower_left_corner.x, y + A_k.lower_left_corner.y);

And that is it. Very similar technique for triangulated figure.

Rectangle selection could be optimized by doing binary search or even more complicated alias method

UPDATE

If your boundary is generic, then the only good way to go is to triangulate your polygon using any good library out there (f.e. Triangle), then select one of the triangles based on area (step 1), then sample uniformly point in the triangle using two random U01 numbers r1 and r2,

P = (1 - sqrt(r1)) * A + (sqrt(r1)*(1 - r2)) * B + (r2*sqrt(r1)) * C

i.e., in pseudocode

r1 = randomU01();
s1 = sqrt(r1);
r2 = randomU01();

x = (1.0-s1)*A.x + s1*(1.0-r2)*B.x + r2*s1*C.x;
y = (1.0-s1)*A.y + s1*(1.0-r2)*B.y + r2*s1*C.y;

return (x,y);