diff options
Diffstat (limited to 'lib/stitches')
| -rw-r--r-- | lib/stitches/ConnectAndSamplePattern.py | 477 | ||||
| -rw-r--r-- | lib/stitches/DebuggingMethods.py | 155 | ||||
| -rw-r--r-- | lib/stitches/LineStringSampling.py | 502 | ||||
| -rw-r--r-- | lib/stitches/PointTransfer.py | 467 | ||||
| -rw-r--r-- | lib/stitches/StitchPattern.py | 223 | ||||
| -rw-r--r-- | lib/stitches/auto_fill.py | 98 | ||||
| -rw-r--r-- | lib/stitches/constants.py | 41 | ||||
| -rw-r--r-- | lib/stitches/fill.py | 64 |
8 files changed, 2005 insertions, 22 deletions
diff --git a/lib/stitches/ConnectAndSamplePattern.py b/lib/stitches/ConnectAndSamplePattern.py new file mode 100644 index 00000000..21a56cd6 --- /dev/null +++ b/lib/stitches/ConnectAndSamplePattern.py @@ -0,0 +1,477 @@ +from shapely.geometry.polygon import LineString, LinearRing +from shapely.geometry import Point, MultiPoint, linestring +from shapely.ops import nearest_points, polygonize +from collections import namedtuple +from depq import DEPQ +import math +from ..stitches import LineStringSampling +from ..stitches import PointTransfer +from ..stitches import constants + +nearest_neighbor_tuple = namedtuple('nearest_neighbor_tuple', ['nearest_point_parent', 'nearest_point_child', 'projected_distance_parent', 'child_node']) + + +# Cuts a closed line so that the new closed line starts at the point with "distance" to the beginning of the old line. +def cut(line, distance): + if distance <= 0.0 or distance >= line.length: + return [LineString(line)] + coords = list(line.coords) + for i, p in enumerate(coords): + if i > 0 and p == coords[0]: + pd = line.length + else: + pd = line.project(Point(p)) + if pd == distance: + if coords[0] == coords[-1]: + return LineString(coords[i:]+coords[1:i+1]) + else: + return LineString(coords[i:]+coords[:i]) + if pd > distance: + cp = line.interpolate(distance) + if coords[0] == coords[-1]: + return LineString([(cp.x, cp.y)] + coords[i:]+coords[1:i]+[(cp.x, cp.y)]) + else: + return LineString([(cp.x, cp.y)] + coords[i:]+coords[:i]) + + +#Takes the offsetted curves organized as tree, connects and samples them. +#Strategy: A connection from parent to child is made where both curves come closest together. +#Input: +#-tree: contains the offsetted curves in a hierachical organized data structure. +#-used_offset: used offset when the offsetted curves were generated +#-stitch_distance: maximum allowed distance between two points after sampling +#-close_point: defines the beginning point for stitching (stitching starts always from the undisplaced curve) +#-offset_by_half: If true the resulting points are interlaced otherwise not. +#Returnvalues: +#-All offsetted curves connected to one line and sampled with points obeying stitch_distance and offset_by_half +#-Tag (origin) of each point to analyze why a point was placed at this position +def connect_raster_tree_nearest_neighbor(tree, used_offset, stitch_distance, close_point, offset_by_half): + + current_coords = tree.val + abs_offset = abs(used_offset) + result_coords = [] + result_coords_origin = [] + + # We cut the current item so that its index 0 is closest to close_point + start_distance = tree.val.project(close_point) + if start_distance > 0: + current_coords = cut(current_coords, start_distance) + tree.val = current_coords + + if not tree.transferred_point_priority_deque.is_empty(): + new_DEPQ = DEPQ(iterable=None, maxlen=None) + for item,priority in tree.transferred_point_priority_deque: + new_DEPQ.insert(item, math.fmod( + priority-start_distance+current_coords.length, current_coords.length)) + tree.transferred_point_priority_deque = new_DEPQ + #print("Gecutted") + + stitching_direction = 1 + # This list should contain a tuple of nearest points between the current geometry + # and the subgeometry, the projected distance along the current geometry, + # and the belonging subtree node + nearest_points_list = [] + + for subnode in tree.children: + point_parent, point_child = nearest_points(current_coords, subnode.val) + proj_distance = current_coords.project(point_parent) + nearest_points_list.append(nearest_neighbor_tuple(nearest_point_parent = point_parent, + nearest_point_child = point_child, + projected_distance_parent = proj_distance, + child_node=subnode)) + nearest_points_list.sort(reverse=False, key=lambda tup: tup.projected_distance_parent) + + if nearest_points_list: + start_distance = min(abs_offset*constants.factor_offset_starting_points, nearest_points_list[0].projected_distance_parent) + end_distance = max(current_coords.length-abs_offset*constants.factor_offset_starting_points, nearest_points_list[-1].projected_distance_parent) + else: + start_distance = abs_offset*constants.factor_offset_starting_points + end_distance = current_coords.length-abs_offset*constants.factor_offset_starting_points + + own_coords, own_coords_origin = LineStringSampling.raster_line_string_with_priority_points(current_coords, start_distance, # We add/subtract an offset to not sample the same point again (avoid double points for start and end) + end_distance, stitch_distance, stitching_direction, tree.transferred_point_priority_deque, abs_offset) + assert(len(own_coords) == len(own_coords_origin)) + own_coords_origin[0] = LineStringSampling.PointSource.ENTER_LEAVING_POINT + own_coords_origin[-1] = LineStringSampling.PointSource.ENTER_LEAVING_POINT + + #tree.val = LineString(own_coords) + #tree.pointsourcelist = own_coords_origin + tree.stitching_direction = stitching_direction + tree.already_rastered = True + + #Next we need to transfer our rastered points to siblings and childs + to_transfer_point_list = [] + to_transfer_point_list_origin = [] + for k in range(1, len(own_coords)-1): #Do not take the first and the last since they are ENTER_LEAVING_POINT points for sure + # if abs(temp[k][0]-5.25) < 0.5 and abs(temp[k][1]-42.9) < 0.5: + # print("HIER gefunden!") + if (not offset_by_half and own_coords_origin[k] == LineStringSampling.PointSource.EDGE_NEEDED): + continue + if own_coords_origin[k] == LineStringSampling.PointSource.ENTER_LEAVING_POINT or own_coords_origin[k] == LineStringSampling.PointSource.FORBIDDEN_POINT: + continue + to_transfer_point_list.append(Point(own_coords[k])) + point_origin = own_coords_origin[k] + to_transfer_point_list_origin.append(point_origin) + + + #since the projection is only in ccw direction towards inner we need to use "-used_offset" for stitching_direction==-1 + PointTransfer.transfer_points_to_surrounding(tree,stitching_direction*used_offset,offset_by_half,stitch_distance, + to_transfer_point_list,to_transfer_point_list_origin,overnext_neighbor=False, + transfer_forbidden_points=False,transfer_to_parent=False,transfer_to_sibling=True,transfer_to_child=True) + + + #We transfer also to the overnext child to get a more straight arrangement of points perpendicular to the stitching lines + if offset_by_half: + PointTransfer.transfer_points_to_surrounding(tree,stitching_direction*used_offset,False,stitch_distance, + to_transfer_point_list,to_transfer_point_list_origin,overnext_neighbor=True, + transfer_forbidden_points=False,transfer_to_parent=False,transfer_to_sibling=True,transfer_to_child=True) + + if not nearest_points_list: + #If there is no child (inner geometry) we can simply take our own rastered coords as result + result_coords = own_coords + result_coords_origin = own_coords_origin + else: + #There are childs so we need to merge their coordinates with our own rastered coords + + #To create a closed ring + own_coords.append(own_coords[0]) + own_coords_origin.append(own_coords_origin[0]) + + + #own_coords does not start with current_coords but has an offset (see call of raster_line_string_with_priority_points) + total_distance = start_distance + current_item_index = 0 + result_coords = [own_coords[0]] + result_coords_origin = [LineStringSampling.PointSource.ENTER_LEAVING_POINT] + for i in range(1, len(own_coords)): + next_distance = math.sqrt((own_coords[i][0]-own_coords[i-1][0])**2 + + (own_coords[i][1]-own_coords[i-1][1])**2) + while (current_item_index < len(nearest_points_list) and + total_distance+next_distance+constants.eps > nearest_points_list[current_item_index].projected_distance_parent): + + item = nearest_points_list[current_item_index] + child_coords, child_coords_origin = connect_raster_tree_nearest_neighbor( + item.child_node, used_offset, stitch_distance, item.nearest_point_child, offset_by_half) + + delta = item.nearest_point_parent.distance(Point(own_coords[i-1])) + if delta > abs_offset*constants.factor_offset_starting_points: + result_coords.append(item.nearest_point_parent.coords[0]) + result_coords_origin.append(LineStringSampling.PointSource.ENTER_LEAVING_POINT) + # reversing avoids crossing when entering and leaving the child segment + result_coords.extend(child_coords[::-1]) + result_coords_origin.extend(child_coords_origin[::-1]) + + + #And here we calculate the point for the leaving + delta = item.nearest_point_parent.distance(Point(own_coords[i])) + if current_item_index < len(nearest_points_list)-1: + delta = min(delta, abs( + nearest_points_list[current_item_index+1].projected_distance_parent-item.projected_distance_parent)) + + if delta > abs_offset*constants.factor_offset_starting_points: + result_coords.append(current_coords.interpolate( + item.projected_distance_parent+abs_offset*constants.factor_offset_starting_points).coords[0]) + result_coords_origin.append(LineStringSampling.PointSource.ENTER_LEAVING_POINT) + + current_item_index += 1 + if i < len(own_coords)-1: + if(Point(result_coords[-1]).distance(Point(own_coords[i])) > abs_offset*constants.factor_offset_remove_points): + result_coords.append(own_coords[i]) + result_coords_origin.append(own_coords_origin[i]) + + # Since current_coords and temp are rastered differently there accumulate errors regarding the current distance. + # Since a projection of each point in temp would be very time consuming we project only every n-th point which resets the accumulated error every n-th point. + if i % 20 == 0: + total_distance = current_coords.project(Point(own_coords[i])) + else: + total_distance += next_distance + + assert(len(result_coords) == len(result_coords_origin)) + return result_coords, result_coords_origin + +#Takes a line and calculates the nearest distance along this line to enter the 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 +#-thresh: The distance between travel_line and next_line needs to below thresh to be a valid point for entering +#Output: +#-tuple - the tuple structure is: (nearest point in travel_line, nearest point in next_line) +def get_nearest_points_closer_than_thresh(travel_line, next_line,thresh): + point_list = list(MultiPoint(travel_line.coords)) + + if point_list[0].distance(next_line) < thresh: + return nearest_points(point_list[0], next_line) + + for i in range(len(point_list)-1): + line_segment = LineString([point_list[i], point_list[i+1]]) + result = nearest_points(line_segment,next_line) + + if result[0].distance(result[1])< thresh: + return result + line_segment = LineString([point_list[-1], point_list[0]]) + result = nearest_points(line_segment,next_line) + + if result[0].distance(result[1])< thresh: + return result + else: + return None + + +#Takes a line and calculates the nearest distance along this line to enter the childs in children_list +#The method calculates the distances along the line and along the reversed line to find the best direction +#which minimizes the overall distance for all childs. +#Input: +#-travel_line: The "parent" line for which the distance should be minimized to enter the childs +#-children_list: contains the childs of travel_line which need to be entered +#-threshold: The distance between travel_line and a child needs to below threshold to be a valid point for entering +#-preferred_direction: Put a bias on the desired travel direction along travel_line. If equals zero no bias is applied. +# preferred_direction=1 means we prefer the direction of travel_line; preferred_direction=-1 means we prefer the opposite direction. +#Output: +#-stitching direction for travel_line +#-list of tuples (one tuple per child). The tuple structure is: ((nearest point in travel_line, nearest point in child), distance along travel_line, belonging child) +def create_nearest_points_list(travel_line, children_list, threshold, threshold_hard,preferred_direction=0): + result_list_in_order = [] + result_list_reversed_order = [] + + travel_line_reversed = LinearRing(travel_line.coords[::-1]) + + weight_in_order = 0 + weight_reversed_order = 0 + for child in children_list: + result = get_nearest_points_closer_than_thresh(travel_line, child.val, threshold) + if result == None: #where holes meet outer borders a distance up to 2*used offset can arise + result = get_nearest_points_closer_than_thresh(travel_line, child.val, threshold_hard) + assert(result != None) + proj = travel_line.project(result[0]) + weight_in_order += proj + result_list_in_order.append(nearest_neighbor_tuple(nearest_point_parent = result[0], + nearest_point_child = result[1], + projected_distance_parent = proj, + child_node = child)) + + result = get_nearest_points_closer_than_thresh(travel_line_reversed, child.val, threshold) + if result == None: #where holes meet outer borders a distance up to 2*used offset can arise + result = get_nearest_points_closer_than_thresh(travel_line_reversed, child.val, threshold_hard) + assert(result != None) + proj = travel_line_reversed.project(result[0]) + weight_reversed_order += proj + result_list_reversed_order.append(nearest_neighbor_tuple(nearest_point_parent = result[0], + nearest_point_child = result[1], + projected_distance_parent = proj, + child_node = child)) + + if preferred_direction == 1: + weight_in_order=min(weight_in_order/2, max(0, weight_in_order-10*threshold)) + if weight_in_order == weight_reversed_order: + return (1, result_list_in_order) + elif preferred_direction == -1: + weight_reversed_order=min(weight_reversed_order/2, max(0, weight_reversed_order-10*threshold)) + if weight_in_order == weight_reversed_order: + return (-1, result_list_reversed_order) + + + if weight_in_order < weight_reversed_order: + return (1, result_list_in_order) + else: + return (-1, result_list_reversed_order) + + +def calculate_replacing_middle_point(line_segment, abs_offset,max_stich_distance): + angles = LineStringSampling.calculate_line_angles(line_segment) + if angles[1] < abs_offset*constants.limiting_angle_straight: + if line_segment.length < max_stich_distance: + return None + else: + return line_segment.interpolate(line_segment.length-max_stich_distance).coords[0] + else: + return line_segment.coords[1] + +#Takes the offsetted curves organized as tree, connects and samples them. +#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 stich afterwards from inner to outer. +#Input: +#-tree: contains the offsetted curves in a hierachical organized data structure. +#-used_offset: used offset when the offsetted curves were generated +#-stitch_distance: maximum allowed distance between two points after sampling +#-close_point: defines the beginning point for stitching (stitching starts always from the undisplaced curve) +#-offset_by_half: If true the resulting points are interlaced otherwise not. +#Returnvalues: +#-All offsetted curves connected to one line and sampled with points obeying stitch_distance and offset_by_half +#-Tag (origin) of each point to analyze why a point was placed at this position +def connect_raster_tree_from_inner_to_outer(tree, used_offset, stitch_distance, close_point, offset_by_half): + + current_coords = tree.val + abs_offset = abs(used_offset) + result_coords = [] + result_coords_origin = [] + + start_distance = tree.val.project(close_point) + # We cut the current path so that its index 0 is closest to close_point + if start_distance > 0: + current_coords = cut(current_coords, start_distance) + tree.val = current_coords + + if not tree.transferred_point_priority_deque.is_empty(): + new_DEPQ = DEPQ(iterable=None, maxlen=None) + for item, priority in tree.transferred_point_priority_deque: + new_DEPQ.insert(item, math.fmod( + priority-start_distance+current_coords.length, current_coords.length)) + tree.transferred_point_priority_deque = new_DEPQ + + #We try to use always the opposite stitching direction with respect to the parent to avoid crossings when entering and leaving the child + parent_stitching_direction = -1 + if tree.parent != None: + parent_stitching_direction = tree.parent.stitching_direction + + #find the nearest point in current_coords and its children and sort it along the stitching direction + stitching_direction, nearest_points_list = create_nearest_points_list(current_coords, tree.children, 1.5*abs_offset,2.05*abs_offset,parent_stitching_direction) + nearest_points_list.sort(reverse=False, key=lambda tup: tup.projected_distance_parent) + + #Have a small offset for the starting and ending to avoid double points at start and end point (since the paths are closed rings) + if nearest_points_list: + start_offset = min(abs_offset*constants.factor_offset_starting_points, nearest_points_list[0].projected_distance_parent) + end_offset = max(current_coords.length-abs_offset*constants.factor_offset_starting_points, nearest_points_list[-1].projected_distance_parent) + else: + start_offset = abs_offset*constants.factor_offset_starting_points + end_offset = current_coords.length-abs_offset*constants.factor_offset_starting_points + + + if stitching_direction == 1: + own_coords, own_coords_origin = LineStringSampling.raster_line_string_with_priority_points(current_coords, start_offset, # We add start_offset to not sample the same point again (avoid double points for start and end) + end_offset, stitch_distance, stitching_direction, tree.transferred_point_priority_deque, abs_offset) + else: + own_coords, own_coords_origin = LineStringSampling.raster_line_string_with_priority_points(current_coords, current_coords.length-start_offset, # We subtract start_offset to not sample the same point again (avoid double points for start and end) + current_coords.length-end_offset, stitch_distance, stitching_direction, tree.transferred_point_priority_deque, abs_offset) + current_coords.coords = current_coords.coords[::-1] + + #Adjust the points origin for start and end (so that they might not be transferred to childs) + #if own_coords_origin[-1] != LineStringSampling.PointSource.HARD_EDGE: + # own_coords_origin[-1] = LineStringSampling.PointSource.ENTER_LEAVING_POINT + #if own_coords_origin[0] != LineStringSampling.PointSource.HARD_EDGE: + # own_coords_origin[0] = LineStringSampling.PointSource.ENTER_LEAVING_POINT + assert(len(own_coords) == len(own_coords_origin)) + + #tree.val = LineString(own_coords) + #tree.pointsourcelist = own_coords_origin + tree.stitching_direction = stitching_direction + tree.already_rastered = True + + + to_transfer_point_list = [] + to_transfer_point_list_origin = [] + for k in range(0, len(own_coords)): #TODO: maybe do not take the first and the last since they are ENTER_LEAVING_POINT points for sure + if (not offset_by_half and own_coords_origin[k] == LineStringSampling.PointSource.EDGE_NEEDED or own_coords_origin[k] == LineStringSampling.PointSource.FORBIDDEN_POINT): + continue + if own_coords_origin[k] == LineStringSampling.PointSource.ENTER_LEAVING_POINT: + continue + to_transfer_point_list.append(Point(own_coords[k])) + to_transfer_point_list_origin.append(own_coords_origin[k]) + + assert(len(to_transfer_point_list) == len(to_transfer_point_list_origin)) + + + #Next we need to transfer our rastered points to siblings and childs + + + #since the projection is only in ccw direction towards inner we need to use "-used_offset" for stitching_direction==-1 + PointTransfer.transfer_points_to_surrounding(tree,stitching_direction*used_offset,offset_by_half,stitch_distance, + to_transfer_point_list,to_transfer_point_list_origin,overnext_neighbor=False, + transfer_forbidden_points=False,transfer_to_parent=False,transfer_to_sibling=True,transfer_to_child=True) + + + #We transfer also to the overnext child to get a more straight arrangement of points perpendicular to the stitching lines + if offset_by_half: + PointTransfer.transfer_points_to_surrounding(tree,stitching_direction*used_offset,False,stitch_distance, + to_transfer_point_list,to_transfer_point_list_origin,overnext_neighbor=True, + transfer_forbidden_points=False,transfer_to_parent=False,transfer_to_sibling=True,transfer_to_child=True) + + if not nearest_points_list: + #If there is no child (inner geometry) we can simply take our own rastered coords as result + result_coords = own_coords + result_coords_origin = own_coords_origin + else: + #There are childs so we need to merge their coordinates with our own rastered coords + + #Create a closed ring for the following code + own_coords.append(own_coords[0]) + own_coords_origin.append(own_coords_origin[0]) + + # own_coords does not start with current_coords but has an offset (see call of raster_line_string_with_priority_points) + total_distance = start_offset + + current_item_index = 0 + result_coords = [own_coords[0]] + result_coords_origin = [own_coords_origin[0]] + + for i in range(1, len(own_coords)): + next_distance = math.sqrt((own_coords[i][0]-own_coords[i-1][0])**2 + + (own_coords[i][1]-own_coords[i-1][1])**2) + while (current_item_index < len(nearest_points_list) and + total_distance+next_distance+constants.eps > nearest_points_list[current_item_index].projected_distance_parent): + #The current and the next point in own_coords enclose the nearest point tuple between this geometry and the child geometry. + #Hence we need to insert the child geometry points here before the next point of own_coords. + item = nearest_points_list[current_item_index] + child_coords, child_coords_origin = connect_raster_tree_from_inner_to_outer( + item.child_node, used_offset, stitch_distance, item.nearest_point_child, offset_by_half) + + #Imagine the nearest point of the child is within a long segment of the parent. Without additonal points + #on the parent side this would cause noticeable deviations. Hence we add here points shortly before and after + #the entering of the child to have only minor deviations to the desired shape. + #Here is the point for the entering: + if(Point(result_coords[-1]).distance(item.nearest_point_parent) > constants.factor_offset_starting_points*abs_offset): + result_coords.append(item.nearest_point_parent.coords[0]) + result_coords_origin.append(LineStringSampling.PointSource.ENTER_LEAVING_POINT) + #if (abs(result_coords[-1][0]-61.7) < 0.2 and abs(result_coords[-1][1]-105.1) < 0.2): + # print("HIIER FOUNDED3") + + #Check whether the number of points of the connecting lines from child to child can be reduced + if len(child_coords) > 1: + point = calculate_replacing_middle_point(LineString([result_coords[-1],child_coords[0],child_coords[1]]),abs_offset,stitch_distance) + #if (abs(result_coords[-1][0]-8.9) < 0.2 and abs(result_coords[-1][1]-8.9) < 0.2): + # print("HIIER FOUNDED3") + if point != None: + #if (abs(point[0]-17.8) < 0.2 and abs(point[1]-17.8) < 0.2): + # print("HIIER FOUNDED3") + result_coords.append(point) + result_coords_origin.append(child_coords_origin[0]) + + result_coords.extend(child_coords[1:]) + result_coords_origin.extend(child_coords_origin[1:]) + else: + result_coords.extend(child_coords) + result_coords_origin.extend(child_coords_origin) + + #And here is the point for the leaving of the child (distance to the own following point should not be too large) + delta = item.nearest_point_parent.distance(Point(own_coords[i])) + if current_item_index < len(nearest_points_list)-1: + delta = min(delta, abs( + nearest_points_list[current_item_index+1].projected_distance_parent-item.projected_distance_parent)) + + if delta > constants.factor_offset_starting_points*abs_offset: + result_coords.append(current_coords.interpolate( + item.projected_distance_parent+2*constants.factor_offset_starting_points*abs_offset).coords[0]) + result_coords_origin.append(LineStringSampling.PointSource.ENTER_LEAVING_POINT) + #check whether this additional point makes the last point of the child unnecessary + point = calculate_replacing_middle_point(LineString([result_coords[-3],result_coords[-2],result_coords[-1]]),abs_offset,stitch_distance) + if point == None: + result_coords.pop(-2) + result_coords_origin.pop(-2) + + #if (abs(result_coords[-1][0]-61.7) < 0.2 and abs(result_coords[-1][1]-105.1) < 0.2): + # print("HIIER FOUNDED3") + + current_item_index += 1 + if i < len(own_coords)-1: + if(Point(result_coords[-1]).distance(Point(own_coords[i])) > abs_offset*constants.factor_offset_remove_points): + result_coords.append(own_coords[i]) + result_coords_origin.append(own_coords_origin[i]) + + # Since current_coords and own_coords are rastered differently there accumulate errors regarding the current distance. + # Since a projection of each point in own_coords would be very time consuming we project only every n-th point which resets the accumulated error every n-th point. + if i % 20 == 0: + total_distance = current_coords.project(Point(own_coords[i])) + else: + total_distance += next_distance + + assert(len(result_coords) == len(result_coords_origin)) + return result_coords, result_coords_origin diff --git a/lib/stitches/DebuggingMethods.py b/lib/stitches/DebuggingMethods.py new file mode 100644 index 00000000..d0f65576 --- /dev/null +++ b/lib/stitches/DebuggingMethods.py @@ -0,0 +1,155 @@ + +import matplotlib.pyplot as plt +from shapely.geometry import Polygon +from shapely.ops import nearest_points, substring, polygonize + +from anytree import PreOrderIter +from shapely.geometry.polygon import orient +#import LineStringSampling as Sampler +import numpy as np +import matplotlib.collections as mcoll +import matplotlib.path as mpath + +# def offset_polygons(polys, offset,joinstyle): +# if polys.geom_type == 'Polygon': +# inners = polys.interiors +# outer = polys.exterior +# polyinners = [] +# for inner in inners: +# inner = inner.parallel_offset(offset,'left', 5, joinstyle, 1) +# polyinners.append(Polygon(inner)) +# outer = outer.parallel_offset(offset,'left', 5, joinstyle, 1) +# return Polygon(outer).difference(MultiPolygon(polyinners)) +# else: +# polyreturns = [] +# for poly in polys: +# inners = poly.interiors +# outer = poly.exterior +# polyinners = [] +# for inner in inners: +# inner = inner.parallel_offset(offset,'left', 5, joinstyle, 1) +# polyinners.append(Polygon(inner)) +# outer = outer.parallel_offset(offset,'left', 5, joinstyle, 1) +# result = Polygon(outer).difference(MultiPolygon(polyinners)) +# polyreturns.append(result) +# return MultiPolygon(polyreturns) + +# For debugging + + +def plot_MultiPolygon(MultiPoly, plt, colorString): + if MultiPoly.is_empty: + return + if MultiPoly.geom_type == 'Polygon': + x2, y2 = MultiPoly.exterior.xy + plt.plot(x2, y2, colorString) + + for inners in MultiPoly.interiors: + x2, y2 = inners.coords.xy + plt.plot(x2, y2, colorString) + else: + for poly in MultiPoly: + x2, y2 = poly.exterior.xy + plt.plot(x2, y2, colorString) + + for inners in poly.interiors: + x2, y2 = inners.coords.xy + plt.plot(x2, y2, colorString) + +# Test whether there are areas which would currently not be stitched but should be stitched + + +def subtractResult(poly, rootPoly, offsetThresh): + poly2 = Polygon(poly) + for node in PreOrderIter(rootPoly): + poly2 = poly2.difference(node.val.buffer(offsetThresh, 5, 3, 3)) + return poly2 + +# Used for debugging - plots all polygon exteriors within an AnyTree which is provided by the root node rootPoly. + + +def drawPoly(rootPoly, colorString): + fig, axs = plt.subplots(1, 1) + axs.axis('equal') + plt.gca().invert_yaxis() + for node in PreOrderIter(rootPoly): + # if(node.id == "hole"): + # node.val = LinearRing(node.val.coords[::-1]) + print("Bounds:") + print(node.val.bounds) + x2, y2 = node.val.coords.xy + plt.plot(x2, y2, colorString) + plt.show(block=True) + + +def drawresult(resultcoords, resultcoords_Origin, colorString): + fig, axs = plt.subplots(1, 1) + axs.axis('equal') + plt.gca().invert_yaxis() + plt.plot(*zip(*resultcoords), colorString) + + colormap = np.array(['r', 'g', 'b', 'c', 'm', 'y', 'k', 'gray', 'm']) + labelmap = np.array(['MUST_USE', 'REGULAR_SPACING', 'INITIAL_RASTERING', 'EDGE_NEEDED', 'NOT_NEEDED', + 'ALREADY_TRANSFERRED', 'ADDITIONAL_TRACKING_POINT_NOT_NEEDED', 'EDGE_RASTERING_ALLOWED', 'EDGE_PREVIOUSLY_SHIFTED']) + + for i in range(0, 8+1): + # if i != Sampler.PointSource.EDGE_NEEDED and i != Sampler.PointSource.INITIAL_RASTERING: + # continue + selection = [] + for j in range(len(resultcoords)): + if i == resultcoords_Origin[j]: + selection.append(resultcoords[j]) + if len(selection) > 0: + plt.scatter(*zip(*selection), c=colormap[i], label=labelmap[i]) + + # plt.scatter(*zip(*resultcoords), + # c=colormap[resultcoords_Origin]) + axs.legend() + plt.show(block=True) + + +# Just for debugging in order to draw the connected line with color gradient + + +def colorline( + x, y, z=None, cmap=plt.get_cmap('copper'), norm=plt.Normalize(0.0, 1.0), + linewidth=3, alpha=1.0): + """ + http://nbviewer.ipython.org/github/dpsanders/matplotlib-examples/blob/master/colorline.ipynb + http://matplotlib.org/examples/pylab_examples/multicolored_line.html + Plot a colored line with coordinates x and y + Optionally specify colors in the array z + Optionally specify a colormap, a norm function and a line width + """ + + # Default colors equally spaced on [0,1]: + if z is None: + z = np.linspace(0.0, 1.0, len(x)) + + # Special case if a single number: + if not hasattr(z, "__iter__"): # to check for numerical input -- this is a hack + z = np.array([z]) + + z = np.asarray(z) + + segments = make_segments(x, y) + lc = mcoll.LineCollection(segments, array=z, cmap=cmap, norm=norm, + linewidth=linewidth, alpha=alpha) + + ax = plt.gca() + ax.add_collection(lc) + + return lc + +# Used by colorline + + +def make_segments(x, y): + """ + Create list of line segments from x and y coordinates, in the correct format + for LineCollection: an array of the form numlines x (points per line) x 2 (x + and y) array + """ + points = np.array([x, y]).T.reshape(-1, 1, 2) + segments = np.concatenate([points[:-1], points[1:]], axis=1) + return segments diff --git a/lib/stitches/LineStringSampling.py b/lib/stitches/LineStringSampling.py new file mode 100644 index 00000000..434c6bbf --- /dev/null +++ b/lib/stitches/LineStringSampling.py @@ -0,0 +1,502 @@ +from sys import path +from shapely.geometry.polygon import LineString +from shapely.geometry import Point +from shapely.ops import substring +import math +import numpy as np +from enum import IntEnum +from ..stitches import constants +from ..stitches import PointTransfer + +#Used to tag the origin of a rastered point +class PointSource(IntEnum): + #MUST_USE = 0 # Legacy + REGULAR_SPACING = 1 # introduced to not exceed maximal stichting distance + #INITIAL_RASTERING = 2 #Legacy + EDGE_NEEDED = 3 # point which must be stitched to avoid to large deviations to the desired path + #NOT_NEEDED = 4 #Legacy + #ALREADY_TRANSFERRED = 5 #Legacy + #ADDITIONAL_TRACKING_POINT_NOT_NEEDED = 6 #Legacy + #EDGE_RASTERING_ALLOWED = 7 #Legacy + #EDGE_PREVIOUSLY_SHIFTED = 8 #Legacy + ENTER_LEAVING_POINT = 9 #Whether this point is used to enter or leave a child + SOFT_EDGE_INTERNAL = 10 #If the angle at a point is <= constants.limiting_angle this point is marked as SOFT_EDGE + HARD_EDGE_INTERNAL = 11 #If the angle at a point is > constants.limiting_angle this point is marked as HARD_EDGE (HARD_EDGES will always be stitched) + PROJECTED_POINT = 12 #If the point was created by a projection (transferred point) of a neighbor it is marked as PROJECTED_POINT + REGULAR_SPACING_INTERNAL = 13 # introduced to not exceed maximal stichting distance + #FORBIDDEN_POINT_INTERNAL=14 #Legacy + SOFT_EDGE = 15 #If the angle at a point is <= constants.limiting_angle this point is marked as SOFT_EDGE + HARD_EDGE = 16 #If the angle at a point is > constants.limiting_angle this point is marked as HARD_EDGE (HARD_EDGES will always be stitched) + FORBIDDEN_POINT=17 #Only relevant for desired interlacing - non-shifted point positions at the next neighbor are marked as forbidden + REPLACED_FORBIDDEN_POINT=18 #If one decides to avoid forbidden points new points to the left and to the right as replacement are created + DIRECT = 19 #Calculated by next neighbor projection + OVERNEXT = 20 #Calculated by overnext neighbor projection + + +# Calculates the angles between adjacent edges at each interior point +#Note that the first and last values in the return array are zero since for the boundary points no angle calculations were possible +def calculate_line_angles(line): + Angles = np.zeros(len(line.coords)) + for i in range(1, len(line.coords)-1): + vec1 = np.array(line.coords[i])-np.array(line.coords[i-1]) + vec2 = np.array(line.coords[i+1])-np.array(line.coords[i]) + vec1length = np.linalg.norm(vec1) + vec2length = np.linalg.norm(vec2) + #if vec1length <= 0: + # print("HIER FEHLER") + + #if vec2length <=0: + # print("HIER FEHLEr") + assert(vec1length >0) + assert(vec2length >0) + scalar_prod=np.dot(vec1, vec2)/(vec1length*vec2length) + scalar_prod = min(max(scalar_prod,-1),1) + #if scalar_prod > 1.0: + # scalar_prod = 1.0 + #elif scalar_prod < -1.0: + # scalar_prod = -1.0 + Angles[i] = math.acos(scalar_prod) + return Angles + +#Rasters a line between start_distance and end_distance. +#Input: +#-line: The line to be rastered +#-start_distance: The distance along the line from which the rastering should start +#-end_distance: The distance along the line until which the rastering should be done +#-maxstitch_distance: The maximum allowed stitch distance +#-stitching_direction: =1 is stitched along line direction, =-1 if stitched in reversed order. Note that +# start_distance > end_distance for stitching_direction = -1 +#-must_use_points_deque: deque with projected points on line from its neighbors. An item of the deque +#is setup as follows: ((projected point on line, LineStringSampling.PointSource), priority=distance along line) +#index of point_origin is the index of the point in the neighboring line +#-abs_offset: used offset between to offsetted curves +#Output: +#-List of tuples with the rastered point coordinates +#-List which defines the point origin for each point according to the PointSource enum. +def raster_line_string_with_priority_points(line, start_distance, end_distance, maxstitch_distance, stitching_direction, must_use_points_deque, abs_offset): + if (abs(end_distance-start_distance) < constants.line_lengh_seen_as_one_point): + return [line.interpolate(start_distance).coords[0]], [PointSource.HARD_EDGE] + + assert (stitching_direction == -1 and start_distance >= end_distance) or ( + stitching_direction == 1 and start_distance <= end_distance) + + deque_points = list(must_use_points_deque) + + linecoords = line.coords + + if start_distance > end_distance: + start_distance, end_distance = line.length - \ + start_distance, line.length-end_distance + linecoords = linecoords[::-1] + for i in range(len(deque_points)): + deque_points[i] = (deque_points[i][0], + line.length-deque_points[i][1]) + else: + deque_points = deque_points[::-1] #Since points with highest priority (=distance along line) are first (descending sorted) + + # Remove all points from the deque which do not fall in the segment [start_distance; end_distance] + while (len(deque_points) > 0 and deque_points[0][1] <= start_distance+min(maxstitch_distance/20, constants.point_spacing_to_be_considered_equal)): + deque_points.pop(0) + while (len(deque_points) > 0 and deque_points[-1][1] >= end_distance-min(maxstitch_distance/20, constants.point_spacing_to_be_considered_equal)): + deque_points.pop() + + +# Ordering in priority queue: +# (point, LineStringSampling.PointSource), priority) + aligned_line = LineString(linecoords) + path_coords = substring(aligned_line, + start_distance, end_distance) + + #aligned line is a line without doubled points. I had the strange situation in which the offset "start_distance" from the line beginning resulted in a starting point which was + # already present in aligned_line causing a doubled point. A double point is not allowed in the following calculations so we need to remove it: + if abs(path_coords.coords[0][0]-path_coords.coords[1][0])<constants.eps and abs(path_coords.coords[0][1]-path_coords.coords[1][1])<constants.eps: + path_coords.coords = path_coords.coords[1:] + if abs(path_coords.coords[-1][0]-path_coords.coords[-2][0])<constants.eps and abs(path_coords.coords[-1][1]-path_coords.coords[-2][1])<constants.eps: + path_coords.coords = path_coords.coords[:-1] + + angles = calculate_line_angles(path_coords) + + current_distance = start_distance + + #Next we merge the line points and the projected (deque) points into one list + merged_point_list = [] + dq_iter = 0 + for point,angle in zip(path_coords.coords,angles): + #if abs(point[0]-40.4) < 0.2 and abs(point[1]-2.3)< 0.2: + # print("GEFUNDEN") + current_distance = start_distance+path_coords.project(Point(point)) + while dq_iter < len(deque_points) and deque_points[dq_iter][1] < current_distance: + #We want to avoid setting points at soft edges close to forbidden points + if deque_points[dq_iter][0].point_source == PointSource.FORBIDDEN_POINT: + #Check whether a previous added point is a soft edge close to the forbidden point + if (merged_point_list[-1][0].point_source == PointSource.SOFT_EDGE_INTERNAL and + abs(merged_point_list[-1][1]-deque_points[dq_iter][1] < abs_offset*constants.factor_offset_forbidden_point)): + item = merged_point_list.pop() + merged_point_list.append((PointTransfer.projected_point_tuple(point=item[0].point, point_source=\ + PointSource.FORBIDDEN_POINT),item[1])) + else: + merged_point_list.append(deque_points[dq_iter]) + dq_iter+=1 + #Check whether the current point is close to a forbidden point + if (dq_iter < len(deque_points) and + deque_points[dq_iter-1][0].point_source == PointSource.FORBIDDEN_POINT and + angle < constants.limiting_angle and + abs(deque_points[dq_iter-1][1]-current_distance) < abs_offset*constants.factor_offset_forbidden_point): + point_source = PointSource.FORBIDDEN_POINT + else: + if angle < constants.limiting_angle: + point_source = PointSource.SOFT_EDGE_INTERNAL + else: + point_source = PointSource.HARD_EDGE_INTERNAL + merged_point_list.append((PointTransfer.projected_point_tuple(point=Point(point), point_source=point_source),current_distance)) + + result_list = [merged_point_list[0]] + + #General idea: Take one point of merged_point_list after another into the current segment until this segment is not simplified to a straight line by shapelys simplify method. + #Then, look at the points within this segment and choose the best fitting one (HARD_EDGE > OVERNEXT projected point > DIRECT projected point) as termination of this segment + # and start point for the next segment (so we do not always take the maximum possible length for a segment) + segment_start_index = 0 + segment_end_index = 1 + forbidden_point_list = [] + while segment_end_index < len(merged_point_list): + #if abs(merged_point_list[segment_end_index-1][0].point.coords[0][0]-67.9) < 0.2 and abs(merged_point_list[segment_end_index-1][0].point.coords[0][1]-161.0)< 0.2: + # print("GEFUNDEN") + + #Collection of points for the current segment + current_point_list = [merged_point_list[segment_start_index][0].point] + + while segment_end_index < len(merged_point_list): + segment_length = merged_point_list[segment_end_index][1]-merged_point_list[segment_start_index][1] + if segment_length > maxstitch_distance+constants.point_spacing_to_be_considered_equal: + new_distance = merged_point_list[segment_start_index][1]+maxstitch_distance + merged_point_list.insert(segment_end_index,(PointTransfer.projected_point_tuple(point=aligned_line.interpolate(new_distance), point_source=\ + PointSource.REGULAR_SPACING_INTERNAL),new_distance)) + if abs(merged_point_list[segment_end_index][0].point.coords[0][0]-12.2) < 0.2 and abs(merged_point_list[segment_end_index][0].point.coords[0][1]-0.9)< 0.2: + print("GEFUNDEN") + segment_end_index+=1 + break + #if abs(merged_point_list[segment_end_index][0].point.coords[0][0]-93.6) < 0.2 and abs(merged_point_list[segment_end_index][0].point.coords[0][1]-122.7)< 0.2: + # print("GEFUNDEN") + + current_point_list.append(merged_point_list[segment_end_index][0].point) + simplified_len = len(LineString(current_point_list).simplify(constants.factor_offset_remove_dense_points*abs_offset,preserve_topology=False).coords) + if simplified_len > 2: #not all points have been simplified - so we need to add it + break + + if merged_point_list[segment_end_index][0].point_source ==PointSource.HARD_EDGE_INTERNAL: + segment_end_index+=1 + break + segment_end_index+=1 + + segment_end_index-=1 + + #Now we choose the best fitting point within this segment + index_overnext = -1 + index_direct = -1 + index_hard_edge = -1 + + iter = segment_start_index+1 + while (iter <= segment_end_index): + if merged_point_list[iter][0].point_source == PointSource.OVERNEXT: + index_overnext = iter + elif merged_point_list[iter][0].point_source == PointSource.DIRECT: + index_direct = iter + elif merged_point_list[iter][0].point_source == PointSource.HARD_EDGE_INTERNAL: + index_hard_edge = iter + iter += 1 + if index_hard_edge != -1: + segment_end_index = index_hard_edge + else: + if index_overnext != -1: + if (index_direct != -1 and index_direct > index_overnext and + (merged_point_list[index_direct][1]-merged_point_list[index_overnext][1]) >= + constants.factor_segment_length_direct_preferred_over_overnext* + (merged_point_list[index_overnext][1]-merged_point_list[segment_start_index][1])): + #We allow to take the direct projected point instead of the overnext projected point if it would result in a + #significant longer segment length + segment_end_index = index_direct + else: + segment_end_index = index_overnext + elif index_direct != -1: + segment_end_index = index_direct + + #Usually OVERNEXT and DIRECT points are close to each other and in some cases both were selected as segment edges + #If they are too close (<abs_offset) we remove one of it + if (((merged_point_list[segment_start_index][0].point_source == PointSource.OVERNEXT and + merged_point_list[segment_end_index][0].point_source == PointSource.DIRECT) or + (merged_point_list[segment_start_index][0].point_source == PointSource.DIRECT and + merged_point_list[segment_end_index][0].point_source == PointSource.OVERNEXT)) and + abs(merged_point_list[segment_end_index][1] - merged_point_list[segment_start_index][1]) < abs_offset): + result_list.pop() + + result_list.append(merged_point_list[segment_end_index]) + #To have a chance to replace all forbidden points afterwards + if merged_point_list[segment_end_index][0].point_source == PointSource.FORBIDDEN_POINT: + forbidden_point_list.append(len(result_list)-1) + + segment_start_index = segment_end_index + segment_end_index+=1 + + return_point_list = [] #[result_list[0][0].point.coords[0]] + return_point_source_list = []#[result_list[0][0].point_source] + + #Currently replacement of forbidden points not satisfying + #result_list = replace_forbidden_points(aligned_line, result_list, forbidden_point_list,abs_offset) + + #Finally we create the final return_point_list and return_point_source_list + for i in range(len(result_list)): + return_point_list.append(result_list[i][0].point.coords[0]) + #if abs(result_list[i][0].point.coords[0][0]-91.7) < 0.2 and abs(result_list[i][0].point.coords[0][1]-106.15)< 0.2: + # print("GEFUNDEN") + if result_list[i][0].point_source == PointSource.HARD_EDGE_INTERNAL: + point_source = PointSource.HARD_EDGE + elif result_list[i][0].point_source == PointSource.SOFT_EDGE_INTERNAL: + point_source = PointSource.SOFT_EDGE + elif result_list[i][0].point_source == PointSource.REGULAR_SPACING_INTERNAL: + point_source = PointSource.REGULAR_SPACING + elif result_list[i][0].point_source == PointSource.FORBIDDEN_POINT: + point_source = PointSource.FORBIDDEN_POINT + else: + point_source = PointSource.PROJECTED_POINT + + return_point_source_list.append(point_source) + + + assert(len(return_point_list) == len(return_point_source_list)) + + #return remove_dense_points(returnpointlist, returnpointsourcelist, maxstitch_distance,abs_offset) + return return_point_list, return_point_source_list + +#Rasters a line between start_distance and end_distance. +#Input: +#-line: The line to be rastered +#-start_distance: The distance along the line from which the rastering should start +#-end_distance: The distance along the line until which the rastering should be done +#-maxstitch_distance: The maximum allowed stitch distance +#-stitching_direction: =1 is stitched along line direction, =-1 if stitched in reversed order. Note that +# start_distance > end_distance for stitching_direction = -1 +#-must_use_points_deque: deque with projected points on line from its neighbors. An item of the deque +#is setup as follows: ((projected point on line, LineStringSampling.PointSource), priority=distance along line) +#index of point_origin is the index of the point in the neighboring line +#-abs_offset: used offset between to offsetted curves +#Output: +#-List of tuples with the rastered point coordinates +#-List which defines the point origin for each point according to the PointSource enum. +def raster_line_string_with_priority_points_graph(line, maxstitch_distance, stitching_direction, must_use_points_deque, abs_offset, offset_by_half): + if (line.length < constants.line_lengh_seen_as_one_point): + return [line.coords[0]], [PointSource.HARD_EDGE] + + deque_points = list(must_use_points_deque) + + linecoords = line.coords + + if stitching_direction==-1: + linecoords = linecoords[::-1] + for i in range(len(deque_points)): + deque_points[i] = (deque_points[i][0], + line.length-deque_points[i][1]) + else: + deque_points = deque_points[::-1] #Since points with highest priority (=distance along line) are first (descending sorted) + +# Ordering in priority queue: +# (point, LineStringSampling.PointSource), priority) + aligned_line = LineString(linecoords) #might be different from line for stitching_direction=-1 + + angles = calculate_line_angles(aligned_line) + #For the first and last point we cannot calculate an angle. Set it to above the limit to make it a hard edge + angles[0] = 1.1*constants.limiting_angle + angles[-1] = 1.1*constants.limiting_angle + + current_distance = 0.0 + + #Next we merge the line points and the projected (deque) points into one list + merged_point_list = [] + dq_iter = 0 + for point,angle in zip(aligned_line.coords,angles): + #if abs(point[0]-52.9) < 0.2 and abs(point[1]-183.4)< 0.2: + # print("GEFUNDEN") + current_distance = aligned_line.project(Point(point)) + while dq_iter < len(deque_points) and deque_points[dq_iter][1] < current_distance: + #We want to avoid setting points at soft edges close to forbidden points + if deque_points[dq_iter][0].point_source == PointSource.FORBIDDEN_POINT: + #Check whether a previous added point is a soft edge close to the forbidden point + if (merged_point_list[-1][0].point_source == PointSource.SOFT_EDGE_INTERNAL and + abs(merged_point_list[-1][1]-deque_points[dq_iter][1] < abs_offset*constants.factor_offset_forbidden_point)): + item = merged_point_list.pop() + merged_point_list.append((PointTransfer.projected_point_tuple(point=item[0].point, point_source=\ + PointSource.FORBIDDEN_POINT),item[1])) + else: + merged_point_list.append(deque_points[dq_iter]) + dq_iter+=1 + #Check whether the current point is close to a forbidden point + if (dq_iter < len(deque_points) and + deque_points[dq_iter-1][0].point_source == PointSource.FORBIDDEN_POINT and + angle < constants.limiting_angle and + abs(deque_points[dq_iter-1][1]-current_distance) < abs_offset*constants.factor_offset_forbidden_point): + point_source = PointSource.FORBIDDEN_POINT + else: + if angle < constants.limiting_angle: + point_source = PointSource.SOFT_EDGE_INTERNAL + else: + point_source = PointSource.HARD_EDGE_INTERNAL + merged_point_list.append((PointTransfer.projected_point_tuple(point=Point(point), point_source=point_source),current_distance)) + + result_list = [merged_point_list[0]] + + #General idea: Take one point of merged_point_list after another into the current segment until this segment is not simplified to a straight line by shapelys simplify method. + #Then, look at the points within this segment and choose the best fitting one (HARD_EDGE > OVERNEXT projected point > DIRECT projected point) as termination of this segment + # and start point for the next segment (so we do not always take the maximum possible length for a segment) + segment_start_index = 0 + segment_end_index = 1 + forbidden_point_list = [] + while segment_end_index < len(merged_point_list): + #if abs(merged_point_list[segment_end_index-1][0].point.coords[0][0]-67.9) < 0.2 and abs(merged_point_list[segment_end_index-1][0].point.coords[0][1]-161.0)< 0.2: + # print("GEFUNDEN") + + #Collection of points for the current segment + current_point_list = [merged_point_list[segment_start_index][0].point] + + while segment_end_index < len(merged_point_list): + segment_length = merged_point_list[segment_end_index][1]-merged_point_list[segment_start_index][1] + if segment_length > maxstitch_distance+constants.point_spacing_to_be_considered_equal: + new_distance = merged_point_list[segment_start_index][1]+maxstitch_distance + merged_point_list.insert(segment_end_index,(PointTransfer.projected_point_tuple(point=aligned_line.interpolate(new_distance), point_source=\ + PointSource.REGULAR_SPACING_INTERNAL),new_distance)) + #if abs(merged_point_list[segment_end_index][0].point.coords[0][0]-12.2) < 0.2 and abs(merged_point_list[segment_end_index][0].point.coords[0][1]-0.9)< 0.2: + # print("GEFUNDEN") + segment_end_index+=1 + break + #if abs(merged_point_list[segment_end_index][0].point.coords[0][0]-34.4) < 0.2 and abs(merged_point_list[segment_end_index][0].point.coords[0][1]-6.2)< 0.2: + # print("GEFUNDEN") + + current_point_list.append(merged_point_list[segment_end_index][0].point) + simplified_len = len(LineString(current_point_list).simplify(constants.factor_offset_remove_dense_points*abs_offset,preserve_topology=False).coords) + if simplified_len > 2: #not all points have been simplified - so we need to add it + break + + if merged_point_list[segment_end_index][0].point_source ==PointSource.HARD_EDGE_INTERNAL: + segment_end_index+=1 + break + segment_end_index+=1 + + segment_end_index-=1 + + #Now we choose the best fitting point within this segment + index_overnext = -1 + index_direct = -1 + index_hard_edge = -1 + + iter = segment_start_index+1 + while (iter <= segment_end_index): + if merged_point_list[iter][0].point_source == PointSource.OVERNEXT: + index_overnext = iter + elif merged_point_list[iter][0].point_source == PointSource.DIRECT: + index_direct = iter + elif merged_point_list[iter][0].point_source == PointSource.HARD_EDGE_INTERNAL: + index_hard_edge = iter + iter += 1 + if index_hard_edge != -1: + segment_end_index = index_hard_edge + else: + if offset_by_half: + index_preferred = index_overnext + index_less_preferred = index_direct + else: + index_preferred = index_direct + index_less_preferred = index_overnext + + if index_preferred != -1: + if (index_less_preferred != -1 and index_less_preferred > index_preferred and + (merged_point_list[index_less_preferred][1]-merged_point_list[index_preferred][1]) >= + constants.factor_segment_length_direct_preferred_over_overnext* + (merged_point_list[index_preferred][1]-merged_point_list[segment_start_index][1])): + #We allow to take the direct projected point instead of the overnext projected point if it would result in a + #significant longer segment length + segment_end_index = index_less_preferred + else: + segment_end_index = index_preferred + elif index_less_preferred != -1: + segment_end_index = index_less_preferred + + #Usually OVERNEXT and DIRECT points are close to each other and in some cases both were selected as segment edges + #If they are too close (<abs_offset) we remove one of it + if (((merged_point_list[segment_start_index][0].point_source == PointSource.OVERNEXT and + merged_point_list[segment_end_index][0].point_source == PointSource.DIRECT) or + (merged_point_list[segment_start_index][0].point_source == PointSource.DIRECT and + merged_point_list[segment_end_index][0].point_source == PointSource.OVERNEXT)) and + abs(merged_point_list[segment_end_index][1] - merged_point_list[segment_start_index][1]) < abs_offset): + result_list.pop() + + result_list.append(merged_point_list[segment_end_index]) + #To have a chance to replace all forbidden points afterwards + if merged_point_list[segment_end_index][0].point_source == PointSource.FORBIDDEN_POINT: + forbidden_point_list.append(len(result_list)-1) + + segment_start_index = segment_end_index + segment_end_index+=1 + + return_point_list = [] #[result_list[0][0].point.coords[0]] + return_point_source_list = []#[result_list[0][0].point_source] + + #Currently replacement of forbidden points not satisfying + result_list = replace_forbidden_points(aligned_line, result_list, forbidden_point_list,abs_offset) + + #Finally we create the final return_point_list and return_point_source_list + for i in range(len(result_list)): + return_point_list.append(result_list[i][0].point.coords[0]) + #if abs(result_list[i][0].point.coords[0][0]-91.7) < 0.2 and abs(result_list[i][0].point.coords[0][1]-106.15)< 0.2: + # print("GEFUNDEN") + if result_list[i][0].point_source == PointSource.HARD_EDGE_INTERNAL: + point_source = PointSource.HARD_EDGE + elif result_list[i][0].point_source == PointSource.SOFT_EDGE_INTERNAL: + point_source = PointSource.SOFT_EDGE + elif result_list[i][0].point_source == PointSource.REGULAR_SPACING_INTERNAL: + point_source = PointSource.REGULAR_SPACING + elif result_list[i][0].point_source == PointSource.FORBIDDEN_POINT: + point_source = PointSource.FORBIDDEN_POINT + else: + point_source = PointSource.PROJECTED_POINT + + return_point_source_list.append(point_source) + + + assert(len(return_point_list) == len(return_point_source_list)) + + #return remove_dense_points(returnpointlist, returnpointsourcelist, maxstitch_distance,abs_offset) + return return_point_list, return_point_source_list + +def replace_forbidden_points(line, result_list, forbidden_point_list_indices, abs_offset): + current_index_shift = 0 #since we add and remove points in the result_list, we need to adjust the indices stored in forbidden_point_list_indices + for index in forbidden_point_list_indices: + #if abs(result_list[index][0].point.coords[0][0]-40.7) < 0.2 and abs(result_list[index][0].point.coords[0][1]-1.3)< 0.2: + # print("GEFUNDEN") + index+=current_index_shift + distance_left = result_list[index][0].point.distance(result_list[index-1][0].point)/2.0 + distance_right = result_list[index][0].point.distance(result_list[(index+1)%len(result_list)][0].point)/2.0 + while distance_left > constants.point_spacing_to_be_considered_equal and distance_right > constants.point_spacing_to_be_considered_equal: + new_point_left_proj = result_list[index][1]-distance_left + if new_point_left_proj < 0: + new_point_left_proj += line.length + new_point_right_proj = result_list[index][1]+distance_right + if new_point_right_proj > line.length: + new_point_right_proj-=line.length + point_left = line.interpolate(new_point_left_proj) + point_right = line.interpolate(new_point_right_proj) + forbidden_point_distance = result_list[index][0].point.distance(LineString([point_left, point_right])) + if forbidden_point_distance < constants.factor_offset_remove_dense_points*abs_offset: + del result_list[index] + result_list.insert(index, (PointTransfer.projected_point_tuple(point=point_right, point_source=\ + PointSource.REPLACED_FORBIDDEN_POINT),new_point_right_proj)) + result_list.insert(index, (PointTransfer.projected_point_tuple(point=point_left, point_source=\ + PointSource.REPLACED_FORBIDDEN_POINT),new_point_left_proj)) + current_index_shift+=1 + break + else: + distance_left/=2.0 + distance_right/=2.0 + return result_list + +if __name__ == "__main__": + line = LineString([(0,0), (1,0), (2,1),(3,0),(4,0)]) + + print(calculate_line_angles(line)*180.0/math.pi) diff --git a/lib/stitches/PointTransfer.py b/lib/stitches/PointTransfer.py new file mode 100644 index 00000000..998282a3 --- /dev/null +++ b/lib/stitches/PointTransfer.py @@ -0,0 +1,467 @@ +from shapely.geometry import Point, MultiPoint +from shapely.geometry.polygon import LineString, LinearRing +from collections import namedtuple +from shapely.ops import nearest_points +import math +from ..stitches import constants +from ..stitches import LineStringSampling + +projected_point_tuple = namedtuple('projected_point_tuple', ['point', 'point_source']) + +#Calculated the nearest interserction point of "bisectorline" with the coordinates of child (child.val). +#It returns the intersection point and its distance along the coordinates of the child or "None, None" if no +#intersection was found. +def calc_transferred_point(bisectorline, child): + result = bisectorline.intersection(child.val) + if result.is_empty: + return None, None + desired_point = Point() + if result.geom_type == 'Point': + desired_point = result + elif result.geom_type == 'LineString': + desired_point = Point(result.coords[0]) + else: + resultlist = list(result) + desired_point = resultlist[0] + if len(resultlist) > 1: + desired_point = nearest_points(result, Point(bisectorline.coords[0]))[0] + + priority = child.val.project(desired_point) + point = desired_point + return point, priority + + +#Takes the current tree item and its rastered points (to_transfer_points) and transfers these points to its parent, siblings and childs +# To do so it calculates the current normal and determines its intersection with the neighbors which gives the transferred points. +#Input: +#-treenode: Tree node whose points stored in "to_transfer_points" shall be transferred to its neighbors. +#-used_offset: The used offset when the curves where offsetted +#-offset_by_half: True if the transferred points shall be interlaced with respect to the points in "to_transfer_points" +#-max_stitching_distance: The maximum allowed stitch distance between two points +#-to_transfer_points: List of points belonging to treenode which shall be transferred - it is assumed that to_transfer_points can be handled as closed ring +#-to_transfer_points_origin: The origin tag of each point in to_transfer_points +#-overnext_neighbor: Transfer the points to the overnext neighbor (gives a more stable interlacing) +#-transfer_forbidden_points: Only allowed for interlacing (offset_by_half): Might be used to transfer points unshifted as forbidden points to the neighbor to avoid a point placing there +#-transfer_to_parent: If True, points will be transferred to the parent +#-transfer_to_sibling: If True, points will be transferred to the siblings +#-transfer_to_child: If True, points will be transferred to the childs +#Output: +#-Fills the attribute "transferred_point_priority_deque" of the siblings and parent in the tree datastructure. An item of the deque +#is setup as follows: ((projected point on line, LineStringSampling.PointSource), priority=distance along line) +#index of point_origin is the index of the point in the neighboring line +def transfer_points_to_surrounding(treenode, used_offset, offset_by_half, max_stitching_distance, to_transfer_points, to_transfer_points_origin=[], + overnext_neighbor = False, transfer_forbidden_points = False, transfer_to_parent=True, transfer_to_sibling=True, transfer_to_child=True): + + assert(len(to_transfer_points)==len(to_transfer_points_origin) or len(to_transfer_points_origin) == 0) + assert((overnext_neighbor and not offset_by_half) or not overnext_neighbor) + assert(not transfer_forbidden_points or transfer_forbidden_points and (offset_by_half or not offset_by_half and overnext_neighbor)) + + if len(to_transfer_points) == 0: + return + + # Get a list of all possible adjacent nodes which will be considered for transferring the points of treenode: + childs_tuple = treenode.children + siblings_tuple = treenode.siblings + # Take only neighbors which have not rastered before + # We need to distinguish between childs (project towards inner) and parent/siblings (project towards outer) + child_list = [] + child_list_forbidden = [] + neighbor_list = [] + neighbor_list_forbidden = [] + + if transfer_to_child: + for child in childs_tuple: + if child.already_rastered == False: + if not overnext_neighbor: + child_list.append(child) + if transfer_forbidden_points: + child_list_forbidden.append(child) + if overnext_neighbor: + for subchild in child.children: + if subchild.already_rastered == False: + child_list.append(subchild) + + if transfer_to_sibling: + for sibling in siblings_tuple: + if sibling.already_rastered == False: + if not overnext_neighbor: + neighbor_list.append(sibling) + if transfer_forbidden_points: + neighbor_list_forbidden.append(sibling) + if overnext_neighbor: + for subchild in sibling.children: + if subchild.already_rastered == False: + neighbor_list.append(subchild) + + if transfer_to_parent and treenode.parent != None: + if treenode.parent.already_rastered == False: + if not overnext_neighbor: + neighbor_list.append(treenode.parent) + if transfer_forbidden_points: + neighbor_list_forbidden.append(treenode.parent) + if overnext_neighbor: + if treenode.parent.parent != None: + if treenode.parent.parent.already_rastered == False: + neighbor_list.append(treenode.parent.parent) + + if not neighbor_list and not child_list: + return + + # Go through all rastered points of treenode and check where they should be transferred to its neighbar + point_list = list(MultiPoint(to_transfer_points)) + point_source_list = to_transfer_points_origin.copy() + + # For a linear ring the last point is the same as the starting point which we delete + # since we do not want to transfer the starting and end point twice + closed_line = LineString(to_transfer_points) + if point_list[0].distance(point_list[-1]) < constants.point_spacing_to_be_considered_equal: + point_list.pop() + if(point_source_list): + point_source_list.pop() + if len(point_list) == 0: + return + else: + # closed line is needed if we offset by half since we need to determine the line + # length including the closing segment + closed_line = LinearRing(to_transfer_points) + + bisectorline_length = abs(used_offset) * \ + constants.transfer_point_distance_factor*(2.0 if overnext_neighbor else 1.0) + + bisectorline_length_forbidden_points = abs(used_offset) * \ + constants.transfer_point_distance_factor + + linesign_child = math.copysign(1, used_offset) + + + i = 0 + currentDistance = 0 + while i < len(point_list): + assert(point_source_list[i] != LineStringSampling.PointSource.ENTER_LEAVING_POINT) + #if abs(point_list[i].coords[0][0]-47) < 0.3 and abs(point_list[i].coords[0][1]-4.5) < 0.3: + # print("HIIIIIIIIIIIERRR") + + # We create a bisecting line through the current point + normalized_vector_prev_x = ( + point_list[i].coords[0][0]-point_list[i-1].coords[0][0]) # makes use of closed shape + normalized_vector_prev_y = ( + point_list[i].coords[0][1]-point_list[i-1].coords[0][1]) + prev_spacing = math.sqrt(normalized_vector_prev_x*normalized_vector_prev_x + + normalized_vector_prev_y*normalized_vector_prev_y) + + normalized_vector_prev_x /= prev_spacing + normalized_vector_prev_y /= prev_spacing + + + normalized_vector_next_x = normalized_vector_next_y = 0 + next_spacing = 0 + while True: + normalized_vector_next_x = ( + point_list[i].coords[0][0]-point_list[(i+1) % len(point_list)].coords[0][0]) + normalized_vector_next_y = ( + point_list[i].coords[0][1]-point_list[(i+1) % len(point_list)].coords[0][1]) + next_spacing = math.sqrt(normalized_vector_next_x*normalized_vector_next_x + + normalized_vector_next_y*normalized_vector_next_y) + if next_spacing < constants.line_lengh_seen_as_one_point: + point_list.pop(i) + if(point_source_list): + point_source_list.pop(i) + currentDistance += next_spacing + continue + + normalized_vector_next_x /= next_spacing + normalized_vector_next_y /= next_spacing + break + + vecx = (normalized_vector_next_x+normalized_vector_prev_x) + vecy = (normalized_vector_next_y+normalized_vector_prev_y) + vec_length = math.sqrt(vecx*vecx+vecy*vecy) + + vecx_forbidden_point = vecx + vecy_forbidden_point = vecy + + # The two sides are (anti)parallel - construct normal vector (bisector) manually: + # If we offset by half we are offseting normal to the next segment + if(vec_length < constants.line_lengh_seen_as_one_point or offset_by_half): + vecx = linesign_child*bisectorline_length*normalized_vector_next_y + vecy = -linesign_child*bisectorline_length*normalized_vector_next_x + + if transfer_forbidden_points: + vecx_forbidden_point = linesign_child*bisectorline_length_forbidden_points*normalized_vector_next_y + vecy_forbidden_point = -linesign_child*bisectorline_length_forbidden_points*normalized_vector_next_x + + else: + vecx *= bisectorline_length/vec_length + vecy *= bisectorline_length/vec_length + + if (vecx*normalized_vector_next_y-vecy * normalized_vector_next_x)*linesign_child < 0: + vecx = -vecx + vecy = -vecy + vecx_forbidden_point = vecx + vecy_forbidden_point = vecy + + assert((vecx*normalized_vector_next_y-vecy * + normalized_vector_next_x)*linesign_child >= 0) + + originPoint = point_list[i] + originPoint_forbidden_point = point_list[i] + if(offset_by_half): + off = currentDistance+next_spacing/2 + if off > closed_line.length: + off -= closed_line.length + originPoint = closed_line.interpolate(off) + + bisectorline_child = LineString([(originPoint.coords[0][0], + originPoint.coords[0][1]), + (originPoint.coords[0][0]+vecx, + originPoint.coords[0][1]+vecy)]) + + bisectorline_neighbor = LineString([(originPoint.coords[0][0], + originPoint.coords[0][1]), + (originPoint.coords[0][0]-vecx, + originPoint.coords[0][1]-vecy)]) + + bisectorline_forbidden_point_child = LineString([(originPoint_forbidden_point.coords[0][0], + originPoint_forbidden_point.coords[0][1]), + (originPoint_forbidden_point.coords[0][0]+vecx_forbidden_point, + originPoint_forbidden_point.coords[0][1]+vecy_forbidden_point)]) + + bisectorline_forbidden_point_neighbor = LineString([(originPoint_forbidden_point.coords[0][0], + originPoint_forbidden_point.coords[0][1]), + (originPoint_forbidden_point.coords[0][0]-vecx_forbidden_point, + originPoint_forbidden_point.coords[0][1]-vecy_forbidden_point)]) + + for child in child_list: + point, priority = calc_transferred_point(bisectorline_child,child) + if point==None: + continue + child.transferred_point_priority_deque.insert(projected_point_tuple(point = point, point_source=LineStringSampling.PointSource.OVERNEXT if overnext_neighbor else LineStringSampling.PointSource.DIRECT), priority) + for child in child_list_forbidden: + point, priority = calc_transferred_point(bisectorline_forbidden_point_child,child) + if point == None: + continue + child.transferred_point_priority_deque.insert(projected_point_tuple(point=point, point_source=LineStringSampling.PointSource.FORBIDDEN_POINT), priority) + + for neighbor in neighbor_list: + point, priority = calc_transferred_point(bisectorline_neighbor,neighbor) + if point==None: + continue + neighbor.transferred_point_priority_deque.insert(projected_point_tuple(point = point, point_source=LineStringSampling.PointSource.OVERNEXT if overnext_neighbor else LineStringSampling.PointSource.DIRECT), priority) + for neighbor in neighbor_list_forbidden: + point, priority = calc_transferred_point(bisectorline_forbidden_point_neighbor,neighbor) + if point == None: + continue + neighbor.transferred_point_priority_deque.insert(projected_point_tuple(point=point, point_source=LineStringSampling.PointSource.FORBIDDEN_POINT), priority) + + i += 1 + currentDistance += next_spacing + + assert(len(point_list) == len(point_source_list)) + +#Calculated the nearest interserction point of "bisectorline" with the coordinates of child. +#It returns the intersection point and its distance along the coordinates of the child or "None, None" if no +#intersection was found. +def calc_transferred_point_graph(bisectorline, edge_geometry): + result = bisectorline.intersection(edge_geometry) + if result.is_empty: + return None, None + desired_point = Point() + if result.geom_type == 'Point': + desired_point = result + elif result.geom_type == 'LineString': + desired_point = Point(result.coords[0]) + else: + resultlist = list(result) + desired_point = resultlist[0] + if len(resultlist) > 1: + desired_point = nearest_points(result, Point(bisectorline.coords[0]))[0] + + priority = edge_geometry.project(desired_point) + point = desired_point + return point, priority + + +#Takes the current tree item and its rastered points (to_transfer_points) and transfers these points to its parent, siblings and childs +# To do so it calculates the current normal and determines its intersection with the neighbors which gives the transferred points. +#Input: +#-treenode: Tree node whose points stored in "to_transfer_points" shall be transferred to its neighbors. +#-used_offset: The used offset when the curves where offsetted +#-offset_by_half: True if the transferred points shall be interlaced with respect to the points in "to_transfer_points" +#-max_stitching_distance: The maximum allowed stitch distance between two points +#-to_transfer_points: List of points belonging to treenode which shall be transferred - it is assumed that to_transfer_points can be handled as closed ring +#-to_transfer_points_origin: The origin tag of each point in to_transfer_points +#-overnext_neighbor: Transfer the points to the overnext neighbor (gives a more stable interlacing) +#-transfer_forbidden_points: Only allowed for interlacing (offset_by_half): Might be used to transfer points unshifted as forbidden points to the neighbor to avoid a point placing there +#-transfer_to_parent: If True, points will be transferred to the parent +#-transfer_to_sibling: If True, points will be transferred to the siblings +#-transfer_to_child: If True, points will be transferred to the childs +#Output: +#-Fills the attribute "transferred_point_priority_deque" of the siblings and parent in the tree datastructure. An item of the deque +#is setup as follows: ((projected point on line, LineStringSampling.PointSource), priority=distance along line) +#index of point_origin is the index of the point in the neighboring line +def transfer_points_to_surrounding_graph(fill_stitch_graph, current_edge, used_offset, offset_by_half, to_transfer_points, + overnext_neighbor = False, transfer_forbidden_points = False, transfer_to_previous=True, transfer_to_next=True): + + assert((overnext_neighbor and not offset_by_half) or not overnext_neighbor) + assert(not transfer_forbidden_points or transfer_forbidden_points and (offset_by_half or not offset_by_half and overnext_neighbor)) + + if len(to_transfer_points) == 0: + return + + + # Take only neighbors which have not rastered before + # We need to distinguish between childs (project towards inner) and parent/siblings (project towards outer) + previous_edge_list = [] + previous_edge_list_forbidden = [] + next_edge_list = [] + next_edge_list_forbidden = [] + + if transfer_to_previous: + previous_neighbors_tuples = current_edge['previous_neighbors'] + for neighbor in previous_neighbors_tuples: + neighbor_edge = fill_stitch_graph[neighbor[0]][neighbor[-1]]['segment'] + if not neighbor_edge['already_rastered']: + if not overnext_neighbor: + previous_edge_list.append(neighbor_edge) + if transfer_forbidden_points: + previous_edge_list_forbidden.append(neighbor_edge) + if overnext_neighbor: + overnext_previous_neighbors_tuples = neighbor_edge['previous_neighbors'] + for overnext_neighbor in overnext_previous_neighbors_tuples: + overnext_neighbor_edge = fill_stitch_graph[overnext_neighbor[0]][overnext_neighbor[-1]]['segment'] + if not overnext_neighbor_edge['already_rastered']: + previous_edge_list.append(overnext_neighbor_edge) + + if transfer_to_next: + next_neighbors_tuples = current_edge['next_neighbors'] + for neighbor in next_neighbors_tuples: + neighbor_edge = fill_stitch_graph[neighbor[0]][neighbor[-1]]['segment'] + if not neighbor_edge['already_rastered']: + if not overnext_neighbor: + next_edge_list.append(neighbor_edge) + if transfer_forbidden_points: + next_edge_list_forbidden.append(neighbor_edge) + if overnext_neighbor: + overnext_next_neighbors_tuples = neighbor_edge['next_neighbors'] + for overnext_neighbor in overnext_next_neighbors_tuples: + overnext_neighbor_edge = fill_stitch_graph[overnext_neighbor[0]][overnext_neighbor[-1]]['segment'] + if not overnext_neighbor_edge['already_rastered']: + next_edge_list.append(overnext_neighbor_edge) + + + if not previous_edge_list and not next_edge_list: + return + + # Go through all rastered points of treenode and check where they should be transferred to its neighbar + point_list = list(MultiPoint(to_transfer_points)) + line = LineString(to_transfer_points) + + bisectorline_length = abs(used_offset) * \ + constants.transfer_point_distance_factor*(2.0 if overnext_neighbor else 1.0) + + bisectorline_length_forbidden_points = abs(used_offset) * \ + constants.transfer_point_distance_factor + + linesign_child = math.copysign(1, used_offset) + + + i = 0 + currentDistance = 0 + while i < len(point_list): + + #if abs(point_list[i].coords[0][0]-47) < 0.3 and abs(point_list[i].coords[0][1]-4.5) < 0.3: + # print("HIIIIIIIIIIIERRR") + + # We create a bisecting line through the current point + normalized_vector_prev_x = ( + point_list[i].coords[0][0]-point_list[i-1].coords[0][0]) # makes use of closed shape + normalized_vector_prev_y = ( + point_list[i].coords[0][1]-point_list[i-1].coords[0][1]) + prev_spacing = math.sqrt(normalized_vector_prev_x*normalized_vector_prev_x + + normalized_vector_prev_y*normalized_vector_prev_y) + + normalized_vector_prev_x /= prev_spacing + normalized_vector_prev_y /= prev_spacing + + + normalized_vector_next_x = normalized_vector_next_y = 0 + next_spacing = 0 + while True: + normalized_vector_next_x = ( + point_list[i].coords[0][0]-point_list[(i+1) % len(point_list)].coords[0][0]) + normalized_vector_next_y = ( + point_list[i].coords[0][1]-point_list[(i+1) % len(point_list)].coords[0][1]) + next_spacing = math.sqrt(normalized_vector_next_x*normalized_vector_next_x + + normalized_vector_next_y*normalized_vector_next_y) + if next_spacing < constants.line_lengh_seen_as_one_point: + point_list.pop(i) + currentDistance += next_spacing + continue + + normalized_vector_next_x /= next_spacing + normalized_vector_next_y /= next_spacing + break + + vecx = (normalized_vector_next_x+normalized_vector_prev_x) + vecy = (normalized_vector_next_y+normalized_vector_prev_y) + vec_length = math.sqrt(vecx*vecx+vecy*vecy) + + vecx_forbidden_point = vecx + vecy_forbidden_point = vecy + + # The two sides are (anti)parallel - construct normal vector (bisector) manually: + # If we offset by half we are offseting normal to the next segment + if(vec_length < constants.line_lengh_seen_as_one_point or offset_by_half): + vecx = linesign_child*bisectorline_length*normalized_vector_next_y + vecy = -linesign_child*bisectorline_length*normalized_vector_next_x + + if transfer_forbidden_points: + vecx_forbidden_point = linesign_child*bisectorline_length_forbidden_points*normalized_vector_next_y + vecy_forbidden_point = -linesign_child*bisectorline_length_forbidden_points*normalized_vector_next_x + + else: + vecx *= bisectorline_length/vec_length + vecy *= bisectorline_length/vec_length + + if (vecx*normalized_vector_next_y-vecy * normalized_vector_next_x)*linesign_child < 0: + vecx = -vecx + vecy = -vecy + vecx_forbidden_point = vecx + vecy_forbidden_point = vecy + + assert((vecx*normalized_vector_next_y-vecy * + normalized_vector_next_x)*linesign_child >= 0) + + originPoint = point_list[i] + originPoint_forbidden_point = point_list[i] + if(offset_by_half): + off = currentDistance+next_spacing/2 + if off > line.length: + break + originPoint = line.interpolate(off) + + bisectorline = LineString([(originPoint.coords[0][0]-vecx, + originPoint.coords[0][1]-vecy), + (originPoint.coords[0][0]+vecx, + originPoint.coords[0][1]+vecy)]) + + bisectorline_forbidden_point = LineString([(originPoint_forbidden_point.coords[0][0]-vecx_forbidden_point, + originPoint_forbidden_point.coords[0][1]-vecy_forbidden_point), + (originPoint_forbidden_point.coords[0][0]+vecx_forbidden_point, + originPoint_forbidden_point.coords[0][1]+vecy_forbidden_point)]) + + + for edge in previous_edge_list+next_edge_list: + point, priority = calc_transferred_point_graph(bisectorline,edge['geometry']) + if point==None: + continue + edge['projected_points'].insert(projected_point_tuple(point = point, point_source=LineStringSampling.PointSource.OVERNEXT if overnext_neighbor else LineStringSampling.PointSource.DIRECT), priority) + for edge_forbidden in previous_edge_list_forbidden+next_edge_list_forbidden: + point, priority = calc_transferred_point_graph(bisectorline_forbidden_point,edge_forbidden['geometry']) + if point == None: + continue + edge_forbidden['projected_points'].insert(projected_point_tuple(point=point, point_source=LineStringSampling.PointSource.FORBIDDEN_POINT), priority) + + + i += 1 + currentDistance += next_spacing diff --git a/lib/stitches/StitchPattern.py b/lib/stitches/StitchPattern.py new file mode 100644 index 00000000..d0a3f7aa --- /dev/null +++ b/lib/stitches/StitchPattern.py @@ -0,0 +1,223 @@ +from shapely.geometry.polygon import LinearRing, LineString +from shapely.geometry import Polygon, MultiLineString +from shapely.ops import polygonize +from shapely.geometry import MultiPolygon +from anytree import AnyNode, PreOrderIter +from shapely.geometry.polygon import orient +from depq import DEPQ +from enum import IntEnum +from ..stitches import ConnectAndSamplePattern +from ..stitches import constants + + + +# Problem: When shapely offsets a LinearRing the start/end point might be handled wrongly since they are only treated as LineString. +# (See e.g. https://i.stack.imgur.com/vVh56.png as a problematic example) +# This method checks first whether the start/end point form a problematic edge with respect to the offset side. If it is not a problematic +# edge we can use the normal offset_routine. Otherwise we need to perform two offsets: +# -offset the ring +# -offset the start/end point + its two neighbors left and right +# Finally both offsets are merged together to get the correct offset of a LinearRing +def offset_linear_ring(ring, offset, side, resolution, join_style, mitre_limit): + coords = ring.coords[:] + # check whether edge at index 0 is concave or convex. Only for concave edges we need to spend additional effort + dx_seg1 = dy_seg1 = 0 + if coords[0] != coords[-1]: + dx_seg1 = coords[0][0]-coords[-1][0] + dy_seg1 = coords[0][1]-coords[-1][1] + else: + dx_seg1 = coords[0][0]-coords[-2][0] + dy_seg1 = coords[0][1]-coords[-2][1] + dx_seg2 = coords[1][0]-coords[0][0] + dy_seg2 = coords[1][1]-coords[0][1] + # use cross product: + crossvalue = dx_seg1*dy_seg2-dy_seg1*dx_seg2 + sidesign = 1 + if side == 'left': + sidesign = -1 + + # We do not need to take care of the joint n-0 since we offset along a concave edge: + if sidesign*offset*crossvalue <= 0: + return ring.parallel_offset(offset, side, resolution, join_style, mitre_limit) + + # We offset along a convex edge so we offset the joint n-0 separately: + if coords[0] != coords[-1]: + coords.append(coords[0]) + offset_ring1 = ring.parallel_offset( + offset, side, resolution, join_style, mitre_limit) + offset_ring2 = LineString((coords[-2], coords[0], coords[1])).parallel_offset( + offset, side, resolution, join_style, mitre_limit) + + # Next we need to merge the results: + if offset_ring1.geom_type == 'LineString': + return LinearRing(offset_ring2.coords[:]+offset_ring1.coords[1:-1]) + else: + # We have more than one resulting LineString for offset of the geometry (ring) = offset_ring1. + # Hence we need to find the LineString which belongs to the offset of element 0 in coords =offset_ring2 + # in order to add offset_ring2 geometry to it: + result_list = [] + thresh = constants.offset_factor_for_adjacent_geometry*abs(offset) + for offsets in offset_ring1: + if(abs(offsets.coords[0][0]-coords[0][0]) < thresh and abs(offsets.coords[0][1]-coords[0][1]) < thresh): + result_list.append(LinearRing( + offset_ring2.coords[:]+offsets.coords[1:-1])) + else: + result_list.append(LinearRing(offsets)) + return MultiLineString(result_list) + + +# Removes all geometries which do not form a "valid" LinearRing (meaning a ring which does not form a straight line) +def take_only_valid_linear_rings(rings): + if(rings.geom_type == 'MultiLineString'): + new_list = [] + for ring in rings: + if len(ring.coords) > 3 or (len(ring.coords) == 3 and ring.coords[0] != ring.coords[-1]): + new_list.append(ring) + if len(new_list) == 1: + return LinearRing(new_list[0]) + else: + return MultiLineString(new_list) + else: + if len(rings.coords) <= 2: + return LinearRing() + elif len(rings.coords) == 3 and rings.coords[0] == rings.coords[-1]: + return LinearRing() + else: + return rings + + +#Since naturally holes have the opposite point ordering than non-holes we make +#all lines within the tree "root" uniform (having all the same ordering direction) +def make_tree_uniform_ccw(root): + for node in PreOrderIter(root): + if(node.id == 'hole'): + node.val.coords = list(node.val.coords)[::-1] + + +#Used to define which stitching strategy shall be used +class StitchingStrategy(IntEnum): + CLOSEST_POINT = 0 + INNER_TO_OUTER = 1 + +# Takes a polygon (which can have holes) as input and creates offsetted versions until the polygon is filled with these smaller offsets. +# These created geometries are afterwards connected to each other and resampled with a maximum stitch_distance. +# The return value is a LineString which should cover the full polygon. +#Input: +#-poly: The shapely polygon which can have holes +#-offset: The used offset for the curves +#-join_style: Join style for the offset - can be round, mitered or bevel (https://shapely.readthedocs.io/en/stable/manual.html#shapely.geometry.JOIN_STYLE) +#For examples look at https://shapely.readthedocs.io/en/stable/_images/parallel_offset.png +#-stitch_distance maximum allowed stitch distance between two points +#-offset_by_half: True if the points shall be interlaced +#-strategy: According to StitchingStrategy you can select between different strategies for the connection between parent and childs +#Output: +#-List of point coordinate tuples +#-Tag (origin) of each point to analyze why a point was placed at this position +def offset_poly(poly, offset, join_style, stitch_distance, offset_by_half, strategy, starting_point): + ordered_poly = orient(poly, -1) + ordered_poly = ordered_poly.simplify( + constants.simplification_threshold, False) + root = AnyNode(id="node", val=ordered_poly.exterior, already_rastered=False, transferred_point_priority_deque=DEPQ( + iterable=None, maxlen=None)) + active_polys = [root] + active_holes = [[]] + + for holes in ordered_poly.interiors: + #print("hole: - is ccw: ", LinearRing(holes).is_ccw) + active_holes[0].append( + AnyNode(id="hole", val=holes, already_rastered=False, transferred_point_priority_deque=DEPQ( + iterable=None, maxlen=None))) + + # counter = 0 + while len(active_polys) > 0: # and counter < 20: + # counter += 1 + # print("New iter") + current_poly = active_polys.pop() + current_holes = active_holes.pop() + poly_inners = [] + + # outer = current_poly.val.parallel_offset(offset,'left', 5, join_style, 10) + outer = offset_linear_ring(current_poly.val, offset, 'left', 5, join_style, 10) + outer = outer.simplify(constants.simplification_threshold, False) + outer = take_only_valid_linear_rings(outer) + + for j in range(len(current_holes)): + # inner = closeLinearRing(current_holes[j].val,offset/2.0).parallel_offset(offset,'left', 5, join_style, 10) + inner = offset_linear_ring( + current_holes[j].val, offset, 'left', 5, join_style, 10) + inner = inner.simplify(constants.simplification_threshold, False) + inner = take_only_valid_linear_rings(inner) + if not inner.is_empty: + poly_inners.append(Polygon(inner)) + if not outer.is_empty: + if len(poly_inners) == 0: + if outer.geom_type == 'LineString': + result = Polygon(outer) + else: + result = MultiPolygon(polygonize(outer)) + else: + if outer.geom_type == 'LineString': + result = Polygon(outer).difference( + MultiPolygon(poly_inners)) + else: + result = MultiPolygon(outer).difference( + MultiPolygon(poly_inners)) + + if not result.is_empty and result.area > offset*offset/10: + result_list = [] + if result.geom_type == 'Polygon': + result_list = [result] + else: + result_list = list(result) + # print("New result_list: ", len(result_list)) + for polygon in result_list: + polygon = orient(polygon, -1) + + if polygon.area < offset*offset/10: + continue + + polygon = polygon.simplify(constants.simplification_threshold, False) + poly_coords = polygon.exterior + # if polygon.exterior.is_ccw: + # hole.coords = list(hole.coords)[::-1] + #poly_coords = polygon.exterior.simplify(constants.simplification_threshold, False) + poly_coords = take_only_valid_linear_rings(poly_coords) + if poly_coords.is_empty: + continue + #print("node: - is ccw: ", LinearRing(poly_coords).is_ccw) + # if(LinearRing(poly_coords).is_ccw): + # print("Fehler!") + node = AnyNode(id="node", parent=current_poly, + val=poly_coords, already_rastered=False, transferred_point_priority_deque=DEPQ( + iterable=None, maxlen=None)) + active_polys.append(node) + hole_node_list = [] + for hole in polygon.interiors: + hole_node = AnyNode( + id="hole", val=hole, already_rastered=False, transferred_point_priority_deque=DEPQ( + iterable=None, maxlen=None)) + for previous_hole in current_holes: + if Polygon(hole).contains(Polygon(previous_hole.val)): + previous_hole.parent = hole_node + hole_node_list.append(hole_node) + active_holes.append(hole_node_list) + 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 previous_hole.parent == None: + previous_hole.parent = current_poly + + + #DebuggingMethods.drawPoly(root, 'r-') + + make_tree_uniform_ccw(root) + # print(RenderTree(root)) + if strategy == StitchingStrategy.CLOSEST_POINT: + connected_line, connected_line_origin = ConnectAndSamplePattern.connect_raster_tree_nearest_neighbor( + root, offset, stitch_distance, starting_point, offset_by_half) + elif strategy == StitchingStrategy.INNER_TO_OUTER: + connected_line, connected_line_origin = ConnectAndSamplePattern.connect_raster_tree_from_inner_to_outer( + root, offset, stitch_distance, starting_point, offset_by_half) + else: + print("Invalid strategy!") + assert(0) + + return connected_line, connected_line_origin diff --git a/lib/stitches/auto_fill.py b/lib/stitches/auto_fill.py index 160d927e..71cfd80f 100644 --- a/lib/stitches/auto_fill.py +++ b/lib/stitches/auto_fill.py @@ -12,14 +12,17 @@ import networkx from shapely import geometry as shgeo from shapely.ops import snap from shapely.strtree import STRtree - +from depq import DEPQ from ..debug import debug from ..stitch_plan import Stitch from ..svg import PIXELS_PER_MM +from ..utils import geometry from ..utils.geometry import Point as InkstitchPoint from ..utils.geometry import line_string_to_point_list -from .fill import intersect_region_with_grating, stitch_row +from .fill import intersect_region_with_grating, intersect_region_with_grating_line, stitch_row from .running_stitch import running_stitch +from .PointTransfer import transfer_points_to_surrounding_graph +from .LineStringSampling import raster_line_string_with_priority_points_graph class PathEdge(object): @@ -49,6 +52,7 @@ class PathEdge(object): @debug.time def auto_fill(shape, + line, angle, row_spacing, end_row_spacing, @@ -58,10 +62,13 @@ def auto_fill(shape, skip_last, starting_point, ending_point=None, - underpath=True): + underpath=True, + offset_by_half=True): + #offset_by_half only relevant for line != None; staggers only relevant for line == None! + fill_stitch_graph = [] try: - fill_stitch_graph = build_fill_stitch_graph(shape, angle, row_spacing, end_row_spacing, starting_point, ending_point) + fill_stitch_graph = build_fill_stitch_graph(shape, line, angle, row_spacing, end_row_spacing, starting_point, ending_point) except ValueError: # Small shapes will cause the graph to fail - min() arg is an empty sequence through insert node return fallback(shape, running_stitch_length) @@ -72,7 +79,7 @@ def auto_fill(shape, travel_graph = build_travel_graph(fill_stitch_graph, shape, angle, underpath) path = find_stitch_path(fill_stitch_graph, travel_graph, starting_point, ending_point) result = path_to_stitches(path, travel_graph, fill_stitch_graph, angle, row_spacing, - max_stitch_length, running_stitch_length, staggers, skip_last) + max_stitch_length, running_stitch_length, staggers, skip_last,line!=None,offset_by_half) return result @@ -106,7 +113,7 @@ def project(shape, coords, outline_index): @debug.time -def build_fill_stitch_graph(shape, angle, row_spacing, end_row_spacing, starting_point=None, ending_point=None): +def build_fill_stitch_graph(shape, line, angle, row_spacing, end_row_spacing, starting_point=None, ending_point=None): """build a graph representation of the grating segments This function builds a specialized graph (as in graph theory) that will @@ -141,18 +148,34 @@ def build_fill_stitch_graph(shape, angle, row_spacing, end_row_spacing, starting debug.add_layer("auto-fill fill stitch") - # Convert the shape into a set of parallel line segments. - rows_of_segments = intersect_region_with_grating(shape, angle, row_spacing, end_row_spacing) - segments = [segment for row in rows_of_segments for segment in row] + if line == None: + # Convert the shape into a set of parallel line segments. + rows_of_segments = intersect_region_with_grating(shape, angle, row_spacing, end_row_spacing) + else: + rows_of_segments = intersect_region_with_grating_line(shape, line, row_spacing, end_row_spacing) + + #segments = [segment for row in rows_of_segments for segment in row] graph = networkx.MultiGraph() - # First, add the grating segments as edges. We'll use the coordinates - # of the endpoints as nodes, which networkx will add automatically. - for segment in segments: - # networkx allows us to label nodes with arbitrary data. We'll - # mark this one as a grating segment. - graph.add_edge(*segment, key="segment", underpath_edges=[]) + + for i in range(len(rows_of_segments)): + for segment in rows_of_segments[i]: + # First, add the grating segments as edges. We'll use the coordinates + # of the endpoints as nodes, which networkx will add automatically. + + # networkx allows us to label nodes with arbitrary data. We'll + # mark this one as a grating segment. + #graph.add_edge(*segment, key="segment", underpath_edges=[]) + previous_neighbors_ = [(seg[0],seg[-1]) for seg in rows_of_segments[i-1] if i > 0] + next_neighbors_ = [(seg[0],seg[-1]) for seg in rows_of_segments[(i+1)% len(rows_of_segments)] if i < len(rows_of_segments)-1] + + graph.add_edge(segment[0],segment[-1], key="segment", underpath_edges=[], + geometry=shgeo.LineString(segment), previous_neighbors = previous_neighbors_, next_neighbors = next_neighbors_, + projected_points=DEPQ(iterable=None, maxlen=None), already_rastered=False) + + +#fill_stitch_graph[start][end]['segment']['underpath_edges'].append(edge) tag_nodes_with_outline_and_projection(graph, shape, graph.nodes()) add_edges_between_outline_nodes(graph, duplicate_every_other=True) @@ -325,7 +348,8 @@ def get_segments(graph): segments = [] for start, end, key, data in graph.edges(keys=True, data=True): if key == 'segment': - segments.append(shgeo.LineString((start, end))) + segments.append(data["geometry"]) + #segments.append(shgeo.LineString((start, end))) return segments @@ -363,7 +387,8 @@ def process_travel_edges(graph, fill_stitch_graph, shape, travel_edges): # segments that _might_ intersect ls. Refining the result is # necessary but the STRTree still saves us a ton of time. if segment.crosses(ls): - start, end = segment.coords + start = segment.coords[0] + end = segment.coords[-1] fill_stitch_graph[start][end]['segment']['underpath_edges'].append(edge) # The weight of a travel edge is the length of the line segment. @@ -614,9 +639,28 @@ def travel(travel_graph, start, end, running_stitch_length, skip_last): # stitch. return stitches[1:] +def stitch_line(stitches, stitching_direction, geometry,projected_points, max_stitch_length,row_spacing,skip_last,offset_by_half): + #print(start_point) + #print(geometry[0]) + #if stitching_direction == -1: + # geometry.coords = geometry.coords[::-1] + stitched_line, stitched_line_origin = raster_line_string_with_priority_points_graph(geometry,max_stitch_length,stitching_direction,projected_points,abs(row_spacing),offset_by_half) + + + stitches.append(Stitch(*stitched_line[0], tags=('fill_row_start',))) + for i in range(1,len(stitched_line)): + stitches.append(Stitch(*stitched_line[i], tags=('fill_row'))) + + if not skip_last: + if stitching_direction==1: + stitches.append(Stitch(*geometry.coords[-1], tags=('fill_row_end',))) + else: + stitches.append(Stitch(*geometry.coords[0], tags=('fill_row_end',))) + @debug.time -def path_to_stitches(path, travel_graph, fill_stitch_graph, angle, row_spacing, max_stitch_length, running_stitch_length, staggers, skip_last): +def path_to_stitches(path, travel_graph, fill_stitch_graph, angle, row_spacing, max_stitch_length, + running_stitch_length, staggers, skip_last, offsetted_line, offset_by_half): path = collapse_sequential_outline_edges(path) stitches = [] @@ -627,7 +671,23 @@ def path_to_stitches(path, travel_graph, fill_stitch_graph, angle, row_spacing, for edge in path: if edge.is_segment(): - stitch_row(stitches, edge[0], edge[1], angle, row_spacing, max_stitch_length, staggers, skip_last) + if offsetted_line: + new_stitches = [] + current_edge = fill_stitch_graph[edge[0]][edge[-1]]['segment'] + path_geometry = current_edge['geometry'] + projected_points = current_edge['projected_points'] + stitching_direction = 1 + if (abs(edge[0][0]-path_geometry.coords[0][0])+abs(edge[0][1]-path_geometry.coords[0][1]) > + abs(edge[0][0]-path_geometry.coords[-1][0])+abs(edge[0][1]-path_geometry.coords[-1][1])): + stitching_direction = -1 + stitch_line(new_stitches, stitching_direction, path_geometry,projected_points, max_stitch_length,row_spacing,skip_last,offset_by_half) + current_edge['already_rastered'] = True + transfer_points_to_surrounding_graph(fill_stitch_graph,current_edge,row_spacing,False,new_stitches,overnext_neighbor=True) + transfer_points_to_surrounding_graph(fill_stitch_graph,current_edge,row_spacing,offset_by_half,new_stitches,overnext_neighbor=False,transfer_forbidden_points=offset_by_half) + + stitches.extend(new_stitches) + else: + stitch_row(stitches, edge[0], edge[1], angle, row_spacing, max_stitch_length, staggers, skip_last) travel_graph.remove_edges_from(fill_stitch_graph[edge[0]][edge[1]]['segment'].get('underpath_edges', [])) else: stitches.extend(travel(travel_graph, edge[0], edge[1], running_stitch_length, skip_last)) diff --git a/lib/stitches/constants.py b/lib/stitches/constants.py new file mode 100644 index 00000000..63746310 --- /dev/null +++ b/lib/stitches/constants.py @@ -0,0 +1,41 @@ +import math + +# Used in the simplify routine of shapely +simplification_threshold = 0.01 + +# If a transferred point is closer than this value to one of its neighbors, it will be checked whether it can be removed +distance_thresh_remove_transferred_point = 0.15 + +# If a line segment is shorter than this threshold it is handled as a single point +line_lengh_seen_as_one_point = 0.05 + +# E.g. to check whether a point is already present in a point list, the point is allowed to be this value in distance apart +point_spacing_to_be_considered_equal = 0.05 + +# Adjacent geometry should have points closer than offset*offset_factor_for_adjacent_geometry to be considered adjacent +offset_factor_for_adjacent_geometry = 1.5 + +# Transfer point distance is used for projecting points from already rastered geometry to adjacent geometry +# (max spacing transfer_point_distance_factor*offset) to get a more regular pattern +transfer_point_distance_factor = 1.5 + +# Used to handle numerical inaccuracies during comparisons +eps = 1E-3 + +factor_offset_starting_points=0.5 #When entering and leaving a child from a parent we introduce an offset of abs_offset*factor_offset_starting_points so + #that entering and leaving points are not lying above each other. + +factor_offset_remove_points=0.5 #if points are closer than abs_offset*factor_offset_remove_points one of it is removed + +fac_offset_edge_shift = 0.25 #if an unshifted relevant edge is closer than abs_offset*fac_offset_edge_shift to the line segment created by the shifted edge, + #the shift is allowed - otherwise the edge must not be shifted. + +limiting_angle = math.pi*15/180.0 #decides whether the point belongs to a hard edge (must use this point during sampling) or soft edge (do not necessarily need to use this point) +limiting_angle_straight = math.pi*0.5/180.0 #angles straighter (smaller) than this are considered as more or less straight (no concrete edges required for path segments having only angles <= this value) + + +factor_offset_remove_dense_points=0.2 #if a point distance to the connected line of its two neighbors is smaller than abs_offset times this factor, this point will be removed if the stitching distance will not be exceeded + +factor_offset_forbidden_point = 1.0 #if a soft edge is closer to a forbidden point than abs_offset*this factor it will be marked as forbidden. + +factor_segment_length_direct_preferred_over_overnext = 0.5 #usually overnext projected points are preferred. If an overnext projected point would create a much smaller segment than a direct projected point we might prefer the direct projected point diff --git a/lib/stitches/fill.py b/lib/stitches/fill.py index 21e35d83..4e1669e9 100644 --- a/lib/stitches/fill.py +++ b/lib/stitches/fill.py @@ -6,12 +6,11 @@ import math import shapely - -from ..stitch_plan import Stitch +from shapely.geometry.linestring import LineString from ..svg import PIXELS_PER_MM from ..utils import Point as InkstitchPoint from ..utils import cache - +from ..stitch_plan import Stitch def legacy_fill(shape, angle, row_spacing, end_row_spacing, max_stitch_length, flip, staggers, skip_last): rows_of_segments = intersect_region_with_grating(shape, angle, row_spacing, end_row_spacing, flip) @@ -89,6 +88,65 @@ def stitch_row(stitches, beg, end, angle, row_spacing, max_stitch_length, stagge if (end - stitches[-1]).length() > 0.1 * PIXELS_PER_MM and not skip_last: stitches.append(end) +def extend_line(line, minx,maxx,miny,maxy): + line = line.simplify(0.01, False) + + upper_left = InkstitchPoint(minx, miny) + lower_right = InkstitchPoint(maxx, maxy) + length = (upper_left - lower_right).length() + + point1 = InkstitchPoint(*line.coords[0]) + point2 = InkstitchPoint(*line.coords[1]) + new_starting_point = point1-(point2-point1).unit()*length + + point3 = InkstitchPoint(*line.coords[-2]) + point4 = InkstitchPoint(*line.coords[-1]) + new_ending_point = point4+(point4-point3).unit()*length + + line = LineString([new_starting_point.as_tuple()]+line.coords[1:-1]+[new_ending_point.as_tuple()]) + + +def intersect_region_with_grating_line(shape, line, row_spacing, end_row_spacing=None, flip=False): + + row_spacing = abs(row_spacing) + (minx, miny, maxx, maxy) = shape.bounds + upper_left = InkstitchPoint(minx, miny) + rows = [] + extend_line(line, minx,maxx,miny,maxy) #extend the line towards the ends to increase probability that all offsetted curves cross the shape + + line_offsetted = line + res = line_offsetted.intersection(shape) + while isinstance(res, (shapely.geometry.GeometryCollection, shapely.geometry.MultiLineString)) or (not res.is_empty and len(res.coords) > 1): + if isinstance(res, (shapely.geometry.GeometryCollection, shapely.geometry.MultiLineString)): + runs = [line_string.coords for line_string in res.geoms if (not line_string.is_empty and len(line_string.coords) > 1)] + else: + runs = [res.coords] + + runs.sort(key=lambda seg: (InkstitchPoint(*seg[0]) - upper_left).length()) + if flip: + runs.reverse() + runs = [tuple(reversed(run)) for run in runs] + + if row_spacing > 0: + rows.append(runs) + else: + rows.insert(0,runs) + line_offsetted = line_offsetted.parallel_offset(row_spacing,'left',5) + if row_spacing < 0: + line_offsetted.coords = line_offsetted.coords[::-1] + line_offsetted = line_offsetted.simplify(0.01, False) + res = line_offsetted.intersection(shape) + if row_spacing > 0 and not isinstance(res, (shapely.geometry.GeometryCollection, shapely.geometry.MultiLineString)): + if (res.is_empty or len(res.coords) == 1): + row_spacing = -row_spacing + #print("Set to right") + line_offsetted = line.parallel_offset(row_spacing,'left',5) + line_offsetted.coords = line_offsetted.coords[::-1] #using negative row spacing leads as a side effect to reversed offsetted lines - here we undo this + line_offsetted = line_offsetted.simplify(0.01, False) + res = line_offsetted.intersection(shape) + + return rows + def intersect_region_with_grating(shape, angle, row_spacing, end_row_spacing=None, flip=False): # the max line length I'll need to intersect the whole shape is the diagonal |