ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/hd/sm_geom.c
(Generate patch)

Comparing ray/src/hd/sm_geom.c (file contents):
Revision 3.7 by gwlarson, Tue Oct 6 18:16:53 1998 UTC vs.
Revision 3.9 by gwlarson, Mon Dec 28 18:07:35 1998 UTC

# Line 27 | Line 27 | FVECT v1,v2;
27   {
28     return(EQUAL(v1[0],v2[0]) && EQUAL(v1[1],v2[1])&& EQUAL(v1[2],v2[2]));
29   }
30 + #if 0
31 + extern FVECT Norm[500];
32 + extern int Ncnt;
33 + #endif
34  
31 –
35   int
36   convex_angle(v0,v1,v2)
37   FVECT v0,v1,v2;
38   {
39      FVECT cp,cp01,cp12,v10,v02;
40      double dp;
41 <    /*
42 <      VSUB(v10,v1,v0);
40 <      VSUB(v02,v0,v2);
41 <      VCROSS(cp,v10,v02);
42 <   */
43 <      /* test sign of (v0Xv1)X(v1Xv2). v1 */
41 >
42 >    /* test sign of (v0Xv1)X(v1Xv2). v1 */
43      VCROSS(cp01,v0,v1);
44      VCROSS(cp12,v1,v2);
45      VCROSS(cp,cp01,cp12);
46          
47      dp = DOT(cp,v1);
48 + #if 0
49 +    VCOPY(Norm[Ncnt++],cp01);
50 +    VCOPY(Norm[Ncnt++],cp12);
51 +    VCOPY(Norm[Ncnt++],cp);
52 + #endif
53      if(ZERO(dp) || dp < 0.0)
54 <       return(FALSE);
54 >      return(FALSE);
55      return(TRUE);
56   }
57  
58   /* calculates the normal of a face contour using Newell's formula. e
59  
60 <               a =  SUMi (yi - yi+1)(zi + zi+1)
60 >               a =  SUMi (yi - yi+1)(zi + zi+1);
61                 b =  SUMi (zi - zi+1)(xi + xi+1)
62                 c =  SUMi (xi - xi+1)(yi + yi+1)
63   */
# Line 70 | Line 74 | int norm;
74    
75    n[1] = (v0[2] - v1[2]) * (v0[0] + v1[0]) +
76             (v1[2] - v2[2]) * (v1[0] + v2[0]) +
77 <           (v2[2] - v0[2]) * (v2[0] + v0[0]);
77 >            (v2[2] - v0[2]) * (v2[0] + v0[0]);
78    
79    n[2] = (v0[1] + v1[1]) * (v0[0] - v1[0]) +
80           (v1[1] + v2[1]) * (v1[0] - v2[0]) +
# Line 78 | Line 82 | int norm;
82  
83    if(!norm)
84       return(0);
81 –
85    
86    mag = normalize(n);
87  
# Line 116 | Line 119 | FVECT p0,p1;
119   }
120  
121  
122 <
122 > double
123   point_on_sphere(ps,p,c)
124   FVECT ps,p,c;
125   {
126 +  double d;
127      VSUB(ps,p,c);    
128 <    normalize(ps);
128 >    d= normalize(ps);
129 >    return(d);
130   }
131  
132  
# Line 179 | Line 184 | intersect_ray_plane(orig,dir,peq,pd,r)
184     double *pd;
185     FVECT r;
186   {
187 <  double t;
187 >  double t,d;
188    int hit;
189      /*
190        Plane is Ax + By + Cz +D = 0:
# Line 190 | Line 195 | intersect_ray_plane(orig,dir,peq,pd,r)
195         line is  l = p1 + (p2-p1)t
196       */
197      /* Solve for t: */
198 <    t =  -(DOT(FP_N(peq),orig) + FP_D(peq))/(DOT(FP_N(peq),dir));
199 <    if(t < 0)
198 >  d = DOT(FP_N(peq),dir);
199 >  if(ZERO(d))
200 >     return(0);
201 >  t =  -(DOT(FP_N(peq),orig) + FP_D(peq))/d;
202 >
203 >  if(t < 0)
204         hit = 0;
205      else
206         hit = 1;
# Line 212 | Line 221 | intersect_ray_oplane(orig,dir,n,pd,r)
221     double *pd;
222     FVECT r;
223   {
224 <  double t;
224 >  double t,d;
225    int hit;
226      /*
227        Plane is Ax + By + Cz +D = 0:
# Line 223 | Line 232 | intersect_ray_oplane(orig,dir,n,pd,r)
232         line is  l = p1 + (p2-p1)t
233       */
234      /* Solve for t: */
235 <    t =  -(DOT(n,orig))/(DOT(n,dir));
235 >    d= DOT(n,dir);
236 >    if(ZERO(d))
237 >       return(0);
238 >    t =  -(DOT(n,orig))/d;
239      if(t < 0)
240         hit = 0;
241      else
# Line 237 | Line 249 | intersect_ray_oplane(orig,dir,n,pd,r)
249    return(hit);
250   }
251  
252 + /* Assumption: know crosses plane:dont need to check for 'on' case */
253 + intersect_edge_coord_plane(v0,v1,w,r)
254 + FVECT v0,v1;
255 + int w;
256 + FVECT r;
257 + {
258 +  FVECT dv;
259 +  int wnext;
260 +  double t;
261  
262 +  VSUB(dv,v1,v0);
263 +  t = -v0[w]/dv[w];
264 +  r[w] = 0.0;
265 +  wnext = (w+1)%3;
266 +  r[wnext] = v0[wnext] + dv[wnext]*t;
267 +  wnext = (w+2)%3;
268 +  r[wnext] = v0[wnext] + dv[wnext]*t;
269 + }
270 +
271   int
272   intersect_edge_plane(e0,e1,peq,pd,r)
273     FVECT e0,e1;
# Line 286 | Line 316 | FVECT p0,p1,p2;
316         that plane and do an in-circle test in the plane
317       */
318      
319 <    /* find the equation of the plane defined by p1-p3 */
319 >    /* find the equation of the plane defined by p0-p2 */
320      tri_plane_equation(p0,p1,p2,&peq,FALSE);
321  
322      /* define a coordinate system on the plane: the x axis is in
# Line 298 | Line 328 | FVECT p0,p1,p2;
328      if(!intersect_vector_plane(p,peq,NULL,np))
329          return(FALSE);
330  
331 <    /* create coordinate system on  plane: p2-p1 defines the x_axis*/
331 >    /* create coordinate system on  plane: p1-p0 defines the x_axis*/
332      VSUB(x_axis,p1,p0);
333      normalize(x_axis);
334      /* The y axis is  */
# Line 309 | Line 339 | FVECT p0,p1,p2;
339      VSUB(p2,p2,p0);
340      VSUB(np,np,p0);
341      
342 <    p1[0] = VLEN(p1);
342 >    p1[0] = DOT(p1,x_axis);
343      p1[1] = 0;
344  
345      d1 = DOT(p2,x_axis);
# Line 338 | Line 368 | int sides[3];
368      /* Find the normal to the triangle ORIGIN:v0:v1 */
369      if(!NTH_BIT(*nset,0))
370      {
371 <        VCROSS(n[0],v1,v0);
371 >        VCROSS(n[0],v0,v1);
372          SET_NTH_BIT(*nset,0);
373      }
374      /* Test the point for sidedness */
# Line 356 | Line 386 | int sides[3];
386      /* Test next edge */
387      if(!NTH_BIT(*nset,1))
388      {
389 <        VCROSS(n[1],v2,v1);
389 >        VCROSS(n[1],v1,v2);
390          SET_NTH_BIT(*nset,1);
391      }
392      /* Test the point for sidedness */
# Line 372 | Line 402 | int sides[3];
402      /* Test next edge */
403      if(!NTH_BIT(*nset,2))
404      {
405 <        VCROSS(n[2],v0,v2);
405 >        VCROSS(n[2],v2,v0);
406          SET_NTH_BIT(*nset,2);
407      }
408      /* Test the point for sidedness */
# Line 398 | Line 428 | FVECT v0,v1,v2,p;
428      double d;
429      FVECT n;  
430  
431 <    VCROSS(n,v1,v0);
431 >    VCROSS(n,v0,v1);
432      /* Test the point for sidedness */
433      d  = DOT(n,p);
434      if(d > 0.0)
435        return(FALSE);
436  
437      /* Test next edge */
438 <    VCROSS(n,v2,v1);
438 >    VCROSS(n,v1,v2);
439      /* Test the point for sidedness */
440      d  = DOT(n,p);
441      if(d > 0.0)
442         return(FALSE);
443  
444      /* Test next edge */
445 <    VCROSS(n,v0,v2);
445 >    VCROSS(n,v2,v0);
446      /* Test the point for sidedness */
447      d  = DOT(n,p);
448      if(d > 0.0)
# Line 506 | Line 536 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
536      if(sides[0][0] == GT_INVALID)
537      {
538        if(!NTH_BIT(nset,0))
539 <        VCROSS(n[0],t1,t0);
539 >        VCROSS(n[0],t0,t1);
540        /* Test the point for sidedness */
541        d  = DOT(n[0],p0);
542        if(d >= 0.0)
# Line 519 | Line 549 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
549      if(sides[0][1] == GT_INVALID)
550      {
551        if(!NTH_BIT(nset,1))
552 <        VCROSS(n[1],t2,t1);
552 >        VCROSS(n[1],t1,t2);
553          /* Test the point for sidedness */
554          d  = DOT(n[1],p0);
555          if(d >= 0.0)
# Line 532 | Line 562 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
562      if(sides[0][2] == GT_INVALID)
563      {
564        if(!NTH_BIT(nset,2))
565 <        VCROSS(n[2],t0,t2);
565 >        VCROSS(n[2],t2,t0);
566        /* Test the point for sidedness */
567        d  = DOT(n[2],p0);
568        if(d >= 0.0)
# Line 548 | Line 578 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
578      if(sides[1][0] == GT_INVALID)
579      {
580        if(!NTH_BIT(nset,0))
581 <        VCROSS(n[0],t1,t0);
581 >        VCROSS(n[0],t0,t1);
582        /* Test the point for sidedness */
583        d  = DOT(n[0],p1);
584        if(d >= 0.0)
# Line 562 | Line 592 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
592      if(sides[1][1] == GT_INVALID)
593      {
594        if(!NTH_BIT(nset,1))
595 <        VCROSS(n[1],t2,t1);
595 >        VCROSS(n[1],t1,t2);
596        /* Test the point for sidedness */
597        d  = DOT(n[1],p1);
598        if(d >= 0.0)
# Line 576 | Line 606 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
606      if(sides[1][2] == GT_INVALID)
607      {
608        if(!NTH_BIT(nset,2))
609 <        VCROSS(n[2],t0,t2);
609 >        VCROSS(n[2],t2,t0);
610        /* Test the point for sidedness */
611        d  = DOT(n[2],p1);
612        if(d >= 0.0)
# Line 592 | Line 622 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
622      if(sides[2][0] == GT_INVALID)
623      {
624        if(!NTH_BIT(nset,0))
625 <        VCROSS(n[0],t1,t0);
625 >        VCROSS(n[0],t0,t1);
626        /* Test the point for sidedness */
627        d  = DOT(n[0],p2);
628        if(d >= 0.0)
# Line 605 | Line 635 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
635      if(sides[2][1] == GT_INVALID)
636      {
637        if(!NTH_BIT(nset,1))
638 <        VCROSS(n[1],t2,t1);
638 >        VCROSS(n[1],t1,t2);
639        /* Test the point for sidedness */
640        d  = DOT(n[1],p2);
641        if(d >= 0.0)
# Line 618 | Line 648 | set_sidedness_tests(t0,t1,t2,p0,p1,p2,test,sides,nset,
648      if(sides[2][2] == GT_INVALID)
649      {
650        if(!NTH_BIT(nset,2))
651 <        VCROSS(n[2],t0,t2);
651 >        VCROSS(n[2],t2,t0);
652        /* Test the point for sidedness */
653        d  = DOT(n[2],p2);
654        if(d >= 0.0)
# Line 637 | Line 667 | FVECT a1,a2,a3,b1,b2,b3;
667    int which,test,n_set[2];
668    int sides[2][3][3],i,j,inext,jnext;
669    int tests[2][3];
670 <  FVECT n[2][3],p,avg[2];
670 >  FVECT n[2][3],p,avg[2],t1,t2,t3;
671  
672 + #ifdef DEBUG
673 +    tri_normal(b1,b2,b3,p,FALSE);
674 +    if(DOT(p,b1) > 0)
675 +      {
676 +       VCOPY(t3,b1);
677 +       VCOPY(t2,b2);
678 +       VCOPY(t1,b3);
679 +      }
680 +    else
681 + #endif
682 +      {
683 +       VCOPY(t1,b1);
684 +       VCOPY(t2,b2);
685 +       VCOPY(t3,b3);
686 +      }
687 +
688    /* Test the vertices of triangle a against the pyramid formed by triangle
689       b and the origin. If any vertex of a is interior to triangle b, or
690       if all 3 vertices of a are ON the edges of b,return TRUE. Remember
691       the results of the edge normal and sidedness tests for later.
692     */
693 < if(vertices_in_stri(a1,a2,a3,b1,b2,b3,&(n_set[0]),n[0],avg[0],sides[1]))
693 > if(vertices_in_stri(a1,a2,a3,t1,t2,t3,&(n_set[0]),n[0],avg[0],sides[1]))
694       return(TRUE);
695    
696 < if(vertices_in_stri(b1,b2,b3,a1,a2,a3,&(n_set[1]),n[1],avg[1],sides[0]))
696 > if(vertices_in_stri(t1,t2,t3,a1,a2,a3,&(n_set[1]),n[1],avg[1],sides[0]))
697       return(TRUE);
698  
699  
700 <  set_sidedness_tests(b1,b2,b3,a1,a2,a3,tests[0],sides[0],n_set[1],n[1]);
700 >  set_sidedness_tests(t1,t2,t3,a1,a2,a3,tests[0],sides[0],n_set[1],n[1]);
701    if(tests[0][0]&tests[0][1]&tests[0][2])
702      return(FALSE);
703  
704 <  set_sidedness_tests(a1,a2,a3,b1,b2,b3,tests[1],sides[1],n_set[0],n[0]);
704 >  set_sidedness_tests(a1,a2,a3,t1,t2,t3,tests[1],sides[1],n_set[0],n[0]);
705    if(tests[1][0]&tests[1][1]&tests[1][2])
706      return(FALSE);
707  
# Line 861 | Line 907 | FVECT a0,a1,b0,b1;
907   }
908  
909  
864 –
910   /* Find the normalized barycentric coordinates of p relative to
911   * triangle v0,v1,v2. Return result in coord
912   */
# Line 879 | Line 924 | double coord[3];
924  
925   }
926  
882 – bary_ith_child(coord,i)
883 – double coord[3];
884 – int i;
885 – {
886 –
887 –  switch(i){
888 –  case 0:
889 –      /* update bary for child */
890 –      coord[0] = 2.0*coord[0]- 1.0;
891 –      coord[1] *= 2.0;
892 –      coord[2] *= 2.0;
893 –      break;
894 –  case 1:
895 –    coord[0] *= 2.0;
896 –    coord[1] = 2.0*coord[1]- 1.0;
897 –    coord[2] *= 2.0;
898 –    break;
899 –  case 2:
900 –    coord[0] *= 2.0;
901 –    coord[1] *= 2.0;
902 –    coord[2] = 2.0*coord[2]- 1.0;
903 –    break;
904 –  case 3:
905 –    coord[0] = 1.0 - 2.0*coord[0];
906 –    coord[1] = 1.0 - 2.0*coord[1];
907 –    coord[2] = 1.0 - 2.0*coord[2];
908 –    break;
909 – #ifdef DEBUG
910 –  default:
911 –    eputs("bary_ith_child():Invalid child\n");
912 –    break;
913 – #endif
914 –  }
915 – }
927  
928  
918 – int
919 – bary_child(coord)
920 – double coord[3];
921 – {
922 –  int i;
929  
924 –  if(coord[0] > 0.5)
925 –  {
926 –      /* update bary for child */
927 –      coord[0] = 2.0*coord[0]- 1.0;
928 –      coord[1] *= 2.0;
929 –      coord[2] *= 2.0;
930 –      return(0);
931 –  }
932 –  else
933 –    if(coord[1] > 0.5)
934 –    {
935 –      coord[0] *= 2.0;
936 –      coord[1] = 2.0*coord[1]- 1.0;
937 –      coord[2] *= 2.0;
938 –      return(1);
939 –    }
940 –    else
941 –      if(coord[2] > 0.5)
942 –      {
943 –        coord[0] *= 2.0;
944 –        coord[1] *= 2.0;
945 –        coord[2] = 2.0*coord[2]- 1.0;
946 –        return(2);
947 –      }
948 –      else
949 –         {
950 –           coord[0] = 1.0 - 2.0*coord[0];
951 –           coord[1] = 1.0 - 2.0*coord[1];
952 –           coord[2] = 1.0 - 2.0*coord[2];
953 –           return(3);
954 –         }
955 – }
956 –
957 – /* Coord was the ith child of the parent: set the coordinate
958 –   relative to the parent
959 – */
930   bary_parent(coord,i)
961 – double coord[3];
962 – int i;
963 – {
964 –
965 –  switch(i) {
966 –  case 0:
967 –    /* update bary for child */
968 –    coord[0] = coord[0]*0.5 + 0.5;
969 –    coord[1] *= 0.5;
970 –    coord[2] *= 0.5;
971 –    break;
972 –  case 1:
973 –    coord[0] *= 0.5;
974 –    coord[1]  = coord[1]*0.5 + 0.5;
975 –    coord[2] *= 0.5;
976 –    break;
977 –    
978 –  case 2:
979 –    coord[0] *= 0.5;
980 –    coord[1] *= 0.5;
981 –    coord[2] = coord[2]*0.5 + 0.5;
982 –    break;
983 –    
984 –  case 3:
985 –    coord[0] = 0.5 - 0.5*coord[0];
986 –    coord[1] = 0.5 - 0.5*coord[1];
987 –    coord[2] = 0.5 - 0.5*coord[2];
988 –    break;
989 – #ifdef DEBUG
990 –  default:
991 –    eputs("bary_parent():Invalid child\n");
992 –    break;
993 – #endif
994 –  }
995 – }
996 –
997 – bary_from_child(coord,child,next)
998 – double coord[3];
999 – int child,next;
1000 – {
1001 – #ifdef DEBUG
1002 –  if(child <0 || child > 3)
1003 –  {
1004 –    eputs("bary_from_child():Invalid child\n");
1005 –    return;
1006 –  }
1007 –  if(next <0 || next > 3)
1008 –  {
1009 –    eputs("bary_from_child():Invalid next\n");
1010 –    return;
1011 –  }
1012 – #endif
1013 –  if(next == child)
1014 –    return;
1015 –
1016 –  switch(child){
1017 –  case 0:
1018 –    switch(next){
1019 –    case 1:
1020 –      coord[0] += 1.0;
1021 –      coord[1] -= 1.0;
1022 –      break;
1023 –    case 2:
1024 –      coord[0] += 1.0;
1025 –      coord[2] -= 1.0;
1026 –      break;
1027 –    case 3:
1028 –      coord[0] *= -1.0;
1029 –      coord[1] = 1 - coord[1];
1030 –      coord[2] = 1 - coord[2];
1031 –      break;
1032 –
1033 –    }
1034 –    break;
1035 –  case 1:
1036 –    switch(next){
1037 –    case 0:
1038 –      coord[0] -= 1.0;
1039 –      coord[1] += 1.0;
1040 –      break;
1041 –    case 2:
1042 –      coord[1] += 1.0;
1043 –      coord[2] -= 1.0;
1044 –      break;
1045 –    case 3:
1046 –      coord[0] = 1 - coord[0];
1047 –      coord[1] *= -1.0;
1048 –      coord[2] = 1 - coord[2];
1049 –      break;
1050 –    }
1051 –    break;
1052 –  case 2:
1053 –    switch(next){
1054 –    case 0:
1055 –      coord[0] -= 1.0;
1056 –      coord[2] += 1.0;
1057 –      break;
1058 –    case 1:
1059 –      coord[1] -= 1.0;
1060 –      coord[2] += 1.0;
1061 –      break;
1062 –    case 3:
1063 –      coord[0] = 1 - coord[0];
1064 –      coord[1] = 1 - coord[1];
1065 –      coord[2] *= -1.0;
1066 –      break;
1067 –    }
1068 –    break;
1069 –  case 3:
1070 –    switch(next){
1071 –    case 0:
1072 –      coord[0] *= -1.0;
1073 –      coord[1] = 1 - coord[1];
1074 –      coord[2] = 1 - coord[2];
1075 –      break;
1076 –    case 1:
1077 –      coord[0] = 1 - coord[0];
1078 –      coord[1] *= -1.0;
1079 –      coord[2] = 1 - coord[2];
1080 –      break;
1081 –    case 2:
1082 –      coord[0] = 1 - coord[0];
1083 –      coord[1] = 1 - coord[1];
1084 –      coord[2] *= -1.0;
1085 –      break;
1086 –    }
1087 –    break;
1088 –  }
1089 – }
1090 –
1091 –
1092 – baryi_parent(coord,i)
931   BCOORD coord[3];
932   int i;
933   {
1096 –
934    switch(i) {
935    case 0:
936      /* update bary for child */
937 <    coord[0] = (coord[0] >> 1) + MAXBCOORD2;
937 >    coord[0] = (coord[0] >> 1) + MAXBCOORD4;
938      coord[1] >>= 1;
939      coord[2] >>= 1;
940      break;
941    case 1:
942      coord[0] >>= 1;
943 <    coord[1]  = (coord[1] >> 1) + MAXBCOORD2;
943 >    coord[1]  = (coord[1] >> 1) + MAXBCOORD4;
944      coord[2] >>= 1;
945      break;
946      
947    case 2:
948      coord[0] >>= 1;
949      coord[1] >>= 1;
950 <    coord[2] = (coord[2] >> 1) + MAXBCOORD2;
950 >    coord[2] = (coord[2] >> 1) + MAXBCOORD4;
951      break;
952      
953    case 3:
954 <    coord[0] = MAXBCOORD2 - (coord[0] >> 1);
955 <    coord[1] = MAXBCOORD2 - (coord[1] >> 1);
956 <    coord[2] = MAXBCOORD2 - (coord[2] >> 1);
954 >    coord[0] = MAXBCOORD4 - (coord[0] >> 1);
955 >    coord[1] = MAXBCOORD4 - (coord[1] >> 1);
956 >    coord[2] = MAXBCOORD4 - (coord[2] >> 1);
957      break;
958   #ifdef DEBUG
959    default:
960 <    eputs("baryi_parent():Invalid child\n");
960 >    eputs("bary_parent():Invalid child\n");
961      break;
962   #endif
963    }
964   }
965  
966 < baryi_from_child(coord,child,next)
966 > bary_from_child(coord,child,next)
967   BCOORD coord[3];
968   int child,next;
969   {
970   #ifdef DEBUG
971    if(child <0 || child > 3)
972    {
973 <    eputs("baryi_from_child():Invalid child\n");
973 >    eputs("bary_from_child():Invalid child\n");
974      return;
975    }
976    if(next <0 || next > 3)
977    {
978 <    eputs("baryi_from_child():Invalid next\n");
978 >    eputs("bary_from_child():Invalid next\n");
979      return;
980    }
981   #endif
# Line 1148 | Line 985 | int child,next;
985    switch(child){
986    case 0:
987        coord[0] = 0;
988 <      coord[1] = MAXBCOORD - coord[1];
989 <      coord[2] = MAXBCOORD - coord[2];
988 >      coord[1] = MAXBCOORD2 - coord[1];
989 >      coord[2] = MAXBCOORD2 - coord[2];
990        break;
991    case 1:
992 <      coord[0] = MAXBCOORD - coord[0];
992 >      coord[0] = MAXBCOORD2 - coord[0];
993        coord[1] = 0;
994 <      coord[2] = MAXBCOORD - coord[2];
994 >      coord[2] = MAXBCOORD2 - coord[2];
995        break;
996    case 2:
997 <      coord[0] = MAXBCOORD - coord[0];
998 <      coord[1] = MAXBCOORD - coord[1];
997 >      coord[0] = MAXBCOORD2 - coord[0];
998 >      coord[1] = MAXBCOORD2 - coord[1];
999        coord[2] = 0;
1000      break;
1001    case 3:
1002      switch(next){
1003      case 0:
1004        coord[0] = 0;
1005 <      coord[1] = MAXBCOORD - coord[1];
1006 <      coord[2] = MAXBCOORD - coord[2];
1005 >      coord[1] = MAXBCOORD2 - coord[1];
1006 >      coord[2] = MAXBCOORD2 - coord[2];
1007        break;
1008      case 1:
1009 <      coord[0] = MAXBCOORD - coord[0];
1009 >      coord[0] = MAXBCOORD2 - coord[0];
1010        coord[1] = 0;
1011 <      coord[2] = MAXBCOORD - coord[2];
1011 >      coord[2] = MAXBCOORD2 - coord[2];
1012        break;
1013      case 2:
1014 <      coord[0] = MAXBCOORD - coord[0];
1015 <      coord[1] = MAXBCOORD - coord[1];
1014 >      coord[0] = MAXBCOORD2 - coord[0];
1015 >      coord[1] = MAXBCOORD2 - coord[1];
1016        coord[2] = 0;
1017        break;
1018      }
# Line 1184 | Line 1021 | int child,next;
1021   }
1022  
1023   int
1024 < baryi_child(coord)
1024 > bary_child(coord)
1025   BCOORD coord[3];
1026   {
1027    int i;
1028  
1029 <  if(coord[0] > MAXBCOORD2)
1029 >  if(coord[0] > MAXBCOORD4)
1030    {
1031        /* update bary for child */
1032 <      coord[0] = (coord[0]<< 1) - MAXBCOORD;
1032 >      coord[0] = (coord[0]<< 1) - MAXBCOORD2;
1033        coord[1] <<= 1;
1034        coord[2] <<= 1;
1035        return(0);
1036    }
1037    else
1038 <    if(coord[1] > MAXBCOORD2)
1038 >    if(coord[1] > MAXBCOORD4)
1039      {
1040        coord[0] <<= 1;
1041 <      coord[1] = (coord[1] << 1) - MAXBCOORD;
1041 >      coord[1] = (coord[1] << 1) - MAXBCOORD2;
1042        coord[2] <<= 1;
1043        return(1);
1044      }
1045      else
1046 <      if(coord[2] > MAXBCOORD2)
1046 >      if(coord[2] > MAXBCOORD4)
1047        {
1048          coord[0] <<= 1;
1049          coord[1] <<= 1;
1050 <        coord[2] = (coord[2] << 1) - MAXBCOORD;
1050 >        coord[2] = (coord[2] << 1) - MAXBCOORD2;
1051          return(2);
1052        }
1053        else
1054           {
1055 <           coord[0] = MAXBCOORD - (coord[0] << 1);
1056 <           coord[1] = MAXBCOORD - (coord[1] << 1);
1057 <           coord[2] = MAXBCOORD - (coord[2] << 1);
1055 >           coord[0] = MAXBCOORD2 - (coord[0] << 1);
1056 >           coord[1] = MAXBCOORD2 - (coord[1] << 1);
1057 >           coord[2] = MAXBCOORD2 - (coord[2] << 1);
1058             return(3);
1059           }
1060   }
1061  
1062   int
1063 < baryi_nth_child(coord,i)
1063 > bary_nth_child(coord,i)
1064   BCOORD coord[3];
1065   int i;
1066   {
# Line 1231 | Line 1068 | int i;
1068    switch(i){
1069    case 0:
1070      /* update bary for child */
1071 <    coord[0] = (coord[0]<< 1) - MAXBCOORD;
1071 >    coord[0] = (coord[0]<< 1) - MAXBCOORD2;
1072      coord[1] <<= 1;
1073      coord[2] <<= 1;
1074      break;
1075    case 1:
1076      coord[0] <<= 1;
1077 <    coord[1] = (coord[1] << 1) - MAXBCOORD;
1077 >    coord[1] = (coord[1] << 1) - MAXBCOORD2;
1078      coord[2] <<= 1;
1079      break;
1080    case 2:
1081      coord[0] <<= 1;
1082      coord[1] <<= 1;
1083 <    coord[2] = (coord[2] << 1) - MAXBCOORD;
1083 >    coord[2] = (coord[2] << 1) - MAXBCOORD2;
1084      break;
1085    case 3:
1086 <    coord[0] = MAXBCOORD - (coord[0] << 1);
1087 <    coord[1] = MAXBCOORD - (coord[1] << 1);
1088 <    coord[2] = MAXBCOORD - (coord[2] << 1);
1086 >    coord[0] = MAXBCOORD2 - (coord[0] << 1);
1087 >    coord[1] = MAXBCOORD2 - (coord[1] << 1);
1088 >    coord[2] = MAXBCOORD2 - (coord[2] << 1);
1089      break;
1090    }
1091   }
1255 –
1256 –
1257 – baryi_children(coord,i,in_tri,rcoord,rin_tri)
1258 – BCOORD coord[3][3];
1259 – int i;
1260 – int in_tri[3];
1261 – BCOORD rcoord[3][3];
1262 – int rin_tri[3];
1263 – {
1264 –  int j;
1265 –
1266 –  for(j=0; j< 3; j++)
1267 –  {
1268 –    if(!in_tri[j])
1269 –    {
1270 –      rin_tri[j]=0;
1271 –      continue;
1272 –    }
1273 –    
1274 –    if(i != 3)
1275 –    {
1276 –      if(coord[j][i] < MAXBCOORD2)
1277 –        {
1278 –          rin_tri[j] = 0;
1279 –          continue;
1280 –        }
1281 –    }
1282 –    else
1283 –      if( !(coord[j][0] <= MAXBCOORD2 && coord[j][1] <= MAXBCOORD2 &&
1284 –            coord[j][2] <= MAXBCOORD2))
1285 –        {
1286 –          rin_tri[j] = 0;
1287 –          continue;
1288 –        }
1289 –      
1290 –    rin_tri[j]=1;
1291 –    VCOPY(rcoord[j],coord[j]);
1292 –    baryi_nth_child(rcoord[j],i);
1293 –  }
1294 –
1295 – }
1296 –
1297 – convert_dtol(b,bi)
1298 – double b[3];
1299 – BCOORD bi[3];
1300 – {
1301 –  int i;
1302 –
1303 –  for(i=0; i < 2;i++)
1304 –  {
1305 –
1306 –    if(b[i] <= 0.0)
1307 –    {
1308 – #ifdef EXTRA_DEBUG
1309 –      if(b[i] < 0.0)
1310 –        printf("under %f\n",b[i]);
1311 – #endif
1312 –      bi[i]= 0;
1313 –    }
1314 –    else
1315 –      if(b[i] >= 1.0)
1316 –      {
1317 – #ifdef EXTRA_DEBUG
1318 –        if(b[i] > 1.0)
1319 –          printf("over %f\n",b[i]);
1320 – #endif
1321 –        bi[i]= MAXBCOORD;
1322 –      }
1323 –      else
1324 –        bi[i] = (BCOORD)(b[i]*MAXBCOORD);
1325 –  }
1326 –  bi[2] = bi[0] +  bi[1];
1327 –  if(bi[2] > MAXBCOORD)
1328 –  {
1329 – #ifdef EXTRA_DEBUG
1330 –      printf("sum over %f\n",b[0]+b[1]);
1331 – #endif
1332 –      bi[2] = 0;
1333 –      bi[1] = MAXBCOORD - bi[0];
1334 –  }
1335 –  else
1336 –    bi[2] = MAXBCOORD - bi[2];
1337 –
1338 – }
1339 –
1340 – /* convert barycentric coordinate b in [-eps,1+eps] to [0,MAXLONG],
1341 –   dir unbounded to [-MAXLONG,MAXLONG]
1342 – */
1343 – bary_dtol(b,db,bi,dbi,t,w)
1344 – double b[3],db[3][3];
1345 – BCOORD bi[3];
1346 – BDIR dbi[3][3];
1347 – TINT t[3];
1348 – int w;
1349 – {
1350 –  int i,id,j,k;
1351 –  double d;
1352 –
1353 –  convert_dtol(b,bi);
1354 –
1355 –  for(j=w; j< 3; j++)
1356 –  {
1357 –    if(t[j] == HUGET)
1358 –    {
1359 –      max_index(db[j],&d);
1360 –      for(i=0; i< 3; i++)
1361 –        dbi[j][i] = (BDIR)(db[j][i]/d*MAXBDIR);
1362 –      break;
1363 –    }
1364 –    else
1365 –    {
1366 –      for(i=0; i< 3; i++)
1367 –          dbi[j][i] = (BDIR)(db[j][i]*MAXBDIR);
1368 –    }
1369 –  }
1370 – }
1371 –
1372 –
1373 – /* convert barycentric coordinate b in [-eps,1+eps] to [0,MAXLONG],
1374 –   dir unbounded to [-MAXLONG,MAXLONG]
1375 – */
1376 – bary_dtol_new(b,db,bi,boi,dbi,t)
1377 – double b[4][3],db[3][3];
1378 – BCOORD bi[3],boi[3][3];
1379 – BDIR dbi[3][3];
1380 – int t[3];
1381 – {
1382 –  int i,id,j,k;
1383 –  double d;
1384 –
1385 –  convert_dtol(b[3],bi);
1386 –
1387 –  for(j=0; j<3;j++)
1388 –  {
1389 –    if(t[j] != 1)
1390 –      continue;
1391 –    convert_dtol(b[j],boi[j]);
1392 –  }
1393 –  for(j=0; j< 3; j++)
1394 –  {
1395 –    k = (j+1)%3;
1396 –    if(t[k]==0)
1397 –      continue;
1398 –    if(t[k] == -1)
1399 –      {
1400 –        max_index(db[j],&d);
1401 –        for(i=0; i< 3; i++)
1402 –          dbi[j][i] = (BDIR)(db[j][i]/d*MAXBDIR);
1403 –        t[k] = 0;
1404 –      }
1405 –    else
1406 –      if(t[j] != 1)
1407 –        for(i=0; i< 3; i++)
1408 –          dbi[j][i] = (BDIR)(db[j][i]*MAXBDIR);
1409 –    else
1410 –      for(i=0; i< 3; i++)
1411 –        dbi[j][i] = boi[k][i] - boi[j][i];
1412 –    
1413 –  }
1414 – }
1415 –
1416 –
1417 – bary_dtolb(b,bi,in_tri)
1418 – double b[3][3];
1419 – BCOORD bi[3][3];
1420 – int in_tri[3];
1421 – {
1422 –  int i,j;
1423 –
1424 –  for(j=0; j<3;j++)
1425 –  {
1426 –    if(!in_tri[j])
1427 –      continue;
1428 –    for(i=0; i < 2;i++)
1429 –    {
1430 –    if(b[j][i] <= 0.0)
1431 –    {
1432 –      bi[j][i]= 0;
1433 –    }
1434 –    else
1435 –      if(b[j][i] >= 1.0)
1436 –      {
1437 –        bi[j][i]= MAXBCOORD;
1438 –      }
1439 –      else
1440 –        bi[j][i] = (BCOORD)(b[j][i]*MAXBCOORD);
1441 –    }
1442 –    bi[j][2] = MAXBCOORD - bi[j][0] - bi[j][1];
1443 –    if(bi[j][2] < 0)
1444 –      {
1445 –        bi[j][2] = 0;
1446 –        bi[j][1] = MAXBCOORD - bi[j][0];
1447 –      }
1448 –  }
1449 – }
1450 –
1451 –
1452 –
1453 –
1454 –
1455 –
1456 –
1457 –
1092  
1093  
1094  

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines