home *** CD-ROM | disk | FTP | other *** search
/ Photo CD Demo 1 / Demo.bin / gems / gemsiii / intell.c < prev    next >
C/C++ Source or Header  |  1992-03-16  |  3KB  |  95 lines
  1. #include    "GraphicsGems.h"
  2.  
  3. /* ----    intell - Intersect a ray with an ellipsoid. -------------------    */
  4. /*                                    */
  5. /*                                    */
  6. /*    Description:                            */
  7. /*        Intell determines the intersection of a ray with an        */
  8. /*        ellipsoid.                            */
  9. /*                                    */
  10. /*    On entry:                            */
  11. /*        raybase = The coordinate defining the base of the        */
  12. /*              intersecting ray.                    */
  13. /*        raycos  = The direction cosines of the above ray.        */
  14. /*        center  = The location of the center of the ellispoid.    */
  15. /*        radius  = The radius of a bounding sphere.            */
  16. /*        Q       = A 4x4 coefficient matrix representing the quadric    */
  17. /*              surface.                        */
  18. /*                                    */
  19. /*    On return:                            */
  20. /*        rin     = The entering distance of the intersection.    */
  21. /*        rout    = The leaving  distance of the intersection.    */
  22. /*                                    */
  23. /*    Returns:  True if the ray intersects the ellipsoid.        */
  24. /*                                    */
  25. /* --------------------------------------------------------------------    */
  26.  
  27.  
  28. boolean    intell    (raybase,raycos,center,radius,Q,rin,rout)
  29.  
  30.     Point3    raybase;        /* Base of the intersection ray    */
  31.     Vector3    raycos;            /* Direction cosines of the ray    */
  32.     Vector3    center;            /* Base point for the quadric    */
  33.     double    radius;            /* Bounding radius        */
  34.     Matrix4    Q;            /* Quadric coefficient matrix    */
  35.     double    *rin;            /* Entering distance        */
  36.     double    *rout;            /* Leaving  distance        */
  37.  
  38. {
  39.     boolean    hit;            /* True if ray intersects quad.    */
  40.     double    qa, qb, qc, qd, qe,    /* Coefficient matrix terms    */
  41.         qf, qg, qh, qi, qj;
  42.     double    k2, k1, k0;        /* Quadric coefficients        */
  43.     double    disc;
  44.     Vector3    D;            /* Ray base to ellipsoid base    */
  45.     double    d;            /* Shortest distance between    */
  46.                     /* ray and ellipsoid        */
  47.  
  48.  
  49. /*    Compute the minimal distance between the ray and the center    */
  50. /*    of the ellipsoid.                        */
  51.  
  52.     V3Sub (¢er,&raybase,&D);
  53.     d   = V3Length (&D);
  54.     V3Normalize (&D);
  55.     d   = d * sqrt (1. - V3Dot(&D,&raycos));
  56.     hit  = (d <= radius);
  57.   
  58.     if  (!hit) return (hit);        /* If outside of sphere    */
  59.  
  60.  
  61.     qa = Q.element[0][0];
  62.     qb = Q.element[0][1] + Q.element[1][0];
  63.     qc = Q.element[0][2] + Q.element[2][0];
  64.     qd = Q.element[0][3] + Q.element[3][0];
  65.     qe = Q.element[1][1];
  66.     qf = Q.element[1][2] + Q.element[2][1];
  67.     qg = Q.element[1][3] + Q.element[3][1];
  68.     qh = Q.element[2][2];
  69.     qi = Q.element[2][3] + Q.element[3][2];
  70.     qj = Q.element[3][3];
  71.  
  72.     k2 = raycos.x * (raycos.x*qa + raycos.y*qb) +
  73.          raycos.y * (raycos.y*qe + raycos.z*qf) +
  74.          raycos.z * (raycos.z*qh + raycos.z*qc);
  75.     k1 = raycos.x * ((raybase.x*2.*qa + raybase.y*qb + raybase.z*qc) + qd) +
  76.          raycos.y * ((raybase.x*qb + raybase.y*2.*qe + raybase.z*qf) + qg) +
  77.          raycos.z * ((raybase.x*qc + raybase.y*qg + raybase.z*2.*qh) + qi);
  78.     k0 = raybase.x * (raybase.x*qa + raybase.y*qb + raybase.z*qc + qd) +
  79.          raybase.y * (               raybase.y*qe + raybase.z*qf + qg) +
  80.          raybase.z * (                              raybase.z*qh + qi) +
  81.          qj;
  82.  
  83.     disc = k1*k1 - 4.*k2*k0;
  84.   
  85.     hit  = (disc >= 0.0);
  86.  
  87.     if  (hit) {                 /* If ray hits quadric    */
  88.         disc  =  sqrt(disc);
  89.         *rin  = (-k1 - disc) / (2.*k2);    /*    entering distance    */
  90.         *rout = (-k1 + disc) / (2.*k2);    /*    leaving distance    */
  91.     }
  92.   
  93.     return (hit);
  94. }
  95.