5

How to determine, using Qhull, which voronoi cells (by index) are "proper" (made of "existing vertices")

I am trying to perform a constrained relaxation using LLoyds algorithm and input generated by scipy.spatial Voronoi (which is a wrapper for Qhull).

In terms of code it looks like:

points = [n  for n in itertools.product(xrange(3),xrange(3))]
vor = Voronoi(points)

vor2 = lloyd(vor) # my relaxation function  - not relevant to the question

The output graph generated by the code looks ok (see below) but the data in the vor structure is not enough to perform the Lloyds relaxation. It is because I should only move the points that are inside the valid voronoi cells ( #4 in the image). The other should be left as they are. Qhull messes the order of points/regions so I can not estimate which region belongs to which point.

Here is the illustration of the problem: 

print vor.vertices
#[[ 0.5  0.5]
# [ 1.5  0.5]
# [ 0.5  1.5]
# [ 1.5  1.5]]

print vor.regions  
# [[], [-1, 0], [-1, 1], [1, -1, 0], [3, -1, 2], [-1, 3], [-1, 2], [3, 2, 0, 1], [2, -1, 0], [3, -1, 1]]

print vor.points
# [[ 0.  0.]
#  [ 0.  1.]
#  [ 0.  2.]
#  [ 1.  0.]
#  [ 1.  1.]
#  [ 1.  2.]
#  [ 2.  0.]
#  [ 2.  1.]
#  [ 2.  2.]]
print vor.point_region
# [1 8 6 3 7 4 2 9 5]

Now i should find somehow that vor.regions[7] is a region that belongs to the point vor.points[4]. How to achieve this?

Voronoi on 3x3 grid

Chris
  • 570
  • 2
  • 9
  • 19

1 Answers1

6

You have a region you know, a point you don't, and you know that vor.point_region[point] == region. For a single region, you can figure out the corresponding point as:

point = np.argwhere(vor.point_region == region)

You can also create a region_point indexing array to figure out multiple points from an array of regions as:

region_point = np.argsort(vor.point_region)
points = region_point[regions-1]
Jaime
  • 65,696
  • 17
  • 124
  • 159
  • Perhaps I did not understand your answer but when I run the command suggested by you: np.argwhere(vor.point_region == 4) the result is: array([[5]]) which is not correct. :( – Chris Jul 14 '13 at 13:21
  • 2
    That's because you should be running `np.argwhere(vor.point_region == 7)`. You know your region (7), it is the point (4) you are after, right? Because to get the region from the point you simply do `vor.point_region[4]`. – Jaime Jul 14 '13 at 15:34