快速检测空间三角形相交算法的代码实现（DevillersGuigue算法）

快速检测空间三⾓形相交算法的代码实现（DevillersGuigue算

法）

Devillers & Guigue算法( 简称Devillers 算法) 通过三⾓形各顶点构成的⾏列式正负的⼏何意义来判断三⾓形中点、线、⾯之间的相对位置关系，从⽽判断两三⾓形是否相交。其基本原理如下：给定空间四个点：a(ax,ay,az),b=(bx,by,bz),c=(cx,cy,cz),d=(dx,dy,dz)，定义⾏列式如下：

[ a, b, c, d] 采⽤右⼿螺旋法则定义了四个空间点的位置关系。[ a, b, c, d] > 0 表⽰ d 在 a、b、c 按逆时针顺序所组成的三⾓形的正法线⽅向( 即上⽅) ; [ a, b, c, d] < 0 表⽰ d 在 △abc的下⽅; [ a, b, c, d] = 0 表⽰四点共⾯。

设两个三⾓形T1和T2，顶点分别为：V10，V11，V12和V20，V21，V22，三⾓形所在的平⾯分别为π1和π2，其法向量分别为N1和N2.算法先判别三⾓形和另⼀个三⾓形所在的平⾯的相互位置关系，提前排除不相交的情况。通过计算[V20,V21,V22,V1i].(i=0,1,2)来判断T1和π2的关系：如果所有的⾏列式的值都不为零且同号，则T1和T2不相交；否则T1和π2相交。相交⼜分为如下⼏种情况：

a)如果所有的⾏列式的值为零，则T1和T2共⾯，转化为共⾯的线段相交问题。

b)如果其中⼀个⾏列式的值为零，⽽其他两个⾏列式同号，则只有⼀个点在平⾯内，测试顶点是否则T2内部，是则相交，否则不相交；

c)否则T1的顶点位于平⾯π2两侧(包含T1的⼀条边在平⾯π2中的情况)。

再按照类似的⽅法对 T 2 和 π 1 作进⼀步的测试。如果通过测试, 则每个三⾓形必有确定的⼀点位于另⼀个三⾓形所在平⾯的⼀侧, ⽽另外两点位于其另⼀侧。算法分别循环置换每个三⾓形的顶点, 以使V10(V20)位于π2(π1)的⼀侧，另两个点位于其另⼀侧；同时对顶点V21，V22(V11,V12)进⾏交换操作，以确保V10(V20)位于π2(π1)的上⽅，即正法线⽅向。经过以上的预排除和置换操作，V10的邻边V10V11，V10V12和V20的邻边V20V21和V20V22与两平⾯的交线L相交于固定形式的点上，分别记为i，j，k，l(i<j,k<l)，如图：

这些点在L上形成的封闭区间为i1=[i,j],i2=[k,l].⾄此，两个三⾓形的相交测试问题转换为封闭区间i1，i2的重叠问题。若重叠则相交，否则不相交。由于交点形式固定，只需满⾜条件k<=j且i<=l即表明区间重叠，条件还可进⼀步缩减为判别式(1)是否成⽴：

[V10,V11,V20,V21]<=0 && [V10,V12,V22,V20]<=0 (1)

以上资料来⾃《计算机应⽤研究》第25卷第10期，邹益胜等的《空间三⾓形快速相交检测算法》下⾯是具体代码实现。

1typedef float float3[3];

3enum TopologicalStructure

5 INTERSECT, NONINTERSECT

6};

8struct Triangle

10 float3 Normal_0;

11 float3 Vertex_1, Vertex_2, Vertex_3;

12};

16 float x, y;

17};

20//三维点拷贝为⼆维点

21static void copy_point( point& p, float3 f )

22{

23 p.x = f[0];

24 p.y = f[1];

25}

27//四点⾏列式

28inline float get_vector4_det( float3 v1, float3 v2, float3 v3, float3 v4 ) 29{

30 float a[3][3];

31 for ( int i = 0; i != 3; ++i )

32 {

33 a[0][i] = v1[i] - v4[i];

34 a[1][i] = v2[i] - v4[i];

35 a[2][i] = v3[i] - v4[i];

36 }

38 return a[0][0] * a[1][1] * a[2][2]

39 + a[0][1] * a[1][2] * a[2][0]

40 + a[0][2] * a[1][0] * a[2][1]

41 - a[0][2] * a[1][1] * a[2][0]

42 - a[0][1] * a[1][0] * a[2][2]

43 - a[0][0] * a[1][2] * a[2][1];

44}

47//利⽤叉积计算点p相对线段p1p2的⽅位

48inline double direction( point p1, point p2, point p ){

49 return ( p.x - p1.x ) * ( p2.y - p1.y ) - ( p2.x - p1.x ) * ( p.y - p1.y ) ; 50}

53//确定与线段p1p2共线的点p是否在线段p1p2上人工呼吸器

54inline int on_segment( point p1, point p2, point p ){

55 double max = p1.x > p2.x ? p1.x : p2.x ;

56 double min = p1.x < p2.x ? p1.x : p2.x ;

57 double max1 = p1.y > p2.y ? p1.y : p2.y ;

58 double min1 = p1.y < p2.y ? p1.y : p2.y ;

球墨铸铁管qiumogg59 if ( p.x >= min && p.x <= max && p.y >= min1 && p.y <= max1 )

60 {

61 return 1;

62 }

63 else

64 {

65 return 0;

66 }

67}

70//判断线段p1p2与线段p3p4是否相交的主函数

71inline int segments_intersert( point p1, point p2, point p3, point p4 ){

72 double d1, d2, d3, d4;

73 d1 = direction( p3, p4, p1 );

74 d2 = direction( p3, p4, p2 );

75 d3 = direction( p1, p2, p3 );

76 d4 = direction( p1, p2, p4 );

77 if ( d1 * d2 < 0 && d3 * d4 < 0 )

78 {

81 else if ( d1 == 0 && on_segment( p3, p4, p1 ) == 1 )

82 {

83 return 1;

84 }

85 else if ( d2 == 0 && on_segment( p3, p4, p2 ) == 1 )

86 {

87 return 1;

88 }

89 else if ( d3 == 0 && on_segment( p1, p2, p3 ) == 1 )

90 {

91 return 1;

92 }

93 else if ( d4 == 0 && on_segment( p1, p2, p4 ) == 1 )

94 {

95 return 1;

96 }

97 return 0;

98}

100

101//判断同⼀平⾯的直线和三⾓形是否相交

102inline bool line_triangle_intersert_inSamePlane( Triangle* tri, float3 f1, float3 f2 ) 103{

104 point p1, p2, p3, p4;

105

106 copy_point( p1, f1 );

107

108 copy_point( p2, f2 );

109

110 copy_point( p3, tri->Vertex_1 );

111

112 copy_point( p4, tri->Vertex_2 );

113

114 if ( segments_intersert( p1, p2, p3, p4 ) )

115 {

116 return true;

117 }

118

119 copy_point( p3, tri->Vertex_2 );

120

121 copy_point( p4, tri->Vertex_3 );

122

123 if ( segments_intersert( p1, p2, p3, p4 ) )

124 {

125 return true;

126 }

127

128 copy_point( p3, tri->Vertex_1 );

129

130 copy_point( p4, tri->Vertex_3 );

131

132 if ( segments_intersert( p1, p2, p3, p4 ) )

133 {

134 return true;

135 }

136

137 return false;

138}

139

140

141inline void get_central_point(float3 ¢ralPoint,Triangle* tri)

142{

143 centralPoint[0]=(tri->Vertex_1[0]+tri->Vertex_2[0]+tri->Vertex_3[0])/3;

144

145 centralPoint[1]=(tri->Vertex_1[1]+tri->Vertex_2[1]+tri->Vertex_3[1])/3;

大灯随动转向146

147 centralPoint[2]=(tri->Vertex_1[2]+tri->Vertex_2[2]+tri->Vertex_3[2])/3;

148}

149

150

151

152

支架搬运车

153//向量之差

154inline void get_vector_diff( float3& aimV, const float3 a, const float3 b )

155{

156 aimV[0] = b[0] - a[0];

157

158 aimV[1] = b[1] - a[1];

159

160 aimV[2] = b[2] - a[2];

161}

162

163//向量内积

164inline float Dot( const float3& v1, const float3& v2 )

165{

166 return v1[0] * v2[0] + v1[1] * v2[1] + v1[2] * v2[2] ;

167}

168

169//重⼼法判断点是否在三⾓形内部

170inline bool is_point_within_triangle( Triangle* tri, float3 point )

171{

172 float3 v0;

173 get_vector_diff( v0, tri->Vertex_1, tri->Vertex_3 );

174 float3 v1;

175 get_vector_diff( v1, tri->Vertex_1, tri->Vertex_2 );

176 float3 v2;

177 get_vector_diff( v2, tri->Vertex_1, point );

178 float dot00 = Dot( v0, v0 ) ;

179 float dot01 = Dot( v0, v1 ) ;

180 float dot02 = Dot( v0, v2 ) ;

181 float dot11 = Dot( v1, v1 ) ;

182 float dot12 = Dot( v1, v2 ) ;

183 float inverDeno = 1 / ( dot00* dot11 - dot01* dot01 ) ;

184 float u = ( dot11* dot02 - dot01* dot12 ) * inverDeno ;

185 if ( u < 0 || u > 1 ) // if u out of range, return directly

186 {

187 return false ;

188 }

189 float v = ( dot00* dot12 - dot01* dot02 ) * inverDeno ;

190 if ( v < 0 || v > 1 ) // if v out of range, return directly

191 {

192 return false ;

193 }

194 return u + v <= 1 ;

195}

196

197

198//判断同⼀平⾯内的三⾓形是否相交

199inline bool triangle_intersert_inSamePlane( Triangle* tri1, Triangle* tri2 )

200{

201 if ( line_triangle_intersert_inSamePlane( tri2, tri1->Vertex_1, tri1->Vertex_2 ) ) 202 {

203 return true;

204 }

205 else if ( line_triangle_intersert_inSamePlane( tri2, tri1->Vertex_2, tri1->Vertex_3 ) ) 206 {

207 return true;

208 }

209 else if ( line_triangle_intersert_inSamePlane( tri2, tri1->Vertex_1, tri1->Vertex_3 ) )

210 {

光纤气体传感器211 return true;

212 }

213 else

214 {

215 float3 centralPoint1,centralPoint2;

216

217 get_central_point(centralPoint1,tri1);

218

219 get_central_point(centralPoint2,tri2);

220

221 if(is_point_within_triangle( tri2, centralPoint1 ) || is_point_within_triangle( tri1, centralPoint2 ) )

222 {

223 return true;

224 }

225

226 return false;

227 }

228}

229

230//Devillers算法主函数

231inline TopologicalStructure judge_triangle_topologicalStructure( Triangle* tri1, Triangle* tri2 )

232{

233 //设tri1所在的平⾯为p1,tri2所在的平⾯为p2

234 float p1_tri2_vertex1 = get_vector4_det( tri1->Vertex_1, tri1->Vertex_2, tri1->Vertex_3, tri2->Vertex_1 ); 235

236 float p1_tri2_vertex2 = get_vector4_det( tri1->Vertex_1, tri1->Vertex_2, tri1->Vertex_3, tri2->Vertex_2 ); 237

238 float p1_tri2_vertex3 = get_vector4_det( tri1->Vertex_1, tri1->Vertex_2, tri1->Vertex_3, tri2->Vertex_3 ); 239

240

241 if ( p1_tri2_vertex1 > 0 && p1_tri2_vertex2 > 0 && p1_tri2_vertex3 > 0 )

242 {

243 return NONINTERSECT;

244 }

245

246 if ( p1_tri2_vertex1 < 0 && p1_tri2_vertex2 < 0 && p1_tri2_vertex3 < 0 )

247 {

248 return NONINTERSECT;

249 }

250

251

252 if ( p1_tri2_vertex1 == 0 && p1_tri2_vertex2 == 0 && p1_tri2_vertex3 == 0 )

253 {

254 if ( triangle_intersert_inSamePlane( tri1, tri2 ) )

255 {

256 return INTERSECT;

257 }

258 else

259 {

260 return NONINTERSECT;

261 }

262 }

263

264

265 if ( p1_tri2_vertex1 == 0 && p1_tri2_vertex2 * p1_tri2_vertex3 > 0 )

266 {

267 if ( is_point_within_triangle( tri1, tri2->Vertex_1 ) )

268 {

269 return INTERSECT;微电脑时间控制器

270 }

271 else

272 {

273 return NONINTERSECT;

本文发布于:2024-09-25 02:27:13，感谢您对本站的认可！

本文链接：https://www.17tex.com/tex/2/191737.html

上一篇：如何根据离散点自动绘制等值线（等高线）之三角形法

下一篇：TriangleCount三角形计数

标签：相交算法线段是否空间重叠测试

留言与评论（共有 0 条评论）