[Numpy-discussion] floating point arithmetic issue

[Pauli Virtanen]

Fri, 30 Jul 2010 19:48:45 +0200, Guillaume Chérel wrote:
> Your solution is really good. It's almost exactly what I was doing, but
> shorter and I didn't know about the mgrid thing.

It's a very brute-force solution and probably won't be very fast with 
many circles or large grids.

> There is an extra optimization step I perform in my code to avoid the
> computation of the many points of the grid (assuming the circles are
> relatively small) which you know only by their x or y coordinate they're
> out of the circle without computing the distance. And that's where I
> need the modulo operator to compute the left-most and top-most points of
> the grid that are inside the circle.

If your circles are quite small, you probably want to clip the "painting" 
to a box not much larger than a single circle:

# untested, something like below
def point_to_index(x, y, pad=0):
    return np.clip(200 * rr / (xmax - xmin) + pad, 0, 200), \
           np.clip(200 * rr / (ymax - ymin) + pad, 0, 200)

i0, j0 = point_to_index(xx - rr, yy - rr, pad=-2)
i1, j1 = point_to_index(xx + rr, yy + rr, pad=2)

box = np.index_exp[i0:j0,i1:j1]
mask[box] |= (grid_x[box] - xx)**2 + (grid_y[box] - yy)**2 < rr**2
# same as: mask[i0:j0,i1:j1] |= (grid_x[i0:j0,i1:j1] ...

-- 
Pauli Virtanen