Hi Pauli, Your solution is really good. It's almost exactly what I was doing, but shorter and I didn't know about the mgrid thing.
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. I'll try with your solution, maybe my optimization isn't that necessary with the mgrid trick. On Fri, Jul 30, 2010 at 7:15 PM, Pauli Virtanen <p...@iki.fi> wrote: > Fri, 30 Jul 2010 14:21:23 +0200, Guillaume Chérel wrote: > [clip] >> As for the details about my problem, I'm trying to compute the total >> surface of overlapping disks. I approximate the surface with a grid and >> count how many points of the grid fall into at least one disk. > > HTH, > > import numpy as np > > # some random circles > x = np.random.randn(200) > y = np.random.randn(200) > radius = np.random.rand(200)*0.1 > > # grid, [-1, 1] x [-1, 1], 200 x 200 points > xmin, xmax = -1, 1 > ymin, ymax = -1, 1 > grid_x, grid_y = np.mgrid[xmin:xmax:200j, ymin:ymax:200j] > > # count > mask = np.zeros((200, 200), bool) > for xx, yy, rr in zip(x, y, radius): > mask |= (grid_x - xx)**2 + (grid_y - yy)**2 < rr**2 > > area = mask.sum() * (xmax-xmin) * (ymax-ymin) / float(mask.size) > > > _______________________________________________ > NumPy-Discussion mailing list > NumPy-Discussion@scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion > _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
