参考:http://www.dyn4j.org/2010/04/gjk-gilbert-johnson-keerthi/
和SAT(分离轴)法一样, GJK可以判断两个凸图形是否重叠. 比起SAT, GJK优在用同一套办法可以处理所有的图形, 而SAT判断某两种不同图形(多边形-多边形/多边形-圆形/圆形-圆形等)都要区别处理.
GJK原理: 如果两个凸图形的闵可夫斯基差包含原点, 那么这两个图形重叠. 所以问题转变成判断一个闵可夫斯基差图形是否包含原点.
闵可夫斯基和就是把两个图形里的点都当作向量, 把两个图形里的向量点都分别相加, 得到一个新的图形.
闵可夫斯基差就是把其中一个图形的向量都取反, 相当于绕原点反转, 然后再把两个图形求闵可夫斯基和.
GJK不会直接求两个图形的闵可夫斯基差, 而是通过下面这个辅助函数, 输入不同的查找方向, 返回闵可夫斯基差上的一个点, 迭代地用这些点组成一个个子图形, 判断子图形是否包含原点.
proc getSupportPoint(s1: Sprite2D, s2: Sprite2D, dir: Vector2D): Vector2D = let invM1 = s1.invMatrix invM2 = s2.invMatrix v1 = s1.matrix * s1.graph.getFarthestPointInDirection(invM1 * dir - invM1 * Vector2D(x: 0, y: 0)) v2 = s2.matrix * s2.graph.getFarthestPointInDirection(invM2 * dir.negate() - invM2 * Vector2D(x: 0, y: 0)) return v1 - v2
上面这个函数会返回闵可夫斯基差图形在某个方向上最远的一个点. 这里先把查找方向变换到物体空间, 由物体自己找到自身在查找方向上的最远点后, 再变换回世界空间.
查找多边形在查找方向上的最远的点的函数如下. 通过把所有顶点投影到查找方向上, 根据投影的大小判断距离.
method getFarthestPointInDirection*(self: Polygon, direction: Vector2D): Vector2D = var curr = self[0] bestProjection: float32 = curr * direction projection: float32 result = curr for i in 1..self.len-1: curr = self[i] projection = curr * direction if projection > bestProjection: bestProjection = projection result = curr
圆形. norm是求单位向量.
method getFarthestPointInDirection*(self: Circle, direction: Vector2D): Vector2D = direction.norm() * self.radius
GJK每次生成的子图形都是单纯形(simplex). 单纯形在2D下是三角形, 在3D下是四面体.
维基上给出的GJK伪代码如下:
function GJK_intersection(shape p, shape q, vector initial_axis): vector A = Support(p, initial_axis) - Support(q, -initial_axis) simplex s = {A} vector D = -A loop: A = Support(p, D) - Support(q, -D) if dot(A, D) < 0: reject s = s ∪ A s, D, contains_origin = NearestSimplex(s) if contains_origin: accept
个人的具体实现如下:
假设闵可夫斯基差的图形是下图的椭圆形.
首先随便用一个初始查找方向, 通过getSupportPoint函数找到单纯形上的第一个点A. 此时, 把向量OA投影到查找方向向量上. 如果投影是负数, 那么原点O一定在闵可夫斯基差外(在下图是过A的红色虚线右侧), 判定为两图形不重叠. 之后每次找到一个单纯形的顶点都如此判断.
然后, 把查找方向取反(或者用AO方向也可), 找到第二个点B. 然后判断原点不在过点B的虚线左侧.
把A和B连线, 求指向AO方向的垂直于AB的一个向量当作查找方向, 得到点C. 现在, 原点可能存在的地方就是三条红色虚线与AB包围的四边形内了.
求垂直于AC的指向AB方向的向量acPerp, 求垂直于BC的指向BA方向的向量bcPerp.
把AO投影到acPerp上, 如果投影为负, 那么原点不在单纯形内. 此时抛弃距原点最远的点B, 用点A和点C当作新的点A‘和点B‘, 重复以上过程.
把OB投影到bcPerp上, 如上判断原点是在BC内侧还是外侧. 如果是在BC外, 则抛弃离原点最远的A, 用B和C作为新的点A和B, 重复以上过程.
如果两个投影都为非负, 那么原点在单纯形ABC内, 判断两个图形重叠.
proc isIntersect*(s1, s2: Sprite2D): bool = var dir = Vector2D(x: 1, y: 0) simplexA, simplexB, simplexC: Vector2D simplexA = getSupportPoint(s1, s2, dir) if simplexA * dir < 0: return false dir = simplexA.negate() simplexB = getSupportPoint(s1, s2, dir) if simplexB * dir < 0: return false let ao = simplexA.negate() ab = simplexB - simplexA dir = tripleProduct(ab, ao, ab) while true: simplexC = getSupportPoint(s1, s2, dir) if simplexC * dir < 0: return false let ba = simplexA - simplexB ac = simplexC - simplexA bc = simplexC - simplexB acPerp = tripleProduct(ac, ba.negate(), ac) bcPerp = tripleProduct(bc, ba, bc) if acPerp * simplexA > 0: simplexB = simplexC dir = acPerp.negate() elif bcPerp * simplexB > 0: simplexA = simplexC dir = bcPerp.negate() else: return true
上面的tripleProduct是连续两次叉乘, A×B×A会得到一个垂直于A的指向B方向的向量.
proc tripleProduct*(v1, v2, v3: Vector2D): Vector2D = result.x = -(v1.x * v2.y - v1.y * v2.x) * v3.y result.y = (v1.x * v2.y - v1.y * v2.x) * v3.x
GJK(Gilbert–Johnson–Keerthi) 判断任意凸图形是否重叠
原文:http://www.cnblogs.com/pngx/p/4395992.html