home *** CD-ROM | disk | FTP | other *** search
/ Stars of Shareware: Raytrace & Morphing / SOS-RAYTRACE.ISO / programm / rad386 / radiosit / src / ray_cyli.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-04-21  |  5.0 KB  |  160 lines
  1. /*
  2. This file contains intersection of ray with Canonical Cylinder, Cone and Ring.
  3. 1)
  4.     Intersection of a ray with canonical Cylinder.
  5.     This canonical cylinder has axis coinciding
  6.     with +ve Z axis with one base on origin with
  7.     the Height H and the Radius R.
  8.  
  9.     The equation :
  10.         x^2 + y^2 - R^2 = 0
  11.  
  12.     "ray_cylinder" computes the intersection and
  13.     returns rayparameter t of intesection
  14.         where t > 0 for valid intersection,
  15.     and surface parameter (u,v) of the intersection point
  16.         where v is the height parameter
  17.         and   u is along the circle in counterclock direction.
  18. 2.
  19.         Ray  Canonical cone intersection.
  20.         Canonical cone has :
  21.                 apex at origin with axis aligned with Z.
  22.                 base diameter and height 1.
  23.                 So the slant angle is 45 degrees.
  24.         In this Cononical form the surface equation :
  25.                 x^2 + y^2 - z^2 = 0
  26.                 i.e. the radius of the crosssection circle  is equal
  27.                      to its z.
  28.         "ray_cone" computes the intersection and returns parameter t of
  29.         intersection
  30.                 where t > 0 for valid intersection,
  31.         and surface parameter (u,v) of the intersection point
  32.                 where v is the height of the canonical cone withing the
  33.                 limit and u is along the crosssectional cirlce in
  34.         counterclockwise direction.
  35. 3.
  36.     Ray Cononical Ring intersection.
  37.     Canonical ring is defined by the ring with
  38.     center at origin and Z-axis as normal.
  39.     Its further defined by its inner radius(R0), outer radius(R1).
  40.     0 <= R0 < R1.
  41. ------
  42. sumant
  43. */
  44. #include <stdio.h>
  45. #include <math.h>
  46. #include "GraphicsGems.h"
  47. #include "data_structure.h"
  48.  
  49. static int check(t,ray_start,ray_direction,lo,hi)
  50. double t;
  51. Point3 *ray_start;
  52. Vector3 *ray_direction;
  53. double lo,hi;
  54. {
  55.     double zposn=ray_start->z + ray_direction->z*t;
  56.     return(zposn >= lo && zposn <= hi);
  57. }
  58.  
  59. double ray_cylinder(ray_start,ray_direction,R,H,tubeflag,pu,pv)
  60. /*
  61.     Tubeflag = 0 for the cylinder whose out surface visible.
  62.            1 for the cylinder whose inside surface visible.
  63.     The intersection algorithm is same except in the surface
  64.         parameterisation.
  65. */
  66. Point3 *ray_start;                                      /* Input */
  67. Vector3 *ray_direction;                                 /* Input */
  68. double R,H;                        /* Input */
  69. int tubeflag;                        /* Input */
  70. double *pu,*pv;                        /* Output */
  71. {
  72.     int quadraticRoots();
  73.     int nroots;
  74.     double t = -1.,roots[2];
  75.     double A,B,C;
  76.     C  = ray_start->x * ray_start->x + ray_start->y * ray_start->y
  77.         - R * R;
  78.     B = 2.*(ray_start->x * ray_direction->x
  79.         + ray_start->y * ray_direction->y);
  80.     A = ray_direction->x * ray_direction->x 
  81.         + ray_direction->y * ray_direction->y;
  82.     if (nroots = quadraticRoots(A,B,C,roots))
  83.         set_t(0.,H);
  84.     if (t > EPSILON){ /* Compute Surface Parameter */
  85.         Point3 p;
  86.         double u;
  87.         point_on_line(*ray_start,*ray_direction,t,p);
  88.         *pv = p.z/H;
  89.         circle_inverse_mapping(R,p.x,p.y,u);
  90. /*
  91.         *pu = (tubeflag)?(1.-u):u;
  92. */
  93.         *pu = u;
  94.         return(t);
  95.     }
  96.     return(-1.);
  97. }
  98.  
  99. double ray_cone(ray_start,ray_direction,distance,cupflag,pu,pv)
  100. Point3 *ray_start;                                      /* Input */
  101. Vector3 *ray_direction;                                 /* Input */
  102. double distance;                                        /* Input */
  103. int cupflag;                                            /* Input */
  104. double *pu,*pv;                                         /* Output */
  105. {
  106.     int quadraticRoots();
  107.     int nroots;
  108.     double t = -1.,roots[2];
  109.     double A,B,C;
  110.     C = ray_start->x * ray_start->x + ray_start->y * ray_start->y
  111.         - ray_start->z * ray_start->z;
  112.     B = 2.*(ray_start->x * ray_direction->x
  113.         + ray_start->y * ray_direction->y
  114.         - ray_start->z * ray_direction->z);
  115.     A = ray_direction->x * ray_direction->x
  116.         + ray_direction->y * ray_direction->y
  117.         - ray_direction->z * ray_direction->z;
  118.     if (nroots = quadraticRoots(A,B,C,roots))
  119.         set_t(distance,1.);
  120.         if (t > EPSILON){ /* Compute Surface Parameter */
  121.                 Point3 p;
  122.                 double u;
  123.                 point_on_line(*ray_start,*ray_direction,t,p);
  124.                 *pv = (p.z-distance)/(1.-distance);
  125.                 circle_inverse_mapping(p.z,p.x,p.y,u);
  126. /*
  127.                 *pu = (cupflag)?(1.-u):u;
  128. */
  129.         *pu = u;
  130.         return(t);
  131.         }
  132.     return(-1.);
  133. }
  134.  
  135. double ray_ring(ray_start,ray_direction,R0,R1,pu,pv)
  136. Point3 *ray_start;                                      /* Input */
  137. Vector3 *ray_direction;                                 /* Input */
  138. double R0,R1;                        /* Input */
  139. double *pu,*pv;                        /* Output */
  140. {
  141.     double t = -1.;
  142.     if (fabs(ray_direction->z) > EPSILON){ /* Possible intersection */
  143.         Point3 p;
  144.         double radius;
  145.         double value = -ray_start->z/ray_direction->z;
  146.         if (value > EPSILON){
  147.             point_on_line(*ray_start,*ray_direction,value,p);
  148.             radius = V3SquaredLength(&p);
  149.             if (radius >= R0*R0 && radius <= R1*R1){
  150.                 t = value;
  151.                 radius = sqrt(radius);
  152.                 *pv = (radius - R0)/(R1-R0);
  153.                 circle_inverse_mapping(radius,p.x,p.y,*pu);
  154.                 return(t);
  155.             }
  156.         }
  157.     } 
  158.     return(-1.);
  159. }
  160.