1 | |
---|
2 | /* |
---|
3 | ----------------------------------------------------------------------------- |
---|
4 | This source file is part of GIMPACT Library. |
---|
5 | |
---|
6 | For the latest info, see http://gimpact.sourceforge.net/ |
---|
7 | |
---|
8 | Copyright (c) 2006 Francisco Leon. C.C. 80087371. |
---|
9 | email: projectileman@yahoo.com |
---|
10 | |
---|
11 | This library is free software; you can redistribute it and/or |
---|
12 | modify it under the terms of EITHER: |
---|
13 | (1) The GNU Lesser General Public License as published by the Free |
---|
14 | Software Foundation; either version 2.1 of the License, or (at |
---|
15 | your option) any later version. The text of the GNU Lesser |
---|
16 | General Public License is included with this library in the |
---|
17 | file GIMPACT-LICENSE-LGPL.TXT. |
---|
18 | (2) The BSD-style license that is included with this library in |
---|
19 | the file GIMPACT-LICENSE-BSD.TXT. |
---|
20 | |
---|
21 | This library is distributed in the hope that it will be useful, |
---|
22 | but WITHOUT ANY WARRANTY; without even the implied warranty of |
---|
23 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the files |
---|
24 | GIMPACT-LICENSE-LGPL.TXT and GIMPACT-LICENSE-BSD.TXT for more details. |
---|
25 | |
---|
26 | ----------------------------------------------------------------------------- |
---|
27 | */ |
---|
28 | |
---|
29 | #include "GIMPACT/gim_trimesh.h" |
---|
30 | |
---|
31 | |
---|
32 | #define FABS(x) (float(fabs(x))) /* implement as is fastest on your machine */ |
---|
33 | |
---|
34 | /* some macros */ |
---|
35 | |
---|
36 | #define CLASSIFY_TRIPOINTS_BY_FACE(v1,v2,v3,faceplane,out_of_face)\ |
---|
37 | { \ |
---|
38 | _distances[0] = DISTANCE_PLANE_POINT(faceplane,v1);\ |
---|
39 | _distances[1] = _distances[0] * DISTANCE_PLANE_POINT(faceplane,v2);\ |
---|
40 | _distances[2] = _distances[0] * DISTANCE_PLANE_POINT(faceplane,v3); \ |
---|
41 | if(_distances[1]>0.0f && _distances[2]>0.0f)\ |
---|
42 | {\ |
---|
43 | out_of_face = 1;\ |
---|
44 | }\ |
---|
45 | else\ |
---|
46 | {\ |
---|
47 | out_of_face = 0;\ |
---|
48 | }\ |
---|
49 | }\ |
---|
50 | |
---|
51 | /* sort so that a<=b */ |
---|
52 | #define SORT(a,b) \ |
---|
53 | if(a>b) \ |
---|
54 | { \ |
---|
55 | float c; \ |
---|
56 | c=a; \ |
---|
57 | a=b; \ |
---|
58 | b=c; \ |
---|
59 | } |
---|
60 | |
---|
61 | |
---|
62 | /* this edge to edge test is based on Franlin Antonio's gem: |
---|
63 | "Faster Line Segment Intersection", in Graphics Gems III, |
---|
64 | pp. 199-202 */ |
---|
65 | #define EDGE_EDGE_TEST(V0,U0,U1) \ |
---|
66 | Bx=U0[i0]-U1[i0]; \ |
---|
67 | By=U0[i1]-U1[i1]; \ |
---|
68 | Cx=V0[i0]-U0[i0]; \ |
---|
69 | Cy=V0[i1]-U0[i1]; \ |
---|
70 | f=Ay*Bx-Ax*By; \ |
---|
71 | d=By*Cx-Bx*Cy; \ |
---|
72 | if((f>0 && d>=0 && d<=f) || (f<0 && d<=0 && d>=f)) \ |
---|
73 | { \ |
---|
74 | e=Ax*Cy-Ay*Cx; \ |
---|
75 | if(f>0) \ |
---|
76 | { \ |
---|
77 | if(e>=0 && e<=f) return 1; \ |
---|
78 | } \ |
---|
79 | else \ |
---|
80 | { \ |
---|
81 | if(e<=0 && e>=f) return 1; \ |
---|
82 | } \ |
---|
83 | } |
---|
84 | |
---|
85 | #define EDGE_AGAINST_TRI_EDGES(V0,V1,U0,U1,U2) \ |
---|
86 | { \ |
---|
87 | float Ax,Ay,Bx,By,Cx,Cy,e,d,f; \ |
---|
88 | Ax=V1[i0]-V0[i0]; \ |
---|
89 | Ay=V1[i1]-V0[i1]; \ |
---|
90 | /* test edge U0,U1 against V0,V1 */ \ |
---|
91 | EDGE_EDGE_TEST(V0,U0,U1); \ |
---|
92 | /* test edge U1,U2 against V0,V1 */ \ |
---|
93 | EDGE_EDGE_TEST(V0,U1,U2); \ |
---|
94 | /* test edge U2,U1 against V0,V1 */ \ |
---|
95 | EDGE_EDGE_TEST(V0,U2,U0); \ |
---|
96 | } |
---|
97 | |
---|
98 | #define POINT_IN_TRI(V0,U0,U1,U2) \ |
---|
99 | { \ |
---|
100 | float a,b,c,d0,d1,d2; \ |
---|
101 | /* is T1 completly inside T2? */ \ |
---|
102 | /* check if V0 is inside tri(U0,U1,U2) */ \ |
---|
103 | a=U1[i1]-U0[i1]; \ |
---|
104 | b=-(U1[i0]-U0[i0]); \ |
---|
105 | c=-a*U0[i0]-b*U0[i1]; \ |
---|
106 | d0=a*V0[i0]+b*V0[i1]+c; \ |
---|
107 | \ |
---|
108 | a=U2[i1]-U1[i1]; \ |
---|
109 | b=-(U2[i0]-U1[i0]); \ |
---|
110 | c=-a*U1[i0]-b*U1[i1]; \ |
---|
111 | d1=a*V0[i0]+b*V0[i1]+c; \ |
---|
112 | \ |
---|
113 | a=U0[i1]-U2[i1]; \ |
---|
114 | b=-(U0[i0]-U2[i0]); \ |
---|
115 | c=-a*U2[i0]-b*U2[i1]; \ |
---|
116 | d2=a*V0[i0]+b*V0[i1]+c; \ |
---|
117 | if(d0*d1>0.0) \ |
---|
118 | { \ |
---|
119 | if(d0*d2>0.0) return 1; \ |
---|
120 | } \ |
---|
121 | } |
---|
122 | |
---|
123 | int coplanar_tri_tri(GIM_TRIANGLE_DATA *tri1, |
---|
124 | GIM_TRIANGLE_DATA *tri2) |
---|
125 | { |
---|
126 | short i0,i1; |
---|
127 | /* first project onto an axis-aligned plane, that maximizes the area */ |
---|
128 | /* of the triangles, compute indices: i0,i1. */ |
---|
129 | PLANE_MINOR_AXES(tri1->m_planes.m_planes[0], i0, i1); |
---|
130 | |
---|
131 | /* test all edges of triangle 1 against the edges of triangle 2 */ |
---|
132 | EDGE_AGAINST_TRI_EDGES(tri1->m_vertices[0],tri1->m_vertices[1],tri2->m_vertices[0],tri2->m_vertices[1],tri2->m_vertices[2]); |
---|
133 | EDGE_AGAINST_TRI_EDGES(tri1->m_vertices[1],tri1->m_vertices[2],tri2->m_vertices[0],tri2->m_vertices[1],tri2->m_vertices[2]); |
---|
134 | EDGE_AGAINST_TRI_EDGES(tri1->m_vertices[2],tri1->m_vertices[0],tri2->m_vertices[0],tri2->m_vertices[1],tri2->m_vertices[2]); |
---|
135 | |
---|
136 | /* finally, test if tri1 is totally contained in tri2 or vice versa */ |
---|
137 | POINT_IN_HULL(tri1->m_vertices[0],(&tri2->m_planes.m_planes[1]),3,i0); |
---|
138 | if(i0==0) return 1; |
---|
139 | |
---|
140 | POINT_IN_HULL(tri2->m_vertices[0],(&tri1->m_planes.m_planes[1]),3,i0); |
---|
141 | if(i0==0) return 1; |
---|
142 | |
---|
143 | return 0; |
---|
144 | } |
---|
145 | |
---|
146 | |
---|
147 | |
---|
148 | #define NEWCOMPUTE_INTERVALS(VV0,VV1,VV2,D0,D1,D2,D0D1,D0D2,A,B,C,X0,X1) \ |
---|
149 | { \ |
---|
150 | if(D0D1>0.0f) \ |
---|
151 | { \ |
---|
152 | /* here we know that D0D2<=0.0 */ \ |
---|
153 | /* that is D0, D1 are on the same side, D2 on the other or on the plane */ \ |
---|
154 | A=VV2; B=(VV0-VV2)*D2; C=(VV1-VV2)*D2; X0=D2-D0; X1=D2-D1; \ |
---|
155 | } \ |
---|
156 | else if(D0D2>0.0f)\ |
---|
157 | { \ |
---|
158 | /* here we know that d0d1<=0.0 */ \ |
---|
159 | A=VV1; B=(VV0-VV1)*D1; C=(VV2-VV1)*D1; X0=D1-D0; X1=D1-D2; \ |
---|
160 | } \ |
---|
161 | else if(D1*D2>0.0f || D0!=0.0f) \ |
---|
162 | { \ |
---|
163 | /* here we know that d0d1<=0.0 or that D0!=0.0 */ \ |
---|
164 | A=VV0; B=(VV1-VV0)*D0; C=(VV2-VV0)*D0; X0=D0-D1; X1=D0-D2; \ |
---|
165 | } \ |
---|
166 | else if(D1!=0.0f) \ |
---|
167 | { \ |
---|
168 | A=VV1; B=(VV0-VV1)*D1; C=(VV2-VV1)*D1; X0=D1-D0; X1=D1-D2; \ |
---|
169 | } \ |
---|
170 | else if(D2!=0.0f) \ |
---|
171 | { \ |
---|
172 | A=VV2; B=(VV0-VV2)*D2; C=(VV1-VV2)*D2; X0=D2-D0; X1=D2-D1; \ |
---|
173 | } \ |
---|
174 | else \ |
---|
175 | { \ |
---|
176 | /* triangles are coplanar */ \ |
---|
177 | return coplanar_tri_tri(tri1,tri2); \ |
---|
178 | } \ |
---|
179 | }\ |
---|
180 | |
---|
181 | |
---|
182 | |
---|
183 | int gim_triangle_triangle_overlap( |
---|
184 | GIM_TRIANGLE_DATA *tri1, |
---|
185 | GIM_TRIANGLE_DATA *tri2) |
---|
186 | { |
---|
187 | vec3f _distances; |
---|
188 | char out_of_face; |
---|
189 | CLASSIFY_TRIPOINTS_BY_FACE(tri1->m_vertices[0],tri1->m_vertices[1],tri1->m_vertices[2],tri2->m_planes.m_planes[0],out_of_face); |
---|
190 | if(out_of_face==1) return 0; |
---|
191 | |
---|
192 | CLASSIFY_TRIPOINTS_BY_FACE(tri2->m_vertices[0],tri2->m_vertices[1],tri2->m_vertices[2],tri1->m_planes.m_planes[0],out_of_face); |
---|
193 | if(out_of_face==1) return 0; |
---|
194 | |
---|
195 | |
---|
196 | float du0=0,du1=0,du2=0,dv0=0,dv1=0,dv2=0; |
---|
197 | float D[3]; |
---|
198 | float isect1[2], isect2[2]; |
---|
199 | float du0du1=0,du0du2=0,dv0dv1=0,dv0dv2=0; |
---|
200 | short index; |
---|
201 | float vp0,vp1,vp2; |
---|
202 | float up0,up1,up2; |
---|
203 | float bb,cc,max; |
---|
204 | |
---|
205 | /* compute direction of intersection line */ |
---|
206 | VEC_CROSS(D,tri1->m_planes.m_planes[0],tri2->m_planes.m_planes[0]); |
---|
207 | |
---|
208 | /* compute and index to the largest component of D */ |
---|
209 | max=(float)FABS(D[0]); |
---|
210 | index=0; |
---|
211 | bb=(float)FABS(D[1]); |
---|
212 | cc=(float)FABS(D[2]); |
---|
213 | if(bb>max) max=bb,index=1; |
---|
214 | if(cc>max) max=cc,index=2; |
---|
215 | |
---|
216 | /* this is the simplified projection onto L*/ |
---|
217 | vp0= tri1->m_vertices[0][index]; |
---|
218 | vp1= tri1->m_vertices[1][index]; |
---|
219 | vp2= tri1->m_vertices[2][index]; |
---|
220 | |
---|
221 | up0= tri2->m_vertices[0][index]; |
---|
222 | up1= tri2->m_vertices[1][index]; |
---|
223 | up2= tri2->m_vertices[2][index]; |
---|
224 | |
---|
225 | /* compute interval for triangle 1 */ |
---|
226 | float a,b,c,x0,x1; |
---|
227 | NEWCOMPUTE_INTERVALS(vp0,vp1,vp2,dv0,dv1,dv2,dv0dv1,dv0dv2,a,b,c,x0,x1); |
---|
228 | |
---|
229 | /* compute interval for triangle 2 */ |
---|
230 | float d,e,f,y0,y1; |
---|
231 | NEWCOMPUTE_INTERVALS(up0,up1,up2,du0,du1,du2,du0du1,du0du2,d,e,f,y0,y1); |
---|
232 | |
---|
233 | float xx,yy,xxyy,tmp; |
---|
234 | xx=x0*x1; |
---|
235 | yy=y0*y1; |
---|
236 | xxyy=xx*yy; |
---|
237 | |
---|
238 | tmp=a*xxyy; |
---|
239 | isect1[0]=tmp+b*x1*yy; |
---|
240 | isect1[1]=tmp+c*x0*yy; |
---|
241 | |
---|
242 | tmp=d*xxyy; |
---|
243 | isect2[0]=tmp+e*xx*y1; |
---|
244 | isect2[1]=tmp+f*xx*y0; |
---|
245 | |
---|
246 | SORT(isect1[0],isect1[1]); |
---|
247 | SORT(isect2[0],isect2[1]); |
---|
248 | |
---|
249 | if(isect1[1]<isect2[0] || isect2[1]<isect1[0]) return 0; |
---|
250 | return 1; |
---|
251 | } |
---|