inkstitch/lib/stitches/contour_fill.py

590 wiersze
21 KiB
Python

from collections import namedtuple
from itertools import chain
import networkx as nx
import numpy as np
import trimesh
from shapely import offset_curve
from shapely.geometry import (GeometryCollection, LineString, MultiPolygon,
Point, Polygon)
from shapely.geometry.polygon import orient
from shapely.ops import nearest_points, polygonize
from shapely.validation import make_valid
from ..stitch_plan import Stitch
from ..utils import DotDict
from ..utils.clamp_path import clamp_path_to_polygon
from ..utils.geometry import (cut, ensure_geometry_collection,
ensure_multi_line_string, reverse_line_string,
roll_linear_ring)
from ..utils.smoothing import smooth_path
from ..utils.threading import check_stop_flag
from .running_stitch import running_stitch
class Tree(nx.DiGraph):
# This lets us do tree.nodes['somenode'].parent instead of the default
# tree.nodes['somenode']['parent'].
node_attr_dict_factory = DotDict
def __init__(self, *args, **kwargs):
self.__node_num = 0
super().__init__(**kwargs)
def generate_node_name(self):
node = self.__node_num
self.__node_num += 1
return node
nearest_neighbor_tuple = namedtuple(
"nearest_neighbor_tuple",
[
"nearest_point_parent",
"nearest_point_child",
"proj_distance_parent",
"child_node",
],
)
def _offset_linear_ring(ring, offset, resolution, join_style, mitre_limit):
ring = Polygon(ring)
result = offset_curve(ring, -offset, resolution, join_style=join_style, mitre_limit=mitre_limit)
result = ensure_multi_line_string(result)
rings = result.simplify(0.01, False)
return _take_only_valid_linear_rings(rings)
def _take_only_valid_linear_rings(rings):
"""
Removes all geometries which do not form a "valid" LinearRing.
A "valid" ring is one that does not form a straight line.
"""
valid_rings = []
for ring in ensure_geometry_collection(rings).geoms:
if len(ring.coords) > 3 or (len(ring.coords) == 3 and ring.coords[0] != ring.coords[-1]):
valid_rings.append(ring)
return GeometryCollection(valid_rings)
def _orient_linear_ring(ring, clockwise=True):
# Unfortunately for us, Inkscape SVGs have an inverted Y coordinate.
# Normally we don't have to care about that, but in this very specific
# case, the meaning of is_ccw is flipped. It actually tests whether
# a ring is clockwise. That makes this logic super-confusing.
if ring.is_ccw != clockwise:
return reverse_line_string(ring)
else:
return ring
def _orient_tree(tree, clockwise=True):
"""
Orient all linear rings in the tree.
Since naturally holes have the opposite point ordering than non-holes we
make all lines within the tree uniform (having all the same ordering
direction)
"""
for node in tree.nodes.values():
node.val = _orient_linear_ring(node.val, clockwise)
def offset_polygon(polygon, offset, join_style, clockwise):
"""
Convert a polygon to a tree of isocontours.
An isocontour is an offset version of the polygon's boundary. For example,
the isocontours of a circle are a set of concentric circles inside the
circle.
This function takes a polygon (which may have holes) as input and creates
isocontours until the polygon is filled completely. The isocontours are
returned as a Tree, with a parent-child relationship indicating that the
parent isocontour contains the child isocontour.
Arguments:
polygon - The shapely Polygon which may have holes
offset - The spacing between isocontours
join_style - Join style used when offsetting the Polygon border to create
isocontours. Can be round, mitered or bevel, as defined by
shapely:
https://shapely.readthedocs.io/en/stable/manual.html#shapely.geometry.JOIN_STYLE
clockwise - If True, isocontour points are in clockwise order; if False, counter-clockwise.
Return Value:
Tree - see above
"""
ordered_polygon = orient(polygon, -1)
tree = Tree()
tree.add_node('root', type='node', parent=None, val=ordered_polygon.exterior)
active_polygons = ['root']
active_holes = [[]]
for hole in ordered_polygon.interiors:
hole_node = tree.generate_node_name()
tree.add_node(hole_node, type="hole", val=hole)
active_holes[0].append(hole_node)
while len(active_polygons) > 0:
check_stop_flag()
current_poly = active_polygons.pop()
current_holes = active_holes.pop()
outer, inners = _offset_polygon_and_holes(tree, current_poly, current_holes, offset, join_style)
polygons = _match_polygons_and_holes(outer, inners)
for polygon in polygons.geoms:
new_polygon, new_holes = _convert_polygon_to_nodes(tree, polygon, parent_polygon=current_poly, child_holes=current_holes)
if new_polygon is not None:
active_polygons.append(new_polygon)
active_holes.append(new_holes)
for previous_hole in current_holes:
# If the previous holes are not
# contained in the new holes they
# have been merged with the
# outer polygon
if not tree.nodes[previous_hole].parent:
tree.nodes[previous_hole].parent = current_poly
tree.add_edge(current_poly, previous_hole)
_orient_tree(tree, clockwise)
return tree
def _offset_polygon_and_holes(tree, poly, holes, offset, join_style):
outer = _offset_linear_ring(
tree.nodes[poly].val,
offset,
resolution=5,
join_style=join_style,
mitre_limit=10,
)
inners = []
for hole in holes:
inner = _offset_linear_ring(
tree.nodes[hole].val,
-offset, # take negative offset for holes
resolution=5,
join_style=join_style,
mitre_limit=10,
)
if not inner.is_empty:
inners.append(Polygon(inner.geoms[0]))
return outer, inners
def _match_polygons_and_holes(outer, inners):
result = MultiPolygon(polygonize(outer.geoms))
if len(inners) > 0:
inners = MultiPolygon(inners)
if not inners.is_valid:
inners = make_valid(MultiPolygon(inners))
result = ensure_geometry_collection(result.difference(inners))
return result
def _convert_polygon_to_nodes(tree, polygon, parent_polygon, child_holes):
if polygon.area < 0.1:
return None, None
polygon = orient(polygon, -1)
valid_rings = _take_only_valid_linear_rings(polygon.exterior)
try:
exterior = valid_rings.geoms[0]
except IndexError:
return None, None
node = tree.generate_node_name()
tree.add_node(node, type='node', parent=parent_polygon, val=exterior)
tree.add_edge(parent_polygon, node)
hole_nodes = []
for hole in polygon.interiors:
hole_node = tree.generate_node_name()
tree.add_node(hole_node, type="hole", val=hole)
for previous_hole in child_holes:
if Polygon(hole).contains(Polygon(tree.nodes[previous_hole].val)):
tree.nodes[previous_hole].parent = hole_node
tree.add_edge(hole_node, previous_hole)
hole_nodes.append(hole_node)
return node, hole_nodes
def _get_nearest_points_closer_than_thresh(travel_line, next_line, threshold):
"""
Find the first point along travel_line that is within threshold of next_line.
Input:
travel_line - The "parent" line for which the distance should
be minimized to enter next_line
next_line - contains the next_line which need to be entered
threshold - The distance between travel_line and next_line needs
to below threshold to be a valid point for entering
Return value:
tuple or None
- the tuple structure is:
(point in travel_line, point in next_line)
- None is returned if there is no point that satisfies the threshold.
"""
# We'll buffer next_line and find the intersection with travel_line.
# Then we'll return the very first point in the intersection,
# matched with a corresponding point on next_line. Fortunately for
# us, intersection of a Polygon with a LineString yields pieces of
# the LineString in the same order as the input LineString.
threshold_area = next_line.buffer(threshold)
portion_within_threshold = travel_line.intersection(threshold_area)
if portion_within_threshold.is_empty:
return None
else:
# Projecting with 0 lets us avoid distinguishing between LineString and
# MultiLineString.
parent_point = Point(portion_within_threshold.interpolate(0))
return nearest_points(parent_point, next_line)
def _create_nearest_points_list(travel_line, tree, children, threshold, threshold_hard):
"""Determine the best place to enter each of parent's children
Arguments:
travel_line - The "parent" line for which the distance should
be minimized to enter each child
children - children of travel_line that need to be entered
threshold - The distance between travel_line and a child should
to be below threshold to be a valid point for entering
threshold_hard - As a last resort, we can accept an entry point
that is this far way
Return value:
list of nearest_neighbor_tuple - indicating where to enter each
respective child
"""
children_nearest_points = []
for child in children:
result = _get_nearest_points_closer_than_thresh(travel_line, tree.nodes[child].val, threshold)
if result is None:
# where holes meet outer borders a distance
# up to 2 * used offset can arise
result = _get_nearest_points_closer_than_thresh(travel_line, tree.nodes[child].val, threshold_hard)
# if we still didn't get a result, ignore this child
# this may lead to oddities, but at least it doesn't fail
if result is None:
continue
proj = travel_line.project(result[0])
children_nearest_points.append(
nearest_neighbor_tuple(
nearest_point_parent=result[0],
nearest_point_child=result[1],
proj_distance_parent=proj,
child_node=child,
)
)
return children_nearest_points
def _find_path_inner_to_outer(tree, node, offset, starting_point, avoid_self_crossing, forward=True):
"""Find a stitch path for this ring and its children.
Strategy: A connection from parent to child is made as fast as possible to
reach the innermost child as fast as possible in order to stitch afterwards
from inner to outer.
This function calls itself recursively to find a stitch path for each child
(and its children).
Arguments:
tree - a Tree of isocontours (as returned by offset_polygon)
offset - offset that was passed to offset_polygon
starting_point - starting point for stitching
avoid_self_crossing - if True, tries to generate a path that does not
cross itself.
forward - if True, this ring will be stitched in its natural direction
(used internally by avoid_self_crossing)
Return value:
LineString -- the stitching path
"""
check_stop_flag()
current_node = tree.nodes[node]
current_ring = current_node.val
if not forward and avoid_self_crossing:
current_ring = reverse_line_string(current_ring)
# reorder the coordinates of this ring so that it starts with
# a point nearest the starting_point
start_distance = current_ring.project(starting_point)
current_ring = roll_linear_ring(current_ring, start_distance)
current_node.val = current_ring
# Find where along this ring to connect to each child.
nearest_points_list = _create_nearest_points_list(
current_ring,
tree,
tree[node],
threshold=1.5 * offset,
threshold_hard=2.05 * offset
)
nearest_points_list.sort(key=lambda tup: tup.proj_distance_parent)
result_coords = []
if not nearest_points_list:
# We have no children, so we're at the center of a spiral. Reversing
# the innermost ring gives a nicer visual appearance.
if not avoid_self_crossing:
current_ring = reverse_line_string(current_ring)
else:
# This is a recursive algorithm. We'll stitch along this ring, pausing
# to jump to each child ring in turn and sew it before continuing on
# this ring. We'll end back where we started.
result_coords.append(current_ring.coords[0])
distance_so_far = 0
for child_connection in nearest_points_list:
# Cut this ring into pieces before and after where this child will connect.
before, after = cut(current_ring, child_connection.proj_distance_parent - distance_so_far)
distance_so_far = child_connection.proj_distance_parent
# Stitch the part leading up to this child.
if before is not None:
result_coords.extend(before.coords)
# Stitch this child. The child will start and end in the same
# place, which should be close to our current location.
child_path = _find_path_inner_to_outer(
tree,
child_connection.child_node,
offset,
child_connection.nearest_point_child,
avoid_self_crossing,
not forward
)
result_coords.extend(child_path.coords)
# Skip ahead a little bit on this ring before resuming. This
# gives a nice spiral pattern, where we spiral out from the
# innermost child.
if after is not None:
skip, after = cut(after, offset)
distance_so_far += offset
current_ring = after
if current_ring is not None:
# skip a little at the end so we don't end exactly where we started.
remaining_length = current_ring.length
if remaining_length > offset:
current_ring, skip = cut(current_ring, current_ring.length - offset)
result_coords.extend(current_ring.coords)
return LineString(result_coords)
def inner_to_outer(tree, polygon, offset,
stitch_length, tolerance, smoothness,
starting_point, avoid_self_crossing,
enable_random, random_sigma, random_seed):
"""Fill a shape with spirals, from innermost to outermost."""
stitch_path = _find_path_inner_to_outer(tree, 'root', offset, starting_point, avoid_self_crossing)
points = [Stitch(*point) for point in stitch_path.coords]
if smoothness > 0:
smoothed = smooth_path(points, smoothness)
points = clamp_path_to_polygon(smoothed, polygon)
stitches = running_stitch(points, stitch_length, tolerance, enable_random, random_sigma, random_seed)
return stitches
def _reorder_linear_ring(ring, start):
distances = ring - start
start_index = np.argmin(np.linalg.norm(distances, axis=1))
return np.roll(ring, -start_index, axis=0)
def _interpolate_linear_rings(ring1, ring2, max_stitch_length, start=None):
"""
Interpolate between two LinearRings
Creates a path from start_point on ring1 and around the rings, ending at a
nearby point on ring2. The path will smoothly transition from ring1 to
ring2 as it travels around the rings.
Inspired by interpolate() from https://github.com/mikedh/pocketing/blob/master/pocketing/polygons.py
Arguments:
ring1 -- LinearRing start point will lie on
ring2 -- LinearRing end point will lie on
max_stitch_length -- maximum stitch length (used to calculate resampling accuracy)
start -- Point on ring1 to start at, as a tuple
Return value: Path interpolated between two LinearRings, as a LineString.
"""
# Resample the two LinearRings so that they are the same number of points
# long. Then take the corresponding points in each ring and interpolate
# between them, gradually going more toward ring2.
#
# This is a little less accurate than the method in interpolate(), but several
# orders of magnitude faster because we're not building and querying a KDTree.
num_points = int(20 * ring1.length / max_stitch_length)
ring1_resampled = trimesh.path.traversal.resample_path(np.array(ring1.coords), count=num_points)
ring2_resampled = trimesh.path.traversal.resample_path(np.array(ring2.coords), count=num_points)
if start is not None:
ring1_resampled = _reorder_linear_ring(ring1_resampled, start)
ring2_resampled = _reorder_linear_ring(ring2_resampled, start)
weights = np.linspace(0.0, 1.0, num_points).reshape((-1, 1))
points = (ring1_resampled * (1.0 - weights)) + (ring2_resampled * weights)
result = LineString(points)
return result.simplify(0.1, False)
def _check_and_prepare_tree_for_valid_spiral(tree):
"""Check whether spiral fill is possible, and tweak if necessary.
Takes a tree consisting of isocontours. If a parent has more than one child
we cannot create a spiral. However, to make the routine more robust, we
allow more than one child if only one of the children has own children. The
other children are removed in this routine then. If the routine returns true,
the tree will have been cleaned up from unwanted children.
If even with these weaker constraints, a spiral is not possible, False is
returned.
"""
def process_node(node):
check_stop_flag()
children = set(tree[node])
if len(children) == 0:
return True
elif len(children) == 1:
child = children.pop()
return process_node(child)
else:
children_with_children = {child for child in children if tree[child]}
if len(children_with_children) > 1:
# Node has multiple children with children, so a perfect spiral is not possible.
# This False value will be returned all the way up the stack.
return False
elif len(children_with_children) == 1:
children_without_children = children - children_with_children
child = children_with_children.pop()
tree.remove_nodes_from(children_without_children)
return process_node(child)
else:
# None of the children has its own children, so we'll just take the longest.
longest = max(children, key=lambda child: tree[child]['val'].length)
shorter_children = children - {longest}
tree.remove_nodes_from(shorter_children)
return process_node(longest)
return process_node('root')
def single_spiral(tree, stitch_length, tolerance, starting_point, enable_random, random_sigma, random_seed):
"""Fill a shape with a single spiral going from outside to center."""
return _spiral_fill(tree, stitch_length, tolerance, starting_point, enable_random, random_sigma, random_seed, _make_spiral)
def double_spiral(tree, stitch_length, tolerance, starting_point, enable_random, random_sigma, random_seed):
"""Fill a shape with a double spiral going from outside to center and back to outside. """
return _spiral_fill(tree, stitch_length, tolerance, starting_point, enable_random, random_sigma, random_seed, _make_fermat_spiral)
def _spiral_fill(tree, stitch_length, tolerance, close_point, enable_random, random_sigma, random_seed, spiral_maker):
starting_point = close_point.coords[0]
rings = _get_spiral_rings(tree)
path = spiral_maker(rings, stitch_length, starting_point)
path = [Stitch(*stitch) for stitch in path]
return running_stitch(path, stitch_length, tolerance, enable_random, random_sigma, random_seed)
def _get_spiral_rings(tree):
rings = []
node = 'root'
while True:
check_stop_flag()
rings.append(tree.nodes[node].val)
children = tree[node]
if len(children) == 0:
break
elif len(children) == 1:
node = list(children)[0]
else:
# We can only really fill a shape with a single spiral if each
# parent has only one child. We'll do our best though, because
# that is probably more helpful to the user than just refusing
# entirely. We'll pick the child that's closest to the center.
parent_center = rings[-1].centroid
node = min(children, key=lambda child: parent_center.distance(tree.nodes[child].val.centroid))
return rings
def _make_fermat_spiral(rings, stitch_length, starting_point):
forward = _make_spiral(rings[::2], stitch_length, starting_point)
back = _make_spiral(rings[1::2], stitch_length, starting_point)
back.reverse()
return chain(forward, back)
def _make_spiral(rings, stitch_length, starting_point):
path = []
spiral_part = None
for ring1, ring2 in zip(rings[:-1], rings[1:]):
check_stop_flag()
spiral_part = _interpolate_linear_rings(ring1, ring2, stitch_length, starting_point)
# skip last to avoid duplicated points
path.extend(spiral_part.coords[:-1])
if spiral_part:
# at the end add last point
path.append(spiral_part.coords[-1])
return path