[comp.graphics] Polygon Representation of a Sphere's Surface

richr@ai.etl.army.mil (Richard Rosenthal) (03/23/89)

I want to have in 3-D space (x, y, z) a polygon near-representation
of the surface of a sphere where each polygon is identical
and regular (if that's the word).

Is an icosahedron the right place to start?

I would like to be able to generate the representation with
increasing numbers of polygons, say first 20, and then
say 80.

Ultimately, I would like to divide up the Earth's surface into
pieces that uniformly tile a sphere's surface with, say,
64,00 tiles (preferable triangles) but maybe upto 233 million tiles!

Any comments, algorithms, or reference suggestions?

/s/ Rich

-- 
Richard Rosenthal                 Internet:  richr@ai.etl.army.mil
Engineer Topographic Labs             UUCP:  ...!ames!ai.etl.army.mil!richr
Ft. Belvoir, VA 22060-5546          BITNET:  richr%ai.etl.army.mil@CUNYVM
+1 202 355 3653                      CSNET:  richr%ai.etl.army.mil@RELAY.CS.NET

foo@titan.rice.edu (Mark Hall) (03/23/89)

In article <270@ai.etl.army.mil> richr@ai.etl.army.mil. (Richard Rosenthal) writes:
>I want to have in 3-D space (x, y, z) a polygon near-representation
>of the surface of a sphere where each polygon is identical
>and regular (if that's the word).
>
>Is an icosahedron the right place to start?
>
>I would like to be able to generate the representation with
>increasing numbers of polygons, say first 20, and then
>say 80.
>Richard Rosenthal                 Internet:  richr@ai.etl.army.mil

  I think you are going to have a problem. As I remember the argument 
 from a class a few years ago: 

  One way to tile a sphere with regular (all sides/angles equal) polyhedra 
 is to take a regular solid centered at the sphere's center. Project 
 the edges of the solid onto the sphere. Works just fine. 

  The kicker is that there are only 5 known regular solids. (Not 
 counting the teapotahedron -  see CACM Feb. 1988 cover)

  If you do discover a tiling with millions of regular shapes, you 
 should be able to run the above process in reverse to create a new
 regular polyhedral solid.  

  Good luck.

 - mark

 

beshers@cs.cs.columbia.edu (Clifford Beshers) (03/23/89)

In article <2908@kalliope.rice.edu> foo@titan.rice.edu (Mark Hall) writes:

   Summary: you may have do discover a new "Platonic" solid
   Lines: 31

   In article <270@ai.etl.army.mil> richr@ai.etl.army.mil. (Richard Rosenthal) writes:
   >I want to have in 3-D space (x, y, z) a polygon near-representation
   >of the surface of a sphere where each polygon is identical
   >and regular (if that's the word).
   >
   >Is an icosahedron the right place to start?
   >
   >I would like to be able to generate the representation with
   >increasing numbers of polygons, say first 20, and then
   >say 80.
   >Richard Rosenthal                 Internet:  richr@ai.etl.army.mil

     I think you are going to have a problem. As I remember the argument 
    from a class a few years ago: 

Not at all.  Start with an icosahedron.  You now have a polygonal
approximation to a sphere where each face is a triangle.
Inscribe a triangle inside each face with its vertices at the
midpoints of the edges, dividing each face into four triangular
regions.  Project the three new points onto the surface of the
sphere.  You now have a a polygonal approximation to a sphere
where each face is a triangle.  Inscribe a triangle...

If you use Gouraud shading, you really don't have to recurse very far
to get something that looks very nice.
--
-----------------------------------------------
Cliff Beshers
Columbia University Computer Science Department
beshers@cs.columbia.edu

flynn@pixel.cps.msu.edu (Patrick J. Flynn) (03/23/89)

In article <2908@kalliope.rice.edu> foo@titan.rice.edu (Mark Hall) writes:
>In article <270@ai.etl.army.mil> richr@ai.etl.army.mil. (Richard Rosenthal) writes:
>>I want to have in 3-D space (x, y, z) a polygon near-representation
>>of the surface of a sphere where each polygon is identical
>>and regular (if that's the word). Is an icosahedron the right place to start?
>>I would like to be able to generate the representation with
>>increasing numbers of polygons, say first 20, and then
>>say 80.
>
>  I think you are going to have a problem. As I remember the argument 
> from a class a few years ago: 
>  One way to tile a sphere with regular (all sides/angles equal) polyhedra 
> is to take a regular solid centered at the sphere's center. Project 
> the edges of the solid onto the sphere. Works just fine. 
>  The kicker is that there are only 5 known regular solids.

Actually, one popular method of constructing a tesselation of the sphere
*starts* with a regular polyhedron (usually the icosahedron) and repeatedly
subdivides each facet into triangles.  The common method is to divide each
edge of a facet into k equal-length segments and construct k**2 equilateral
triangles from the original.  Each of these triangles can then be subdivided
to the desired resolution.

Ballard and Brown, "Computer Vision", p. 492 has a recipe for constructing
the icosahedron.

Also see

C. Brown, "Fast Display of well-tesselated surfaces", Computers and Graphics,
v. 4, p. 77-85, 1979.
--
Pat Flynn, CS, Mich. State U. | "What kind of chump do you take me for?" -Nick
flynn@cps.msu.edu             | "First Class!" -Rocky Rococo
(517) 353-4638                |

djs@mit-amt (David J. Sturman) (03/24/89)

In article <BESHERS.89Mar22205747@cs.cs.columbia.edu> beshers@cs.cs.columbia.edu
 (Clifford Beshers) writes:
 >In article <2908@kalliope.rice.edu> foo@titan.rice.edu (Mark Hall) writes:
 >
 >   In article <270@ai.etl.army.mil> richr@ai.etl.army.mil. (Richard Rosenthal)
 writes:
 >   >I want to have in 3-D space (x, y, z) a polygon near-representation
 >   >of the surface of a sphere where each polygon is identical
 >   >and regular (if that's the word).
 >   >
 >   >Is an icosahedron the right place to start?
 >   > ...
 >
 >Not at all.  Start with an icosahedron.  You now have a polygonal
 >approximation to a sphere where each face is a triangle.
 >Inscribe a triangle inside each face with its vertices at the
 >midpoints of the edges, dividing each face into four triangular
 >regions.  Project the three new points onto the surface of the
 >sphere.  You now have a a polygonal approximation to a sphere
 >where each face is a triangle.  Inscribe a triangle...
 >

 It seems to me that this will give you a polyhedral approximation to the
 sphere but NOT a regular polyhedron.  I suspect that the projection of
 the new points to the the surface of the sphere will unevenly
 deform the new triangles, since one vertex of each outer triangle is
 already on the sphere while all three vertices of the inner triangle
 are moved by the projection.  Thus, the interior triangle will no longer be
 isomorphic to the outer triangles.

 I have not worked out the geometry, and am working on intuition at this
 point.

 David Sturman
 MIT Media Lab
 djs@media-lab.media.mit.edu

kevinwu@lionsgate.Sun.COM (Kevin Wu) (03/25/89)

In article <3661@mit-amt> djs@media-lab.media.mit.edu (David J. Sturman) writes:
> It seems to me that this will give you a polyhedral approximation to the
> sphere but NOT a regular polyhedron.

David is correct.  You can use group theory to prove that there exists no
regular polyhedron with more than 20 facets.  I think that the original
poster really wanted a subjectively good polygonal approximation to a sphere,
not necessarily a regular polyhedron.  Subdividing the facets of a regular
icosahedron as described earlier gives a polyhedron that looks pretty close
to a sphere, especially with smooth shading.


--
Kevin Wu
Sun Microsystems, Inc.			Internet:  kevinwu@sun.com
Mail Stop 8-01				UUCP:      ...!sun!kevinwu

leech@alanine.cs.unc.edu (Jonathan Leech) (03/25/89)

In article <270@ai.etl.army.mil> richr@ai.etl.army.mil. (Richard Rosenthal) writes:
>I want to have in 3-D space (x, y, z) a polygon near-representation
>of the surface of a sphere where each polygon is identical
>and regular (if that's the word).
>...
>Any comments, algorithms, or reference suggestions?

    This request comes up (seemingly) once a month or so.  For a
change, here's the code to do it instead of a description.  Note that
the triangles generated are all equilateral, but not all the same
size. This program generates tesselations with

    #triangles(n'th approximation) = 8 * 4^(n-1),

since the starting point is an octahedron. Starting with other
platonic solids would change the coefficient, of course.

    The program retains generated triangles in core. It could be
easily modified to print out results on the fly if memory is limited.

    Code modifications will be required to generate your preferred
database format - see the comments.

P.S. Eugene, or whoever is keeping the commonly asked questions list,
how about adding this?

/*
 * sphere - generate a polygon mesh approximating a sphere by
 *  recursive subdivision. First approximation is an octahedron;
 *  each level of refinement increases the number of polygons by
 *  a factor of 4.
 * Level 3 (128 polygons) is a good tradeoff if gouraud
 *  shading is used to render the database.
 *
 * Usage: sphere [level]
 *	level is an integer >= 1 setting the recursion level (default 1).
 *
 * Notes:
 *
 *  The triangles are generated with vertices in clockwise order as
 *  viewed from the outside in a right-handed coordinate system.
 *  To reverse the order, compile with COUNTERCLOCKWISE defined.
 *
 *  Shared vertices are not retained, so numerical errors may produce
 *  cracks between polygons at high subdivision levels.
 *
 *  The subroutines print_object() and print_triangle() should
 *  be changed to generate whatever the desired database format is.
 *  If UNC is defined, a PPHIGS text archive is generated.
 *
 * Jon Leech 3/24/89
 */
#include <stdio.h>
#include <math.h>

typedef struct {
    double  x, y, z;
} point;

typedef struct {
    point     pt[3];	/* Vertices of triangle */
    double    area;	/* Unused; might be used for adaptive subdivision */
} triangle;

typedef struct {
    int       npoly;	/* # of polygons in object */
    triangle *poly;	/* Polygons in no particular order */
} object;

/* Six equidistant points lying on the unit sphere */
#define XPLUS {  1,  0,  0 }	/*  X */
#define XMIN  { -1,  0,  0 }	/* -X */
#define YPLUS {  0,  1,  0 }	/*  Y */
#define YMIN  {  0, -1,  0 }	/* -Y */
#define ZPLUS {  0,  0,  1 }	/*  Z */
#define ZMIN  {  0,  0, -1 }	/* -Z */

/* Vertices of a unit octahedron */
triangle octahedron[] = {
    { XPLUS, ZPLUS, YPLUS }, 0.0,
    { YPLUS, ZPLUS, XMIN  }, 0.0,
    { XMIN , ZPLUS, YMIN  }, 0.0,
    { YMIN , ZPLUS, XPLUS }, 0.0,
    { XPLUS, YPLUS, ZMIN  }, 0.0,
    { YPLUS, XMIN , ZMIN  }, 0.0,
    { XMIN , YMIN , ZMIN  }, 0.0,
    { YMIN , XPLUS, ZMIN  }, 0.0
};

/* An octahedron */
object oct = {
    sizeof(octahedron) / sizeof(octahedron[0]),
    &octahedron[0]
};

/* Forward declarations */
point *normalize(/* point *p */);
point *midpoint(/* point *a, point *b */);
void print_object(/* object *obj, int level */);
void print_triangle(/* triangle *t */);
double sqr(/* double x */);
double area_of_triangle(/* triangle *t */);

extern char *malloc(/* unsigned */);

main(ac, av)
int ac;
char *av[];
{
    object *old, *new;
    int     i, level, maxlevels = 1;

    if (ac > 1)
	if ((maxlevels = atoi(av[1])) < 1) {
	    fprintf(stderr, "%s: # of levels must be >= 1\n", av[0]);
	    exit(1);
	}

#ifdef COUNTERCLOCKWISE
    /* Reverse order of points in each triangle */
    for (i = 0; i < oct.npoly; i++) {
	point tmp;
		      tmp = oct.poly[i].pt[0];
	oct.poly[i].pt[0] = oct.poly[i].pt[2];
	oct.poly[i].pt[2] = tmp;
    }
#endif

    old = &oct;

    /* Subdivide each starting triangle (maxlevels - 1) times */
    for (level = 1; level < maxlevels; level++) {
	/* Allocate a new object */
	new = (object *)malloc(sizeof(object));
	if (new == NULL) {
	    fprintf(stderr, "%s: Out of memory on subdivision level %d\n",
		av[0], level);
	    exit(1);
	}
	new->npoly = old->npoly * 4;

	/* Allocate 4* the number of points in the current approximation */
	new->poly  = (triangle *)malloc(new->npoly * sizeof(triangle));
	if (new->poly == NULL) {
	    fprintf(stderr, "%s: Out of memory on subdivision level %d\n",
		av[0], level);
	    exit(1);
	}

	/* Subdivide each polygon in the old approximation and normalize
	 *  the new points thus generated to lie on the surface of the unit
	 *  sphere.
	 * Each input triangle with vertices labelled [0,1,2] as shown
	 *  below will be turned into four new triangles:
	 *
	 *			Make new points
	 *			    a = (0+2)/2
	 *			    b = (0+1)/2
	 *			    c = (1+2)/2
	 *	  1
	 *	 /\		Normalize a, b, c
	 *	/  \
	 *    b/____\ c		Construct new triangles
	 *    /\    /\		    [0,b,a]
	 *   /	\  /  \		    [b,1,c]
	 *  /____\/____\	    [a,b,c]
	 * 0	  a	2	    [a,c,2]
	 */
	for (i = 0; i < old->npoly; i++) {
	    triangle
		 *oldt = &old->poly[i],
		 *newt = &new->poly[i*4];
	    point a, b, c;

	    a = *normalize(midpoint(&oldt->pt[0], &oldt->pt[2]));
	    b = *normalize(midpoint(&oldt->pt[0], &oldt->pt[1]));
	    c = *normalize(midpoint(&oldt->pt[1], &oldt->pt[2]));

	    newt->pt[0] = oldt->pt[0];
	    newt->pt[1] = b;
	    newt->pt[2] = a;
	    newt++;

	    newt->pt[0] = b;
	    newt->pt[1] = oldt->pt[1];
	    newt->pt[2] = c;
	    newt++;

	    newt->pt[0] = a;
	    newt->pt[1] = b;
	    newt->pt[2] = c;
	    newt++;

	    newt->pt[0] = a;
	    newt->pt[1] = c;
	    newt->pt[2] = oldt->pt[2];
	}

	if (level > 1) {
	    free(old->poly);
	    free(old);
	}

	/* Continue subdividing new triangles */
	old = new;
    }

    /* Print out resulting approximation */
    print_object(old, maxlevels);
}

/* Normalize a point p */
point *normalize(p)
point *p;
{
    static point r;
    double mag;

    r = *p;
    mag = r.x * r.x + r.y * r.y + r.z * r.z;
    if (mag != 0.0) {
	mag = 1.0 / sqrt(mag);
	r.x *= mag;
	r.y *= mag;
	r.z *= mag;
    }

    return &r;
}

/* Return the average of two points */
point *midpoint(a, b)
point *a, *b;
{
    static point r;

    r.x = (a->x + b->x) * 0.5;
    r.y = (a->y + b->y) * 0.5;
    r.z = (a->z + b->z) * 0.5;

    return &r;
}

/* Write out all polygons in an object */
void print_object(obj, level)
object *obj;
int level;
{
    int i;

#ifdef UNC  /* Generate PPHIGS archive format */
    int dx, dy, dz;

    printf("structure sphere%d posted {\n", level);
    printf("\tcolor polygon {\n");
    printf("\t\t200 100  50   0     50 100 200   0\n");
    printf("\t};\n");

    switch (level) {
	case 1:
	    dx = -2000; dy =  2000; dz = 0;
	    break;
	case 2:
	    dx =  2000; dy =  2000; dz = 0;
	    break;
	case 3:
	    dx = -2000; dy = -2000; dz = 0;
	    break;
	case 4:
	    dx =  2000; dy = -2000; dz = 0;
	    break;
	case 5:
	    dx =     0; dy =	 0; dz = 0;
	    break;
	default:
	    dx = dy = dz = 0;
	    break;
    }

    printf("\tmatrix Pre scale 1000 1000 1000;\n");
    printf("\tmatrix Pre translate %d %d %d ;\n", dx, dy, dz);
#endif	/* UNC */

    /* Spit out coordinates for each triangle */
    for (i = 0; i < obj->npoly; i++)
	print_triangle(&obj->poly[i]);

#ifdef UNC
    printf("};\n");
#endif	/* UNC */
}

/* Write a single polygon */
void print_triangle(t)
triangle *t;
{
    int i;

#ifdef UNC
    printf("\tpolygon 3 {\n");
    for (i = 0; i < 3; i++)
	printf("\t\t%g\t%g\t%g %g\t%g\t%g ;\n",
	    t->pt[i].x, t->pt[i].y, t->pt[i].z,     /* Point */
	    t->pt[i].x, t->pt[i].y, t->pt[i].z);    /* Normal */
    printf("\t};\n");
#else
    /* Insert code for your output format here.
     * Triangle vertices are in t->pt[0..2].{x,y,z}
     */
    static int errflag = 0;
    fprintf(stderr, "Error! You must modify print_triangle() to generate triangles in\n");
    fprintf(stderr, "\tthe desired format!\n");
    exit(1);
#endif
}

double sqr(x) double x; { return x * x ; }

/* Returns area of triangle (currently unused) */
double area_of_triangle(t)
triangle *t;
{
    double a;

    a  = sqr(t->pt[0].y * (t->pt[1].z - t->pt[2].z) -
	     t->pt[1].y * (t->pt[0].z - t->pt[2].z) +
	     t->pt[2].y * (t->pt[0].z - t->pt[1].z));
    a += sqr(t->pt[0].z * (t->pt[1].x - t->pt[2].x) -
	     t->pt[1].z * (t->pt[0].x - t->pt[2].x) +
	     t->pt[2].z * (t->pt[0].x - t->pt[1].x));
    a += sqr(t->pt[0].x * (t->pt[1].y - t->pt[2].y) -
	     t->pt[1].x * (t->pt[0].y - t->pt[2].y) +
	     t->pt[2].x * (t->pt[0].y - t->pt[1].y));

    return 0.5 * sqrt(a);
}
--
    Jon Leech (leech@cs.unc.edu)    __@/
    SUSHIDO: the Way of the Tuna

jdchrist@watcgl.waterloo.edu (Dan Christensen) (03/27/89)

In article <7409@thorin.cs.unc.edu> leech@alanine.UUCP (Jonathan Leech) writes:
>    This request comes up (seemingly) once a month or so.  For a
>change, here's the code to do it instead of a description.  Note that
>the triangles generated are all equilateral, but not all the same
>size. This program generates tesselations with
>
>--
>    Jon Leech (leech@cs.unc.edu)    __@/
>    SUSHIDO: the Way of the Tuna

Are you sure that the triangles are all equilateral?  I wrote a program a
while ago that seems to do exactly what yours does and they some of the
triangles that were generated looked isoceles.  In fact, after several
iterations, the triangles became extremely thin and so the polygonization
was not very efficient.  Starting with a 20 sided regular figure (I forget
its name) helps.

----
Dan Christensen, Computer Graphics Lab,	         jdchrist@watcgl.uwaterloo.ca
University of Waterloo, Waterloo, Ont.	         jdchrist@watcgl.waterloo.edu

leech@zeta.cs.unc.edu (Jonathan Leech) (03/27/89)

In article <8836@watcgl.waterloo.edu> jdchrist@watcgl.waterloo.edu (Dan Christensen) writes:
>Are you sure that the triangles are all equilateral?  I wrote a program a
>while ago that seems to do exactly what yours does and they some of the
>triangles that were generated looked isoceles.

    Oops, yeah, you're correct. They looked equilateral :-)
--
    Jon Leech (leech@cs.unc.edu)    __@/
    ``Those what cannot remedy the past can pretend to repeal it."
	- Attributed to Santa Ana by Owl