The Math Forum



Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.


Math Forum » Discussions » Software » comp.soft-sys.matlab

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Topic: Voronoi diagram on a sphere
Replies: 11   Last Post: Apr 14, 2013 3:44 AM

Advanced Search

Back to Topic List Back to Topic List Jump to Tree View Jump to Tree View   Messages: [ Previous | Next ]
Abel Brown

Posts: 62
Registered: 10/4/07
Re: Voronoi diagram on a sphere
Posted: Nov 6, 2009 10:10 AM
  Click to see the message monospaced in plain text Plain Text   Click to reply to this topic Reply

"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!!



Point your RSS reader here for a feed of the latest messages in this topic.

[Privacy Policy] [Terms of Use]

© The Math Forum at NCTM 1994-2018. All Rights Reserved.