home *** CD-ROM | disk | FTP | other *** search
/ rtsi.com / 2014.01.www.rtsi.com.tar / www.rtsi.com / OS9 / OSK / GRAPHICS / rayshade.lzh / triangle.c < prev    next >
Text File  |  1990-05-08  |  7KB  |  271 lines
  1. /*
  2.  * triangle.c
  3.  *
  4.  * Copyright (C) 1989, Craig E. Kolb
  5.  *
  6.  * This software may be freely copied, modified, and redistributed,
  7.  * provided that this copyright notice is preserved on all copies.
  8.  *
  9.  * There is no warranty or other guarantee of fitness for this software,
  10.  * it is provided solely .  Bug reports or fixes may be sent
  11.  * to the author, who may or may not act on them as he desires.
  12.  *
  13.  * You may not include this software in a program or other software product
  14.  * without supplying the source, or without informing the end-user that the
  15.  * source is available for no extra charge.
  16.  *
  17.  * If you modify this software, you should include a notice giving the
  18.  * name of the person performing the modification, the date of modification,
  19.  * and the reason for such modification.
  20.  *
  21.  * $Id: triangle.c,v 3.0.1.2 90/04/04 18:55:45 craig Exp $
  22.  *
  23.  * $Log:    triangle.c,v $
  24.  * Revision 3.0.1.2  90/04/04  18:55:45  craig
  25.  * patch5: Test for edge-on intersection is more robust.
  26.  * patch5: Lint removal.
  27.  * 
  28.  * Revision 3.0.1.1  89/12/06  16:33:32  craig
  29.  * patch2: Added calls to new error/warning routines.
  30.  * 
  31.  * Revision 3.0  89/10/27  02:06:07  craig
  32.  * Baseline for first official release.
  33.  * 
  34.  */
  35. #include <stdio.h>
  36. #include <math.h>
  37. #include "constants.h"
  38. #include "typedefs.h"
  39. #include "funcdefs.h"
  40.  
  41. #define within(x, a)        (((a) <= 0 && (x) >= (a) && (x) <= 0) || \
  42.                     ((a) > 0 && (x) >= 0 && (x) <= (a)))
  43. /*
  44.  * Create and return reference to a triangle.
  45.  */
  46. Object *
  47. maktri(type, surf, p1, p2, p3, n1, n2, n3)
  48. int    type;
  49. char    *surf;
  50. Vector    *p1, *p2, *p3, *n1, *n2, *n3;
  51. {
  52.     Triangle *triangle;
  53.     Primitive *prim;
  54.     Vector vc1, vc2, vc3, ptmp, anorm;
  55.     Object *newobj;
  56.     double indexval;
  57.  
  58.     prim = mallocprim();
  59.     triangle = (Triangle *)Malloc(sizeof(Triangle));
  60.     prim->objpnt.p_triangle = triangle;
  61.     newobj = new_object(NULL, (char)type, (char *)prim, (Trans *)NULL);
  62.     prim->surf = find_surface(surf);
  63.  
  64.     if (type == PHONGTRI) {
  65.         prim->type = PHONGTRI;
  66.         (void)normalize(n1);
  67.         (void)normalize(n2);
  68.         (void)normalize(n3);
  69.         triangle->vnorm = (Vector *)Malloc(3 * sizeof(Vector));
  70.         triangle->vnorm[0] = *n1;
  71.         triangle->vnorm[1] = *n2;
  72.         triangle->vnorm[2] = *n3;
  73.         triangle->b = (double *)Malloc(3 * sizeof(double));
  74.     }
  75.     else
  76.         prim->type = TRIANGLE;
  77.  
  78.     vecsub(*p2, *p1, &vc1);
  79.     vecsub(*p3, *p2, &vc2);
  80.     vecsub(*p1, *p3, &vc3);
  81.  
  82.     /* Find plane normal. */
  83.     rawcrossp(&triangle->nrm, &vc1, &vc2);
  84.     ptmp = triangle->nrm;
  85.     if (normalize(&ptmp) == 0.) {
  86.         yywarning("Degenerate triangle.\n");
  87.         free((char *)prim);
  88.         free((char *)triangle);
  89.         free((char *)newobj);
  90.         return (Object *)0;
  91.     }
  92.  
  93.     triangle->d = dotp(&triangle->nrm, p1);
  94.  
  95.     triangle->p1 = *p1;
  96.     triangle->p2 = *p2;
  97.     triangle->p3 = *p3;
  98.  
  99.     triangle->e1 = vc1;
  100.     triangle->e2 = vc2;
  101.     triangle->e3 = vc3;
  102.  
  103.     if (type == PHONGTRI && dotp(&triangle->vnorm[0], &ptmp) < 0.) {
  104.         /*
  105.          * Reverse direction of surface normal on Phong
  106.          * triangle if the surface normal points "away"
  107.          * from the first vertex normal.
  108.          */
  109.         scalar_prod(-1., triangle->nrm, &triangle->nrm);
  110.         triangle->d = -triangle->d;
  111.         scalar_prod(-1., triangle->e1, &triangle->e1);
  112.         scalar_prod(-1., triangle->e2, &triangle->e2);
  113.         scalar_prod(-1., triangle->e3, &triangle->e3);
  114.     }
  115.     /*
  116.      * Find "dominant" part of normal vector.
  117.      */
  118.     anorm.x = abs(triangle->nrm.x);
  119.     anorm.y = abs(triangle->nrm.y);
  120.     anorm.z = abs(triangle->nrm.z);
  121.     indexval = max(anorm.y, anorm.z);
  122.     indexval = max(anorm.x, indexval);
  123.     if (indexval == anorm.x)
  124.         triangle->index = XNORMAL;
  125.     else if (indexval == anorm.y)
  126.         triangle->index = YNORMAL;
  127.     else
  128.         triangle->index = ZNORMAL;
  129.  
  130.     return newobj;
  131. }
  132.  
  133. /*
  134.  * Intersect ray with triangle.  See Snyder & Barr for details on
  135.  * how this works.
  136.  */
  137. double
  138. inttri(pos, ray, obj)
  139. register Vector           *pos, *ray;
  140. Primitive       *obj;
  141. {
  142.     register Triangle     *triangle;
  143.     double qi1, qi2, s, k, b1, b2, b3;
  144.     extern unsigned long primtests[];
  145.  
  146.     primtests[obj->type]++;
  147.  
  148.     triangle = obj->objpnt.p_triangle;
  149.     /*
  150.      * Plane intersection.
  151.      */
  152.     k = dotp(&triangle->nrm, ray);
  153.     if (equal(k, 0.))
  154.         return 0.;
  155.     s = (triangle->d - dotp(&triangle->nrm, pos)) / k;
  156.     if (s <= 0.)
  157.         return 0.;
  158.  
  159.     if (triangle->index == XNORMAL) {
  160.         qi1 = pos->y + s * ray->y;
  161.         qi2 = pos->z + s * ray->z;
  162.         b1 = triangle->e2.y * (qi2 - triangle->p2.z) -
  163.                 triangle->e2.z * (qi1 - triangle->p2.y);
  164.         if (!within(b1, triangle->nrm.x))
  165.             return 0.;
  166.         b2 = triangle->e3.y * (qi2 - triangle->p3.z) -
  167.                 triangle->e3.z * (qi1 - triangle->p3.y);
  168.         if (!within(b2, triangle->nrm.x))
  169.             return 0.;
  170.         b3 = triangle->e1.y * (qi2 - triangle->p1.z) -
  171.                 triangle->e1.z * (qi1 - triangle->p1.y);
  172.         if (!within(b3, triangle->nrm.x))
  173.             return 0.;
  174.     } else if (triangle->index == YNORMAL) {
  175.         qi1 = pos->x + s * ray->x;
  176.         qi2 = pos->z + s * ray->z;
  177.         b1 = triangle->e2.z * (qi1 - triangle->p2.x) -
  178.             triangle->e2.x * (qi2 - triangle->p2.z);
  179.         if (!within(b1, triangle->nrm.y))
  180.             return 0.;
  181.         b2 = triangle->e3.z * (qi1 - triangle->p3.x) -
  182.             triangle->e3.x * (qi2 - triangle->p3.z);
  183.         if (!within(b2, triangle->nrm.y))
  184.             return 0.;
  185.         b3 = triangle->e1.z * (qi1 - triangle->p1.x) -
  186.             triangle->e1.x * (qi2 - triangle->p1.z);
  187.         if (!within(b3, triangle->nrm.y))
  188.             return 0.;
  189.     } else {
  190.         qi1 = pos->x + s * ray->x;
  191.         qi2 = pos->y + s * ray->y;
  192.         b1 = triangle->e2.x * (qi2 - triangle->p2.y) -
  193.                 triangle->e2.y * (qi1 - triangle->p2.x);
  194.         if (!within(b1, triangle->nrm.z))
  195.             return 0.;
  196.  
  197.         b2 = triangle->e3.x * (qi2 - triangle->p3.y) -
  198.                 triangle->e3.y * (qi1 - triangle->p3.x);
  199.         if (!within(b2, triangle->nrm.z))
  200.             return 0.;
  201.  
  202.         b3 = triangle->e1.x * (qi2 - triangle->p1.y) -
  203.                 triangle->e1.y * (qi1 - triangle->p1.x);
  204.         if (!within(b3, triangle->nrm.z))
  205.             return 0.;
  206.     }
  207.     /*
  208.      * Take abs value if there was an intersection.
  209.      */
  210.     if (obj->type == PHONGTRI) {
  211.         triangle->b[0] = abs(b1);
  212.         triangle->b[1] = abs(b2);
  213.         triangle->b[2] = abs(b3);
  214.     }
  215.     return s;
  216. }
  217.  
  218. /*ARGSUSED*/
  219. nrmtri(pos, obj, nrm)
  220. Vector           *pos, *nrm;
  221. Primitive       *obj;
  222. {
  223.     Triangle *tri;
  224.  
  225.     /*
  226.      * Normals will be normalized later...
  227.      */
  228.     if (obj->type == TRIANGLE) {
  229.         *nrm = obj->objpnt.p_triangle->nrm;
  230.     } else {
  231.         /*
  232.          * Interpolate normals of Phong-shaded triangles.
  233.          */
  234.         tri = obj->objpnt.p_triangle;
  235.         nrm->x = tri->b[0]*tri->vnorm[0].x+tri->b[1]*tri->vnorm[1].x+
  236.             tri->b[2]*tri->vnorm[2].x;
  237.         nrm->y = tri->b[0]*tri->vnorm[0].y+tri->b[1]*tri->vnorm[1].y+
  238.             tri->b[2]*tri->vnorm[2].y;
  239.         nrm->z = tri->b[0]*tri->vnorm[0].z+tri->b[1]*tri->vnorm[1].z+
  240.             tri->b[2]*tri->vnorm[2].z;
  241.     }
  242. }
  243.  
  244. triextent(o, bounds)
  245. Primitive *o;
  246. double bounds[2][3];
  247. {
  248.     Triangle *tri;
  249.  
  250.     tri = o->objpnt.p_triangle;
  251.  
  252.     bounds[LOW][X] = bounds[HIGH][X] = tri->p1.x;
  253.     bounds[LOW][Y] = bounds[HIGH][Y] = tri->p1.y;
  254.     bounds[LOW][Z] = bounds[HIGH][Z] = tri->p1.z;
  255.  
  256.     if (tri->p2.x < bounds[LOW][X]) bounds[LOW][X] = tri->p2.x;
  257.     if (tri->p2.x > bounds[HIGH][X]) bounds[HIGH][X] = tri->p2.x;
  258.     if (tri->p3.x < bounds[LOW][X]) bounds[LOW][X] = tri->p3.x;
  259.     if (tri->p3.x > bounds[HIGH][X]) bounds[HIGH][X] = tri->p3.x;
  260.  
  261.     if (tri->p2.y < bounds[LOW][Y]) bounds[LOW][Y] = tri->p2.y;
  262.     if (tri->p2.y > bounds[HIGH][Y]) bounds[HIGH][Y] = tri->p2.y;
  263.     if (tri->p3.y < bounds[LOW][Y]) bounds[LOW][Y] = tri->p3.y;
  264.     if (tri->p3.y > bounds[HIGH][Y]) bounds[HIGH][Y] = tri->p3.y;
  265.  
  266.     if (tri->p2.z < bounds[LOW][Z]) bounds[LOW][Z] = tri->p2.z;
  267.     if (tri->p2.z > bounds[HIGH][Z]) bounds[HIGH][Z] = tri->p2.z;
  268.     if (tri->p3.z < bounds[LOW][Z]) bounds[LOW][Z] = tri->p3.z;
  269.     if (tri->p3.z > bounds[HIGH][Z]) bounds[HIGH][Z] = tri->p3.z;
  270. }
  271.