home *** CD-ROM | disk | FTP | other *** search
/ APDL Public Domain 1 / APDL_PD1A.iso / raytrace / dbw_render / source / c / HIT < prev    next >
Encoding:
Text File  |  1992-10-23  |  10.0 KB  |  378 lines
  1. /************************************************************************
  2.  *                                                                      *
  3.  *                    Copyright (c) 1987, David B. Wecker               *
  4.  *                            All Rights Reserved                       *
  5.  *                                                                      *
  6.  * This file is part of DBW_Render                                      *
  7.  *                                                                      *
  8.  * DBW_Render is distributed in the hope that it will be useful, but    *
  9.  * WITHOUT ANY WARRANTY. No author or distributor accepts               *
  10.  * responsibility to anyone for the consequences of using it or for     *
  11.  * whether it serves any particular purpose or works at all, unless     *
  12.  * he says so in writing. Refer to the DBW_Render General Public        *
  13.  * License for full details.                                            *
  14.  *                                                                      *
  15.  * Everyone is granted permission to copy, modify and redistribute      *
  16.  * DBW_Render, but only under the conditions described in the           *
  17.  * DBW_Render General Public License. A copy of this license is         *
  18.  * supposed to have been given to you along with DBW_Render so you      *
  19.  * can know your rights and responsibilities. It should be in a file    *
  20.  * named COPYING. Among other things, the copyright notice and this     *
  21.  * notice must be preserved on all copies.                              *
  22.  ************************************************************************
  23.  
  24. #define MODULE_HIT
  25. #include "ray.h"
  26.  
  27. void findnormal(np,p,n)
  28. node    *np;
  29. vector  p,
  30. n;
  31. {
  32.      vector ripple,fuzzy,bevelx,bevely,bevelz,tmp;
  33.      float t;
  34.      int w,j,it;
  35.  
  36.      switch (np->kind) 
  37.      {
  38.      case SPHERE   :
  39.           SPHERENORMAL(sptr(np)->center,p,n);
  40.           break;
  41.      case TRIANGLE :
  42.           PLANENORMAL(tptr(np)->ve,tptr(np)->vp,n);
  43.           break;
  44.      case QUAD     :
  45.           PLANENORMAL(qptr(np)->ve,qptr(np)->vp,n);
  46.           break;
  47.      case RING     :
  48.           PLANENORMAL(rptr(np)->ve,rptr(np)->vp,n);
  49.           break;
  50.      }
  51.  
  52.      if (np->attr.tex != 0) 
  53.      { /* don't bother checking if there's no texture */
  54.           if (np->attr.tex == 4) 
  55.           { /* fiddle with normal for all ripples */
  56.                for (w = 0; w < numwaves; w++) 
  57.                {
  58.                     calcripple(p,w,ripple);  /* calculate the wave perturbation */
  59.                     vecsum(ripple,n,n);  /* add ripple bend to normal */
  60.                }
  61.                normalize(n);  /* make sure it's a unit vector after all that */
  62.           }
  63.  
  64.           /* fiddle with normal for one ripple */
  65.           else if (np->attr.tex >= 10 && np->attr.tex <= 19) 
  66.           {
  67.                calcripple(p,np->attr.tex - 10,ripple);  /* calc wave perturbation */
  68.                vecsum(ripple,n,n);  /* add ripple bend to normal */
  69.                normalize(n);  /* make sure it's still a unit vector */
  70.           }
  71.  
  72.           /* check for rough pebbly surface */
  73.           else if (np->attr.tex >= 90 && np->attr.tex <= 99) 
  74.           {
  75.                j = np->attr.tex - 90;  /* select the desired pebble finish */
  76.                vecscale(pebble[j].zoom,p,tmp);
  77.                fuzzy[0] = turbulence(tmp);
  78.                vecscale(1.5,tmp,tmp);
  79.                fuzzy[1] = turbulence(tmp);
  80.                vecscale(1.5,tmp,tmp);
  81.                fuzzy[2] = turbulence(tmp);
  82.                if (fuzzy[0] > 0.7)
  83.                     fuzzy[0] = 0.7;
  84.                if (fuzzy[1] > 0.7)
  85.                     fuzzy[1] = 0.7;
  86.                if (fuzzy[2] > 0.7)
  87.                     fuzzy[2] = 0.7;
  88.                if (rnd() < 0.5)
  89.                     fuzzy[0] = -fuzzy[0];
  90.                if (rnd() < 0.5)
  91.                     fuzzy[1] = -fuzzy[1];
  92.                if (rnd() < 0.5)
  93.                     fuzzy[2] = -fuzzy[2];
  94.                vecscale(pebble[j].scale,fuzzy,fuzzy);
  95.                vecsum(n,fuzzy,n);
  96.                normalize(n);
  97.           }
  98.      }
  99.  
  100.      /*---------------------------------------------------------------------*/
  101.      /* Add any general normal perturbations                                */
  102.  
  103.      if (np->attr.fuz > 0.0) 
  104.      {
  105.           /* Perturb the normal randomly to produce fuzzy surfaces */
  106.           fuzzy[0] = rnd();  /* 0..1 */
  107.           fuzzy[1] = rnd();
  108.           fuzzy[2] = rnd();
  109.           if (rnd() < 0.5)
  110.                fuzzy[0] = -fuzzy[0];
  111.           if (rnd() < 0.5)
  112.                fuzzy[1] = -fuzzy[1];
  113.           if (rnd() < 0.5)
  114.                fuzzy[2] = -fuzzy[2];
  115.  
  116.           /* 'fuzzy' is now approximately a random unit vector */
  117.           vecscale(rnd() * np->attr.fuz,fuzzy,fuzzy);
  118.           vecsum(fuzzy,n,n); /* vector addition of fuzz compunent to true normal */
  119.           normalize(n);  /* Make sure it's still a unit vector */
  120.      }
  121. }
  122.  
  123. int hitcylinder(cp,eye,d,p,t)
  124. cylinder *cp;  /* the cylinder      */
  125. vector eye;    /* source ray origin */
  126. vector d;      /* source ray        */
  127. vector p;
  128. float  *t;
  129. {
  130.      float aa,bb,cc,radical,a2,b2,c2,dist1,dist2;
  131.      vector otop,oeye,p1,p2,opoint1,opoint2;
  132.  
  133.      /*
  134.       *   translate bottom of cylinder to 0,0,0
  135.       *   get translated eyepoint and cylinder top
  136.       */
  137.  
  138.      VECSUB(eye,    cp->bottom, oeye);
  139.      VECSUB(cp->top,cp->bottom, otop);  /* ? */
  140.  
  141.  
  142.      a2 = cp->a*cp->a;
  143.      b2 = cp->b*cp->b;
  144.      c2 = cp->c*cp->c;
  145.  
  146.      if (a2 < SMALL) return( FALSE );
  147.      if (b2 < SMALL) return( FALSE );
  148.      if (c2 < SMALL) return( FALSE );
  149.  
  150.  
  151.      aa  = d[0]*d[0] / a2;
  152.      aa += d[2]*d[2] / b2;
  153.      aa -= d[1]*d[1] / c2;
  154.  
  155.      bb  = d[0]*oeye[0] / a2;
  156.      bb += d[2]*oeye[2] / b2;
  157.      bb -= d[1]*oeye[1] / c2;
  158.  
  159.      cc  = oeye[0]*oeye[0] / a2;
  160.      cc += oeye[2]*oeye[2] / b2;
  161.      cc -= oeye[1]*oeye[1] / c2;
  162.  
  163. /*   
  164.  *   descriminate < 0,  ray misses cylinder
  165.  *   descriminate == 0, ray grazes cylinder
  166.  */
  167.  
  168.      if ((radical = (bb*bb) - (4.0*aa*cc)) < 0.0)
  169.           return( FALSE );
  170.  
  171.      radical = sqrt(radical);
  172.  
  173.      aa = 2.0 * aa;
  174.  
  175.      if ( aa < SMALL ) 
  176.           return( FALSE );  
  177.  
  178.      /* the roots of the quadratic */
  179.  
  180.      dist1 = (-bb + radical) / aa;
  181.      dist2 = (-bb - radical) / aa;
  182.  
  183.  
  184.      VECSCALE(dist1,d,p1);
  185.      VECSUM(p1,oeye,p1);
  186.      VECSCALE(dist2,d,p2);
  187.      VECSUM(p2,oeye,p2);
  188.      VECSUB(oeye,p1,opoint1);
  189.      VECSUB(oeye,p2,opoint2);
  190.      dist1 = NORM(opoint1);
  191.      dist2 = NORM(opoint2);
  192.  
  193.      if (dist1 < dist2) 
  194.      { 
  195.           VECCOPY(p1,p); 
  196.           *t = dist1; 
  197.      }
  198.      else
  199.      { 
  200.           VECCOPY(p2,p); 
  201.           *t = dist2; 
  202.      }
  203.  
  204.      return( TRUE );
  205. }
  206.  
  207. int hitsphere(center,radius,eye,d,p,t)
  208. vector  center;     /* where it is  */
  209. vector  eye;        /* where we are */
  210. vector  d;
  211. vector  p;
  212. float   radius;     /* how big it is */
  213. float   *t;
  214. {
  215.      float   r_r,d_r,t2,radical;
  216.      vector  r;
  217.  
  218.      VECSUB(center,eye,r);
  219.      r_r = DOT(r,r);
  220.      d_r = DOT(d,r);
  221.  
  222.      if ((radical = (d_r*d_r) + (radius*radius) - r_r) < 0.0) 
  223.           return( FALSE );
  224.  
  225.      radical = sqrt(radical);
  226.  
  227.      if (d_r < radical) 
  228.      { 
  229.           *t = d_r + radical; 
  230.           t2 = d_r - radical; 
  231.      }
  232.      else
  233.      { 
  234.           *t = d_r - radical; 
  235.           t2 = d_r + radical; 
  236.      }
  237.  
  238.      if (fabs(*t) < SMALL) 
  239.           *t = t2;
  240.  
  241.      if (*t <= 0) 
  242.           return( FALSE );
  243.  
  244.      VECSCALE((*t),d,p);
  245.      VECSUM(p,eye,p);
  246.      return( TRUE );
  247. }
  248.  
  249. int hitplane(p,ve,vp,eye,d,sfs,inter)
  250. vector  p,
  251.   ve,
  252.   vp,
  253.   eye,
  254.   d,
  255.   sfs,
  256.   inter;
  257. {
  258.      vector h;
  259.      float det;
  260.  
  261.      det  = ve[0] * ((vp[1] * d[2]) - (vp[2] * d[1]));
  262.      det -= vp[0] * ((ve[1] * d[2]) - (ve[2] * d[1]));
  263.      det +=  d[0] * ((ve[1] * vp[2])- (ve[2] * vp[1]));
  264.  
  265.      if (det == 0.0)
  266.           return FALSE;
  267.  
  268.      vecsub(eye,p,h);
  269.  
  270.      sfs[2]  = h[0] * ((ve[1] * vp[2]) - (ve[2] * vp[1]));
  271.      sfs[2] -= h[1] * ((ve[0] * vp[2]) - (ve[2] * vp[0]));
  272.      sfs[2] += h[2] * ((ve[0] * vp[1]) - (ve[1] * vp[0]));
  273.      sfs[2] /= det;
  274.      sfs[2]  = -sfs[2];
  275.  
  276.      if (sfs[2] < SMALL)
  277.           return FALSE;
  278.  
  279.      sfs[0]  = h[0] * ((vp[1] * d[2]) - (vp[2] * d[1]));
  280.      sfs[0] -= h[1] * ((vp[0] * d[2]) - (vp[2] * d[0]));
  281.      sfs[0] += h[2] * ((vp[0] * d[1]) - (vp[1] * d[0]));
  282.      sfs[0] /= det;
  283.  
  284.      sfs[1]  = h[0] * ((ve[2] * d[1]) - (ve[1] * d[2]));
  285.      sfs[1] += h[1] * ((ve[0] * d[2]) - (ve[2] * d[0]));
  286.      sfs[1] -= h[2] * ((ve[0] * d[1]) - (ve[1] * d[0]));
  287.      sfs[1] /= det;
  288.  
  289.      VECCOPY(d,inter);
  290.      VECSCALE(sfs[2],inter,inter);
  291.      VECSUM(inter,eye,inter);
  292.  
  293.      return TRUE;
  294. }
  295.  
  296. int hittriangle(tp,eye,d,p,t)
  297. triangle    *tp;
  298. vector      eye,
  299. d,
  300. p;
  301. float       *t;
  302. {
  303.      vector sfs;
  304.      int hit;
  305.  
  306.      hit = hitplane(tp->position,tp->ve,tp->vp,eye,d,sfs,p) &&
  307.        sfs[0] > 0.0  && sfs[1] > 0.0 &&  sfs[0] + sfs[1] <= 1.0;
  308.      *t = sfs[2];
  309.  
  310.      return hit;
  311. }
  312.  
  313. int hitquad(qp,eye,d,p,t)
  314. quad        *qp;
  315. vector      eye,
  316. d,
  317. p;
  318. float       *t;
  319. {
  320.      vector sfs;
  321.      int hit = FALSE;
  322.  
  323.      if (hitplane(qp->position,qp->ve,qp->vp,eye,d,sfs,p))
  324.           hit = sfs[0] > 0.0 && sfs[1] > 0.0 && sfs[0] <= 1.0 && sfs[1] <= 1.0;
  325.      *t = sfs[2];
  326.      return hit;
  327. }
  328.  
  329. int hitring(rp,eye,d,p,t)
  330. ring        *rp;
  331. vector      eye,
  332. d,
  333. p;
  334. float       *t;
  335. {
  336.      vector  sfs;
  337.      float   r;
  338.      int     hit;
  339.  
  340.      hit = hitplane(rp->position,rp->ve,rp->vp,eye,d,sfs,p) &&
  341.        (r = sfs[0] * sfs[0] + sfs[1] * sfs[1]) <= rp->maxrad &&
  342.        r >= rp->minrad;
  343.      *t = sfs[2];
  344.      return( hit );
  345. }
  346.  
  347. void shell(v,v1,v2,n)
  348. float       v[];
  349. vector      v1[];
  350. node        *v2[];
  351. int         n;
  352. {
  353.      int     gap,i,j;
  354.      node    *temp2;
  355.      vector  temp1;
  356.      float   temp;
  357.  
  358.      sorts++;  /* Keep statistics */
  359.      sort_size += n;
  360.  
  361.      for (gap = n/2; gap > 0; gap /=2)
  362.           for (i = gap; i < n; i++)
  363.                for (j = i-gap; j>=0 && v[j] > v[j+gap]; j -= gap) 
  364.                {
  365.                     temp = v[j];
  366.                     v[j] = v[j+gap];
  367.                     v[j+gap] = temp;
  368.  
  369.                     VECCOPY(v1[j],temp1);
  370.                     VECCOPY(v1[j+gap],v1[j]);
  371.                     VECCOPY(temp1,v1[j+gap]);
  372.  
  373.                     temp2 = v2[j];
  374.                     v2[j] = v2[j+gap];
  375.                     v2[j+gap] = temp2;
  376.                }
  377. }
  378.