home *** CD-ROM | disk | FTP | other *** search
/ RISC DISC 3 / RISC_DISC_3.iso / resources / etexts / gems / gemsi / boundsphere.c < prev    next >
C/C++ Source or Header  |  1992-09-08  |  4KB  |  134 lines
  1. /*
  2. An Efficient Bounding Sphere
  3. by Jack Ritter
  4. from "Graphics Gems", Academic Press, 1990
  5. */
  6.  
  7. /* Routine to calculate tight bounding sphere over    */
  8. /* a set of points in 3D */
  9. /* This contains the routine find_bounding_sphere(), */
  10. /* the struct definition, and the globals used for parameters. */
  11. /* The abs() of all coordinates must be < BIGNUMBER */
  12. /* Code written by Jack Ritter and Lyle Rains. */
  13.  
  14. #include <stdio.h>
  15. #include <math.h>
  16. #include "GraphicsGems.h"
  17.  
  18. #define BIGNUMBER 100000000.0          /* hundred million */
  19.  
  20. /* GLOBALS. These are used as input and output parameters. */
  21.  
  22. struct Point3Struct caller_p,cen;
  23. double rad;
  24. int NUM_POINTS;
  25.  
  26. /* Call with no parameters. Caller must set NUM_POINTS */
  27. /* before calling. */
  28. /* Caller must supply the routine GET_iTH_POINT(i), which loads his */
  29. /* ith point into the global struct caller_p. (0 <= i < NUM_POINTS). */
  30. /* The calling order of the points is irrelevant. */
  31. /* The final bounding sphere will be put into the globals */
  32. /* cen and rad. */
  33.  
  34.  
  35. find_bounding_sphere()
  36. {
  37. register int i;
  38. register double dx,dy,dz;
  39. register double rad_sq,xspan,yspan,zspan,maxspan;
  40. double old_to_p,old_to_p_sq,old_to_new;
  41. struct Point3Struct xmin,xmax,ymin,ymax,zmin,zmax,dia1,dia2;
  42.  
  43.  
  44. /* FIRST PASS: find 6 minima/maxima points */
  45. xmin.x=ymin.y=zmin.z= BIGNUMBER; /* initialize for min/max compare */
  46. xmax.x=ymax.y=zmax.z= -BIGNUMBER;
  47. for (i=0;i<NUM_POINTS;i++)
  48.     {
  49.     GET_iTH_POINT(i); /* load global struct caller_p with */
  50.                     /* his ith point. */
  51.     if (caller_p.x<xmin.x)
  52.         xmin = caller_p; /* New xminimum point */
  53.     if (caller_p.x>xmax.x)
  54.         xmax = caller_p;
  55.     if (caller_p.y<ymin.y)
  56.         ymin = caller_p;
  57.     if (caller_p.y>ymax.y)
  58.         ymax = caller_p;
  59.     if (caller_p.z<zmin.z)
  60.         zmin = caller_p;
  61.     if (caller_p.z>zmax.z)
  62.         zmax = caller_p;
  63.     }
  64. /* Set xspan = distance between the 2 points xmin & xmax (squared) */
  65. dx = xmax.x - xmin.x;
  66. dy = xmax.y - xmin.y;
  67. dz = xmax.z - xmin.z;
  68. xspan = dx*dx + dy*dy + dz*dz;
  69.  
  70. /* Same for y & z spans */
  71. dx = ymax.x - ymin.x;
  72. dy = ymax.y - ymin.y;
  73. dz = ymax.z - ymin.z;
  74. yspan = dx*dx + dy*dy + dz*dz;
  75.  
  76. dx = zmax.x - zmin.x;
  77. dy = zmax.y - zmin.y;
  78. dz = zmax.z - zmin.z;
  79. zspan = dx*dx + dy*dy + dz*dz;
  80.  
  81. /* Set points dia1 & dia2 to the maximally separated pair */
  82. dia1 = xmin; dia2 = xmax; /* assume xspan biggest */
  83. maxspan = xspan;
  84. if (yspan>maxspan)
  85.     {
  86.     maxspan = yspan;
  87.     dia1 = ymin; dia2 = ymax;
  88.     }
  89. if (zspan>maxspan)
  90.     {
  91.     dia1 = zmin; dia2 = zmax;
  92.     }
  93.  
  94.  
  95. /* dia1,dia2 is a diameter of initial sphere */
  96. /* calc initial center */
  97. cen.x = (dia1.x+dia2.x)/2.0;
  98. cen.y = (dia1.y+dia2.y)/2.0;
  99. cen.z = (dia1.z+dia2.z)/2.0;
  100. /* calculate initial radius**2 and radius */
  101. dx = dia2.x-cen.x; /* x component of radius vector */
  102. dy = dia2.y-cen.y; /* y component of radius vector */
  103. dz = dia2.z-cen.z; /* z component of radius vector */
  104. rad_sq = dx*dx + dy*dy + dz*dz;
  105. rad = sqrt(rad_sq);
  106.  
  107. /* SECOND PASS: increment current sphere */
  108.  
  109. for (i=0;i<NUM_POINTS;i++)
  110.     {
  111.     GET_iTH_POINT(i); /* load global struct caller_p  */
  112.                     /* with his ith point. */
  113.     dx = caller_p.x-cen.x;
  114.     dy = caller_p.y-cen.y;
  115.     dz = caller_p.z-cen.z;
  116.     old_to_p_sq = dx*dx + dy*dy + dz*dz;
  117.     if (old_to_p_sq > rad_sq)     /* do r**2 test first */
  118.         {     /* this point is outside of current sphere */
  119.         old_to_p = sqrt(old_to_p_sq);
  120.         /* calc radius of new sphere */
  121.         rad = (rad + old_to_p) / 2.0;
  122.         rad_sq = rad*rad;     /* for next r**2 compare */
  123.         old_to_new = old_to_p - rad;
  124.         /* calc center of new sphere */
  125.         cen.x = (rad*cen.x + old_to_new*caller_p.x) / old_to_p;
  126.         cen.y = (rad*cen.y + old_to_new*caller_p.y) / old_to_p;
  127.         cen.z = (rad*cen.z + old_to_new*caller_p.z) / old_to_p;
  128.         /* Suppress if desired */
  129.         printf("\n New sphere: cen,rad = %f %f %f   %f",
  130.             cen.x,cen.y,cen.z, rad);
  131.         }    
  132.     }
  133. }              /* end of find_bounding_sphere()  */
  134.