
Re: Voronoi diagram on a sphere
Posted:
Nov 6, 2009 10:10 AM


"Damian Sheehy" <Damian.Sheehy@mathworks.com> wrote in message <hd1a6r$rq0$1@fred.mathworks.com>... > Hi Abel, > > If you have release in R2009a or later, the computational geometry tools > make this task relatively easy. > The Voronoi diagram is the dual of the Delaunay triangulation; the > triangulation lying on the on the surface of the sphere. > > % Create the points on the surface of the sphere > theth_a = 2*pi*rand(100,1); > phi = pi*rand(100,1); > x = cos(theth_a).*sin(phi); > y = sin(theth_a).*sin(phi); > z = cos(phi); > > % Create a Delaunay triangulation of the point set. > % This is a solid triangulation composed of tetrahedra. > dt = DelaunayTri(x,y,z); > > % Extract the surface (boundary) triangulation, in this instance it is > actually the convex hull > % so you could use the convexHull method also. > ch = dt.freeBoundary(); > > % Create a Triangulation Representation from this surface triangulation > % We will use it to compute the location of the voronoi vertices > (circumcenter of the triangles), > % and the dual edges from the triangle neighbor information. > > tr = TriRep(ch,x,y,z); > numt = size(tr,1); > T = (1:numt)'; > neigh = tr.neighbors(); > cc = tr.circumcenters(); > xcc = cc(:,1); > ycc = cc(:,2); > zcc = cc(:,3); > idx1 = T < neigh(:,1); > idx2 = T < neigh(:,2); > idx3 = T < neigh(:,3); > neigh = [T(idx1) neigh(idx1,1); T(idx2) neigh(idx2,2); T(idx3) > neigh(idx3,3)]'; > plot3(xcc(neigh), ycc(neigh), zcc(neigh),'r'); > axis equal; > > This is just a 3D extension of the following example; > http://www.mathworks.com/products/demos/shipping/matlab/demoDelaunayTri.html#24 > > Best regards, > > Damian
Damian
Thank you very much! This is exactly what I needed!
Ultimately, I need to calculate the area of each of these polygon (which is no problem). Typically, using [V,C] = voronoin(X), the index of the i'th polygon
% faces polygonIndx = C{i}
Then to get the vertex information just index into V
vX = V(polygonIndx,1); vY = V(polygonIndx,2);
Is there a way, using your example, to get the vertices for each polygon?
Should I just computed the convexhull of each of the original points using the Voronoi points computed?
Thanks!!

