Forster-Overfelt版Greiner-Horman多边形裁剪算法的伪代码有什么问题?
问题
我正在尝试理解并实现Greiner-Horman多边形裁剪算法的Forster-Overfelt版本。我读过其他关于这个算法的帖子,但是我还是无法让它正常工作。
我知道有些地方出错了,因为即使是一个简单的例子,它也会产生两个多边形的错误交集,而这个例子没有退化情况:
subjpoly = [(0,0),(6,0),(6,6),(0,6),(0,0)]
clippoly = [(1,4),(3,8),(5,4),(5,10),(1,10),(1,4)]
这产生了一个交集:
[ [(5.0, 6.0), (4.0, 6.0), (5, 4), (5.0, 6.0)],
[(1.0, 6.0), (2.0, 6.0), (4.0, 6.0)] ]
可视化效果如下:
所以这个问题不是关于某一段代码或语言语法,而是关于理解算法并将其转化为伪代码。谢谢大家的帮助!
算法
这里介绍的算法直接基于Forster-Oberfelt论文第4.2节的内容,应该与之相似,但显然我遗漏了什么,导致结果不正确。
第一部分:
首先,循环遍历被裁剪的多边形和裁剪的多边形,标记每个顶点的位置为“在内”、“在外”或“在另一多边形上”。
for s in subj.iter():
s.loc = testLocation(s, clip)
for c in clip.iter():
c.loc = testLocation(c, subj)
第二部分:
接下来,循环遍历被裁剪多边形的交点。
for s in subj.iter():
if s.intersect:
子部分2.1:
处理每个被裁剪多边形的交点,要么取消它们的交点标记,要么标记为进入或退出,并对邻近的交点做同样的处理。注意:文章中解释的算法只说明了如何标记主要的被裁剪多边形,但没有说明如何标记邻近的交点,所以在这里我假设两者使用相同的规则进行标记。
mark(s)
mark(s.neighbour)
标记处理规则定义如下:
def mark(s):
curlocs = (s.prev.loc,s.next.loc)
neighlocs = (s.neighbour.prev.loc,s.neighbour.next.loc)
# in in
if curlocs == ("in","in"):
if neighlocs == ("in","in")\
or neighlocs == ("out","out")\
or neighlocs == ("on","on"):
s.intersect = False
else:
s.entry = True
# out out
elif curlocs == ("out","out"):
if neighlocs == ("in","in")\
or neighlocs == ("out","out")\
or neighlocs == ("on","on"):
s.intersect = False
else:
s.entry = False
# on on
elif curlocs == ("on","on"):
if neighlocs == ("in","in")\
or neighlocs == ("out","out")\
or neighlocs == ("on","on"):
s.intersect = False
else:
# label opposite of neighbour
# NOTE: this is not specified in the article,
# but one cannot take the opposite of the neighbour's entry flag
# if the neighbour hasn't been marked yet,
# thus the decision to mark the neighbour first
mark(s.neighbour)
s.entry = not s.neighbour
# partial exit
elif curlocs == ("in","on")\
or curlocs == ("on","out"):
s.entry = False
# partial entry
elif curlocs == ("on","in")\
or curlocs == ("out","on"):
s.entry = True
# normal exit
elif curlocs == ("in","out"):
s.entry = False
# normal entry
elif curlocs == ("out","in"):
s.entry = True
子部分2.2:
最后,确保当前交点和邻近交点没有相同的进入或退出标记;如果有,就禁用它们的交点标记并更改位置标记。
if s.entry and s.neighbour.entry:
s.intersect = False
s.neighbour.intersect = False
s.loc = "in"
elif not s.entry and not s.neighbour.entry:
s.intersect = False
s.neighbour.intersect = False
s.loc = "out"
附加问题
一个附加问题是如何让这个算法支持并集和交集操作,因为原始的Greiner算法通过简单地反转初始的进入/退出标记来支持并集,但这个Forster算法并没有使用这样的标记?
1 个回答
再说一句关于并集和交集的事。主要的意思是,并集操作和交集操作的方向是相反的。也就是说,如果在进行交集操作时需要沿着多边形向后移动,那么在进行并集操作时就要向前移动,反之亦然。
接下来讲讲算法:首先,我们先看看这个算法的大致框架。我这里的算法每次交集操作只会生成一个多边形,所以你需要调整一下,才能生成多个多边形。
'''
The following is an adaptation of the above Greiner-Hormann* algorithm to deal
with degenerate cases. The adaptation was briefly described by Liu et al.**
*Greiner, G. and Hormann K., Efficient Clipping of Arbitrary Polygons, ACM
Trans. on Graphics, 17(2), 1998, pp.71-83
**Liu, Y. K., Wang X. Q., Bao S. Z., Gombosi M., and Zalik B, An Algorithm for
Polygon Clipping and for Determining Polygon Intersections and Unions, Comp. &
Geo, 33, pp. 589-598, 2007
'''
def clip(subject, constraint):
subject, constraint = inside_outside(subject, constraint) #label vertices as inside or outside
subject, constraint = poly_inters(subject, constraint) #find intersections
subject, constraint = label(subject, constraint) #label intersections and entry or exit and possibly remove
flag = True #loop flag
#set our current location to the first point in subject
current = subject.first
#loop through our polygon until we have found the first intersection
while flag:
current = current.next
#Either an intersection has been found or no intersections found
if current.intersect or current.pt == subject.first.pt:
flag = False
#reset our flag for the new loop
flag = True
#If a point lies outside of c and there was an intersection clip s
if current.intersect:
append(clipped, current.pt)
While flag:
#Entry
if current.en:
clipped, current = forward(clipped, current)
#Exit
else:
clipped, current = backward(clipped, current)
#Check to see if we have completed a polygon
if current.pt == clipped.first.pt:
#check if the polygon intersect at a point
if clipped.num_vertices is not 1:
#remove the last vertex because it is also the first
remove(clipped, clipped.last)
#we have created our polygon so we can exit
flag = .FALSE.
#change to the neighbor polygon since we are at a new intersection
current = current.neighbor
#Check if one polygon is contained wholly within the other
elif contained(subject, constraint):
clipped = subject
elif contained(subject, constraint):
clipped = constraint
return clipped
现在我们可以讨论标记的问题。下面的代码是一个循环,用来标记交点是内部还是外部。它并不包含判断内部或外部的逻辑,只是操作的顺序。
#label intersections as entry or exit
def label(poly1, poly2):
#cycle through first polygon and label intersections as en or ex
current = poly1.first
for i in range(0,poly1.num_vertices):
if current.intersect:
current = intersect_cases(current)
#Make sure current is still an intersection
if current.isect:
current.neighbor = intersect_cases(current.neighbor)
#if the intersection is en/en or ex/ex
if current.en == current.neighbor.en:
current = remove_inter(current)
current = current.next #move to the next point
return poly1, poly2
最后,我们要处理各种情况的标记。
#deal with the cases
#on/on, on/out, on/in, out/on, out/out, out/in, in/on, in/out, in/in
def intersect_cases(current):
neighbor = current.neighbor
#on/on
if current.prev.intersect and current.next.intersect:
#Determine what to do based on the neighbor
#en tag is the opposite of the neighbor's en tag
if neighbor.prev.intersect and neighbor.next.intersect:
current = remove_inter(current)
current.en = True
neighbor.en = True
elif neighbor.prev.intersect and not neighbor.next.en:
current.en = False
elif neighbor.prev.intersect and neighbor.next.en:
current.en = True
elif not neighbor.prev.en and neighbor.next.intersect:
current.en = False
elif not (neighbor.prev.en or neighbor.next.en):
current = remove_inter(current)
current.en = True
neighbor.en = False
elif not neighbor.prev.en and neighbor.next.en:
current.en = False
elif neighbor.prev.en and neighbor.next.isect:
current.en = True
elif neighbor.prev.en and not neighbor.next.en:
current.en = True
elif neighbor.prev.en and neighbor.next.en:
current = remove_inter(current)
current.en = False
neighbor.en = True
#on/out
elif current.prev.intersect and not current.next.en:
current.en = False
#on/in
elif current.prev.intersect and current.next.en:
current.en = True
#out/on
elif not current.prev.en and current.next.intersect:
current.en = True
#out/out
elif not (current.prev%en or current.next.en):
if neighbor.prev%intersect and neighbor.next.intersect:
current = remove_inter(current)
neighbor.en = True
elif neighbor.prev.en == neighbor.next.en:
current = remove_inter(current)
else:
if neighbor.prev.en and not neighbor.next.en:
current.en = True
else:
current.en = False
#out/in
elif not current.prev.en and current.next.en:
current.en = True
#in/on
elif current.prev.en and current.next.intersect:
current.en = False
#in/out
elif current.prev.en and not current.next.en:
current.en = False
#in/in
elif current.prev.en and current.next.en:
if neighbor.prev.intersect and neighbor.next.intersect:
current = remove_inter(current)
neighbor.en = False
elif neighbor.prev.en == neighbor.next.en:
current = remove_inter(current)
else:
if neighbor.prev.en and not neighbor.next.en:
current.en = True
else:
current.en = False
return current
上面的代码还没有经过测试,也不是为了效率而写的,而是为了让人更容易阅读和理解。