# This file is part of pymadcad, distributed under license LGPL v3
''' This module provides the chamfer and bevel functions, and the associated tools.
`multicut`, `chamfer`, and `bevel` are built in the same cutting algorithm. It cuts the mesh faces by propagation from the given edges. The cutting planes are determined by an offset vector from the original primitive (point or edge).
Most of the time you don't need to set the offset yourself. It can be automatically calculated by several methods, depending on the shape you want to get. Those methods are called `cutters` and are executed using `planeoffsets`_.
'''
from .mathutils import *
from .mesh import *
from . import generation as gt
from . import hashing
from . import text
from . import settings
from .triangulation import triangulation_outline
from .blending import blenditer, match_length
from numbers import Number
from functools import singledispatch
__all__ = [ 'chamfer', 'bevel', 'multicut',
'mesh_cut', 'web_cut', 'planeoffsets',
'tangentend', 'tangentcorner', 'tangentjunction',
'cutter_width', 'cutter_distance', 'cutter_depth', 'cutter_radius',
]
# ---- common cut methods ----
[docs]@singledispatch
def multicut(mesh, indices, cutter):
''' Cut a Mesh/Web/Wire around the given edges/points, using the given cutter '''
raise TypeError('wrong argument type: {}'.format(type(mesh)))
[docs]@singledispatch
def chamfer(mesh, indices, cutter):
''' Chamfer a Mesh/Web/Wire around the given edges/points, using the given cutter '''
raise TypeError('wrong argument type: {}'.format(type(mesh)))
[docs]@singledispatch
def bevel(mesh, indices, cutter):
''' Bevel a Mesh/Web/Wire around the given edges/points, using the given cutter '''
raise TypeError('wrong argument type: {}'.format(type(mesh)))
# ---- cut methods -----
[docs]def cutter_width(width, fn1, fn2):
''' Plane offset for a cut based on the width of the bevel '''
n = normalize(fn1+fn2)
s = dot(fn1,n)
return -width/2 * sqrt(1/s**2 - 1) * n
[docs]def cutter_distance(depth, fn1, fn2):
"""Plane offset for a cut based on the distance along the side faces"""
n = normalize(fn1 +fn2)
cos_b = dot(fn1, n)
cos_a = sqrt(1-cos_b**2)
return -depth * n * cos_a
[docs]def cutter_depth(dist, fn1, fn2):
''' Plane offset for a cut based on the distance to the cutted edge '''
return -dist * cross(normalize(cross(fn1,fn2)), fn1-fn2)
[docs]def cutter_radius(depth, fn1, fn2):
''' Plane offset for a cut based on the angle between faces '''
n = normalize(fn1 + fn2)
s = dot(fn1,n)
return -depth * (1/s - s) * n
# ----- mesh operations ------
[docs]def planeoffsets(mesh, edges, cutter):
''' Compute the offsets for cutting planes using the given method
cutter is a tuple or a function
- function(fn1,fn2) -> offset
fn1, fn2 are the adjacents face normals
offset is the distance from segment to plane times the normal to the plane
- ('method', distance)
the method is the string name of the method (a function named 'cutter_'+method existing in this module)
distance depends on the method and is the numeric parameter of the method
'''
if isinstance(cutter, dict):
return cutter
else:
cutter = interpretcutter(cutter)
# get adjacent faces' normals to lines
offsets = {e: [None,None] for e in edges}
for f in mesh.faces:
for e in ((f[0],f[1]), (f[1],f[2]), (f[2],f[0])):
if e in offsets: offsets[e][0] = mesh.facenormal(f)
elif (e[1],e[0]) in offsets: offsets[(e[1],e[0])][1] = mesh.facenormal(f)
# compute offsets (displacement on depth)
for e, adjacents in offsets.items():
offset = cutter(*adjacents)
if dot(cross(*adjacents), mesh.points[e[0]]-mesh.points[e[1]]) > 0:
offset = -offset
offsets[e] = offset
return offsets
def interpretcutter(cutter):
if isinstance(cutter, tuple):
func = globals()['cutter_'+cutter[0]]
arg = cutter[1]
return lambda fn1,fn2: func(arg, fn1, fn2)
elif callable(cutter):
return cutter
else:
raise TypeError("cutter must be a callable or a tuple (name, param)")
[docs]def mesh_cut(mesh, start, cutplane, stops, conn, prec, removal, cutghost=True):
''' Propagation cut for an edge
:start: the edge or point to start propagation from
:cutplane: the plane cutting the faces. Its normal must be oriented toward the propagation area.
:stops: the planes stopping the propagation. Their normal must be oriented toward the propagation area.
:removal: the set in wich the function will put the indices of faces inside
:cutghost: whether the function should propagate on faces already marked for removal (previously or during the propagation)
'''
pts = mesh.points
ptmgr = hashing.PointSet(prec, manage=pts)
stops = list(filter(lambda e:e, stops))
# find the intersections axis between the planes
axis = [intersection_plane_plane(cutplane, stop) for stop in stops]
# prepare propagation
if isinstance(start, tuple):
front = [start, (start[1],start[0])]
elif isinstance(start, Number):
front = [e for e in conn if start in e]
else:
raise TypeError('wrong start: {}, cut can only start from a point or an edge'.format(start))
intersections = set()
seen = set()
while front:
# propagation affairs
frontedge = front.pop()
if frontedge not in conn: continue
fi = conn[frontedge]
if fi in seen: continue
f = mesh.faces[fi]
# do not process empty faces (its a topology singularity), but propagate
if faceheight(mesh, fi) <= prec:
seen.add(fi)
for j in range(3):
front.append((f[j],f[j-1]))
continue
# find the intersection of the triangle with the common axis to the cutplane and the stop plane
fpts = mesh.facepoints(fi)
p = None
for a in axis:
# if the intersection is on a triangle's edge, p is None
p = intersection_axis_face(a, fpts, prec)
if p: break
if p:
# mark cutplane change
# empty faces can be generated, but they are necessary to keep the connectivity working
pi = ptmgr.add(p)
l = len(mesh.faces)
mesh.faces[fi] = (f[0], f[1], pi)
mesh.faces.append((f[1], f[2], pi))
mesh.faces.append((f[2], f[0], pi))
mesh.tracks.append(mesh.tracks[fi])
mesh.tracks.append(mesh.tracks[fi])
registerface(mesh, conn, fi)
registerface(mesh, conn, l)
registerface(mesh, conn, l+1)
front.append(frontedge)
continue
# mark this face as processed
seen.add(fi)
# point side for propagation
# True if the point is over or on the cutplane
goodside = [False]*3
for j,pi in enumerate(f):
goodside[j] = dot(pts[pi]-cutplane[0], cutplane[1]) >= -prec
# abort conditions
if not any(goodside):
continue
outside = any( all( dot(pts[pi]-stop[0], stop[1]) <= -prec
for pi in f)
for stop in stops)
if outside:
continue
# True if the point is over or on for all stop planes
goodx = [False]*3
for j,pi in enumerate(f):
goodx[j] = all(dot(pts[pi]-stop[0], stop[1]) >= -prec for stop in stops)
# intersections of triangle's edges with the plane
cut = [None]*3
for j in range(3):
cut[j] = intersection_edge_plane((pts[f[j]], pts[f[j-2]]), cutplane, prec)
# don't intersect if the 2 points are at the corner of the triangle
for j in range(3):
if cut[j-1] and cut[j] and distance(cut[j-1], cut[j]) <= prec:
cut[j] = None
if all(goodside) and any(goodx):
removal.add(fi)
if fi not in removal or cutghost: # NOTE: this is a ugly trick to make the current multicut implementation work with the same function for corners and edges
# cut the face
for j in range(3):
if cut[j] and cut[j-1]:
# cut only if the intersection segment is in the delimited area
if any( dot(cut[j-1]-stop[0], stop[1]) <= prec
and dot(cut[j] -stop[0], stop[1]) <= prec
for stop in stops):
continue
# cut the face (create face even for non kept side, necessary for propagation)
p1 = ptmgr.add(cut[j])
p2 = ptmgr.add(cut[j-1])
l = len(mesh.faces)
mesh.faces[fi] = (p1, f[j-2], f[j-1])
mesh.faces.append((p1, f[j-1], p2))
mesh.faces.append((p1, p2, f[j-0]))
mesh.tracks.append(mesh.tracks[fi])
mesh.tracks.append(mesh.tracks[fi])
registerface(mesh, conn, fi)
registerface(mesh, conn, l)
registerface(mesh, conn, l+1)
seen.update((fi, l, l+1))
# mark to remove the faces outside
if dot(pts[f[j]]-cutplane[0], cutplane[1]) < 0:
removal.add(fi)
removal.add(l)
intersections.add((p2,p1))
else:
removal.add(l+1)
removal.discard(fi)
intersections.add((p1,p2))
break
# propagate
for j in range(3):
front.append((f[j],f[j-1]))
return intersections
def propagate(mesh, edges, conn, prec):
front = list(edges)
edges = set(edges)
seen = set()
count = 0
while front:
edge = front.pop()
if edge not in conn: continue
fi = conn[edge]
if fi in seen: continue
seen.add(fi)
count += 1
if faceheight(mesh, fi) <= prec: continue
a,b,c = arrangeface(mesh.faces[fi], edge[0])
if c in edge: continue
if (b,c) not in edges: front.append((c,b))
if (c,a) not in edges: front.append((a,c))
return seen
@multicut.register(Mesh)
def mesh_multicut(mesh, edges, cutter, conn=None, prec=None, removal=None):
''' General purpose edge cutting function.
cut the given edges and the crossing corners, resolving the interference issues.
cutter can either be a argument for planeoffsets or more directly a dict `{edgekey: offset vector}` for each passed edge.
edges must be unoriented.
'''
# perpare working objects
if conn is None: conn = connef(mesh.faces)
if prec is None: prec = mesh.precision()
if removal is None:
removal = set()
final = True
else:
final = False
if isinstance(edges, Web): edges = edges.edges
edges = [edgekey(*e) for e in edges]
offsets = planeoffsets(mesh, edges, cutter)
normals = mesh.vertexnormals()
juncconn = connpp(edges)
pts = mesh.points
debmesh = mesh
separators = {} # plans separateurs pour chaque arete
for i, prox in juncconn.items():
if len(prox) != 2: continue
o0 = offsets[edgekey(i,prox[0])]
o1 = offsets[edgekey(i,prox[1])]
axis = intersection_plane_plane((pts[prox[0]]+o0, o0), (pts[prox[1]]+o1, o1), pts[i]) # fail when the two planes have the same orientation
if axis:
n = normalize(cross(pts[i]-axis[0], axis[1]))
if dot(n, pts[prox[0]]-pts[i]) < 0: n = -n
separators[(i,prox[0])] = (axis[0], n)
separators[(i,prox[1])] = (axis[0],-n)
## display cut planes in the mesh (debug)
#grp = len(debmesh.groups)
#debmesh.groups.append(None)
#m = len(debmesh.points)
#debmesh.points.append(pts[i])
#debmesh.points.append(axis[0]+axis[1])
#debmesh.points.append(axis[0]-axis[1])
#debmesh.faces.append((m,m+1,m+2))
#debmesh.tracks.append(grp)
else:
separators[(i,prox[0])] = None
separators[(i,prox[1])] = None
# fix the missing separators in case of straight lines
for e,sep in separators.items():
i,j = e
while not sep:
prox = juncconn[i]
if len(prox) != 2: raise Exception('unable to compute all separators')
j,i = i, prox[0] if prox[0] != j else prox[1]
separators[e] = sep = separators[(i,j)]
corners = {} # offset pour chaque coin
# pour chaque jonction
for junc, prox in juncconn.items():
if len(prox) < 3: continue
# calculer le plan de distance max demandé par les aretes adjacentes
offset = offsets.get(junc, 0.)
# pour chaque arete adjacente (dans chaque sens)
for b in prox:
for c in prox:
a = junc
if b == c or b == a or c == a: continue
# calculer la distance
ab = normalize(pts[b]-pts[a])
ac = normalize(pts[c]-pts[a])
n = cross(ab, ac)
nb = cross(n, ab)
nc = cross(ac, n)
ob = offsets[edgekey(a,b)]
oc = offsets[edgekey(a,c)]
ib = dot(ob,ob) / dot(ob, nb)
ic = dot(oc,oc) / dot(oc, nc)
s = dot(nb,ac)
if not s: continue
abl = (ib*dot(ab,ac) + ic) / s
offset = dot(abl*ab + ib*nb, normals[a])
# assignation du plan de coupe
corners[junc] = offset*normals[junc]
plane = (pts[junc] + corners[junc], -normals[junc])
for p in prox:
separators[(junc,p)] = plane
outlines = {}
# couper les aretes
for edge, offset in offsets.items():
outlines[edge] = mesh_cut(mesh, edge,
(pts[edge[0]]+offset, -normalize(offset)),
(separators.get(edge), separators.get((edge[1], edge[0]))),
conn, prec, removal, True)
# couper les sommets
for corner,offset in corners.items():
outlines[corner] = mesh_cut(mesh, corner,
(pts[corner]+offset, -normalize(offset)),
(),
conn, prec, removal, False)
if final:
frontier = []
for edges in outlines.values():
frontier.extend(edges)
finalize(mesh, outlines, propagate(mesh, frontier, conn, prec), prec)
#finalize(mesh, outlines, removal, prec) # NOTE this is faster but also leads to not removing all faces
return outlines
def finalize(mesh, outlines, removal, prec):
''' Function finalizing the mesh for multicut
removing faces and simplifying outlines
'''
# simplify cuts
merges = {}
for e,cuts in outlines.items():
reindex = line_simplification(Web(mesh.points, cuts), prec)
lines = []
for a,b in cuts:
c,d = reindex.get(a,a), reindex.get(b,b)
if c != d: lines.append((c,d))
outlines[e] = lines
merges.update(reindex)
# apply merges, but keeping the empty faces, to keep the removal list valid
for i,f in enumerate(mesh.faces):
mesh.faces[i] = (
merges.get(f[0],f[0]),
merges.get(f[1],f[1]),
merges.get(f[2],f[2]),
)
# delete inside faces and empty ones
removefaces(mesh, lambda fi: fi in removal or faceheight(mesh,fi) <= prec)
def arrangeface(f, p):
if p == f[1]: return f[1],f[2],f[0]
elif p == f[2]: return f[2],f[0],f[1]
else: return f
def registerface(mesh, conn, fi):
f = mesh.faces[fi]
if f[0] == f[1] or f[1] == f[2] or f[2] == f[0]:
return
for e in ((f[0],f[1]), (f[1],f[2]), (f[2],f[0])):
conn[e] = fi
def unregisterface(mesh, conn, fi):
f = mesh.faces[fi]
for e in ((f[0],f[1]), (f[1],f[2]), (f[2],f[0])):
if e in conn: del conn[e]
def segmentsdict(line):
segments = {}
for i in range(len(line)-1):
segments[(line[i],line[i+1])] = i
return segments
def intersection_plane_plane(p0, p1, neigh=None):
''' Return the intersection axis between two planes '''
if not neigh: neigh = p0[0]
d = cross(p0[1], p1[1])
if length(d) <= NUMPREC: return None
return vec3(
inverse(transpose(dmat3(p0[1], p1[1], d)))
* dvec3(dot(p0[0],p0[1]), dot(p1[0],p1[1]), dot(neigh,d))
), normalize(d)
def intersection_edge_plane(edge, axis, prec):
''' Return the intersection point of an edge with a plane, or None if it doesn't exist
In case of an intersection at an extremity of the edge, return that edge point
'''
a, b = edge
da = dot(axis[0]-a, axis[1])
db = dot(axis[0]-b, axis[1])
if abs(da) <= prec: return a
if abs(db) <= prec: return b
if da * db > 0: return None # the segment doesn't reach the plane
edgedir = b-a
proj = dot(edgedir, axis[1])
if abs(proj) <= prec: return None # the segment is parallel to the plane
return a + da * edgedir / proj
def intersection_axis_face(axis, face, prec):
''' Intersection between an axis and a triangle
In case of intersection with an edge of the triangle, return None
'''
n = cross(face[1]-face[0], face[2]-face[0])
unp = dot(n, axis[1])
if abs(unp) <= prec: return None
p = axis[0] + dot(face[0]-axis[0], n) / unp * axis[1]
for i in range(3):
if dot(p - face[i], normalize(cross(n,face[i]-face[i-1]))) <= prec:
return None
return p
def removefaces(mesh, crit):
''' Remove faces whose indices are present in faces, (for huge amount, prefer pass faces as a set) '''
newfaces = typedlist(dtype=uvec3)
newtracks = typedlist(dtype='I')
for i in range(len(mesh.faces)):
if not crit(i):
newfaces.append(mesh.faces[i])
newtracks.append(mesh.tracks[i])
mesh.faces = newfaces
mesh.tracks = newtracks
def removeedges(mesh, crit):
''' Remove faces whose indices are present in faces, (for huge amount, prefer pass faces as a set) '''
newedges = typedlist(dtype=uvec2)
newtracks = typedlist(dtype='I')
for i in range(len(mesh.edges)):
if not crit(i):
newedges.append(mesh.edges[i])
newtracks.append(mesh.tracks[i])
mesh.edges = newedges
mesh.tracks = newtracks
#def facesurf(mesh, fi):
#o,x,y = mesh.facepoints(fi)
#return length(cross(x-o, y-o))
#def faceangle(mesh, fi):
#o,x,y = mesh.facepoints(fi)
#return length(cross(normalize(x-o), normalize(y-o)))
def faceheight(mesh, fi):
f = mesh.facepoints(fi)
m = inf
for i in range(3):
a,b = f[i]-f[i-1], f[i]-f[i-2]
if not dot(b,b): return 0
m = min(m, dot(a,a) - dot(a,b)**2/dot(b,b))
return sqrt(max(0,m))
def faceheight(mesh, fi):
f = mesh.facepoints(fi)
m = inf
for i in range(3):
l = length(f[i-2]-f[i])
if not l: return 0
h = length(cross(f[i-2]-f[i], f[i-1]-f[i])) / l
if h < m: m = h
return m
@chamfer.register(Mesh)
def mesh_chamfer(mesh, edges, cutter):
''' Create a chamfer on the given suite of points, create faces are planes.
cutter is described in function planeoffsets()
'''
# cut faces
segments = mesh_multicut(mesh, edges, cutter)
# create junctions
group = len(mesh.groups)
mesh.groups.append(None)
pts = mesh.points
# edges for corners surfaces
corners = {}
def cornerreg(c, edge):
if c not in corners: corners[c] = []
corners[c].append(edge)
for e,s in segments.items():
if not s: continue
# corner cuts
if isinstance(e, Number):
if e not in corners: corners[e] = []
corners[e].extend(s)
# edge cuts
else:
# assemble the intersection segments in a single outline
frags = suites(s)
if len(frags) == 1:
lp = frags[0]
cornerreg(e[0], (lp[0],lp[-1]))
elif len(frags) == 2:
lp = frags[0] + frags[1]
# get the sides in the edge direction
l,r = frags
if dot(pts[l[-1]]-pts[l[0]], pts[e[1]]-pts[e[0]]) < 0 or dot(pts[r[-1]]-pts[r[0]], pts[e[1]]-pts[e[0]]) > 0:
l,r = r,l
cornerreg(e[0], (l[0], r[-1]))
cornerreg(e[1], (r[0], l[-1]))
else:
raise Exception('a cutted edge has more than 2 cutted sides')
# triangulate cutted edge
faces = triangulation_outline(Wire(mesh.points, lp)).faces
mesh.faces.extend(faces)
mesh.tracks.extend([group]*len(faces))
for c,s in corners.items():
# triangulate cutted corner
faces = triangulation_outline(Wire(mesh.points, suites(s)[0])).faces
mesh.faces.extend(faces)
mesh.tracks.extend([group]*len(faces))
@bevel.register(Mesh)
def mesh_bevel(mesh, edges, cutter, resolution=None):
''' Create a chamfer on the given suite of points, create faces are planes.
cutter is described in function planeoffsets()
'''
if isinstance(edges, Web): edges = edges.edges
edges = [edgekey(*e) for e in edges]
conn = connef(mesh.faces)
enormals = {}
for e in edges:
enormals[e] = (
mesh.facenormal(conn[e]),
mesh.facenormal(conn[(e[1],e[0])]),
)
# cut faces
segments = mesh_multicut(mesh, edges, cutter, conn)
mesh.check()
conn = connef(mesh.faces)
pts = mesh.points
# edges for corners surfaces
corners = {}
def cornerreg(c, edge):
if c not in corners: corners[c] = []
corners[c].append(edge)
# create junctions
normals = {} # face normal at each point
junctions = []
ends = []
for e,s in segments.items():
if not s: continue
# corner cuts
if isinstance(e, Number):
if e not in corners: corners[e] = []
corners[e].extend(s)
# edge cuts
else:
# assemble the intersection segments in a single outline
frags = suites(s)
# one loop: bevel extremity
if len(frags) == 1:
lp = frags[0]
# find a separation to put a junction between 2 sides
side = noproject(pts[lp[0]]-pts[lp[-1]], normalize(pts[e[1]] - pts[e[0]]))
for i in range(1,len(lp)):
if dot(pts[lp[i-1]]-pts[e[0]], side) * dot(pts[lp[i]]-pts[e[0]], side) < 0:
l,r = lp[:i], lp[i:]
break
cornerreg(e[0], (lp[0],lp[-1]))
if (lp[i], lp[i-1]) in conn:
ends.append((lp[i-1], lp[i]))
# two parts: cutted edge
elif len(frags) == 2:
l,r = frags
if dot(pts[l[-1]]-pts[l[0]], pts[e[1]]-pts[e[0]]) < 0 or dot(pts[r[-1]]-pts[r[0]], pts[e[1]]-pts[e[0]]) > 0:
l,r = r,l
cornerreg(e[0], (l[0], r[-1]))
cornerreg(e[1], (r[0], l[-1]))
else:
raise Exception('a cutted edge has more than 2 cutted sides')
# match the curves
r.reverse()
match = list(match_length(Wire(pts,l), Wire(pts,r)))
# prepare tangents
ll = Wire(pts,l).length()
x = 0.
for i in range(len(match)):
if i: x += distance(pts[match[i-1][0]], pts[match[i][0]]) / ll
# get the tangents from the triangle between the cuts and the cutting edge
l,r = match[i]
o = interpol1(pts[e[1]], pts[e[0]], x)
plane = cross(pts[r]-o, pts[l]-o)
normals[l] = enormals[e][1] + normals.get(l, 0)
normals[r] = enormals[e][0] + normals.get(r, 0)
junctions.append(match)
# normals neighbooring corners
for e,s in corners.items():
for a,b in s:
if (b,a) in conn:
n = mesh.facenormal(conn[(b,a)])
normals[b] = normals.get(b, 0) + n
normals[a] = normals.get(a, 0) + n
# normals neighbooring cut ends
for a,b in ends:
n = mesh.facenormal(conn[(b,a)])
normals[a] = noproject(normals[a], n)
normals[b] = noproject(normals[b], n)
# normalize all contributions
for k,n in normals.items():
normals[k] = normalize(n)
# compute resolution based on adjacent faces to edges
div = 0
c = 0
for match in junctions:
for a,b in match:
div += settings.curve_resolution(
arclength(a,b,normals[a],normals[b]),
anglebt(normals[a], normals[b]),
resolution)
c += 1
div //= c
new = Mesh()
# round cutted corner
for e,s in corners.items():
if len(s) < 3: continue
# assemble the intersection segments in a single outline
lp = suites(s)[0]
if lp[-1] == lp[0]: lp.pop() # remove the redundant point if there is (not all the time due to some BUG)
new += tangentcorner(pts, lp, normals, div)
# fill gap between existing surface and round extremities
for edge in ends:
n = mesh.facenormal(conn[(edge[1],edge[0])])
new += tangentend(pts, edge, normals, div)
# round cutted edges
for match in junctions:
# put a junction
new += tangentjunction(pts, match, normals, div)
new.groups = [None]
new.tracks = typedlist.full(0, len(new.faces), 'I')
mesh += new
mesh.mergeclose()
# ---- web operations -----
from .nprint import nprint
[docs]def web_cut(web, start, cutplane, conn, prec, removal):
''' Propagation cut for a point
:start: the edge or point to start propagation from
:cutplane: the plane cutting the faces. Its normal must be oriented toward the propagation area.
'''
intersections = []
# propagate from point to point
seen = set()
front = list(conn[start])
while front:
ei = front.pop()
if ei in seen: continue
seen.add(ei)
e = web.edges[ei]
# look intersection with the cutplane
ep = web.edgepoints(e)
side = [dot(p-cutplane[0], cutplane[1]) for p in ep]
if side[0]*side[1] <= prec:
p = intersection_edge_plane(ep, cutplane, prec)
pi = len(web.points)
web.points.append(p)
intersections.append(pi)
if side[0] <= prec:
if distance2(p,ep[0]) < prec**2:
pi = e[0]
web.edges[ei] = (pi, e[1])
web.edges.append((e[0], pi))
web.tracks.append(web.tracks[ei])
conn.discard(e[0], ei)
conn.add(e[0], len(web.edges)-1)
else:
if distance2(p,ep[1]) < prec**2:
pi = e[1]
web.edges[ei] = (e[0], pi)
web.edges.append((pi, e[1]))
web.tracks.append(web.tracks[ei])
conn.discard(e[1], ei)
conn.add(e[1], len(web.edges)-1)
else:
# propagate
for p in e:
front.extend(conn[p])
removal.update(seen)
return intersections
@multicut.register(Web)
def web_multicut(web, points, cutter, conn=None, prec=None, removal=None):
if conn is None: conn = connpe(web.edges)
if prec is None: prec = web.precision()
if isinstance(points, Wire): points = points.indices
removal = set()
cutter = interpretcutter(cutter)
pts = web.points
intersections = {}
for pi in points:
dirs = []
for ni in conn[pi]:
e = web.edges[ni]
dirs.append(normalize( pts[pi] - pts[e[e[0] == pi]] ))
if not dirs: # the cut has removed everything connected
continue
normal = normalize(sum(dirs))
c = sum(dot(normal,d) for d in dirs) / len(dirs)
s = sqrt(max(0, 1-c**2))
x,_,_ = dirbase(normal)
offset = cutter(c*x + s*normal, -c*x + s*normal)
if dot(offset, normal) > 0: offset = -offset
plane = (pts[pi]+offset, -normalize(offset))
intersections[pi] = web_cut(web, pi, plane, conn, prec, removal)
removeedges(web, lambda i: i in removal)
return intersections
@chamfer.register(Web)
def web_chamfer(web, points, cutter):
holes = web_multicut(web, points, cutter)
g = len(web.groups)
web.groups.append(None)
conn = connpe(web.edges)
def link(e):
if web.edges[next(conn[e[0]])][0] == e[0]: return (e[1],e[0])
else: return e
for cuts in holes.values():
if len(cuts) == 2:
web.edges.append(link(tuple(cuts)))
web.tracks.append(g)
else:
pi = len(web.points)
web.points.append(sum(web.points[c] for c in cuts) / len(cuts))
web.edges.extend( link((c,pi)) for c in cuts )
web.tracks.extend( [g] * len(cuts) )
@bevel.register(Web)
def web_bevel(obj, points, cutter, resolution=None):
holes = web_multicut(obj, points, cutter)
conn = connpe(obj.edges)
pts = obj.points
div = 6
g = len(obj.groups)
obj.groups.append(None)
def tangentat(p):
''' Find the tangent axis to the given cut point,
the third element is whether the original edge direction is opposite to the given tangent
'''
e = obj.edges[next(conn[p])]
if e[0] == p: return pts[p], normalize(pts[p] - pts[e[1]]), True
else: return pts[p], normalize(pts[p] - pts[e[0]]), False
for cuts in holes.values():
if len(cuts) == 2:
# place an arc between the two ends
s0 = tangentat(cuts[0])
s1 = tangentat(cuts[1])
corner = web(tangentarc(s0[:2], s1[:2], resolution))
if s0[2]: corner = corner.flip()
obj += corner
else:
# find an approximate sphere center
tangents = [tangentat(i) for i in cuts]
normal = normalize(sum(t[1] for t in tangents))
med = sum(pts[i] for i in cuts) / len(cuts)
center = sum(pts[i] + unproject(med-pts[i], normal) for i in cuts)
radius = sum(distance(pts[i], center) for i in cuts)
# place summit
pi = len(pts)
top = center + normal*radius
pts.append(top)
# place arcs
for t in tangents:
u = (top, noproject(t[0]-top, normal))
corner = web(tangentarc(t[:2], u, resolution))
if t[2]: corner = corner.flip()
obj += corner
# ---- wire operations -----
@multicut.register(Wire)
def wire_multicut(wire: Wire, points, cutter):
if isinstance(points, Wire):
points = points.indices
if not isinstance(points, set):
points = set(points)
prec = wire.precision()
cutter = interpretcutter(cutter)
if not wire.tracks:
wire.tracks = typedlist.full(0, len(wire.indices), "I")
g = len(wire.groups)
wire.groups.append(None)
cuts = []
closed = wire.indices[0] == wire.indices[-1]
for index, origin in enumerate(wire.indices):
if origin not in points:
continue
# get point location in the wire
if not closed:
if origin == wire.indices[0] or origin == wire.indices[-1]:
raise MeshError("a chamfer/bevel cannot have only one side")
# compute cut plane
ref = wire[index-1] - wire[index]
p0, p1 = wire[index:index+2]
# handle special case of closed wire
prior_p = wire[index-1]
if prior_p == p0:
prior_p = wire[index-2]
t0, t1 = normalize(p0 - prior_p), normalize(p1 - p0)
axis = cross(t0, t1)
offset = cutter(normalize(cross(axis, t0)), normalize(cross(axis, t1)))
if dot(offset, t0) > 0:
offset = -offset
cutplane = (p0 + offset, -normalize(offset))
def find_intersection(it):
for ii0, ii1 in it:
if ii0 == -1: # handle closed wires
ii0 = -2
p0, p1 = wire[ii0], wire[ii1]
ps = intersection_edge_plane((p0, p1), cutplane, prec)
if ps:
if isfinite(ps):
return ii1, ps
raise ValueError("no intersection found")
# propagate backward
back_range = range(index, -closed, -1)
back_it = zip(map(lambda x: x-1, back_range), back_range)
iiback, ps = find_intersection(back_it)
cuts.append((iiback, iiback + 1)) # point interval including start and excluding end
# propagate forward
wlen = len(wire.indices)
forward_range = range(index, wlen-1)
forward_it = zip(forward_range, map(lambda x: x+1, forward_range))
iiforward, pe = find_intersection(forward_it)
# update wire
m = len(wire.points)
wire.points.append(ps)
wire.points.append(pe)
wire.indices[iiback: iiforward] = [m, m + 1]
wire.tracks[iiback: iiforward] = [g, wire.tracks[iiforward-1]]
if origin == wire.indices[-1]:
wire.indices[-1] = wire.indices[0]
wire.tracks[-1] = wire.tracks[0]
return cuts
@chamfer.register(Wire)
def wire_chamfer(wire, points, cutter):
wire_multicut(wire, points, cutter)
@bevel.register(Wire)
def wire_bevel(wire, points, cutter, resolution=None):
closed = wire.indices[0] == wire.indices[-1]
cuts = wire_multicut(wire, points, cutter)
g = len(wire.groups) - 1
wire.groups[g] = None
cuts.reverse()
for ii0, ii1 in cuts:
p0 = wire[ii0]
p1 = wire[ii1]
tpii0 = ii0 - 1 # indices index if prior point
tpii1 = ii1 + 1 # indices index if prior point
if tpii0 < 0:
tpii0 = -2
t0 = normalize(p0 - wire[tpii0])
t1 = normalize(p1 - wire[tpii1])
old_l = len(wire.points)
wire.points.extend( tangentarc((p0,t0), (p1,t1), resolution))
wire.indices[ii0:ii0+2] = range(old_l, len(wire.points))
wire.tracks[ii0:ii0+1] = [g] * (len(wire.points)-old_l-1)
if closed and ii0==0:
wire.indices[-1] == wire.indices[0]
# --- mesh tangent generation functions ----
def tangentarc(s0, s1, resolution=None):
''' Create a tangent curve to both given axis.
return the curve as a list of points
'''
p0, t0 = s0
p1, t1 = s1
z = cross(t0,t1)
factor = arclength(p0, p1, cross(z,t0), cross(t1,z))
s0 = (p0, factor*t0)
s1 = (p1, factor*t1)
div = settings.curve_resolution(factor, anglebt(t0,-t1), resolution)
return [ interpol2(s0, s1, j/(div+1)) for j in range(div+2) ]
[docs]def tangentend(points, edge, normals, div):
''' Join a tangent surface resulting of `tangentcorner` or `tangentjunction` to a straight edge e
`normals` is the same dict as for tangentcorner and tangentjunction
'''
l, r = edge
pl, pr = points[l], points[r]
dist = arclength(points[l], points[r], normals[l], normals[r])
e = points[r]-points[l]
tl = normalize(noproject( e, normals[l])) * dist
tr = normalize(noproject(-e, normals[r])) * dist
def infos():
for i in range(div+2):
x = i/(div+1)
yield interpol2((pl,tl), (pr,tr), x), interpol1(pl,pr, x)
return blenditer(infos(), 0, interpol1)
[docs]def tangentcorner(pts, lp, normals, div):
''' Create a rounded surface tangent to the loop given
`normals` is a dict {point: normal}
'''
new = Mesh()
# compute the common points to all triangular regions
n = normalize(sum(normals[i] for i in lp))
med = sum(pts[i] for i in lp) / len(lp)
center = sum(pts[i] + unproject(med-pts[i], normals[i]) for i in lp) / len(lp)
radius = sum(distance(pts[i], center) for i in lp) / len(lp)
# place a central point
p = len(pts)
pts.append(center + n*radius)
normals[p] = n
# triangulate
for i in range(len(lp)):
a,b,c = lp[i-1], lp[i], p
new += gt.icosurface(
(pts[a],pts[b],pts[c]),
(normals[a],normals[b],normals[c]),
resolution=('div',div))
return new
[docs]def tangentjunction(points, match, normals, div):
''' Create a surface joining the given couples of points, tangent to the two sides
`normals` is a dict {point: normal}
'''
def infos():
# estimate normals
for l,r in match:
dist = arclength(points[l], points[r], normals[l], normals[r])
e = points[r]-points[l]
tl = normalize(noproject( e, normals[l])) * dist
tr = normalize(noproject(-e, normals[r])) * dist
yield (points[r], tr), (points[l], tl)
return blenditer(infos(), div, interpol2)