Skip to content

Instantly share code, notes, and snippets.

@kevinjardine

kevinjardine/tiler.py

Created July 24, 2011 12:52
Show Gist options
  • Star 2 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save kevinjardine/1102587 to your computer and use it in GitHub Desktop.
Save kevinjardine/1102587 to your computer and use it in GitHub Desktop.
Download ZIP
A universal tiler and polygon decompositor
Raw
tiler.py
# To use this, you will want to change the location of the saved data
# at least and possibly the functions called below (at the end of the file).
# You first need to install the Python Shapely library for manipulating polygons:
# http://trac.gispython.org/lab/wiki/Shapely
from math import pi, sin, cos, atan2
import shapely.geometry
import sys, copy
def position_opposite(large_order,small_order):
sym_lookup = {}
s0 = large_order/small_order
#opposites only make sense if s0 is even
if ((s0 % 2) == 0):
for i in range(0,s0):
if i % 2 == 0:
sym_lookup[i] = i
else:
sym_lookup[i] = (i+(s0/2)-1) % s0
else:
sym_lookup = position_sequential(large_order,small_order)
return sym_lookup
def position_6(large_order,small_order):
#return [3, 2, 1, 0, 5, 4]
return [4,5,0,1,2,3]
def position_sequential(large_order,small_order):
sym_lookup = {}
s0 = large_order/small_order
for i in range(0,s0):
sym_lookup[i] = i
return sym_lookup
def position_alternating(large_order,small_order):
sym_lookup = {}
s0 = large_order/small_order
for i in range(0,s0/2):
sym_lookup[i] = 2*i
for i in range(0,s0/2):
sym_lookup[i+s0/2] = 2*i+1
return sym_lookup
class Tiling(object):
"""Stores and manipulates a tile patch"""
def __init__(self):
self.pc = 0
self.polygons = []
self.ec = 0
self.edges = []
self.vc = 0
self.vertices = []
self.spc = 0
self.shapely_polygons = []
self.stack = []
def add_to_shapely_polygons(self,pe):
poly = []
for e in pe:
poly.append(e[0])
self.shapely_polygons.append(shapely.geometry.Polygon(poly))
self.spc += 1
def get_vertex_list(self,v0,e0,p0):
error = False
e = self.edges[e0]
v = self.vertices[v0]
poly = [[v.x,v.y]]
v1 = e.v1
if v1 == v0:
v1 = e.v0
i = 0
# cycle around the vertices until we return to v0
while v1 != v0:
if i > len(self.vertices):
error = True
break
v = self.vertices[v1]
poly.append([v.x,v.y])
e0 = self.get_other_edge2(v1,e0,p0)
v1 = self.get_other_vertex(v1,e0)
i += 1
if error:
print "cannot compute vertex list"
sys.exit()
else:
return poly
def get_edge_list(self,v0,e0,p0):
error = False
e = self.edges[e0]
poly = [e0]
v1 = e.v1
if v0 == v1:
v1 = e.v0
i = 0
# cycle around the vertices until we return to v0
while v1 != v0:
e0 = self.get_other_edge2(v1,e0,p0)
if e0 == -1 or i > len(self.vertices):
error = True
break
v1 = self.get_other_vertex(v1,e0)
poly.append(e0)
i += 1
if error:
print "cannot compute edge list"
sys.exit()
else:
return poly
def convert_to_shapely_polygon(self,p,r=True):
e0 = self.polygons[p].edges[0]
e = self.edges[e0]
v0 = e.v0
poly = self.get_vertex_list(v0,e0,p)
if poly == False or len(poly) < 3:
#sanity checks
return False
else:
return shapely.geometry.Polygon(poly)
def add_to_shapely_polygons2(self,poly):
self.shapely_polygons.append(shapely.geometry.Polygon(poly))
self.spc += 1
def polygon_to_svg(self,p):
if self.shapely_polygons[p] == False:
return ''
points = []
for x,y in list(self.shapely_polygons[p].exterior.coords):
points.append(str(20*x+500)+","+str(20*y+500))
r = '<polygon style="stroke:#000000;stroke-width:1;fill:#F0F0F0" points="'
r += " ".join(points)
r += '"/>'+"\n"
return r
def labels_to_svg(self):
r = ''
for v in self.vertices:
r += '<text x="'+str(20*v.x+500)
r += '" y="'+str(20*v.y+500)
r += '">('+str(v.num)+")</text>\n"
for e in self.edges:
vx = (0.5)*(self.vertices[e.v0].x+self.vertices[e.v1].x)
vy = (0.5)*(self.vertices[e.v0].y+self.vertices[e.v1].y)
r += '<text x="'+str(20*vx+500)
r += '" y="'+str(20*vy+500)
r += '">'+str(e.num)+"</text>\n"
return r
# TODO: change the data structures to make this more efficient to compute?
def get_other_edge2(self,v0,e0,p0):
for e in self.polygons[p0].edges:
if e != e0 and (self.edges[e].v0 == v0 or self.edges[e].v1 == v0):
return e
return -1
def add_to_vertex_sweep(self,v,e,p):
vertex = self.vertices[v]
if v in []:
print "@@@************ adding to sweeps", v, e, p
print "current sweep_e", vertex.sweep_e
print "current sweep_p", vertex.sweep_p
e_prev = -1
p_prev = -1
e2 = -1
p2 = -1
if len(vertex.sweep_e) <= 1:
# the sweep direction is not yet established
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
if e not in vertex.sweep_e:
vertex.sweep_e.append(e)
else:
# need to add only in the sweep direction
if vertex.sweep_e[-1] == e:
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
elif vertex.sweep_e[0] == e:
if p not in vertex.sweep_p:
vertex.sweep_p = [p] + vertex.sweep_p
else:
# try the polygon method
if vertex.sweep_p[-1] == p:
if e not in vertex.sweep_e:
vertex.sweep_e.append(e)
elif vertex.sweep_p[0] == p:
if e not in vertex.sweep_e:
vertex.sweep_e = [e] + vertex.sweep_e
else:
#try the other edge method
#TODO: should the get_other_edge code be here?
e2 = self.get_other_edge2(v,e,p)
if vertex.sweep_e[-1] == e2:
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
elif vertex.sweep_e[0] == e2:
if p not in vertex.sweep_p:
vertex.sweep_p = [p] + vertex.sweep_p
else:
# try the other polygon method
p2 = self.get_other_polygon(e,p)
if vertex.sweep_p[-1] == p2:
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
if e not in vertex.sweep_e:
vertex.sweep_e.append(e)
elif vertex.sweep_p[0] == p2:
if p not in vertex.sweep_p:
vertex.sweep_p = [p] + vertex.sweep_p
if e not in vertex.sweep_e:
vertex.sweep_e = [e] + vertex.sweep_e
else:
# try the previous polygon method
p_prev = self.get_previous_polygon(v,e,p)
if vertex.sweep_p[-1] == p_prev:
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
if e not in vertex.sweep_e:
vertex.sweep_e.append(e)
elif vertex.sweep_p[0] == p_prev:
if p not in vertex.sweep_p:
vertex.sweep_p = [p] + vertex.sweep_p
if e not in vertex.sweep_e:
vertex.sweep_e = [e] + vertex.sweep_e
else:
# try the previous edge method
e_prev = self.get_other_edge2(v,e,p_prev)
p_prev = self.get_previous_polygon(v,e,p)
if vertex.sweep_p[-1] == p_prev:
if p not in vertex.sweep_p:
vertex.sweep_p.append(p)
if e not in vertex.sweep_e:
vertex.sweep_e.append(e)
elif vertex.sweep_p[0] == p_prev:
if p not in vertex.sweep_p:
vertex.sweep_p = [p] + vertex.sweep_p
if e not in vertex.sweep_e:
vertex.sweep_e = [e] + vertex.sweep_e
else:
pass
def add_to_sweeps(self,edge_ids,e0,p0):
for e in edge_ids:
if e != e0 and e != -1:
self.add_to_vertex_sweep(self.edges[e].v0,e,p0)
self.add_to_vertex_sweep(self.edges[e].v1,e,p0)
def get_other_polygon(self,e,p):
plist = self.polygons_for([e])
if len(plist) < 2:
#print "warning: cannot find both polygons for", e, p
return plist[0]
if plist[0] == p:
return plist[1]
else:
return plist[0]
def get_previous_polygon(self,v,e0,p):
elist = self.polygons[p].edges
for e in elist:
if e != e0 and (self.edges[e].v0 == v or self.edges[e].v1 == v):
return self.get_other_polygon(e,p)
return -1
def get_previous_vertex(self,v0,e0):
for p in self.polygons:
elist = p.edges
if e0 in elist:
for e in elist:
if e != e0:
if self.edges[e].v0 == v0:
return self.edges[e].v1
elif self.edges[e].v1 == v0:
return self.edges[e].v0
return -1
def create_polygon(self,e0,p,poly_edges):
#TODO: try to simplify this
#In particular, some of the new_vertex calls are redundant
edge_ids = self.get_edge_ids(poly_edges)
self.add_to_sweeps(edge_ids,e0,p)
vl = []
for i in range(0,len(edge_ids)):
e = edge_ids[i]
if e == -1:
# need a new edge
vfi = self.new_vertex(poly_edges[i][0])
vsi = self.new_vertex(poly_edges[i][1])
e = self.new_edge(vfi,vsi)
else:
vfi = self.edges[e].v0
vsi = self.edges[e].v1
if e != e0:
self.add_to_polygon(p,e)
self.add_to_vertex_sweep(vfi,e,p)
self.add_to_vertex_sweep(vsi,e,p)
for pe in poly_edges:
vi = self.new_vertex(pe[0])
if vi not in vl:
vl.append(vi)
return vl
def get_edge_ids(self,poly_edges):
edge_ids = []
for pe in poly_edges:
e = self.get_edge_from_poly_edge(pe)
edge_ids.append(e)
return edge_ids
def add_to_polygon(self,p,e):
self.polygons[p].edges.append(e)
def remove_from_polygon(self,p,e):
v0 = self.edges[e].v0
v1 = self.edges[e].v1
self.polygons[p].edges.remove(e)
#search vertex list looking for duplicates
# TODO: find something more efficient?
def new_vertex(self,v):
for vi in range(0,self.vc):
if not different_vertex(v,self.vertices[vi]):
return vi
v = Vertex(v[0],v[1])
v.num = self.vc
self.vertices.append(v)
self.vc += 1
return v.num
def new_edge(self,v0,v1):
for ei in range(0,self.ec):
if not different_edge(v0,v1,self.edges[ei]):
return ei
e = Edge(v0,v1)
e.num = self.ec
self.edges.append(e)
self.ec += 1
return e.num
def new_polygon(self,e0,ptype,make_internal=True):
p = Polygon(e0,ptype,make_internal)
p.num = self.pc
self.polygons.append(p)
self.pc += 1
return p.num
def get_edge_from_poly_edge(self,pe):
for ei in range(0,self.ec):
e = self.edges[ei]
v0 = [self.vertices[e.v0].x, self.vertices[e.v0].y]
v1 = [self.vertices[e.v1].x, self.vertices[e.v1].y]
t0 = equal_vertex(v0,pe[0]) and equal_vertex(v1,pe[1])
t1 = equal_vertex(v0,pe[1]) and equal_vertex(v1,pe[0])
if t0 or t1:
return ei
return -1
def start_polygon(self,n):
v0 = self.new_vertex([-0.5,0.0])
v1 = self.new_vertex([0.5,0.0])
e = self.new_edge(v1,v0)
ab,external,le = self.get_leading_edge2(v0,0)
self.complete_at_vertex(v0,[n],external,1,[],[])
def star_around(self,a):
v0 = self.new_vertex([-0.5,0.0])
v1 = self.new_vertex([0.5,0.0])
e = self.new_edge(v0,v1)
n = int(2*pi/a + 0.5)
for i in range(0,n):
p = compute_rhombus(v1,e,a)
e = get_other_edge(v1,e,p)
def edges_for(self,v):
ev = []
for e in self.edges:
if (e.v0 == v) or (e.v1 == v):
ev.append(e.num)
return ev
# elist is a list of edges
def polygons_for(self,elist):
pe = []
for p in self.polygons:
for e in elist:
if e in p.edges:
pe.append(p.num)
break
return pe
def compute_circle_gap_at_vertex(self,v):
at = 0
vertex = self.vertices[v]
for p in vertex.sweep_p:
g = vertex.sweep_lu[p]
if g > 0:
ng = (pi / g) * (g-2)
#ng = d[g]['p']
else:
ng = -g
at += ng
return 2*pi - at
def get_poly_edge_from_edge(self,e):
v0 = self.vertices[e[0]]
v1 = self.vertices[e[1]]
return [[v0.x,v0.y],[v1.x,v1.y]]
def validate_sweep(self,v,species,rspecies):
if len(species) == 0:
# an empty species turns off validation
return True
# search for a solution to the sweep
sweep = self.get_sweep(v)
# if the sweep is complete, backstep to make sure that this is a legal vertex
if abs(self.compute_circle_gap_at_vertex(v)) < 0.0001:
sweep = sweep[1:]
#print "in validate_sweep, checking sweep", convert_sweep(sweep)
for cspecies in [species, rspecies]:
for i in range(0,len(cspecies)):
#print "against vertex type", convert_sweep(cspecies[i])
pl = get_plist_from_sweep(sweep,cspecies[i])
if len(pl) != 0:
return True
return False
def find_overlapping_polygons(self,pl):
overlaps = []
p0 = []
for e in pl:
p0.append(e[0])
poly = shapely.geometry.Polygon(p0)
i = 0
for p in self.shapely_polygons:
try:
if poly.intersects(p):
pintersect = poly.intersection(p)
if (type(pintersect) is shapely.geometry.polygon.Polygon) and pintersect.area > 0.01:
#print "overlap detected", i, p,poly,pintersect,pintersect.area
overlaps.append(i)
except:
print "overlap intersect function error"
return []
i += 1
return overlaps
def punch(self,v,e,n):
external = [self.edges[e].v0,self.edges[e].v1]
vraw = [self.vertices[v].x,self.vertices[v].y]
e0 = self.get_poly_edge_from_edge(external)
parity = -1
if n > 0:
r = compute_polygon(vraw,e0,n,parity)
else:
r = compute_rhombus(vraw,e0,-n,parity)
# only do something if the new polygon overlaps exactly one
# other polygon
overlaps = self.find_overlapping_polygons(r[0])
if len(overlaps) == 0:
print "Error: no overlap in punch"
sys.exit()
elif len(overlaps) > 1:
print "Error: multiple overlaps in punch"
sys.exit()
else :
ps = t.polygons_for([e])
if len(ps) > 0:
p0 = ps[0]
else:
print "Error in punch: cannot find existing polygon for edge ", e
sys.exit()
poly_edges = r[0]
root_edge = r[1]
a = r[2]
e0 = self.get_edge_from_poly_edge(root_edge)
p = self.new_polygon(e0,n)
vl = self.create_polygon(e0,p,poly_edges)
# add the sweep angles to the vertices
# TODO: this code is wrong for punches
if n > 0:
for vi in vl:
vertex = self.vertices[vi]
vertex.sweep_lu[p] = n
vertex.sweep_plu[p] = n
else:
# rhombuses have two angles, so the sweep code is a bit more complex
if v == vl[0] or v == vl[2]:
ang = [a,pi-a,a,pi-a]
else:
ang = [pi-a,a,pi-a,a]
for i in range(0,len(vl)):
vertex = self.vertices[vl[i]]
vertex.sweep_lu[p] = -ang[i]
if ang[i] > (pi/2):
vertex.sweep_plu[p] = int(2*pi/(pi-ang[i])+0.5)
else:
vertex.sweep_plu[p] = -int((2*pi/ang[i])+0.5)
self.add_to_shapely_polygons(poly_edges)
esp = self.polygons[p].edges
esp0 = copy.deepcopy(self.polygons[p0].edges)
for e0 in esp:
if e0 in esp0:
self.remove_from_polygon(p0,e0)
else:
self.polygons[p0].edges.append(e0)
self.shapely_polygons[p0] = self.convert_to_shapely_polygon(p0)
return True
def split(self,v0,e0,p0):
#take the edges from the end of e0 to the beginning of e1
# copy it and translate it to the beginning of e0
e1 = self.find_equal_slope_edge(v0,e0,p0)
edges = self.get_edge_list(v0,e0,p0)
if edges == False:
return False
es = []
for e in edges:
es.append(e)
if e == e1:
break
v = self.vertices[v0]
v0x = v.x
v0y = v.y
v1 = self.get_other_vertex(v0,e0)
v = self.vertices[v1]
v1x = v.x
v1y = v.y
vx = v0x - v1x
vy = v0y - v1y
es2 = copy.deepcopy(es)
for e1 in es[1:-1]:
e = self.edges[e1]
pe1 = self.get_poly_edge_from_edge([e.v0,e.v1])
nv0 = self.new_vertex([pe1[0][0]+vx,pe1[0][1]+vy])
nv1 = self.new_vertex([pe1[1][0]+vx,pe1[1][1]+vy])
e = self.new_edge(nv0,nv1)
if e not in es2:
es2.append(e)
self.add_to_polygon(p0,e)
# make a new polygon out of the edges
# sets ptype to 1 to signify a worm
p = self.new_polygon(es2[0],1)
for e in es2[1:]:
self.add_to_polygon(p,e)
s = self.convert_to_shapely_polygon(p)
self.shapely_polygons.append(s)
self.spc += 1
#remove edges from original polygon
for e in es:
self.remove_from_polygon(p0,e)
self.shapely_polygons[p0] = self.convert_to_shapely_polygon(p0)
def clip_all(self,p0):
"""clip worm completely into rhombs"""
n = len(self.polygons[p0].edges)
while n > 4:
self.clip(p0)
n = len(self.polygons[p0].edges)
def clip(self,p0):
num_edges = len(self.edges)
es = copy.deepcopy(self.polygons[p0].edges)
v0 = self.edges[es[0]].v0
vs = self.get_vertex_list(v0,es[0],p0)
n = len(es)/2
v0 = self.edges[es[1]].v0
v1 = self.edges[es[n+1]].v0
e3 = self.new_edge(v0,v1)
if e3 < num_edges:
#sanity check
#oops, try the other vertices
v0 = self.edges[es[1]].v1
v1 = self.edges[es[n+1]].v1
e3 = self.new_edge(v0,v1)
# make a new worm out of the edges
p = self.new_polygon(es[0],1)
self.remove_from_polygon(p0,es[0])
self.add_to_polygon(p,es[1])
self.remove_from_polygon(p0,es[1])
self.add_to_polygon(p,e3)
self.add_to_polygon(p,es[n+1])
self.remove_from_polygon(p0,es[n+1])
self.polygons[p0].edges = [e3] + self.polygons[p0].edges
s = self.convert_to_shapely_polygon(p)
self.shapely_polygons.append(s)
self.shapely_polygons[p0] = self.convert_to_shapely_polygon(p0)
def find_equal_slope_edge(self,v0,e0,p0):
es = self.get_edge_list(v0,e0,p0)
vs = self.get_vertex_list(v0,e0,p0)
if vs == False or len(vs) < 3:
# sanity check
return None
v0 = vs[0]
v1 = vs[1]
vax = v1[0]
vay = v1[1]
dx0 = vax - v0[0]
dy0 = vay - v0[1]
i = 1
for v in vs[2:]:
dx = vax - v[0]
dy = vay - v[1]
if abs(dx) < 0.0001 and abs(dx0) < 0.0001:
return es[i]
elif abs(dx-dx0) < 0.0001 and abs(dy-dy0) < 0.0001:
return es[i]
vax = v[0]
vay = v[1]
i += 1
dx = v0[0] - vax
dy = v0[1] - vay
if abs(dx) < 0.0001 and abs(dx0) < 0.0001:
return es[i]
elif abs(dx-dx0) < 0.0001 and abs(dy-dy0) < 0.0001:
return es[i]
return None
def rhombic_triangle(self,m):
sa = 1.0/(3*m)
self.start_polygon(-pi*sa)
stn = -1
for i in range(2,m):
stn = stn + i - 1
stn2 = 2*stn
a = -pi*i*sa
self.extend_edge(stn+i-1,stn2+1,[a])
for n in range(0,i-2):
self.extend_edge(stn+n,stn2+2*(n+1),[a])
self.extend_edge(stn+2*i-1,stn2+2*(i-1),[a])
i = m
stn = stn + i - 1
stn2 = 2*stn
self.extend_edge(stn+i-1,stn2+1,[3])
for n in range(0,i-2):
self.extend_edge(stn+n,stn2+2*(n+1),[3])
self.extend_edge(stn+2*i-1,stn2+2*(i-1),[3])
def fan(self,m):
sa = 1.0/(m)
self.start_polygon(-pi*sa)
stn = -1
for i in range(2,m):
stn = stn + i - 1
stn2 = 2*stn
a = -pi*i*sa
self.extend_edge(stn+i-1,stn2+1,[a])
for n in range(0,i-2):
self.extend_edge(stn+n,stn2+2*(n+1),[a])
self.extend_edge(stn+2*i-1,stn2+2*(i-1),[a])
def universal_tiler_tu(self,n,decomp_polygon):
self.rhombic_triangle(n)
tri = copy.deepcopy(self)
v0 = (n*(n-1) + 4*n)/2
e0 = (2*n*(n-1) + 6*n)/2 -1
inc = len(self.edges) - n
for i in range(0,5):
self.merge(v0,e0+i*inc,tri,v0,e0-1)
incv = v0 - n
for i in range(0,5):
if i == 0:
extra = 0
else:
extra = 1
self.merge(v0-1+i*incv+extra, inc-1+i*inc+extra,tri,v0-1,inc-2)
self.merge(0, inc-1+5*inc,tri,v0-1,inc-2)
self.merge(0,2*(inc-1+5*inc),decomp_polygon,0,0)
self.merge(v0-1,14*n+7*n*(n-1)-2,decomp_polygon,0,0)
def get_decompose_data(self,v0,e0,p0):
es0 = self.get_edge_list(v0,e0,p0)
if len(es0) < 2:
#sanity check
print "unable to get decompose data"
sys.exit()
es0a = es0[1]
v0a = self.get_other_vertex(v0,es0a)
es0b = es0[-1]
v0b = self.get_other_vertex(v0,es0b)
return [[v0a,es0a,-1],[v0b,es0b,1]]
# given adjacent edges, find a common vertex
def get_common_vertex(self,e0,e1):
e0 = self.edges[e0]
e1 = self.edges[e1]
if e0.v0 == e1.v0:
return e0.v0
elif e0.v0 == e1.v1:
return e0.v0
elif e0.v1 == e1.v0:
return e0.v1
elif e0.v1 == e1.v1:
return e0.v1
else:
return False
def decompose_zonagon(self,p0,as_worms=True):
n = len(self.polygons[p0].edges)
if n % 2 == 1:
#sanity check, this is not a zonagon!
return
for i in range(0,(n/2)-2):
e0 = self.polygons[p0].edges[0]
v0 = self.edges[e0].v0
es = self.get_edge_list(v0,e0,p0)
va = self.get_common_vertex(es[0],es[1])
d0 = self.get_decompose_data(va, es[0],p0)
self.split(d0[0][0],d0[0][1],p0)
p1 = self.pc - 1
p = self.polygons[p1]
if as_worms:
p.ptype = 1
else:
self.clip_all(p1)
def get_sym_edge_index(self,i,symclasses):
for j in range(0,len(symclasses)):
if i == symclasses[j]:
return j
return False
def get_split_data(self,v0,e0,p0,p1,do_extra_inc=False):
es0 = self.get_edge_list(v0,e0,p0)
n = len(es0)
extra_inc = 0
if n % 2 == 0:
m = (n/2)-1
if do_extra_inc:
extra_inc = 1
else:
m = n/2
es1 = self.polygons[p1].edges
found = False
# seek along the p0 edges starting with e0
# until we find an edge that is in p1
for i in range(1,n):
if es0[i] in es1:
found = True
break
if not found:
#sanity check
print "Error: could not find edge in get_split_data"
sys.exit()
i0 = i
# now try the same thing but in the other direction
for i in range(1,n):
if es0[n-i] in es1:
found = True
break
if not found:
#sanity check
print "Error: could not find edge in get_split_data"
sys.exit()
i1 = i
# we have found our split edges in both directions,
# so return the results
if i0 > 0:
i0b = i0 - 1
else:
i0b = n - 1
if i0 < n - 1:
i0a = i0 + 1
else:
i0a = 0
if i1 > 0:
i1a = i1 - 1
else:
i1a = n - 1
if i1 < n - 1:
i1b = i1 + 1
else:
i1b = 0
va = self.get_common_vertex(es0[i0a],es0[i0])
r0 = [va,es0[i0],m-i0+1,es0[i0a]-es0[i0]]
vb = self.get_common_vertex(es0[n-i1a-1],es0[n-i1-1])
r1 = [vb,es0[n-i1],m-i1+1+extra_inc,es0[i1]-es0[i1b]]
return [r0,r1]
def decompose_polygon(self,p0,n,small_order,position_function=position_opposite,eliminate_sym_classes=True):
# By default, this function tries to arrange the
# small order regular polygons in rings consisting of rhombs
# and two regular polygons of opposite symmetry classes.
# If the large order polygon has an even number of sides,
# this ensures that the polygon inside each ring is
# a zonagon and so can optionally be decomposed into rhombs.
if n == 0:
return
large_order = len(self.polygons[p0].edges)
s0 = large_order/small_order
sym_lookup = position_function(large_order,small_order)
if small_order % 2 == 0:
if eliminate_sym_classes and (large_order == n*small_order):
m1 = n - 1
else:
m1 = n
else:
if large_order == n*small_order:
m1 = n - 1
else:
m1 = n
print sym_lookup
for i in range(0,m1):
vs,es,symclasses = self.report(p0,large_order,small_order)
qi0 = self.get_sym_edge_index(sym_lookup[i],symclasses)
e0 = self.polygons[p0].edges[qi0]
v0 = self.edges[e0].v0
if qi0 == 0:
qi1 = len(es)-1
else:
qi1 = qi0 - 1
va = self.get_common_vertex(es[qi0],es[qi1])
self.punch(va, es[qi0],small_order)
p1 = self.pc - 1
r = self.get_split_data(va,es[qi0],p1,p0,eliminate_sym_classes)
for k in range(0,2):
s,f,m1,inc = r[k]
for p in range(0,m1):
self.split(s+p*inc,f+p*inc,p0)
self.report(p0,large_order,small_order)
def clip_all_worms(self):
# clip the worms into rhombs
num_polygons = self.pc
for p in range(1,num_polygons):
poly = self.polygons[p]
if len(poly.edges) > 4 and poly.ptype == 1:
self.clip_all(p)
# a specialized function that creates and decomposes a regular polygon
# with 6n sides where n is odd.
def decompose6_odd(self,n):
#self.start_polygon(6*n)
m = n/2
self.punch(0,0,n)
p0 = self.pc - 1
self.punch(3*n+1,3*n,n)
p1 = p0 + 1
d0 = self.get_decompose_data(0,0,p0)
d1 = self.get_decompose_data(3*n+1,3*n,p1)
d = [d0[0],d1[0],d1[1],d0[1]]
for i in range(0,4):
s,f,inc = d[i]
for p in range(0,m):
self.split(s+p*inc,f+p*inc,0)
self.report(0,6*n,n)
e0 = self.polygons[0].edges[0]
v0 = self.edges[e0].v0
es = self.get_edge_list(v0,e0,0)
va = self.get_common_vertex(es[n-1],es[n])
vb = self.get_common_vertex(es[3*n-1],es[3*n])
self.punch(va,es[n-1],n)
p0 = self.pc - 1
self.punch(vb,es[3*n-1],n)
p1 = p0 + 1
d0 = self.get_decompose_data(va,es[n-1],p0)
d1 = self.get_decompose_data(vb,es[3*n-1],p1)
d = [d0[1],d1[1],d0[0],d1[0]]
for i in range(0,4):
s,f,inc = d[i]
for p in range(0,m):
self.split(s+p*inc,f+p*inc,0)
self.report(0,6*n,n)
e0 = self.polygons[0].edges[0]
v0 = self.edges[e0].v0
es = self.get_edge_list(v0,e0,0)
va = self.get_common_vertex(es[n+m],es[n+m+1])
self.punch(va, es[n+m],n)
p0 = self.pc - 1
d0 = self.get_decompose_data(va, es[n+m],p0)
d = [d0[0],d0[1]]
for i in range(0,2):
s,f,inc = d[i]
for p in range(0,m):
self.split(s+p*inc,f+p*inc,0)
def complete_at_vertex(self,v,alist,external,parity,species,rspecies):
success = True
vraw = [self.vertices[v].x,self.vertices[v].y]
e0 = self.get_poly_edge_from_edge(external)
results = []
new_polygons = []
for n in alist:
if self.insufficient_gap(v,n):
success = False
break
if n > 0:
r = compute_polygon(vraw,e0,n,parity)
else:
r = compute_rhombus(vraw,e0,-n,parity)
if overlaps_polygons(r[0], self.shapely_polygons+new_polygons):
if v not in []:
success = False
break
e0 = get_other_edge(vraw,e0,r[0])
if len(e0) <= 1:
success = True
break
#so far, so good
results.append([n,r])
# add to temporary shapely list
spoly = []
for e in r[0]:
spoly.append(e[0])
new_polygons.append(shapely.geometry.Polygon(spoly))
if success:
# it is possible to add all the polygons in the list without overlap,
# so now we will attempt to add them to the tiling
#save state
state = [self.vc,self.ec,self.pc]
for poly in results:
n = poly[0]
r = poly[1]
poly_edges = r[0]
root_edge = r[1]
a = r[2]
e0 = self.get_edge_from_poly_edge(root_edge)
p = self.new_polygon(e0,n)
vl = self.create_polygon(e0,p,poly_edges)
# now attempt to add the sweep angles to the vertices
if n > 0:
for vi in vl:
vertex = self.vertices[vi]
vertex.sweep_lu[p] = n
vertex.sweep_plu[p] = n
if not self.validate_sweep(vi,species,rspecies):
success = False
vb = vi
break
else:
# rhombuses have two angles, so the sweep code is a bit more complex
if v == vl[0] or v == vl[2]:
ang = [a,pi-a,a,pi-a]
else:
ang = [pi-a,a,pi-a,a]
for i in range(0,len(vl)):
vertex = self.vertices[vl[i]]
vertex.sweep_lu[p] = -ang[i]
if ang[i] > (pi/2):
vertex.sweep_plu[p] = int(2*pi/(pi-ang[i])+0.5)
else:
vertex.sweep_plu[p] = -int((2*pi/ang[i])+0.5)
if not self.validate_sweep(vl[i],species,rspecies):
success = False
vb = vl[i]
break
if success:
self.add_to_shapely_polygons(poly_edges)
else:
# backtrack
sb = self.get_sweep(vb)
self.vc,self.ec,self.pc = state
self.vertices = self.vertices[:self.vc]
self.edges = self.edges[:self.ec]
self.polygons = self.polygons[:self.pc]
self.shapely_polygons = self.shapely_polygons[:self.pc]
self.clean_sweeps()
break
return success
def get_external_edges(self,v):
ev = []
for e in self.edges:
if (e.v0 == v) or (e.v1 == v):
if not e.internal:
ev.append(e.num)
if len(ev) < 3:
return ev
# We have too many edges so need to determine which
# ones are really external. If an edge belongs to
# two polygons, it is internal.
ep = {}
for p in self.polygons:
for ei in ev:
if ei in p.edges:
if ei not in ep:
ep[ei] = []
ep[ei].append(p.num)
ee = []
for ei in ev:
if ei in ep:
if len(ep[ei]) == 2:
# TODO: make the next line work with backtracking
#edges[ei].internal = True
pass
else:
ee.append(ei)
return ee
def get_leading_edge2(self,v,gap):
ee = self.get_external_edges(v)
if len(ee) == 0:
return [0,[],-1]
elif len(ee) == 1:
ne = self.edges[ee[0]]
nv0 = self.vertices[ne.v0]
nv1 = self.vertices[ne.v1]
return [0,[ne.v1,ne.v0],ne.num]
e0 = self.edges[ee[0]]
e1 = self.edges[ee[1]]
if e0.v0 == v:
v00 = self.vertices[e0.v1]
v01 = self.vertices[v]
else:
v00 = self.vertices[e0.v0]
v01 = self.vertices[v]
if e1.v0 == v:
v11 = self.vertices[e1.v1]
v10 = self.vertices[v]
else:
v11 = self.vertices[e1.v0]
v10 = self.vertices[v]
if abs(gap-pi) <= 0.001:
# need to resolve ambiguity
v2 = self.get_previous_vertex(v,ee[0])
vx = self.vertices[v2].x
vy = self.vertices[v2].y
ab = angle_between([v00.x-v01.x,v00.y-v01.y],[vx-v10.x,vy-v10.y])
if ab > 0:
return [ab,[v00.num,v],e0.num]
else:
return [ab,[v,v00.num],e0.num]
else:
ab = angle_between([v00.x-v01.x,v00.y-v01.y],[v11.x-v10.x,v11.y-v10.y])
if abs(abs(ab) - gap) < 0.001:
if ab > 0:
return [ab,[v,v00.num],e0.num]
else:
return [ab,[v00.num,v],e0.num]
else:
if ab < 0:
return [ab,[v,v00.num],e0.num]
else:
return [ab,[v00.num,v],e0.num]
def insufficient_gap(self,v,n):
if n < 0:
ang = -n
else:
ang = (pi / n) * (n-2)
#ang = d[n]['p']
gap = self.compute_circle_gap_at_vertex(v)
return ang - gap > 0.001
# make sure that the proposed tiling is legal
# at both ends of an edge
def verify_external(self,v,e):
if self.edges[e].internal:
return False
v2 = self.get_other_vertex(v,e)
gap2 = self.compute_circle_gap_at_vertex(v2)
return gap2 > 0
def dump_polygon(self,p):
print p, str(len(self.polygons[p].edges)) + "-gon"
for e in self.polygons[p].edges:
v0 = self.vertices[edges[e].v0]
v1 = self.vertices[edges[e].v1]
print "Edge[",e,"]:",self.edges[e].internal,self.edges[e].v0,[v0.x,v0.y],self.edges[e].v1,[v1.x,v1.y]
def get_other_vertex(self,v,e):
edge = self.edges[e]
if edge.v0 != v:
return edge.v0
else:
return edge.v1
def get_sweep(self,v):
sweep = []
for p in self.vertices[v].sweep_p:
sweep.append(self.vertices[v].sweep_lu[p])
return sweep
def get_psweep(self,v):
sweep = []
for p in self.vertices[v].sweep_p:
sweep.append(self.vertices[v].sweep_plu[p])
return sweep
def extend_edge(self,v,e,plist):
external = [self.edges[e].v0,self.edges[e].v1]
if not self.complete_at_vertex(v,plist,external,1,[],[]):
print "unable to complete vertex",v,"in extend_edge"
def add_to_vertices(self,sweep,plist):
for vi in range(0,self.vc):
self.add_to_vertex(vi,sweep,plist)
def clean_sweeps(self):
for v in self.vertices:
new_e = []
vertex = self.vertices[v.num]
for e in vertex.sweep_e:
if e < self.ec:
new_e.append(e)
vertex.sweep_e = new_e
new_p = []
for p in vertex.sweep_p:
if p < self.pc:
new_p.append(p)
else:
if p in vertex.sweep_lu:
del vertex.sweep_lu[p]
vertex.sweep_p = new_p
def clean_sweeps2(self,current_vc,todo):
for v in self.vertices:
new_e = []
vertex = self.vertices[v.num]
for e in vertex.sweep_e:
if e < self.ec:
new_e.append(e)
vertex.sweep_e = new_e
new_p = []
for p in vertex.sweep_p:
if p < self.pc:
new_p.append(p)
else:
if p in vertex.sweep_lu:
del vertex.sweep_lu[p]
if len(vertex.sweep_p) != len(new_p):
todo.append(v.num)
vertex.sweep_p = new_p
todo2 = []
for i in todo:
if i < vc and i not in todo2:
todo2.append(i)
return todo2
def insert_slope_class(self,n,dx,dy,ec):
ec2 = []
found = False
for [dx1,dy1,nums] in ec:
if abs(dx1) < 0.0001 and abs(dx) < 0.0001:
ec2.append([dx1,dy1,nums + [n]])
found = True
elif abs(dx) >= 0.0001 and abs(dx1) >= 0.0001:
if abs((dy/dx)-(dy1/dx1)) < 0.0001:
ec2.append([dx1,dy1,nums + [n]])
found = True
else:
ec2.append([dx1,dy1,nums])
else:
ec2.append([dx1,dy1,nums])
if not found:
ec2.append([dx,dy,[n]])
return ec2
def get_angle_data(self,large_order,small_order,v0,vs):
if large_order % 2 == 0:
r = large_order/(2*pi)
else:
r = large_order/pi
n = large_order/small_order
v0 = self.vertices[v0]
v0 = [v0.x,v0.y]
va = v0
ec = []
angles = []
symclasses = []
i = 0
for v in vs[1:]:
dx = v[0] - va[0]
dy = v[1] - va[1]
a = r*atan2(dy,dx)
angles.append(a)
if a > 0:
s = int(a+0.5) % n
else:
s = int(a-0.5) % n
symclasses.append(s)
ec = self.insert_slope_class(i,dx,dy,ec)
va = v
i += 1
dx = v0[0] - va[0]
dy = v0[1] - va[1]
a = r*atan2(dy,dx)
angles.append(a)
if a > 0:
s = int(a+0.5) % n
else:
s = int(a-0.5) % n
symclasses.append(s)
ec = self.insert_slope_class(i,dx,dy,ec)
return [angles, symclasses,ec]
def get_data(self,p0,large_order,small_order):
e0 = self.polygons[p0].edges[0]
v0 = self.edges[e0].v0
es = self.get_edge_list(v0,e0,p0)
vs = self.get_vertex_list(v0,e0,p0)
[angles, symclasses,ec] = self.get_angle_data(large_order,small_order,v0,vs)
return [vs,es,symclasses]
def report(self,p0,large_order,small_order):
"""A version of get_data which prints out the information as well"""
print
print "****** Polygon report ******"
print "Polygon",p0,"has",len(self.polygons[p0].edges),"edges"
e0 = self.polygons[p0].edges[0]
v0 = self.edges[e0].v0
es = self.get_edge_list(v0,e0,p0)
print "edge list",es
vs = self.get_vertex_list(v0,e0,p0)
print "vertex list",vs
[angles, symclasses,ec] = self.get_angle_data(large_order,small_order,v0,vs)
print "with angles",angles
print "with symclasses",symclasses
print "and",len(ec),"slope classes",ec
return [vs,es,symclasses]
def merge(self,v0,e0,t1,v1,e1):
theta = self.get_merge_rotate(v0,e0,t1,v1,e1)
vxi,vyi = self.get_merge_translate(v0,e0,t1,v1,e1)
for p in t1.polygons:
self.merge_polygon(v0,t1,p,vxi,vyi,theta)
def get_merge_translate(self,v0,e0,t1,v1,e1):
v0x = self.vertices[v0].x
v0y = self.vertices[v0].y
v1x = t1.vertices[v1].x
v1y = t1.vertices[v1].y
return (v0x-v1x), (v0y-v1y)
def get_merge_rotate(self,v0,e0,t1,v1,e1):
vxi,vyi = self.get_merge_translate(v0,e0,t1,v1,e1)
# get the other vertices
v0other = self.get_other_vertex(v0,e0)
v1other = t1.get_other_vertex(v1,e1)
#translate the other vertices (mapping v0 onto (0,0))
v0x = self.vertices[v0other].x - self.vertices[v0].x
v0y = self.vertices[v0other].y - self.vertices[v0].y
v1x = t1.vertices[v1other].x + vxi - self.vertices[v0].x
v1y = t1.vertices[v1other].y + vyi - self.vertices[v0].y
#return the angle between the other vertices
return angle_between([v1x,v1y],[v0x,v0y])
def merge_polygon(self,v0,t1,p,vxi,vyi,theta):
# get the vertices and map them into the new coordinates
# we first use vxi, vyi to translate t2 coordinates into self
# coordinates, then rotate the new polygons around
# v0x, voy by theta radians
v0x = self.vertices[v0].x
v0y = self.vertices[v0].y
edges = []
for ei in p.edges:
v0 = t1.vertices[t1.edges[ei].v0]
v1 = t1.vertices[t1.edges[ei].v1]
x0,y0 = rotate([v0.x+vxi-v0x,v0.y+vyi-v0y],theta)
x1,y1 = rotate([v1.x+vxi-v0x,v1.y+vyi-v0y],theta)
v0n = self.new_vertex([x0+v0x,y0+v0y])
v1n = self.new_vertex([x1+v0x,y1+v0y])
edges.append(self.new_edge(v0n,v1n))
pedges = []
for x,y in list(t1.shapely_polygons[p.num].exterior.coords):
x0,y0 = rotate([x+vxi-v0x,y+vyi-v0y],theta)
pedges.append([x0+v0x,y0+v0y])
# we have the edges, so create the new polygon
p0 = self.new_polygon(edges[0],p.ptype,False)
for ei in edges[1:]:
self.add_to_polygon(p0,ei)
self.add_to_shapely_polygons2(pedges)
def toSVG(self,include_labels=False):
xml = """<?xml version="1.0" encoding="UTF-8" standalone="no"?>
<svg width="1500" height="1500"
xmlns:dc="http://purl.org/dc/elements/1.1/"
xmlns:cc="http://web.resource.org/cc/"
xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"
xmlns:svg="http://www.w3.org/2000/svg"
xmlns="http://www.w3.org/2000/svg">
>
<defs>
<style type="text/css">
<![CDATA[
text {
font-family: Arial, sans-serif;
font-size:2;
background-color:#ffff00;
fill:#ff0000;
}
]]>
</style>
</defs>
"""
for p in range(0,self.pc):
xml += self.polygon_to_svg(p)
if include_labels:
xml += self.labels_to_svg()
xml += "</svg>"
return xml
class Vertex(object):
"""Vertex class with public x, y, num and sweep attributes """
def __init__(self,x,y):
self.x = x
self.y = y
self.num = 0
self.sweep = []
self.sweep_e = []
self.sweep_p = []
self.sweep_lu = {}
self.sweep_plu = {}
self.parity = 0
class Edge(object):
"""Edge class with public v0, v1, num and internal attributes """
def __init__(self,v0,v1,internal=False):
self.v0 = v0
self.v1 = v1
self.num = 0
self.internal = False
class Polygon(object):
"""Polygon class with num, edges list and polygon/rhombus type """
def __init__(self,e0,ptype,make_internal=True):
self.edges = [e0]
self.ptype = ptype
self.num = 0
class LineSegment(object):
"""line segment is an object with four values x1, y1, x2,y2"""
def __init__(self,x1,y1,x2,y2):
self.x1 = x1
self.y1 = y1
self.x2 = x2
self.y2 = y2
def compute_polygon(vs,e0,n,parity=1):
if n < 0:
return compute_rhombus(vs,e0,n,parity)
if parity == 1:
v0 = e0[1]
v1 = e0[0]
else:
v0 = e0[0]
v1 = e0[1]
a = pi/n
#compute the edges
r = 1/(2*sin(a))
r_internal = r*cos(a)
# compute edge angle
ea = atan2(v1[1]-v0[1],v1[0]-v0[0])
phi = a - ea
x = 0.5
y = -r_internal
xc, yc = rotate([x,y],ea)
centre = [xc+v0[0], yc+v0[1]]
poly_edges = get_polygon_edges(r,centre,phi,n)
return [poly_edges,e0,a]
#compute a new rhombus given an edge and an angle
# and add it to the end of the polygon list
# rewrite and add correct sweep bit plus fix angles
def compute_rhombus(vs,e0,a,parity=1):
if parity == 1:
v0 = e0[0]
v1 = e0[1]
else:
v0 = e0[1]
v1 = e0[0]
# the following trick generates the correct rhombus shape based
# upon the vertex ordering (vertex vs should always have angle a)
if equal_vertex(vs,v0):
flip = True
a2 = a
else:
flip = True
a2 = pi - a
# compute edge angle
ea = atan2(v1[1]-v0[1],v1[0]-v0[0])
h = cos(0.5*(pi-a2))
w = sin(0.5*(pi-a2))
vlist = [[-w,0],[0,-h],[w,0],[0,h]]
rot = []
for v in vlist:
x0 = v[0]+w
y0 = v[1]
r = rotate([x0,y0],ea+(a2/2))
x1 = r[0]+v0[0]
y1 = r[1]+v0[1]
rot.append([x1,y1])
if flip:
pe = [[rot[1],rot[0]],[rot[2],rot[1]],[rot[3],rot[2]],[rot[0],rot[3]]]
else:
pe = [[rot[0],rot[1]],[rot[1],rot[2]],[rot[2],rot[3]],[rot[3],rot[0]]]
return [pe,e0,a]
def overlaps_polygons(pl,sp,ignore=False):
p0 = []
for e in pl:
p0.append(e[0])
poly = shapely.geometry.Polygon(p0)
i = 0
for p in sp:
try:
if poly.intersects(p):
pintersect = poly.intersection(p)
if (type(pintersect) is shapely.geometry.polygon.Polygon) and pintersect.area > 0.01:
return True
except:
print "overlap intersect function error"
return True
i += 1
return False
# compute a new regular polygon given an edge and an order
# and add it to the end of the polygon list
def rotate(v,theta):
xc = cos(theta) * v[0] - sin(theta) * v[1]
yc = sin(theta) * v[0] + cos(theta) * v[1]
return xc, yc
def get_polygon_edges(r,centre,phi,n):
pe = []
p0 = [centre[0]+r*sin(phi),centre[1]+ r*cos(phi)]
for i in range(1,n):
p1 = [centre[0]+r*sin(phi+(2*i*pi/n)), centre[1]+ r*cos(phi+(2*i*pi/n))]
pe.append([p0,p1])
p0 = p1
pe.append([p1,pe[0][0]])
return pe
def different_vertex(v0,v1):
return (abs(v0[0] - v1.x) > 0.0001) or (abs(v0[1] - v1.y) > 0.0001)
def equal_vertex(v0,v1):
return (abs(v0[0] - v1[0]) < 0.0001) and (abs(v0[1] - v1[1]) < 0.0001)
def different_edge(v0,v1,e):
return (v0 != e.v0 or v1 != e.v1) and (v0 != e.v1 or v1 != e.v0)
def angle_between(v0,v1):
return atan2(v0[0]*v1[1] - v0[1]*v1[0], v0[0]*v1[0]+v0[1]*v1[1])
# edges are equal if the vertices are equal in either order
def equal_poly_edges(e,et):
t0 = equal_vertex(e[0],et[0]) and equal_vertex(e[1],et[1])
t1 = equal_vertex(e[0],et[1]) and equal_vertex(e[1],et[0])
return t0 or t1
def convert_sweep(s):
s2 = []
for a in s:
if a > 0:
s2.append("p"+str(a))
else:
s2.append("r"+str(-180*a/pi))
return s2
def rh(n,m=2):
if n > 0:
return -m*pi/n
else:
n = -n
return (-(n-m)*pi/n)
def get_plist_from_sweep(sweep,sp):
lsp = len(sp)
for i in range(0,lsp):
sp0 = i
sp1 = i
spi = i
plist = []
success = True
for a in sweep:
if abs(sp[spi] - a) < 0.0001:
spi += 1
if spi >= lsp:
spi = 0
sp1 = spi
else:
success = False
break
if success:
break
if success:
if sp1 <= sp0:
return sp[sp1:sp0]
else:
return sp[sp1:]+sp[:sp0]
else:
return []
def equal_sweep(s0,s1):
if len(s0) == len(s1):
for i in range(0,len(s0)):
if abs(s0[i]-s1[i]) > 0.0001:
return False
return True
return False
# v - raw vertex
# et - raw root edge
# poly_edges - list of raw polygon edges
def get_other_edge(v,et,poly_edges):
for e in poly_edges:
if not equal_poly_edges(e,et) and (equal_vertex(e[0],v) or equal_vertex(e[1],v)):
return e
return []
### You will want to change the code after this point
# hendecagon tiling example
# You can change n below to generate a translational unit
# for a tiling containing any n-gon
print "Computing parts for translational unit."
print
print "Decomposing 6n polygon."
n = 11
t = Tiling()
t.start_polygon(6*n)
#t.decompose6_odd(n)
t.decompose_polygon(0,6,n,position_6)
#t.decompose_zonagon(0)
t.clip_all_worms()
decomp = copy.deepcopy(t)
print
print "Computing translational unit. This may take some time."
t = Tiling()
t.universal_tiler_tu(n,decomp)
print t.pc, "polygons"
xml = t.toSVG()
# You need to change where the results are saved.
fd = open("d:/tess/output.svg","w")
fd.write(xml)
fd.close()
print "Tiling output completed."
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment