From fa9745d8058657e92078ebeaffd76b8d2dcc782f Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Thu, 9 May 2024 18:24:24 -0700 Subject: [PATCH 01/14] check in wip code for first hobby splines --- solid/extrude_along_path.py | 156 +++++++--- solid/solidpython.py | 287 ++++++++++------- solid/splines.py | 605 +++++++++++++++++++++++++++++------- solid/test/test_splines.py | 166 ++++++++-- 4 files changed, 899 insertions(+), 315 deletions(-) diff --git a/solid/extrude_along_path.py b/solid/extrude_along_path.py index bdd94760..d88b939e 100644 --- a/solid/extrude_along_path.py +++ b/solid/extrude_along_path.py @@ -10,38 +10,43 @@ FacetIndices = Tuple[int, int, int] Point3Transform = Callable[[Point3, Optional[float], Optional[float]], Point3] + # ========================== # = Extrusion along a path = # ========================== -def extrude_along_path( shape_pts:Points, - path_pts:Points, - scales:Sequence[Union[Vector2, float, Tuple2]] = None, - rotations: Sequence[float] = None, - transforms: Sequence[Point3Transform] = None, - connect_ends = False, - cap_ends = True) -> OpenSCADObject: - ''' +def extrude_along_path( + shape_pts: Points, + path_pts: Points, + scales: Sequence[Union[Vector2, float, Tuple2]] = None, + rotations: Sequence[float] = None, + transforms: Sequence[Callable] = None, + connect_ends=False, + cap_ends=True, + transform_args: List = None, + transform_kwargs: Dict = None, +) -> OpenSCADObject: + """ Extrude the curve defined by shape_pts along path_pts. -- For predictable results, shape_pts must be planar, convex, and lie in the XY plane centered around the origin. *Some* nonconvexity (e.g, star shapes) and nonplanarity will generally work fine - + -- len(scales) should equal len(path_pts). No-op if not supplied - Each entry may be a single number for uniform scaling, or a pair of + Each entry may be a single number for uniform scaling, or a pair of numbers (or Point2) for differential X/Y scaling If not supplied, no scaling will occur. - + -- len(rotations) should equal 1 or len(path_pts). No-op if not supplied. Each point in shape_pts will be rotated by rotations[i] degrees at each point in path_pts. Or, if only one rotation is supplied, the shape will be rotated smoothly over rotations[0] degrees in the course of the extrusion - + -- len(transforms) should be 1 or be equal to len(path_pts). No-op if not supplied. - Each entry should be have the signature: + Each entry should be have the signature: def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 where path_norm is in [0,1] and expresses progress through the extrusion and loop_norm is in [0,1] and express progress through a single loop of the extrusion - + -- if connect_ends is True, the first and last loops of the extrusion will be joined, which is useful for toroidal geometries. Overrides cap_ends @@ -49,11 +54,10 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 will be connected to the centroid of that loop. For planar, convex shapes, this works nicely. If shape is less planar or convex, some self-intersection may happen. Not applied if connect_ends is True - ''' - + """ - polyhedron_pts:Points= [] - facet_indices:List[Tuple[int, int, int]] = [] + polyhedron_pts: Points = [] + facet_indices: List[Tuple[int, int, int]] = [] # Make sure we've got Euclid Point3's for all elements shape_pts = euclidify(shape_pts, Point3) @@ -66,7 +70,7 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 tangent_path_points: List[Point3] = [] # If first & last points are the same, let's close the shape - first_last_equal = ((path_pts[0] - path_pts[-1]).magnitude_squared() < EPSILON) + first_last_equal = (path_pts[0] - path_pts[-1]).magnitude_squared() < EPSILON if first_last_equal: connect_ends = True path_pts = path_pts[:][:-1] @@ -77,11 +81,14 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 first = Point3(*(path_pts[0] - (path_pts[1] - path_pts[0]))) last = Point3(*(path_pts[-1] - (path_pts[-2] - path_pts[-1]))) tangent_path_points = [first] + path_pts + [last] - tangents = [tangent_path_points[i+2] - tangent_path_points[i] for i in range(len(path_pts))] + tangents = [ + tangent_path_points[i + 2] - tangent_path_points[i] + for i in range(len(path_pts)) + ] for which_loop in range(len(path_pts)): # path_normal is 0 at the first path_pts and 1 at the last - path_normal = which_loop/ (len(path_pts) - 1) + path_normal = which_loop / (len(path_pts) - 1) path_pt = path_pts[which_loop] tangent = tangents[which_loop] @@ -89,21 +96,53 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 rotate_degrees = None if rotations: - rotate_degrees = rotations[which_loop] if len(rotations) > 1 else rotations[0] * path_normal + rotate_degrees = ( + rotations[which_loop] + if len(rotations) > 1 + else rotations[0] * path_normal + ) transform_func = None if transforms: - transform_func = transforms[which_loop] if len(transforms) > 1 else transforms[0] + transform_func = ( + transforms[which_loop] if len(transforms) > 1 else transforms[0] + ) + # Default to no *args or **kwargs + this_transform_func_args = None + if transform_args: + this_transform_func_args = ( + transform_args[which_loop] + if len(transform_args) > 1 + else transform_args[0] + ) + this_transform_func_kwargs = None + if transform_kwargs: + this_transform_func_kwargs = ( + transform_kwargs[which_loop] + if len(transform_kwargs) > 1 + else transform_kwargs[0] + ) this_loop = shape_pts[:] - this_loop = _scale_loop(this_loop, scale) + this_loop = _scale_loop( + this_loop, + scale, + ) this_loop = _rotate_loop(this_loop, rotate_degrees) - this_loop = _transform_loop(this_loop, transform_func, path_normal) - - this_loop = transform_to_point(this_loop, dest_point=path_pt, dest_normal=tangent, src_up=src_up) + this_loop = _transform_loop( + this_loop, + transform_func, + path_normal, + this_transform_func_args, + this_transform_func_kwargs, + ) + + this_loop = transform_to_point( + this_loop, dest_point=path_pt, dest_normal=tangent, src_up=src_up + ) loop_start_index = which_loop * shape_pt_count - if (which_loop < len(path_pts) - 1): + if which_loop < len(path_pts) - 1: loop_facets = _loop_facet_indices(loop_start_index, shape_pt_count) facet_indices += loop_facets @@ -117,46 +156,61 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 elif cap_ends: # OpenSCAD's polyhedron will automatically triangulate faces as needed. - # So just include all points at each end of the tube - last_loop_start_index = len(polyhedron_pts) - shape_pt_count + # So just include all points at each end of the tube + last_loop_start_index = len(polyhedron_pts) - shape_pt_count start_loop_indices = list(reversed(range(shape_pt_count))) - end_loop_indices = list(range(last_loop_start_index, last_loop_start_index + shape_pt_count)) + end_loop_indices = list( + range(last_loop_start_index, last_loop_start_index + shape_pt_count) + ) facet_indices.append(start_loop_indices) facet_indices.append(end_loop_indices) - return polyhedron(points=euc_to_arr(polyhedron_pts), faces=facet_indices) # type: ignore + return polyhedron(points=euc_to_arr(polyhedron_pts), faces=facet_indices) # type: ignore + -def _loop_facet_indices(loop_start_index:int, loop_pt_count:int, next_loop_start_index=None) -> List[FacetIndices]: +def _loop_facet_indices( + loop_start_index: int, loop_pt_count: int, next_loop_start_index=None +) -> List[FacetIndices]: facet_indices: List[FacetIndices] = [] # nlsi == next_loop_start_index if next_loop_start_index == None: next_loop_start_index = loop_start_index + loop_pt_count - loop_indices = list(range(loop_start_index, loop_pt_count + loop_start_index)) + [loop_start_index] - next_loop_indices = list(range(next_loop_start_index, loop_pt_count + next_loop_start_index )) + [next_loop_start_index] + loop_indices = list(range(loop_start_index, loop_pt_count + loop_start_index)) + [ + loop_start_index + ] + next_loop_indices = list( + range(next_loop_start_index, loop_pt_count + next_loop_start_index) + ) + [next_loop_start_index] for i, (a, b) in enumerate(zip(loop_indices[:-1], loop_indices[1:])): - c, d = next_loop_indices[i: i+2] + c, d = next_loop_indices[i : i + 2] # OpenSCAD's polyhedron will accept quads and do its own triangulation with them, - # so we could just append (a,b,d,c). + # so we could just append (a,b,d,c). # However, this lets OpenSCAD (Or CGAL?) do its own triangulation, leading # to some strange outcomes. Prefer to do our own triangulation. # c--d # |\ | # | \| - # a--b + # a--b # facet_indices.append((a,b,d,c)) - facet_indices.append((a,b,c)) - facet_indices.append((b,d,c)) + facet_indices.append((a, b, c)) + facet_indices.append((b, d, c)) return facet_indices -def _rotate_loop(points:Sequence[Point3], rotation_degrees:float=None) -> List[Point3]: + +def _rotate_loop( + points: Sequence[Point3], rotation_degrees: float = None +) -> List[Point3]: if rotation_degrees is None: return points - up = Vector3(0,0,1) + up = Vector3(0, 0, 1) rads = radians(rotation_degrees) return [p.rotate_around(up, rads) for p in points] -def _scale_loop(points:Sequence[Point3], scale:Union[float, Point2, Tuple2]=None) -> List[Point3]: + +def _scale_loop( + points: Sequence[Point3], scale: Union[float, Point2, Tuple2] = None +) -> List[Point3]: if scale is None: return points @@ -164,7 +218,14 @@ def _scale_loop(points:Sequence[Point3], scale:Union[float, Point2, Tuple2]=None scale = [scale] * 2 return [Point3(point.x * scale[0], point.y * scale[1], point.z) for point in points] -def _transform_loop(points:Sequence[Point3], transform_func:Point3Transform = None, path_normal:float = None) -> List[Point3]: + +def _transform_loop( + points: Sequence[Point3], + transform_func: Callable = None, + path_normal: float = None, + *args, + **kwargs, +) -> List[Point3]: # transform_func is a function that takes a point and optionally two floats, # a `path_normal`, in [0,1] that indicates where this loop is in a path extrusion, # and `loop_normal` in [0,1] that indicates where this point is in a list of points @@ -174,8 +235,9 @@ def _transform_loop(points:Sequence[Point3], transform_func:Point3Transform = No result = [] for i, p in enumerate(points): # i goes from 0 to 1 across points - loop_normal = i/(len(points) -1) - new_p = transform_func(p, path_normal, loop_normal) + loop_normal = i / (len(points) - 1) + # print(f"args:\t{args}") + # print(f"kwargs:\t{kwargs}") + new_p = transform_func(p, path_normal, loop_normal, *args, **kwargs) result.append(new_p) return result - diff --git a/solid/solidpython.py b/solid/solidpython.py index 4fc935f2..f7440352 100755 --- a/solid/solidpython.py +++ b/solid/solidpython.py @@ -27,10 +27,10 @@ import re PathStr = Union[Path, str] -AnimFunc = Callable[[Optional[float]], 'OpenSCADObject'] +AnimFunc = Callable[[Optional[float]], "OpenSCADObject"] # These are features added to SolidPython but NOT in OpenSCAD. # Mark them for special treatment -non_rendered_classes = ['hole', 'part'] +non_rendered_classes = ["hole", "part", "stem", "nose"] # Words reserved in Python but not OpenSCAD # Re: https://github.com/SolidCode/SolidPython/issues/99 @@ -54,10 +54,10 @@ def __init__(self, name: str, params: dict): self.is_part_root = False self.traits: Dict[str, Dict[str, float]] = {} - def add_trait(self, trait_name:str, trait_data:Dict[str, float]): + def add_trait(self, trait_name: str, trait_data: Dict[str, float]): self.traits[trait_name] = trait_data - def get_trait(self, trait_name:str) -> Optional[Dict[str, float]]: + def get_trait(self, trait_name: str) -> Optional[Dict[str, float]]: return self.traits.get(trait_name) def set_hole(self, is_hole: bool = True) -> "OpenSCADObject": @@ -68,7 +68,9 @@ def set_part_root(self, is_root: bool = True) -> "OpenSCADObject": self.is_part_root = is_root return self - def find_hole_children(self, path: List["OpenSCADObject"] = None) -> List["OpenSCADObject"]: + def find_hole_children( + self, path: List["OpenSCADObject"] = None + ) -> List["OpenSCADObject"]: """ Because we don't force a copy every time we re-use a node (e.g a = cylinder(2, 6); b = right(10) (a) @@ -101,22 +103,24 @@ def set_modifier(self, m: str) -> "OpenSCADObject": Used to add one of the 4 single-character modifiers: #(debug) !(root) %(background) or *(disable) """ - string_vals = {'disable': '*', - 'debug': '#', - 'background': '%', - 'root': '!', - '*': '*', - '#': '#', - '%': '%', - '!': '!'} - - self.modifier = string_vals.get(m.lower(), '') + string_vals = { + "disable": "*", + "debug": "#", + "background": "%", + "root": "!", + "*": "*", + "#": "#", + "%": "%", + "!": "!", + } + + self.modifier = string_vals.get(m.lower(), "") return self def _render(self, render_holes: bool = False) -> str: """ NOTE: In general, you won't want to call this method. For most purposes, - you really want scad_render(), + you really want scad_render(), Calling obj._render won't include necessary 'use' or 'include' statements """ # First, render all children @@ -166,8 +170,8 @@ def _render_str_no_children(self) -> str: # OpenSCAD doesn't have a 'segments' argument, but it does # have '$fn'. Swap one for the other - if 'segments' in self.params: - self.params['$fn'] = self.params.pop('segments') + if "segments" in self.params: + self.params["$fn"] = self.params.pop("segments") valid_keys = self.params.keys() @@ -230,16 +234,18 @@ def _render_hole_children(self) -> str: # with union in the hole segment of the compiled tree. # And if you figure out a better way to explain this, # please, please do... because I think this works, but I - # also think my rationale is shaky and imprecise. + # also think my rationale is shaky and imprecise. # -ETJ 19 Feb 2013 s = s.replace("intersection", "union") s = s.replace("difference", "union") return s - def add(self, child: Union["OpenSCADObject", Sequence["OpenSCADObject"]]) -> "OpenSCADObject": + def add( + self, child: Union["OpenSCADObject", Sequence["OpenSCADObject"]] + ) -> "OpenSCADObject": """ - if child is a single object, assume it's an OpenSCADObjects and + if child is a single object, assume it's an OpenSCADObjects and add it to self.children if child is a list, assume its members are all OpenSCADObjects and @@ -264,8 +270,8 @@ def set_parent(self, parent: "OpenSCADObject"): self.parent = parent def add_param(self, k: str, v: float) -> "OpenSCADObject": - if k == '$fn': - k = 'segments' + if k == "$fn": + k = "segments" self.params[k] = v return self @@ -281,8 +287,8 @@ def copy(self) -> "OpenSCADObject": # Python can't handle an '$fn' argument, while openSCAD only wants # '$fn'. Swap back and forth as needed; the final renderer will # sort this out. - if '$fn' in self.params: - self.params['segments'] = self.params.pop('$fn') + if "$fn" in self.params: + self.params["segments"] = self.params.pop("$fn") other = type(self)(**self.params) other.set_modifier(self.modifier) @@ -345,12 +351,9 @@ def _repr_png_(self) -> Optional[bytes]: tmp.write(scad_text) tmp.close() tmp_png.close() - subprocess.Popen([ - "openscad", - "--preview", - "-o", tmp_png.name, - tmp.name - ]).communicate() + subprocess.Popen( + ["openscad", "--preview", "-o", tmp_png.name, tmp.name] + ).communicate() with open(tmp_png.name, "rb") as png: png_data = png.read() @@ -368,11 +371,13 @@ class IncludedOpenSCADObject(OpenSCADObject): to the scad file it's included from. """ - def __init__(self, name, params, include_file_path, use_not_include=False, **kwargs): + def __init__( + self, name, params, include_file_path, use_not_include=False, **kwargs + ): self.include_file_path = self._get_include_path(include_file_path) - use_str = 'use' if use_not_include else 'include' - self.include_string = f'{use_str} <{self.include_file_path}>\n' + use_str = "use" if use_not_include else "include" + self.include_string = f"{use_str} <{self.include_file_path}>\n" # Just pass any extra arguments straight on to OpenSCAD; it'll accept # them @@ -393,26 +398,31 @@ def _get_include_path(self, include_file_path): return os.path.abspath(whole_path) # No loadable SCAD file was found in sys.path. Raise an error - raise ValueError(f"Unable to find included SCAD file: {include_file_path} in sys.path") + raise ValueError( + f"Unable to find included SCAD file: {include_file_path} in sys.path" + ) # ========================================= # = Rendering Python code to OpenSCAD code= # ========================================= -def _find_include_strings(obj: Union[IncludedOpenSCADObject, OpenSCADObject]) -> Set[str]: +def _find_include_strings( + obj: Union[IncludedOpenSCADObject, OpenSCADObject] +) -> Set[str]: include_strings = set() if isinstance(obj, IncludedOpenSCADObject): include_strings.add(obj.include_string) for child in obj.children: include_strings.update(_find_include_strings(child)) - # We also accept IncludedOpenSCADObject instances as parameters to functions, + # We also accept IncludedOpenSCADObject instances as parameters to functions, # so search in obj.params as well for param in obj.params.values(): if isinstance(param, OpenSCADObject): include_strings.update(_find_include_strings(param)) return include_strings -def scad_render(scad_object: OpenSCADObject, file_header: str = '') -> str: + +def scad_render(scad_object: OpenSCADObject, file_header: str = "") -> str: # Make this object the root of the tree root = scad_object @@ -421,18 +431,21 @@ def scad_render(scad_object: OpenSCADObject, file_header: str = '') -> str: include_strings = _find_include_strings(root) # and render the string - includes = ''.join(include_strings) + "\n" + includes = "".join(include_strings) + "\n" scad_body = root._render() - if file_header and not file_header.endswith('\n'): - file_header += '\n' + if file_header and not file_header.endswith("\n"): + file_header += "\n" return file_header + includes + scad_body -def scad_render_animated(func_to_animate: AnimFunc, - steps: int =20, - back_and_forth: bool=True, - file_header: str='') -> str: + +def scad_render_animated( + func_to_animate: AnimFunc, + steps: int = 20, + back_and_forth: bool = True, + file_header: str = "", +) -> str: # func_to_animate takes a single float argument, _time in [0, 1), and # returns an OpenSCADObject instance. # @@ -465,7 +478,7 @@ def scad_render_animated(func_to_animate: AnimFunc, scad_obj = func_to_animate(_time=0) # type: ignore include_strings = _find_include_strings(scad_obj) # and render the string - includes = ''.join(include_strings) + "\n" + includes = "".join(include_strings) + "\n" rendered_string = file_header + includes @@ -486,28 +499,36 @@ def scad_render_animated(func_to_animate: AnimFunc, scad_obj = func_to_animate(_time=eval_time) # type: ignore scad_str = indent(scad_obj._render()) - rendered_string += f"if ($t >= {time} && $t < {end_time}){{" \ - f" {scad_str}\n" \ - f"}}\n" + rendered_string += ( + f"if ($t >= {time} && $t < {end_time}){{" f" {scad_str}\n" f"}}\n" + ) return rendered_string -def scad_render_animated_file(func_to_animate:AnimFunc, - steps: int=20, - back_and_forth: bool=True, - filepath: Optional[str]=None, - out_dir: PathStr=None, - file_header: str='', - include_orig_code: bool=True) -> str: - rendered_string = scad_render_animated(func_to_animate, steps, - back_and_forth, file_header) - return _write_code_to_file(rendered_string, filepath, out_dir=out_dir, - include_orig_code=include_orig_code) - -def scad_render_to_file(scad_object: OpenSCADObject, - filepath: PathStr=None, - out_dir: PathStr=None, - file_header: str='', - include_orig_code: bool=True) -> str: + +def scad_render_animated_file( + func_to_animate: AnimFunc, + steps: int = 20, + back_and_forth: bool = True, + filepath: Optional[str] = None, + out_dir: PathStr = None, + file_header: str = "", + include_orig_code: bool = True, +) -> str: + rendered_string = scad_render_animated( + func_to_animate, steps, back_and_forth, file_header + ) + return _write_code_to_file( + rendered_string, filepath, out_dir=out_dir, include_orig_code=include_orig_code + ) + + +def scad_render_to_file( + scad_object: OpenSCADObject, + filepath: PathStr = None, + out_dir: PathStr = None, + file_header: str = "", + include_orig_code: bool = True, +) -> str: header = file_header if include_orig_code: version = _get_version() @@ -517,20 +538,23 @@ def scad_render_to_file(scad_object: OpenSCADObject, rendered_string = scad_render(scad_object, header) return _write_code_to_file(rendered_string, filepath, out_dir, include_orig_code) -def _write_code_to_file(rendered_string: str, - filepath: PathStr=None, - out_dir: PathStr=None, - include_orig_code: bool=True) -> str: + +def _write_code_to_file( + rendered_string: str, + filepath: PathStr = None, + out_dir: PathStr = None, + include_orig_code: bool = True, +) -> str: try: calling_file = Path(calling_module(stack_depth=3).__file__).absolute() # Output path is determined four ways: # -- If filepath is supplied, use filepath - # -- If no filepath is supplied but an out_dir is supplied, + # -- If no filepath is supplied but an out_dir is supplied, # give the calling file a .scad suffix and put it in out_dir # -- If neither filepath nor out_dir are supplied, give the new # file a .scad suffix and put it next to the calling file - # -- If no path info is supplied and we can't find a calling file - # (i.e, this is being called from an interactive terminal), + # -- If no path info is supplied and we can't find a calling file + # (i.e, this is being called from an interactive terminal), # write a file to Path.cwd() / 'solid.scad' out_path = Path() if filepath: @@ -539,10 +563,10 @@ def _write_code_to_file(rendered_string: str, odp = Path(out_dir) if not odp.exists(): odp.mkdir() - out_path = odp / calling_file.with_suffix('.scad').name + out_path = odp / calling_file.with_suffix(".scad").name else: - out_path = calling_file.with_suffix('.scad') - + out_path = calling_file.with_suffix(".scad") + if include_orig_code: rendered_string += sp_code_in_scad_comment(calling_file) except AttributeError as e: @@ -557,32 +581,34 @@ def _write_code_to_file(rendered_string: str, odp = Path(out_dir) if out_dir else Path.cwd() if not odp.exists(): odp.mkdir() - out_path = odp / 'solid.scad' + out_path = odp / "solid.scad" out_path.write_text(rendered_string) return out_path.absolute().as_posix() + def _get_version() -> str: """ Returns SolidPython version Returns '' if no version can be found """ - version = '' + version = "" try: # if SolidPython is installed use `pkg_resources` - version = pkg_resources.get_distribution('solidpython').version + version = pkg_resources.get_distribution("solidpython").version except pkg_resources.DistributionNotFound: # if the running SolidPython is not the one installed via pip, # try to read it from the project setup file version_pattern = re.compile(r"version = ['\"]([^'\"]*)['\"]") - version_file_path = Path(__file__).parent.parent / 'pyproject.toml' + version_file_path = Path(__file__).parent.parent / "pyproject.toml" if version_file_path.exists(): version_match = version_pattern.search(version_file_path.read_text()) if version_match: version = version_match.group(1) return version + def sp_code_in_scad_comment(calling_file: PathStr) -> str: """ Once a SCAD file has been created, it's difficult to reconstruct @@ -597,14 +623,16 @@ def sp_code_in_scad_comment(calling_file: PathStr) -> str: # to create a given file; That would future-proof any given SP-created # code because it would point to the relevant dependencies as well as # the actual code - pyopenscad_str = (f"\n" - f"/***********************************************\n" - f"********* SolidPython code: **********\n" - f"************************************************\n" - f" \n" - f"{pyopenscad_str} \n" - f" \n" - f"************************************************/\n") + pyopenscad_str = ( + f"\n" + f"/***********************************************\n" + f"********* SolidPython code: **********\n" + f"************************************************\n" + f" \n" + f"{pyopenscad_str} \n" + f" \n" + f"************************************************/\n" + ) return pyopenscad_str @@ -621,20 +649,21 @@ def parse_scad_callables(filename: str) -> List[dict]: args = [] kwargs = [] - #for some reason solidpython needs to treat all openscad arguments as if - #they where optional. I don't know why, but at least to pass the tests - #it's neccessary to handle it like this !?!?! + # for some reason solidpython needs to treat all openscad arguments as if + # they where optional. I don't know why, but at least to pass the tests + # it's neccessary to handle it like this !?!?! for p in c.parameters: kwargs.append(p.name) - #if p.optional: + # if p.optional: # kwargs.append(p.name) - #else: + # else: # args.append(p.name) - callables.append({'name': c.name, 'args': args, 'kwargs': kwargs}) + callables.append({"name": c.name, "args": args, "kwargs": kwargs}) return callables + def calling_module(stack_depth: int = 2) -> ModuleType: """ Returns the module *2* back in the frame stack. That means: @@ -642,7 +671,7 @@ def calling_module(stack_depth: int = 2) -> ModuleType: for module A. This means that we have to know exactly how far back in the stack - our desired module is; if code in module B calls another function in + our desired module is; if code in module B calls another function in module B, we have to increase the stack_depth argument to account for this. @@ -657,13 +686,16 @@ def calling_module(stack_depth: int = 2) -> ModuleType: import __main__ as calling_mod # type: ignore return calling_mod -def new_openscad_class_str(class_name: str, - args: Sequence[str] = None, - kwargs: Sequence[str] = None, - include_file_path: Optional[str] = None, - use_not_include: bool = True) -> str: - args_str = '' - args_pairs = '' + +def new_openscad_class_str( + class_name: str, + args: Sequence[str] = None, + kwargs: Sequence[str] = None, + include_file_path: Optional[str] = None, + use_not_include: bool = True, +) -> str: + args_str = "" + args_pairs = "" args = args or [] kwargs = kwargs or [] @@ -675,7 +707,7 @@ def new_openscad_class_str(class_name: str, args = map(_subbed_keyword, args) # type: ignore for arg in args: - args_str += ', ' + arg + args_str += ", " + arg args_pairs += f"'{arg}':{arg}, " # kwargs have a default value defined in their SCAD versions. We don't @@ -683,7 +715,7 @@ def new_openscad_class_str(class_name: str, # that one is defined. kwargs = map(_subbed_keyword, kwargs) # type: ignore for kwarg in kwargs: - args_str += f', {kwarg}=None' + args_str += f", {kwarg}=None" args_pairs += f"'{kwarg}':{kwarg}, " if include_file_path: @@ -694,21 +726,26 @@ def new_openscad_class_str(class_name: str, # NOTE the explicit import of 'solid' below. This is a fix for: # https://github.com/SolidCode/SolidPython/issues/20 -ETJ 16 Jan 2014 - result = (f"import solid\n" - f"class {class_name}(solid.IncludedOpenSCADObject):\n" - f" def __init__(self{args_str}, **kwargs):\n" - f" solid.IncludedOpenSCADObject.__init__(self, '{class_name}', {{{args_pairs} }}, include_file_path='{include_file_str}', use_not_include={use_not_include}, **kwargs )\n" - f" \n" - f"\n") + result = ( + f"import solid\n" + f"class {class_name}(solid.IncludedOpenSCADObject):\n" + f" def __init__(self{args_str}, **kwargs):\n" + f" solid.IncludedOpenSCADObject.__init__(self, '{class_name}', {{{args_pairs} }}, include_file_path='{include_file_str}', use_not_include={use_not_include}, **kwargs )\n" + f" \n" + f"\n" + ) else: - result = (f"class {class_name}(OpenSCADObject):\n" - f" def __init__(self{args_str}):\n" - f" OpenSCADObject.__init__(self, '{class_name}', {{{args_pairs }}})\n" - f" \n" - f"\n") + result = ( + f"class {class_name}(OpenSCADObject):\n" + f" def __init__(self{args_str}):\n" + f" OpenSCADObject.__init__(self, '{class_name}', {{{args_pairs }}})\n" + f" \n" + f"\n" + ) return result + def _subbed_keyword(keyword: str) -> str: """ Append an underscore to any python reserved word. @@ -730,18 +767,24 @@ def _subbed_keyword(keyword: str) -> str: new_key = "__" + keyword[1:] if new_key != keyword: - print(f"\nFound OpenSCAD code that's not compatible with Python. \n" - f"Imported OpenSCAD code using `{keyword}` \n" - f"can be accessed with `{new_key}` in SolidPython\n") + print( + f"\nFound OpenSCAD code that's not compatible with Python. \n" + f"Imported OpenSCAD code using `{keyword}` \n" + f"can be accessed with `{new_key}` in SolidPython\n" + ) return new_key + def _unsubbed_keyword(subbed_keyword: str) -> str: """ Remove trailing underscore for already-subbed python reserved words. Remove prepending underscore if remaining identifier starts with a digit. No-op for all other strings: e.g. 'or_' => 'or', 'other_' => 'other_' """ - if subbed_keyword.endswith("_") and subbed_keyword[:-1] in PYTHON_ONLY_RESERVED_WORDS: + if ( + subbed_keyword.endswith("_") + and subbed_keyword[:-1] in PYTHON_ONLY_RESERVED_WORDS + ): return subbed_keyword[:-1] elif subbed_keyword.startswith("__"): @@ -755,18 +798,21 @@ def _unsubbed_keyword(subbed_keyword: str) -> str: return subbed_keyword + # now that we have the base class defined, we can do a circular import from . import objects + def py2openscad(o: Union[bool, float, str, Iterable]) -> str: if type(o) == bool: return str(o).lower() if type(o) == float: return f"{o:.10f}" # type: ignore if type(o) == str: - return f'\"{o}\"' # type: ignore + return f'"{o}"' # type: ignore if type(o).__name__ == "ndarray": import numpy # type: ignore + return numpy.array2string(o, separator=",", threshold=1000000000) if isinstance(o, IncludedOpenSCADObject): return o._render()[1:-1] @@ -782,5 +828,6 @@ def py2openscad(o: Union[bool, float, str, Iterable]) -> str: return s return str(o) + def indent(s: str) -> str: return s.replace("\n", "\n\t") diff --git a/solid/splines.py b/solid/splines.py index ef4310f4..60ad98d0 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -1,9 +1,21 @@ #! /usr/bin/env python -from math import pow - -from solid import union, circle, cylinder, polygon, color, OpenSCADObject, translate, linear_extrude, polyhedron +from math import copysign, pow, cos, sin, pi, atan2, acos +from itertools import pairwise +from pudb import set_trace + +from solid import ( + union, + circle, + cylinder, + polygon, + color, + OpenSCADObject, + translate, + linear_extrude, + polyhedron, +) from solid.utils import bounding_box, right, Red, Tuple3, euclidify -from euclid3 import Vector2, Vector3, Point2, Point3 +from euclid3 import Vector2, Vector3, Point2, Point3, Matrix4 from typing import Sequence, Tuple, Union, List, cast @@ -20,27 +32,31 @@ Vec23 = Union[Vector2, Vector3] FourPoints = Tuple[Point23Input, Point23Input, Point23Input, Point23Input] -SEGMENTS = 48 +SEGMENTS = 359 +DEBUG = True DEFAULT_SUBDIVISIONS = 10 DEFAULT_EXTRUDE_HEIGHT = 1 + # ======================= # = CATMULL-ROM SPLINES = # ======================= -def catmull_rom_polygon(points: Sequence[Point23Input], - subdivisions: int = DEFAULT_SUBDIVISIONS, - extrude_height: float = DEFAULT_EXTRUDE_HEIGHT, - show_controls: bool =False, - center: bool=True) -> OpenSCADObject: +def catmull_rom_polygon( + points: Sequence[Point23Input], + subdivisions: int = DEFAULT_SUBDIVISIONS, + extrude_height: float = DEFAULT_EXTRUDE_HEIGHT, + show_controls: bool = False, + center: bool = True, +) -> OpenSCADObject: """ - Return a closed OpenSCAD polygon object through all of `points`, - extruded to `extrude_height`. If `show_controls` is True, return red + Return a closed OpenSCAD polygon object through all of `points`, + extruded to `extrude_height`. If `show_controls` is True, return red cylinders at each of the specified control points; this makes it easier to move determine which points should move to get a desired shape. NOTE: if `extrude_height` is 0, this function returns a 2D `polygon()` - object, which OpenSCAD can only combine with other 2D objects + object, which OpenSCAD can only combine with other 2D objects (e.g. `square`, `circle`, but not `cube` or `cylinder`). If `extrude_height` is nonzero, the object returned will be 3D and only combine with 3D objects. """ @@ -53,20 +69,23 @@ def catmull_rom_polygon(points: Sequence[Point23Input], shape += control_points(points, extrude_height, center) return shape -def catmull_rom_points( points: Sequence[Point23Input], - subdivisions:int = DEFAULT_SUBDIVISIONS, - close_loop: bool=False, - start_tangent: Vec23 = None, - end_tangent: Vec23 = None) -> List[Point3]: + +def catmull_rom_points( + points: Sequence[Point23Input], + subdivisions: int = DEFAULT_SUBDIVISIONS, + close_loop: bool = False, + start_tangent: Vec23 = None, + end_tangent: Vec23 = None, +) -> List[Point3]: """ - Return a smooth set of points through `points`, with `subdivisions` points - between each pair of control points. - - If `close_loop` is False, `start_tangent` and `end_tangent` can specify - tangents at the open ends of the returned curve. If not supplied, tangents + Return a smooth set of points through `points`, with `subdivisions` points + between each pair of control points. + + If `close_loop` is False, `start_tangent` and `end_tangent` can specify + tangents at the open ends of the returned curve. If not supplied, tangents will be colinear with first and last supplied segments - Credit due: Largely taken from C# code at: + Credit due: Largely taken from C# code at: https://www.habrador.com/tutorials/interpolation/1-catmull-rom-splines/ retrieved 20190712 """ @@ -77,37 +96,46 @@ def catmull_rom_points( points: Sequence[Point23Input], points_list = list([euclidify(p, Point3) for p in points]) if close_loop: - cat_points = euclidify([points_list[-1]] + points_list + points_list[0:2], Point3) + cat_points = euclidify( + [points_list[-1]] + points_list + points_list[0:2], Point3 + ) else: # Use supplied tangents or just continue the ends of the supplied points start_tangent = start_tangent or (points_list[1] - points_list[0]) start_tangent = euclidify(start_tangent, Vector3) end_tangent = end_tangent or (points_list[-2] - points_list[-1]) end_tangent = euclidify(end_tangent, Vector3) - cat_points = [points_list[0]+ start_tangent] + points_list + [points_list[-1] + end_tangent] + cat_points = ( + [points_list[0] + start_tangent] + + points_list + + [points_list[-1] + end_tangent] + ) last_point_range = len(cat_points) - 3 if close_loop else len(cat_points) - 3 for i in range(0, last_point_range): include_last = True if i == last_point_range - 1 else False - controls = cat_points[i:i+4] + controls = cat_points[i : i + 4] # If we're closing a loop, controls needs to wrap around the end of the array points_needed = 4 - len(controls) if points_needed > 0: controls += cat_points[0:points_needed] controls_tuple = cast(FourPoints, controls) - catmull_points += _catmull_rom_segment(controls_tuple, subdivisions, include_last) + catmull_points += _catmull_rom_segment( + controls_tuple, subdivisions, include_last + ) return catmull_points -def _catmull_rom_segment(controls: FourPoints, - subdivisions: int, - include_last=False) -> List[Point3]: + +def _catmull_rom_segment( + controls: FourPoints, subdivisions: int, include_last=False +) -> List[Point3]: """ Returns `subdivisions` Points between the 2nd & 3rd elements of `controls`, on a quadratic curve that passes through all 4 control points. If `include_last` is True, return `subdivisions` + 1 points, the last being - controls[2]. + controls[2]. No reason to call this unless you're trying to do something very specific """ @@ -121,18 +149,21 @@ def _catmull_rom_segment(controls: FourPoints, p0, p1, p2, p3 = [euclidify(p, Point3) for p in controls] a = 2 * p1 b = p2 - p0 - c = 2* p0 - 5*p1 + 4*p2 - p3 - d = -p0 + 3*p1 - 3*p2 + p3 + c = 2 * p0 - 5 * p1 + 4 * p2 - p3 + d = -p0 + 3 * p1 - 3 * p2 + p3 for i in range(num_points): - t = i/subdivisions + t = i / subdivisions pos = 0.5 * (a + (b * t) + (c * t * t) + (d * t * t * t)) positions.append(Point3(*pos)) return positions -def catmull_rom_patch_points(patch:Tuple[PointInputs, PointInputs], - subdivisions:int = DEFAULT_SUBDIVISIONS, - index_start:int = 0) -> Tuple[List[Point3], List[FaceTrio]]: + +def catmull_rom_patch_points( + patch: Tuple[PointInputs, PointInputs], + subdivisions: int = DEFAULT_SUBDIVISIONS, + index_start: int = 0, +) -> Tuple[List[Point3], List[FaceTrio]]: verts: List[Point3] = [] faces: List[FaceTrio] = [] @@ -142,9 +173,11 @@ def catmull_rom_patch_points(patch:Tuple[PointInputs, PointInputs], strip_length = len(cm_points_a) for i in range(subdivisions + 1): - frac = i/subdivisions - verts += list([affine_combination(a,b, frac) for a,b in zip(cm_points_a, cm_points_b)]) - a_start = i*strip_length + index_start + frac = i / subdivisions + verts += list( + [affine_combination(a, b, frac) for a, b in zip(cm_points_a, cm_points_b)] + ) + a_start = i * strip_length + index_start b_start = a_start + strip_length # This connects the verts we just created to the verts we'll make on the # next loop. So don't calculate for the last loop @@ -153,18 +186,26 @@ def catmull_rom_patch_points(patch:Tuple[PointInputs, PointInputs], return verts, faces -def catmull_rom_patch(patch:Tuple[PointInputs, PointInputs], subdivisions:int = DEFAULT_SUBDIVISIONS) -> OpenSCADObject: + +def catmull_rom_patch( + patch: Tuple[PointInputs, PointInputs], subdivisions: int = DEFAULT_SUBDIVISIONS +) -> OpenSCADObject: faces, vertices = catmull_rom_patch_points(patch, subdivisions) return polyhedron(faces, vertices) -def catmull_rom_prism( control_curves:Sequence[PointInputs], - subdivisions:int = DEFAULT_SUBDIVISIONS, - closed_ring:bool = True, - add_caps:bool = True, - smooth_edges: bool = False ) -> polyhedron: + +def catmull_rom_prism( + control_curves: Sequence[PointInputs], + subdivisions: int = DEFAULT_SUBDIVISIONS, + closed_ring: bool = True, + add_caps: bool = True, + smooth_edges: bool = False, +) -> polyhedron: if smooth_edges: - return catmull_rom_prism_smooth_edges(control_curves, subdivisions, closed_ring, add_caps) + return catmull_rom_prism_smooth_edges( + control_curves, subdivisions, closed_ring, add_caps + ) verts: List[Point3] = [] faces: List[FaceTrio] = [] @@ -172,8 +213,8 @@ def catmull_rom_prism( control_curves:Sequence[PointInputs], curves = list([euclidify(c) for c in control_curves]) if closed_ring: curves.append(curves[0]) - - curve_length = (len(curves[0]) -1) * subdivisions + 1 + + curve_length = (len(curves[0]) - 1) * subdivisions + 1 for i, (a, b) in enumerate(zip(curves[:-1], curves[1:])): index_start = len(verts) - curve_length first_new_vert = curve_length @@ -181,7 +222,9 @@ def catmull_rom_prism( control_curves:Sequence[PointInputs], index_start = 0 first_new_vert = 0 - new_verts, new_faces = catmull_rom_patch_points((a,b), subdivisions=subdivisions, index_start=index_start) + new_verts, new_faces = catmull_rom_patch_points( + (a, b), subdivisions=subdivisions, index_start=index_start + ) # new_faces describes all the triangles in the patch we just computed, # but new_verts shares its first curve_length vertices with the last @@ -191,7 +234,7 @@ def catmull_rom_prism( control_curves:Sequence[PointInputs], if closed_ring and add_caps: bot_indices = range(0, len(verts), curve_length) - top_indices = range(curve_length-1, len(verts), curve_length) + top_indices = range(curve_length - 1, len(verts), curve_length) bot_centroid, bot_faces = centroid_endcap(verts, bot_indices) verts.append(bot_centroid) @@ -201,14 +244,17 @@ def catmull_rom_prism( control_curves:Sequence[PointInputs], top_centroid, top_faces = centroid_endcap(verts, top_indices, invert=True) verts.append(top_centroid) faces += top_faces - + p = polyhedron(faces=faces, points=verts, convexity=3) return p -def catmull_rom_prism_smooth_edges( control_curves:Sequence[PointInputs], - subdivisions:int = DEFAULT_SUBDIVISIONS, - closed_ring:bool = True, - add_caps:bool = True ) -> polyhedron: + +def catmull_rom_prism_smooth_edges( + control_curves: Sequence[PointInputs], + subdivisions: int = DEFAULT_SUBDIVISIONS, + closed_ring: bool = True, + add_caps: bool = True, +) -> polyhedron: verts: List[Point3] = [] faces: List[FaceTrio] = [] @@ -217,11 +263,15 @@ def catmull_rom_prism_smooth_edges( control_curves:Sequence[PointInputs], curves = list([euclidify(c) for c in control_curves]) - expanded_curves = [catmull_rom_points(c, subdivisions, close_loop=False) for c in curves] + expanded_curves = [ + catmull_rom_points(c, subdivisions, close_loop=False) for c in curves + ] expanded_length = len(expanded_curves[0]) for i in range(expanded_length): contour_controls = [c[i] for c in expanded_curves] - contour = catmull_rom_points(contour_controls, subdivisions, close_loop=closed_ring) + contour = catmull_rom_points( + contour_controls, subdivisions, close_loop=closed_ring + ) verts += contour contour_length = len(contour) @@ -233,9 +283,11 @@ def catmull_rom_prism_smooth_edges( control_curves:Sequence[PointInputs], # are pointed outwards for the test cases I ran. I think if control # curves were specified clockwise rather than counter-clockwise, all # of the faces would be pointed inwards - new_faces = face_strip_list(b_start, a_start, length=contour_length, close_loop=closed_ring) + new_faces = face_strip_list( + b_start, a_start, length=contour_length, close_loop=closed_ring + ) faces += new_faces - + if closed_ring and add_caps: bot_indices = range(0, contour_length) top_indices = range(len(verts) - contour_length, len(verts)) @@ -247,32 +299,36 @@ def catmull_rom_prism_smooth_edges( control_curves:Sequence[PointInputs], # top endcap; otherwise both endcaps would point to the same centroid point top_centroid, top_faces = centroid_endcap(verts, top_indices, invert=True) verts.append(top_centroid) - faces += top_faces + faces += top_faces p = polyhedron(faces=faces, points=verts, convexity=3) return p + # ================== # = BEZIER SPLINES = # ================== -# Ported from William A. Adams' Bezier OpenSCAD code at: +# Ported from William A. Adams' Bezier OpenSCAD code at: # https://www.thingiverse.com/thing:8443 -def bezier_polygon( controls: FourPoints, - subdivisions:int = DEFAULT_SUBDIVISIONS, - extrude_height:float = DEFAULT_EXTRUDE_HEIGHT, - show_controls: bool = False, - center: bool = True) -> OpenSCADObject: - ''' + +def bezier_polygon( + controls: FourPoints, + subdivisions: int = DEFAULT_SUBDIVISIONS, + extrude_height: float = DEFAULT_EXTRUDE_HEIGHT, + show_controls: bool = False, + center: bool = True, +) -> OpenSCADObject: + """ Return an OpenSCAD object representing a closed quadratic Bezier curve. - If extrude_height == 0, return a 2D `polygon()` object. - If extrude_height > 0, return a 3D extrusion of specified height. + If extrude_height == 0, return a 2D `polygon()` object. + If extrude_height > 0, return a 3D extrusion of specified height. Note that OpenSCAD won't render 2D & 3D objects together correctly, so pick one and use that. - ''' + """ points = bezier_points(controls, subdivisions) # OpenSCAD can'ts handle Point3s in creating a polygon. Convert them to Point2s - # Note that this prevents us from making polygons outside of the XY plane, + # Note that this prevents us from making polygons outside of the XY plane, # even though a polygon could reasonably be in some other plane while remaining 2D points = list((Point2(p.x, p.y) for p in points)) shape: OpenSCADObject = polygon(points) @@ -280,23 +336,28 @@ def bezier_polygon( controls: FourPoints, shape = linear_extrude(extrude_height, center=center)(shape) if show_controls: - control_objs = control_points(controls, extrude_height=extrude_height, center=center) + control_objs = control_points( + controls, extrude_height=extrude_height, center=center + ) shape += control_objs - + return shape -def bezier_points(controls: FourPoints, - subdivisions: int = DEFAULT_SUBDIVISIONS, - include_last: bool = True) -> List[Point3]: + +def bezier_points( + controls: FourPoints, + subdivisions: int = DEFAULT_SUBDIVISIONS, + include_last: bool = True, +) -> List[Point3]: """ Returns a list of `subdivisions` (+ 1, if `include_last` is True) points - on the cubic bezier curve defined by `controls`. The curve passes through + on the cubic bezier curve defined by `controls`. The curve passes through controls[0] and controls[3] If `include_last` is True, the last point returned will be controls[3]; if False, (useful for linking several curves together), controls[3] won't be included - Ported from William A. Adams' Bezier OpenSCAD code at: + Ported from William A. Adams' Bezier OpenSCAD code at: https://www.thingiverse.com/thing:8443 """ # TODO: enable a smooth curve through arbitrarily many points, as described at: @@ -305,51 +366,353 @@ def bezier_points(controls: FourPoints, points: List[Point3] = [] last_elt = 1 if include_last else 0 for i in range(subdivisions + last_elt): - u = i/subdivisions + u = i / subdivisions points.append(_point_along_bez4(*controls, u)) return points -def _point_along_bez4(p0: Point23Input, p1: Point23Input, p2: Point23Input, p3: Point23Input, u:float) -> Point3: + +def _point_along_bez4( + p0: Point23Input, p1: Point23Input, p2: Point23Input, p3: Point23Input, u: float +) -> Point3: p0 = euclidify(p0) p1 = euclidify(p1) p2 = euclidify(p2) p3 = euclidify(p3) - x = _bez03(u)*p0.x + _bez13(u)*p1.x + _bez23(u)*p2.x + _bez33(u)*p3.x - y = _bez03(u)*p0.y + _bez13(u)*p1.y + _bez23(u)*p2.y + _bez33(u)*p3.y - z = _bez03(u)*p0.z + _bez13(u)*p1.z + _bez23(u)*p2.z + _bez33(u)*p3.z + x = _bez03(u) * p0.x + _bez13(u) * p1.x + _bez23(u) * p2.x + _bez33(u) * p3.x + y = _bez03(u) * p0.y + _bez13(u) * p1.y + _bez23(u) * p2.y + _bez33(u) * p3.y + z = _bez03(u) * p0.z + _bez13(u) * p1.z + _bez23(u) * p2.z + _bez33(u) * p3.z return Point3(x, y, z) -def _bez03(u:float) -> float: - return pow((1-u), 3) -def _bez13(u:float) -> float: - return 3*u*(pow((1-u),2)) +def _bez03(u: float) -> float: + return pow((1 - u), 3) + + +def _bez13(u: float) -> float: + return 3 * u * (pow((1 - u), 2)) + + +def _bez23(u: float) -> float: + return 3 * (pow(u, 2)) * (1 - u) + -def _bez23(u:float) -> float: - return 3*(pow(u,2))*(1-u) +def _bez33(u: float) -> float: + return pow(u, 3) -def _bez33(u:float) -> float: - return pow(u,3) # ================ # = HOBBY CURVES = # ================ +# Based on this post +# https://www.jakelow.com/blog/hobby-curves +# Comments mostly transferred to respective positions +# This code implements Hobby's algorithm for fitting a cubic Bézier curve onto +# a sequence of points in two dimensions. + +# Hobby's algorithm was devised by John D. Hobby and Donald Knuth for use in +# the METAFONT program. Their original paper on the algorithm, titled “Smooth, +# Easy to Compute Interpolating Splines”, was published in 'Discrete and +# Computational Geometry' vol. 1 in 1986. A copy can be found at this URL: +# https://link.springer.com/content/pdf/10.1007/BF02187690.pdf + + +# Based on the paper "Typographers, programmers and mathematicians, or the case +# of an æsthetically pleasing interpolation" by Bogusław Jackowski, published in +# TUGboat, Volume 34 (2013), No. 2, and available at the following URL: +# https://tug.org/TUGboat/tb34-2/tb107jackowski.pdf +def hobby_points( + points: List[Point3], omega: float = 0, close_loop: bool = False +) -> List[Point2]: + """Hobby's algorithm: given a set of points, fit a Bézier spline to them. + + The chosen splines tend to have pleasing, rounded shapes. + + Parameters: + - points: an array of points as [x, y] pairs + - omega: a number between 0 and 1 (inclusive); controls how much curl + there will be at the endpoints of the curve + + Returns: an array of points as [x, y] pairs that define a Bézier spline. + + The output array will have 3n - 2 points where n is the number of input points. + The output will contain every point in the input (these become knots in the + Bézier spline), interspersed with pairs of new points which define positions + of handle points on the spline. All points are in the same coordinate space. + """ + # n is defined such that the points can be numbered P[0] ... P[n]. + # such that there are a total of n-1 points + if close_loop: + points.append(points[0]) + points.insert(0, points[-2]) + n = len(points) - 1 + if n < 2: + raise ValueError( + "Solving for Hobby points only possible for two or more points." + ) + # chords[i] is the vector from P[i] to P[i+1]. + # d[i] is the length of the ith chord. + # Chords from first to last point + chords = [ + Vector3(*(next_point - point).xyz) for point, next_point in pairwise(points) + ] + # if close_loop, prepend chord from last point to first + # and append chord from first to last + # or maybe the other way around? + # chords.insert(0, Vector3(*(points[-1] - points[0]).xyz)) + # chords.append(Vector3(*(points[0] - points[-1]).xyz)) + # d = [point.distance(next_point) for point, next_point in pairwise(points)] + + d = [chord.magnitude() for chord in chords] + if min(d) <= 0: + # cannot support successive points being the same + raise ValueError("No chord can be zero-length.") + + # gamma[i] is the signed turning angle at P[i], i.e. the angle between + # the chords from P[i-1] to P[i] and from P[i] to P[i+1]. + # gamma[0] is undefined; gamma[n] is artificially defined to be zero + gamma = [ + # vec3_angle_between(chord, next_chord) + # for chord, next_chord in pairwise(chords) + chord.angle(next_chord) + for chord, next_chord in pairwise(chords) + ] + # gamma.insert(0, copysign(1, gamma[0])) + gamma.insert(0, gamma[-1]) + gamma.append(gamma[1]) + + normals = [ + vec3_perp_norm(chord, next_chord) for chord, next_chord in pairwise(chords) + ] + normals.insert(0, normals[0]) + normals.append(normals[-1]) + + # Set up the system of linear equations (Jackowski, formula 38). + # We're representing this system as a tridiagonal matrix, because + # we can solve such a system in O(n) time using the Thomas algorithm. + # + # Here, A, B, and C are the matrix diagonals and D is the right-hand side. + # See Wikipedia for a more detailed explanation: + # https://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm + A = [0] + B = [2.0 + omega] + C = [2.0 * omega + 1] + D = [-1.0 * C[0] * gamma[1]] + + for i in range(1, n): + A.append(1.0 / d[i - 1]) + B.append((2.0 * d[i - 1] + 2 * d[i]) / (d[i - 1] * d[i])) + C.append(1.0 / d[i]) + D.append( + (-1.0 * (2 * gamma[i] * d[i] + gamma[i + 1] * d[i - 1])) / (d[i - 1] * d[i]) + ) + A.append(2.0 * omega + 1) + B.append(2.0 + omega) + C.append(0) + D.append(0) + + if DEBUG: + alpha, c_prime, d_prime = thomas(A, B, C, D) + else: + alpha = thomas(A, B, C, D) + + # Use alpha (the chord angle) and gamma (the turning angle of the chord + # polyline) to solve for beta at each point (beta is like alpha, but for + # the chord and handle vector arriving at P[i] rather than leaving from it). + beta = [] + for i in range(n - 1): + beta.append(-1 * gamma[i + 1] - alpha[i + 1]) + beta.append(-1 * alpha[n]) + + # Now we have the angles between the handle vector and the chord + # both arrriving at and leaving from each point, we can solve for the + # positions of the handle (control) points themselves. + c0 = [] + c1 = [] + for i in range(n): + # Compute the magnitudes of the handle vectors at this point + a = (rho(alpha[i], beta[i]) * d[i]) / 3.0 + b = (rho(beta[i], alpha[i]) * d[i]) / 3.0 + + # Use the magnitutes, chords, and turning angles to find + # the positions of the control points in the global coord space. + c0_vec = chords[i].rotate_around(normals[i], alpha[i]) + c1_vec = chords[i].rotate_around(normals[i + 1], -1 * beta[i]) + + c0_new_point = points[i] + c0_vec.normalize() * a + c1_new_point = points[i + 1] - c1_vec.normalize() * b + c0.append(c0_new_point) + c1.append(c1_new_point) + + # c0_first_point = (c0[0] / a) * b + # c0_first_point = points[0] + ( + # chords[0].rotate_around(normals[0], alpha[n - 1]).normalize() + # * ((rho(alpha[0], beta[0]) * d[0]) / 3.0) + # ) + # # CORRECT + # c0_last_point = points[-2] + ( + # chords[-1].rotate_around(normals[n], alpha[n - 1]).normalize() + # # * a + # * ((rho(alpha[-1], beta[-1]) * d[-1]) / 3.0) + # ) + # c1_last_point = points[0] - ( + # chords[-1].rotate_around(normals[n - 1], pi / 3 - alpha[0]).normalize() + # * ((rho(beta[0], alpha[0]) * d[0]) / 3.0) + # # * ((rho(beta[-1], beta[-1]) * d[0]) / 3.0) + # ) + # Finally gather up and return the spline points (knots & control points) + # as a List of Point2s + # res_controls = [(points[i], c0[i], c1[i]) for i in range(0, n - 1)] + res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 1)] + # c0 at points[0] is None as there is no leading control + # c1 at points[-1] is None as there is no leading control + # res_controls[0] = (points[1], c0[-1], c1[0]) + # res_controls[-1] = (c0[0], c1[-1]) + res_controls.insert( + 0, + (points[0], c0[-1], c1[1]), + ) + # res_controls.append( + # (points[n - 2], c0[0], c1[-1]), + # ) + # res_controls.append((points[n - 1],)) + # og_angle = + res_controls.append((points[n - 1],)) + if DEBUG: + return [ + [Point3(*p.xyz) for p in flatten(res_controls)], + c0, + c1, + chords, + d, + gamma, + alpha, + c_prime, + d_prime, + beta, + A, + B, + C, + D, + normals, + ] + return [Point3(*p.xyz) for p in flatten(res_controls)] + + +# Once the angles alpha and beta have been computed for each knot (which +# determine the direction from the knot to each of its neighboring handles), +# this function is used to compute the lengths of the vectors from the knot to +# those handles. Combining the length and angle together lets us solve for the +# handle positions. +# +# The exact choice of function is somewhat arbitrary. The aim is to return +# handle lengths that produce a Bézier curve which is a good approximation of a +# circular arc for points near the knot. +# +# Hobby and Knuth both proposed multiple candidate functions. This code uses +# the function from Jackowski formula 28, due to its simplicity. For other +# choices see Jackowski, section 5. +def rho(alpha: float, beta: float) -> float: + """ "Velocity function" to compute the length of handles for a bezier spline.""" + c = 2.0 / 3 + return 2.0 / (1 + c * cos(beta) + (1 - c) * cos(alpha)) + + +def thomas(A: Sequence, B: Sequence, C: Sequence, D: Sequence) -> Sequence: + """Thomas algorithm: Solve a system of equations encoded in a tridiagonal matrix. + + A, B, and C are diagonals of the matrix. B is the main diagonal. + D is the vector on the right-hand-side of the equation. + + Both B and D will have n elements. The arrays A and C will have + length n as well, but each has one fewer element than B (the values + A[0] and C[n-1] are undefined).""" + + # Note: n is defined here so that B[n] is valid, i.e. we are solving + # a system of n+1 equations. + n = len(B) - 1 + + # Step 1: forward sweep to eliminate A[I] from each Equation + # allocate lists for modified C and D coefficients, c_prime, d_prime + c_prime = [C[0] / B[0]] + d_prime = [D[0] / B[0]] + + for i in range(1, n + 1): + denom = B[i] - c_prime[i - 1] * A[i] + c_prime.append(C[i] / denom) + d_prime.append((D[i] - d_prime[i - 1] * A[i]) / denom) + + # Step 2: back substitution to solve for x + X = [] + X.append(d_prime[n]) + # Work back to front, solving for x[i] using x[-1] + for i in range(n - 1, -1, -1): + X.append(d_prime[i] - c_prime[i] * X[-1]) + if DEBUG: + return [X[::-1], c_prime, d_prime] + else: + return X[::-1] + # =========== # = HELPERS = # =========== -def control_points(points: Sequence[Point23], extrude_height:float=0, center:bool=True, points_color:Tuple3=Red) -> OpenSCADObject: +def vec3_perp_norm(v: Vector3, other: Vector3) -> Vector3: + """Return the normalized cross product of vectors V and OTHER. + + Used to fetch a normalized vector perpendicular to + both V and OTHER. + """ + return v.cross(other).normalized() + + +def vec2_angle_between(v: Vector3, other: Vector3) -> float: + """Signed angle between the xy components of V and OTHER.""" + return atan2(other.y * v.x - other.x * v.y, v.x * other.x + v.y * other.y) + + +def vec3_angle_between(v: Vector3, other: Vector3) -> float: + """Signed angle between the xy components of V and OTHER.""" + dot_product = v.dot(other) + cross_product = v.cross(other) + return atan2(cross_product.magnitude(), dot_product) + + +def vec3_rotate(v: Vector3, theta: float) -> Vector3: + """Rotate a vector V about its origin by angle THETA.""" + rx = Matrix4.new_rotatex(theta) + ry = Matrix4.new_rotatey(theta) + rz = Matrix4.new_rotatez(theta) + return v * rz * ry * rx + + +def vec2_rotate(v: Vector3, theta: float) -> Vector3: + """Rotate the x and y components of V by THETA while preserving V.z""" + ca = cos(theta) + sa = sin(theta) + return Vector3(*[v.x * ca - v.y * sa, v.x * sa + v.y * ca, v.z]) + + +def flatten(matr: Sequence[Sequence], keep_none_values: bool = True) -> List: + return [item for row in matr for item in row if item or keep_none_values] + + +def control_points( + points: Sequence[Point23], + extrude_height: float = 0, + center: bool = True, + points_color: Tuple3 = Red, +) -> OpenSCADObject: """ Return a list of red cylinders/circles (depending on `extrude_height`) at - a supplied set of 2D points. Useful for visualizing and tweaking a curve's + a supplied set of 2D points. Useful for visualizing and tweaking a curve's control points """ # Figure out how big the circles/cylinders should be based on the spread of points min_bb, max_bb = bounding_box(points) outline_w = max_bb[0] - min_bb[0] outline_h = max_bb[1] - min_bb[1] - r = min(outline_w, outline_h) / 20 # + r = min(min(outline_w, outline_h) / 20, 0.25) # if extrude_height == 0: c = circle(r=r) else: @@ -358,7 +721,10 @@ def control_points(points: Sequence[Point23], extrude_height:float=0, center:boo controls = color(points_color)([translate((p.x, p.y, 0))(c) for p in points]) return controls -def face_strip_list(a_start:int, b_start:int, length:int, close_loop:bool=False) -> List[FaceTrio]: + +def face_strip_list( + a_start: int, b_start: int, length: int, close_loop: bool = False +) -> List[FaceTrio]: # If a_start is the index of the vertex at one end of a row of points in a surface, # and b_start is the index of the vertex at the same end of the next row of points, # return a list of lists of indices describing faces for the whole row: @@ -373,15 +739,16 @@ def face_strip_list(a_start:int, b_start:int, length:int, close_loop:bool=False loop = length - 1 for a, b in zip(range(a_start, a_start + loop), range(b_start, b_start + loop)): - faces.append((a, b+1, b)) - faces.append((a, a+1, b+1)) + faces.append((a, b + 1, b)) + faces.append((a, a + 1, b + 1)) if close_loop: - faces.append((a+loop, b, b+loop)) - faces.append((a+loop, a, b)) + faces.append((a + loop, b, b + loop)) + faces.append((a + loop, a, b)) return faces -def fan_endcap_list(cap_points:int=3, index_start:int=0) -> List[FaceTrio]: - ''' + +def fan_endcap_list(cap_points: int = 3, index_start: int = 0) -> List[FaceTrio]: + """ Return a face-triangles list for the endpoint of a tube with cap_points points We construct a fan of triangles all starting at point index_start and going to each point in turn. @@ -401,45 +768,49 @@ def fan_endcap_list(cap_points:int=3, index_start:int=0) -> List[FaceTrio]: 3 returns: [(0,1,2), (0,2,3), (0,3,4), (0,4,5)] - ''' + """ faces: List[FaceTrio] = [] for i in range(index_start + 1, index_start + cap_points - 1): - faces.append((index_start, i, i+1)) + faces.append((index_start, i, i + 1)) return faces -def centroid_endcap(tube_points:Sequence[Point3], indices:Sequence[int], invert:bool = False) -> Tuple[Point3, List[FaceTrio]]: + +def centroid_endcap( + tube_points: Sequence[Point3], indices: Sequence[int], invert: bool = False +) -> Tuple[Point3, List[FaceTrio]]: # tube_points: all points in a polyhedron tube # indices: the indexes of the points at the desired end of the tube - # invert: if True, invert the order of the generated faces. One endcap in + # invert: if True, invert the order of the generated faces. One endcap in # each pair should be inverted # # Return all the triangle information needed to make an endcap polyhedron # - # This is sufficient for some moderately concave polygonal endcaps, + # This is sufficient for some moderately concave polygonal endcaps, # (a star shape, say), but wouldn't be enough for more irregularly convex - # polygons (anyplace where a segment from the centroid to a point on the + # polygons (anyplace where a segment from the centroid to a point on the # polygon crosses an edge of the polygon) faces: List[FaceTrio] = [] center = centroid([tube_points[i] for i in indices]) centroid_index = len(tube_points) - - for a,b in zip(indices[:-1], indices[1:]): + + for a, b in zip(indices[:-1], indices[1:]): faces.append((centroid_index, a, b)) faces.append((centroid_index, indices[-1], indices[0])) if invert: - faces = list((reversed(f) for f in faces)) # type: ignore + faces = list((reversed(f) for f in faces)) # type: ignore return (center, faces) -def centroid(points:Sequence[Point23]) -> Point23: - total = Point3(0,0,0) + +def centroid(points: Sequence[Point23]) -> Point23: + total = Point3(0, 0, 0) for p in points: total += p total /= len(points) return total -def affine_combination(a:Point23, b:Point23, fraction:float) -> Point23: - # Return a Point[23] between a & b, where fraction==0 => a, fraction==1 => b - return (1-fraction) * a + fraction*b +def affine_combination(a: Point23, b: Point23, fraction: float) -> Point23: + # Return a Point[23] between a & b, where fraction==0 => a, fraction==1 => b + return (1 - fraction) * a + fraction * b diff --git a/solid/test/test_splines.py b/solid/test/test_splines.py index 10641044..c7129812 100755 --- a/solid/test/test_splines.py +++ b/solid/test/test_splines.py @@ -4,27 +4,39 @@ from solid.test.ExpandedTestCase import DiffOutput from solid import * from solid.utils import euclidify -from solid.splines import catmull_rom_points, catmull_rom_prism, bezier_points, bezier_polygon +from solid.splines import ( + catmull_rom_points, + catmull_rom_prism, + bezier_points, + bezier_polygon, + hobby_points, +) from euclid3 import Point2, Point3, Vector2, Vector3 from math import pi SEGMENTS = 8 + class TestSplines(DiffOutput): def setUp(self): self.points = [ - Point3(0,0), - Point3(1,1), - Point3(2,1), + Point3(0, 0), + Point3(1, 1), + Point3(2, 1), + ] + self.points_raw = [ + (0, 0), + (1, 1), + (2, 1), ] - self.points_raw = [ (0,0), (1,1), (2,1), ] self.bezier_controls = [ - Point3(0,0), - Point3(1,1), - Point3(2,1), - Point3(2,-1), + Point3(0, 0), + Point3(1, 1), + Point3(2, 1), + Point3(2, -1), ] - self.bezier_controls_raw = [ (0,0), (1,1), (2,1), (2,-1) ] + self.bezier_controls_raw = [(0, 0), (1, 1), (2, 1), (2, -1)] + self.hobby_omega = 0.69 self.subdivisions = 2 def assertPointsListsEqual(self, a, b): @@ -32,24 +44,46 @@ def assertPointsListsEqual(self, a, b): self.assertEqual(str_list(a), str_list(b)) def test_catmull_rom_points(self): - expected = [Point3(0.00, 0.00), Point3(0.38, 0.44), Point3(1.00, 1.00), Point3(1.62, 1.06), Point3(2.00, 1.00)] - actual = catmull_rom_points(self.points, subdivisions=self.subdivisions, close_loop=False) + expected = [ + Point3(0.00, 0.00), + Point3(0.38, 0.44), + Point3(1.00, 1.00), + Point3(1.62, 1.06), + Point3(2.00, 1.00), + ] + actual = catmull_rom_points( + self.points, subdivisions=self.subdivisions, close_loop=False + ) self.assertPointsListsEqual(expected, actual) # TODO: verify we always have the right number of points for a given call - # verify that `close_loop` always behaves correctly + # verify that `close_loop` always behaves correctly # verify that start_tangent and end_tangent behavior is correct def test_catmull_rom_points_raw(self): # Verify that we can use raw sequences of floats as inputs (e.g [(1,2), (3.2,4)]) - # rather than sequences of Point2s - expected = [Point3(0.00, 0.00), Point3(0.38, 0.44), Point3(1.00, 1.00), Point3(1.62, 1.06), Point3(2.00, 1.00)] - actual = catmull_rom_points(self.points_raw, subdivisions=self.subdivisions, close_loop=False) - self.assertPointsListsEqual(expected, actual) + # rather than sequences of Point2s + expected = [ + Point3(0.00, 0.00), + Point3(0.38, 0.44), + Point3(1.00, 1.00), + Point3(1.62, 1.06), + Point3(2.00, 1.00), + ] + actual = catmull_rom_points( + self.points_raw, subdivisions=self.subdivisions, close_loop=False + ) + self.assertPointsListsEqual(expected, actual) def test_catmull_rom_points_3d(self): - points = [Point3(-1,-1,0), Point3(0,0,1), Point3(1,1,0)] - expected = [Point3(-1.00, -1.00, 0.00), Point3(-0.62, -0.62, 0.50), Point3(0.00, 0.00, 1.00), Point3(0.62, 0.62, 0.50), Point3(1.00, 1.00, 0.00)] + points = [Point3(-1, -1, 0), Point3(0, 0, 1), Point3(1, 1, 0)] + expected = [ + Point3(-1.00, -1.00, 0.00), + Point3(-0.62, -0.62, 0.50), + Point3(0.00, 0.00, 1.00), + Point3(0.62, 0.62, 0.50), + Point3(1.00, 1.00, 0.00), + ] actual = catmull_rom_points(points, subdivisions=2) self.assertPointsListsEqual(expected, actual) @@ -63,34 +97,104 @@ def test_bezier_points_raw(self): # rather than sequences of Point2s expected = [Point3(0.00, 0.00), Point3(1.38, 0.62), Point3(2.00, -1.00)] actual = bezier_points(self.bezier_controls_raw, subdivisions=self.subdivisions) - self.assertPointsListsEqual(expected, actual) + self.assertPointsListsEqual(expected, actual) def test_bezier_points_3d(self): - # verify that we get a valid bezier curve back even when its control points + # verify that we get a valid bezier curve back even when its control points # are outside the XY plane and aren't coplanar - controls_3d = [Point3(-2,-1, 0), Point3(-0.5, -0.5, 1), Point3(0.5, 0.5, 1), Point3(2,1,0)] + controls_3d = [ + Point3(-2, -1, 0), + Point3(-0.5, -0.5, 1), + Point3(0.5, 0.5, 1), + Point3(2, 1, 0), + ] actual = bezier_points(controls_3d, subdivisions=self.subdivisions) - expected = [Point3(-2.00, -1.00, 0.00),Point3(0.00, 0.00, 0.75), Point3(2.00, 1.00, 0.00)] + expected = [ + Point3(-2.00, -1.00, 0.00), + Point3(0.00, 0.00, 0.75), + Point3(2.00, 1.00, 0.00), + ] + self.assertPointsListsEqual(expected, actual) + + def test_hobby_points(self): + expected = [ + Point3(0.00, 0.00, 0.00), + Point3(0.18, 0.46, 0.00), + Point3(0.53, 0.84, 0.00), + Point3(1.00, 1.00, 0.00), + Point3(1.32, 1.11, 0.00), + Point3(1.67, 1.09, 0.00), + Point3(2.00, 1.00, 0.00), + Point3(2.94, 0.74, 0.00), + Point3(2.92, -0.59, 0.00), + Point3(2.00, -1.00, 0.00), + ] + + actual = hobby_points(self.bezier_controls, self.hobby_omega) + self.assertPointsListsEqual(expected, actual) + + def test_hobby_points_raw(self): + expected = [ + Point3(0.00, 0.00, 0.00), + Point3(0.18, 0.46, 0.00), + Point3(0.53, 0.84, 0.00), + Point3(1.00, 1.00, 0.00), + Point3(1.32, 1.11, 0.00), + Point3(1.67, 1.09, 0.00), + Point3(2.00, 1.00, 0.00), + Point3(2.94, 0.74, 0.00), + Point3(2.92, -0.59, 0.00), + Point3(2.00, -1.00, 0.00), + ] + actual = hobby_points(self.bezier_controls, self.hobby_omega) + self.assertPointsListsEqual(expected, actual) + + def test_hobby_points_3d(self): + controls_3d = [ + Point3(-2, -1, 0), + Point3(-0.5, -0.5, 1), + Point3(0.5, 0.5, 1), + Point3(2, 1, 0), + ] + expected = expected = [ + Point3(-2.00, -1.00, 0.00), + Point3(-1.56, -1.06, 0.49), + Point3(-1.01, -0.90, 0.88), + Point3(-0.50, -0.50, 1.00), + Point3(-0.13, -0.22, 1.08), + Point3(0.16, 0.17, 1.01), + Point3(0.50, 0.50, 1.00), + Point3(1.00, 0.98, 0.99), + Point3(1.60, 1.16, 0.55), + Point3(2.00, 1.00, 0.00), + ] + actual = hobby_points(controls_3d, self.hobby_omega) self.assertPointsListsEqual(expected, actual) def test_catmull_rom_prism(self): sides = 3 - UP = Vector3(0,0,1) + UP = Vector3(0, 0, 1) control_points = [[10, 10, 0], [10, 10, 5], [8, 8, 15]] cat_tube = [] - angle_step = 2*pi/sides + angle_step = 2 * pi / sides for i in range(sides): - rotated_controls = list((euclidify(p, Point3).rotate_around(UP, angle_step*i) for p in control_points)) + rotated_controls = list( + ( + euclidify(p, Point3).rotate_around(UP, angle_step * i) + for p in control_points + ) + ) cat_tube.append(rotated_controls) - poly = catmull_rom_prism(cat_tube, self.subdivisions, closed_ring=True, add_caps=True) - actual = (len(poly.params['points']), len(poly.params['faces'])) + poly = catmull_rom_prism( + cat_tube, self.subdivisions, closed_ring=True, add_caps=True + ) + actual = (len(poly.params["points"]), len(poly.params["faces"])) expected = (37, 62) self.assertEqual(expected, actual) - -if __name__ == '__main__': - unittest.main() \ No newline at end of file +if __name__ == "__main__": + unittest.main() From 7438d2cb9481ba1fc5afe1564bc9a6cceb033953 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Thu, 9 May 2024 20:32:17 -0700 Subject: [PATCH 02/14] fix optional closed loop in 3d --- solid/splines.py | 39 ++++++++++++++++++++++++++++----------- 1 file changed, 28 insertions(+), 11 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 60ad98d0..e96beebd 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -442,9 +442,9 @@ def hobby_points( """ # n is defined such that the points can be numbered P[0] ... P[n]. # such that there are a total of n-1 points - if close_loop: - points.append(points[0]) - points.insert(0, points[-2]) + # if close_loop: + points.append(points[0]) + points.insert(0, points[-2]) n = len(points) - 1 if n < 2: raise ValueError( @@ -563,21 +563,38 @@ def hobby_points( # Finally gather up and return the spline points (knots & control points) # as a List of Point2s # res_controls = [(points[i], c0[i], c1[i]) for i in range(0, n - 1)] - res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 1)] - # c0 at points[0] is None as there is no leading control - # c1 at points[-1] is None as there is no leading control + res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 2)] # res_controls[0] = (points[1], c0[-1], c1[0]) # res_controls[-1] = (c0[0], c1[-1]) - res_controls.insert( - 0, - (points[0], c0[-1], c1[1]), - ) + if close_loop: + # Insert the curve from the last original point (points[0]) + # to the first original point, + # res_controls.insert(0, (points[1], c0[1], c1[2])) + res_controls.insert(0, (points[0], c0[-1], c1[1])) + # res_controls[1] = ( + # points[1], + # c0[-1], + # c1[-1], + # ) + # res_controls.insert(0, (c0[-1], c1[-1], points[0])) + res_controls.append( + (points[n - 2], c0[n - 2], c1[n - 2]), + ) + # res_controls.append( + # (points[0], c0[n - 1], c1[-1]), + # ) + res_controls.append((points[0],)) + else: + res_controls.insert(0, (points[0], c0[0], c1[1])) + res_controls.append((points[n - 2],)) + + set_trace() + # res_controls.insert(0, (points[1], c0[0], c1[-1])) # res_controls.append( # (points[n - 2], c0[0], c1[-1]), # ) # res_controls.append((points[n - 1],)) # og_angle = - res_controls.append((points[n - 1],)) if DEBUG: return [ [Point3(*p.xyz) for p in flatten(res_controls)], From 2629678099863ec7a4245d54ddeb77177d21b05b Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Thu, 9 May 2024 21:56:57 -0700 Subject: [PATCH 03/14] set correct unit test values --- solid/test/test_splines.py | 41 ++++++++++++++------------------------ 1 file changed, 15 insertions(+), 26 deletions(-) diff --git a/solid/test/test_splines.py b/solid/test/test_splines.py index c7129812..ecb9e884 100755 --- a/solid/test/test_splines.py +++ b/solid/test/test_splines.py @@ -119,34 +119,25 @@ def test_bezier_points_3d(self): def test_hobby_points(self): expected = [ Point3(0.00, 0.00, 0.00), - Point3(0.18, 0.46, 0.00), - Point3(0.53, 0.84, 0.00), + Point3(0.07, 0.50, 0.00), + Point3(0.52, 0.82, 0.00), Point3(1.00, 1.00, 0.00), - Point3(1.32, 1.11, 0.00), - Point3(1.67, 1.09, 0.00), - Point3(2.00, 1.00, 0.00), - Point3(2.94, 0.74, 0.00), - Point3(2.92, -0.59, 0.00), + Point3(1.33, 1.12, 0.00), + Point3(1.71, 1.19, 0.00), Point3(2.00, -1.00, 0.00), ] - actual = hobby_points(self.bezier_controls, self.hobby_omega) + actual = hobby_points(self.bezier_controls, self.hobby_omega, close_loop=False) self.assertPointsListsEqual(expected, actual) def test_hobby_points_raw(self): expected = [ Point3(0.00, 0.00, 0.00), - Point3(0.18, 0.46, 0.00), - Point3(0.53, 0.84, 0.00), - Point3(1.00, 1.00, 0.00), - Point3(1.32, 1.11, 0.00), - Point3(1.67, 1.09, 0.00), + Point3(-0.10, 0.54, 0.00), + Point3(0.43, 0.91, 0.00), Point3(2.00, 1.00, 0.00), - Point3(2.94, 0.74, 0.00), - Point3(2.92, -0.59, 0.00), - Point3(2.00, -1.00, 0.00), ] - actual = hobby_points(self.bezier_controls, self.hobby_omega) + actual = hobby_points(self.points_raw, self.hobby_omega, close_loop=False) self.assertPointsListsEqual(expected, actual) def test_hobby_points_3d(self): @@ -156,19 +147,17 @@ def test_hobby_points_3d(self): Point3(0.5, 0.5, 1), Point3(2, 1, 0), ] - expected = expected = [ + expected = [ Point3(-2.00, -1.00, 0.00), - Point3(-1.56, -1.06, 0.49), - Point3(-1.01, -0.90, 0.88), + Point3(-1.75, -1.03, 0.62), + Point3(-1.04, -0.90, 0.86), Point3(-0.50, -0.50, 1.00), - Point3(-0.13, -0.22, 1.08), - Point3(0.16, 0.17, 1.01), - Point3(0.50, 0.50, 1.00), - Point3(1.00, 0.98, 0.99), - Point3(1.60, 1.16, 0.55), + Point3(-0.12, -0.22, 1.10), + Point3(0.11, 0.24, 1.13), Point3(2.00, 1.00, 0.00), ] - actual = hobby_points(controls_3d, self.hobby_omega) + + actual = hobby_points(controls_3d, self.hobby_omega, close_loop=False) self.assertPointsListsEqual(expected, actual) def test_catmull_rom_prism(self): From 7184d550a102e6b355fdf429d0a4686340ee04b2 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Thu, 9 May 2024 21:57:20 -0700 Subject: [PATCH 04/14] fix OBOB in hobby results, cleanup --- solid/splines.py | 130 ++++++++--------------------------------------- 1 file changed, 20 insertions(+), 110 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index e96beebd..44ac0e74 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -442,9 +442,14 @@ def hobby_points( """ # n is defined such that the points can be numbered P[0] ... P[n]. # such that there are a total of n-1 points - # if close_loop: + # Ensure points are 3d + if not isinstance(points[0], Point3): + points = [Point3(*point) for point in points] + + # Track edge control points to smooth things out points.append(points[0]) points.insert(0, points[-2]) + n = len(points) - 1 if n < 2: raise ValueError( @@ -456,13 +461,6 @@ def hobby_points( chords = [ Vector3(*(next_point - point).xyz) for point, next_point in pairwise(points) ] - # if close_loop, prepend chord from last point to first - # and append chord from first to last - # or maybe the other way around? - # chords.insert(0, Vector3(*(points[-1] - points[0]).xyz)) - # chords.append(Vector3(*(points[0] - points[-1]).xyz)) - # d = [point.distance(next_point) for point, next_point in pairwise(points)] - d = [chord.magnitude() for chord in chords] if min(d) <= 0: # cannot support successive points being the same @@ -470,14 +468,7 @@ def hobby_points( # gamma[i] is the signed turning angle at P[i], i.e. the angle between # the chords from P[i-1] to P[i] and from P[i] to P[i+1]. - # gamma[0] is undefined; gamma[n] is artificially defined to be zero - gamma = [ - # vec3_angle_between(chord, next_chord) - # for chord, next_chord in pairwise(chords) - chord.angle(next_chord) - for chord, next_chord in pairwise(chords) - ] - # gamma.insert(0, copysign(1, gamma[0])) + gamma = [chord.angle(next_chord) for chord, next_chord in pairwise(chords)] gamma.insert(0, gamma[-1]) gamma.append(gamma[1]) @@ -511,10 +502,7 @@ def hobby_points( C.append(0) D.append(0) - if DEBUG: - alpha, c_prime, d_prime = thomas(A, B, C, D) - else: - alpha = thomas(A, B, C, D) + alpha = thomas(A, B, C, D) # Use alpha (the chord angle) and gamma (the turning angle of the chord # polyline) to solve for beta at each point (beta is like alpha, but for @@ -544,75 +532,24 @@ def hobby_points( c0.append(c0_new_point) c1.append(c1_new_point) - # c0_first_point = (c0[0] / a) * b - # c0_first_point = points[0] + ( - # chords[0].rotate_around(normals[0], alpha[n - 1]).normalize() - # * ((rho(alpha[0], beta[0]) * d[0]) / 3.0) - # ) - # # CORRECT - # c0_last_point = points[-2] + ( - # chords[-1].rotate_around(normals[n], alpha[n - 1]).normalize() - # # * a - # * ((rho(alpha[-1], beta[-1]) * d[-1]) / 3.0) - # ) - # c1_last_point = points[0] - ( - # chords[-1].rotate_around(normals[n - 1], pi / 3 - alpha[0]).normalize() - # * ((rho(beta[0], alpha[0]) * d[0]) / 3.0) - # # * ((rho(beta[-1], beta[-1]) * d[0]) / 3.0) - # ) # Finally gather up and return the spline points (knots & control points) # as a List of Point2s - # res_controls = [(points[i], c0[i], c1[i]) for i in range(0, n - 1)] res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 2)] - # res_controls[0] = (points[1], c0[-1], c1[0]) - # res_controls[-1] = (c0[0], c1[-1]) if close_loop: # Insert the curve from the last original point (points[0]) # to the first original point, - # res_controls.insert(0, (points[1], c0[1], c1[2])) - res_controls.insert(0, (points[0], c0[-1], c1[1])) - # res_controls[1] = ( - # points[1], - # c0[-1], - # c1[-1], - # ) - # res_controls.insert(0, (c0[-1], c1[-1], points[0])) + res_controls.insert(0, (points[1], c0[1], c1[1])) res_controls.append( (points[n - 2], c0[n - 2], c1[n - 2]), ) - # res_controls.append( - # (points[0], c0[n - 1], c1[-1]), - # ) - res_controls.append((points[0],)) + res_controls.append( + (points[n - 1], c0[n - 1], c1[0]), + ) + res_controls.append((points[1],)) else: - res_controls.insert(0, (points[0], c0[0], c1[1])) - res_controls.append((points[n - 2],)) - - set_trace() - # res_controls.insert(0, (points[1], c0[0], c1[-1])) - # res_controls.append( - # (points[n - 2], c0[0], c1[-1]), - # ) - # res_controls.append((points[n - 1],)) - # og_angle = - if DEBUG: - return [ - [Point3(*p.xyz) for p in flatten(res_controls)], - c0, - c1, - chords, - d, - gamma, - alpha, - c_prime, - d_prime, - beta, - A, - B, - C, - D, - normals, - ] + # Set the first control point group's c1 to the + res_controls.insert(0, (points[1], c0[1], c1[1])) + res_controls.append((points[n - 1],)) return [Point3(*p.xyz) for p in flatten(res_controls)] @@ -665,10 +602,10 @@ def thomas(A: Sequence, B: Sequence, C: Sequence, D: Sequence) -> Sequence: # Work back to front, solving for x[i] using x[-1] for i in range(n - 1, -1, -1): X.append(d_prime[i] - c_prime[i] * X[-1]) - if DEBUG: - return [X[::-1], c_prime, d_prime] - else: - return X[::-1] + # if DEBUG: + # return [X[::-1], c_prime, d_prime] + # else: + return X[::-1] # =========== @@ -683,33 +620,6 @@ def vec3_perp_norm(v: Vector3, other: Vector3) -> Vector3: return v.cross(other).normalized() -def vec2_angle_between(v: Vector3, other: Vector3) -> float: - """Signed angle between the xy components of V and OTHER.""" - return atan2(other.y * v.x - other.x * v.y, v.x * other.x + v.y * other.y) - - -def vec3_angle_between(v: Vector3, other: Vector3) -> float: - """Signed angle between the xy components of V and OTHER.""" - dot_product = v.dot(other) - cross_product = v.cross(other) - return atan2(cross_product.magnitude(), dot_product) - - -def vec3_rotate(v: Vector3, theta: float) -> Vector3: - """Rotate a vector V about its origin by angle THETA.""" - rx = Matrix4.new_rotatex(theta) - ry = Matrix4.new_rotatey(theta) - rz = Matrix4.new_rotatez(theta) - return v * rz * ry * rx - - -def vec2_rotate(v: Vector3, theta: float) -> Vector3: - """Rotate the x and y components of V by THETA while preserving V.z""" - ca = cos(theta) - sa = sin(theta) - return Vector3(*[v.x * ca - v.y * sa, v.x * sa + v.y * ca, v.z]) - - def flatten(matr: Sequence[Sequence], keep_none_values: bool = True) -> List: return [item for row in matr for item in row if item or keep_none_values] From f12617833db900b73fcd1e970100c495fe34875f Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Thu, 9 May 2024 22:27:17 -0700 Subject: [PATCH 05/14] remove old debug output --- solid/splines.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 44ac0e74..14c5ce15 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -602,9 +602,6 @@ def thomas(A: Sequence, B: Sequence, C: Sequence, D: Sequence) -> Sequence: # Work back to front, solving for x[i] using x[-1] for i in range(n - 1, -1, -1): X.append(d_prime[i] - c_prime[i] * X[-1]) - # if DEBUG: - # return [X[::-1], c_prime, d_prime] - # else: return X[::-1] From 8af9491251347e746cebbd195e058e65271ebdf3 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 00:27:14 -0700 Subject: [PATCH 06/14] initial hobby splines examples --- solid/examples/hobby_splines_example.py | 56 +++++++++++++++++++++++++ 1 file changed, 56 insertions(+) create mode 100644 solid/examples/hobby_splines_example.py diff --git a/solid/examples/hobby_splines_example.py b/solid/examples/hobby_splines_example.py new file mode 100644 index 00000000..3478fd6d --- /dev/null +++ b/solid/examples/hobby_splines_example.py @@ -0,0 +1,56 @@ +import sys +import os +from solid import polygon, scad_render_to_file +from solid.splines import SEGMENTS, bezier_points, control_points +from solid.examples.path_extrude_example import circle_points +from euclid3 import Vector3, Point3 + +from solid import splines +from solid.utils import right, extrude_along_path + +SCALE = 30 +OMEGA = 1.0 +SUBDIVISIONS = 12 +SEGMENTS = 12 + + +def assembly(): + knot_points = [ + Point3(0, 0, -1), + Point3(1, 3, -2), + Point3(2, 1, 1), + Point3(3, -1, 2), + Point3(0, -4, 0), + ] + knot_points = [kp * SCALE for kp in knot_points] + path = splines.hobby_points(knot_points, OMEGA, close_loop=True) + path_open = splines.hobby_points(knot_points, OMEGA, close_loop=False) + bzp = [] # Closed Hobby spline control points used to make beziers + bzp_open = [] # Open Hobby spline control for same thing + for i in range(0, len(path) - 3, 3): + # Controls take the form of: Knot, control vec, control vec, knot + controls = path[i : i + 4] + controls_open = path_open[i : i + 4] + bzp += bezier_points(controls, subdivisions=SUBDIVISIONS) + if len(controls_open) == 4: + # PATH_OPEN may have fewer segments than PATH so assume last valid + # group of 4 points is the final segment of the open curve + bzp_open += bezier_points(controls_open, subdivisions=SUBDIVISIONS) + + assembly = polygon(bzp) + right(5 * SCALE)( + extrude_along_path(circle_points(0.3 * SCALE, SEGMENTS), bzp) + + right(5 * SCALE)( + extrude_along_path(circle_points(0.3 * SCALE, SEGMENTS), bzp_open) + + right(5 * SCALE)(polygon(bzp_open)) + ) + ) + return assembly + + +if __name__ == "__main__": + out_dir = sys.argv[1] if len(sys.argv) > 1 else os.curdir + + a = assembly() + + out_path = scad_render_to_file(a, out_dir=out_dir, include_orig_code=True) + print(f"{__file__}: SCAD file written to: \n{out_path}") From eb7e8fde51896079b1db0feb45d08808a7abe8a7 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 00:44:47 -0700 Subject: [PATCH 07/14] clarify return type and function purpose --- solid/splines.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 14c5ce15..2fb4ca1f 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -423,7 +423,7 @@ def _bez33(u: float) -> float: # https://tug.org/TUGboat/tb34-2/tb107jackowski.pdf def hobby_points( points: List[Point3], omega: float = 0, close_loop: bool = False -) -> List[Point2]: +) -> List[Point3]: """Hobby's algorithm: given a set of points, fit a Bézier spline to them. The chosen splines tend to have pleasing, rounded shapes. @@ -612,7 +612,7 @@ def vec3_perp_norm(v: Vector3, other: Vector3) -> Vector3: """Return the normalized cross product of vectors V and OTHER. Used to fetch a normalized vector perpendicular to - both V and OTHER. + both V and OTHER about which to rotate when forming smooth curves. """ return v.cross(other).normalized() From 6267314b970ac5a2fd50a44ec8b2a62198934429 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 10:43:28 -0700 Subject: [PATCH 08/14] cleanup open vs closed loop results --- solid/splines.py | 10 +++------- 1 file changed, 3 insertions(+), 7 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 2fb4ca1f..2f0ee1c1 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -534,20 +534,16 @@ def hobby_points( # Finally gather up and return the spline points (knots & control points) # as a List of Point2s - res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 2)] + res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 1)] if close_loop: - # Insert the curve from the last original point (points[0]) - # to the first original point, + # Insert the curve from the last input point + # to the first input point, res_controls.insert(0, (points[1], c0[1], c1[1])) - res_controls.append( - (points[n - 2], c0[n - 2], c1[n - 2]), - ) res_controls.append( (points[n - 1], c0[n - 1], c1[0]), ) res_controls.append((points[1],)) else: - # Set the first control point group's c1 to the res_controls.insert(0, (points[1], c0[1], c1[1])) res_controls.append((points[n - 1],)) return [Point3(*p.xyz) for p in flatten(res_controls)] From 3511b4765bcce022ed396699f91acdb56b8b32f5 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 11:23:04 -0700 Subject: [PATCH 09/14] curvier example, fix pytest expected path --- solid/examples/hobby_splines_example.py | 6 +++--- solid/test/test_splines.py | 20 +++++++++++++++++--- 2 files changed, 20 insertions(+), 6 deletions(-) diff --git a/solid/examples/hobby_splines_example.py b/solid/examples/hobby_splines_example.py index 3478fd6d..2b17b0fd 100644 --- a/solid/examples/hobby_splines_example.py +++ b/solid/examples/hobby_splines_example.py @@ -17,10 +17,10 @@ def assembly(): knot_points = [ Point3(0, 0, -1), - Point3(1, 3, -2), - Point3(2, 1, 1), + Point3(1, -2, -2), + Point3(2, 0, 1), Point3(3, -1, 2), - Point3(0, -4, 0), + Point3(3, 1, -2), ] knot_points = [kp * SCALE for kp in knot_points] path = splines.hobby_points(knot_points, OMEGA, close_loop=True) diff --git a/solid/test/test_splines.py b/solid/test/test_splines.py index ecb9e884..dd275976 100755 --- a/solid/test/test_splines.py +++ b/solid/test/test_splines.py @@ -124,6 +124,9 @@ def test_hobby_points(self): Point3(1.00, 1.00, 0.00), Point3(1.33, 1.12, 0.00), Point3(1.71, 1.19, 0.00), + Point3(2.00, 1.00, 0.00), + Point3(2.63, 0.60, 0.00), + Point3(2.35, -0.28, 0.00), Point3(2.00, -1.00, 0.00), ] @@ -133,11 +136,19 @@ def test_hobby_points(self): def test_hobby_points_raw(self): expected = [ Point3(0.00, 0.00, 0.00), - Point3(-0.10, 0.54, 0.00), - Point3(0.43, 0.91, 0.00), + Point3(0.07, 0.50, 0.00), + Point3(0.52, 0.82, 0.00), + Point3(1.00, 1.00, 0.00), + Point3(1.33, 1.12, 0.00), + Point3(1.71, 1.19, 0.00), Point3(2.00, 1.00, 0.00), + Point3(2.63, 0.60, 0.00), + Point3(2.35, -0.28, 0.00), + Point3(2.00, -1.00, 0.00), ] - actual = hobby_points(self.points_raw, self.hobby_omega, close_loop=False) + actual = hobby_points( + self.bezier_controls_raw, self.hobby_omega, close_loop=False + ) self.assertPointsListsEqual(expected, actual) def test_hobby_points_3d(self): @@ -154,6 +165,9 @@ def test_hobby_points_3d(self): Point3(-0.50, -0.50, 1.00), Point3(-0.12, -0.22, 1.10), Point3(0.11, 0.24, 1.13), + Point3(0.50, 0.50, 1.00), + Point3(1.01, 0.84, 0.83), + Point3(1.55, 0.88, 0.44), Point3(2.00, 1.00, 0.00), ] From 68e4166d9ecdd1737cda30163f20b5538cff1833 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 11:28:26 -0700 Subject: [PATCH 10/14] remove unwanted pudb import --- solid/splines.py | 2 -- 1 file changed, 2 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 2f0ee1c1..3464509e 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -1,8 +1,6 @@ #! /usr/bin/env python from math import copysign, pow, cos, sin, pi, atan2, acos from itertools import pairwise -from pudb import set_trace - from solid import ( union, circle, From d77180ad33f4af8141f740463a61a6532d307610 Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 11:51:24 -0700 Subject: [PATCH 11/14] undo formatting-only applied in development --- solid/extrude_along_path.py | 156 +++++++++++------------------------- 1 file changed, 47 insertions(+), 109 deletions(-) diff --git a/solid/extrude_along_path.py b/solid/extrude_along_path.py index d88b939e..bdd94760 100644 --- a/solid/extrude_along_path.py +++ b/solid/extrude_along_path.py @@ -10,43 +10,38 @@ FacetIndices = Tuple[int, int, int] Point3Transform = Callable[[Point3, Optional[float], Optional[float]], Point3] - # ========================== # = Extrusion along a path = # ========================== -def extrude_along_path( - shape_pts: Points, - path_pts: Points, - scales: Sequence[Union[Vector2, float, Tuple2]] = None, - rotations: Sequence[float] = None, - transforms: Sequence[Callable] = None, - connect_ends=False, - cap_ends=True, - transform_args: List = None, - transform_kwargs: Dict = None, -) -> OpenSCADObject: - """ +def extrude_along_path( shape_pts:Points, + path_pts:Points, + scales:Sequence[Union[Vector2, float, Tuple2]] = None, + rotations: Sequence[float] = None, + transforms: Sequence[Point3Transform] = None, + connect_ends = False, + cap_ends = True) -> OpenSCADObject: + ''' Extrude the curve defined by shape_pts along path_pts. -- For predictable results, shape_pts must be planar, convex, and lie in the XY plane centered around the origin. *Some* nonconvexity (e.g, star shapes) and nonplanarity will generally work fine - + -- len(scales) should equal len(path_pts). No-op if not supplied - Each entry may be a single number for uniform scaling, or a pair of + Each entry may be a single number for uniform scaling, or a pair of numbers (or Point2) for differential X/Y scaling If not supplied, no scaling will occur. - + -- len(rotations) should equal 1 or len(path_pts). No-op if not supplied. Each point in shape_pts will be rotated by rotations[i] degrees at each point in path_pts. Or, if only one rotation is supplied, the shape will be rotated smoothly over rotations[0] degrees in the course of the extrusion - + -- len(transforms) should be 1 or be equal to len(path_pts). No-op if not supplied. - Each entry should be have the signature: + Each entry should be have the signature: def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 where path_norm is in [0,1] and expresses progress through the extrusion and loop_norm is in [0,1] and express progress through a single loop of the extrusion - + -- if connect_ends is True, the first and last loops of the extrusion will be joined, which is useful for toroidal geometries. Overrides cap_ends @@ -54,10 +49,11 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 will be connected to the centroid of that loop. For planar, convex shapes, this works nicely. If shape is less planar or convex, some self-intersection may happen. Not applied if connect_ends is True - """ + ''' + - polyhedron_pts: Points = [] - facet_indices: List[Tuple[int, int, int]] = [] + polyhedron_pts:Points= [] + facet_indices:List[Tuple[int, int, int]] = [] # Make sure we've got Euclid Point3's for all elements shape_pts = euclidify(shape_pts, Point3) @@ -70,7 +66,7 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 tangent_path_points: List[Point3] = [] # If first & last points are the same, let's close the shape - first_last_equal = (path_pts[0] - path_pts[-1]).magnitude_squared() < EPSILON + first_last_equal = ((path_pts[0] - path_pts[-1]).magnitude_squared() < EPSILON) if first_last_equal: connect_ends = True path_pts = path_pts[:][:-1] @@ -81,14 +77,11 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 first = Point3(*(path_pts[0] - (path_pts[1] - path_pts[0]))) last = Point3(*(path_pts[-1] - (path_pts[-2] - path_pts[-1]))) tangent_path_points = [first] + path_pts + [last] - tangents = [ - tangent_path_points[i + 2] - tangent_path_points[i] - for i in range(len(path_pts)) - ] + tangents = [tangent_path_points[i+2] - tangent_path_points[i] for i in range(len(path_pts))] for which_loop in range(len(path_pts)): # path_normal is 0 at the first path_pts and 1 at the last - path_normal = which_loop / (len(path_pts) - 1) + path_normal = which_loop/ (len(path_pts) - 1) path_pt = path_pts[which_loop] tangent = tangents[which_loop] @@ -96,53 +89,21 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 rotate_degrees = None if rotations: - rotate_degrees = ( - rotations[which_loop] - if len(rotations) > 1 - else rotations[0] * path_normal - ) + rotate_degrees = rotations[which_loop] if len(rotations) > 1 else rotations[0] * path_normal transform_func = None if transforms: - transform_func = ( - transforms[which_loop] if len(transforms) > 1 else transforms[0] - ) - # Default to no *args or **kwargs - this_transform_func_args = None - if transform_args: - this_transform_func_args = ( - transform_args[which_loop] - if len(transform_args) > 1 - else transform_args[0] - ) - this_transform_func_kwargs = None - if transform_kwargs: - this_transform_func_kwargs = ( - transform_kwargs[which_loop] - if len(transform_kwargs) > 1 - else transform_kwargs[0] - ) + transform_func = transforms[which_loop] if len(transforms) > 1 else transforms[0] this_loop = shape_pts[:] - this_loop = _scale_loop( - this_loop, - scale, - ) + this_loop = _scale_loop(this_loop, scale) this_loop = _rotate_loop(this_loop, rotate_degrees) - this_loop = _transform_loop( - this_loop, - transform_func, - path_normal, - this_transform_func_args, - this_transform_func_kwargs, - ) - - this_loop = transform_to_point( - this_loop, dest_point=path_pt, dest_normal=tangent, src_up=src_up - ) + this_loop = _transform_loop(this_loop, transform_func, path_normal) + + this_loop = transform_to_point(this_loop, dest_point=path_pt, dest_normal=tangent, src_up=src_up) loop_start_index = which_loop * shape_pt_count - if which_loop < len(path_pts) - 1: + if (which_loop < len(path_pts) - 1): loop_facets = _loop_facet_indices(loop_start_index, shape_pt_count) facet_indices += loop_facets @@ -156,61 +117,46 @@ def transform_func(p:Point3, path_norm:float, loop_norm:float): Point3 elif cap_ends: # OpenSCAD's polyhedron will automatically triangulate faces as needed. - # So just include all points at each end of the tube - last_loop_start_index = len(polyhedron_pts) - shape_pt_count + # So just include all points at each end of the tube + last_loop_start_index = len(polyhedron_pts) - shape_pt_count start_loop_indices = list(reversed(range(shape_pt_count))) - end_loop_indices = list( - range(last_loop_start_index, last_loop_start_index + shape_pt_count) - ) + end_loop_indices = list(range(last_loop_start_index, last_loop_start_index + shape_pt_count)) facet_indices.append(start_loop_indices) facet_indices.append(end_loop_indices) - return polyhedron(points=euc_to_arr(polyhedron_pts), faces=facet_indices) # type: ignore - + return polyhedron(points=euc_to_arr(polyhedron_pts), faces=facet_indices) # type: ignore -def _loop_facet_indices( - loop_start_index: int, loop_pt_count: int, next_loop_start_index=None -) -> List[FacetIndices]: +def _loop_facet_indices(loop_start_index:int, loop_pt_count:int, next_loop_start_index=None) -> List[FacetIndices]: facet_indices: List[FacetIndices] = [] # nlsi == next_loop_start_index if next_loop_start_index == None: next_loop_start_index = loop_start_index + loop_pt_count - loop_indices = list(range(loop_start_index, loop_pt_count + loop_start_index)) + [ - loop_start_index - ] - next_loop_indices = list( - range(next_loop_start_index, loop_pt_count + next_loop_start_index) - ) + [next_loop_start_index] + loop_indices = list(range(loop_start_index, loop_pt_count + loop_start_index)) + [loop_start_index] + next_loop_indices = list(range(next_loop_start_index, loop_pt_count + next_loop_start_index )) + [next_loop_start_index] for i, (a, b) in enumerate(zip(loop_indices[:-1], loop_indices[1:])): - c, d = next_loop_indices[i : i + 2] + c, d = next_loop_indices[i: i+2] # OpenSCAD's polyhedron will accept quads and do its own triangulation with them, - # so we could just append (a,b,d,c). + # so we could just append (a,b,d,c). # However, this lets OpenSCAD (Or CGAL?) do its own triangulation, leading # to some strange outcomes. Prefer to do our own triangulation. # c--d # |\ | # | \| - # a--b + # a--b # facet_indices.append((a,b,d,c)) - facet_indices.append((a, b, c)) - facet_indices.append((b, d, c)) + facet_indices.append((a,b,c)) + facet_indices.append((b,d,c)) return facet_indices - -def _rotate_loop( - points: Sequence[Point3], rotation_degrees: float = None -) -> List[Point3]: +def _rotate_loop(points:Sequence[Point3], rotation_degrees:float=None) -> List[Point3]: if rotation_degrees is None: return points - up = Vector3(0, 0, 1) + up = Vector3(0,0,1) rads = radians(rotation_degrees) return [p.rotate_around(up, rads) for p in points] - -def _scale_loop( - points: Sequence[Point3], scale: Union[float, Point2, Tuple2] = None -) -> List[Point3]: +def _scale_loop(points:Sequence[Point3], scale:Union[float, Point2, Tuple2]=None) -> List[Point3]: if scale is None: return points @@ -218,14 +164,7 @@ def _scale_loop( scale = [scale] * 2 return [Point3(point.x * scale[0], point.y * scale[1], point.z) for point in points] - -def _transform_loop( - points: Sequence[Point3], - transform_func: Callable = None, - path_normal: float = None, - *args, - **kwargs, -) -> List[Point3]: +def _transform_loop(points:Sequence[Point3], transform_func:Point3Transform = None, path_normal:float = None) -> List[Point3]: # transform_func is a function that takes a point and optionally two floats, # a `path_normal`, in [0,1] that indicates where this loop is in a path extrusion, # and `loop_normal` in [0,1] that indicates where this point is in a list of points @@ -235,9 +174,8 @@ def _transform_loop( result = [] for i, p in enumerate(points): # i goes from 0 to 1 across points - loop_normal = i / (len(points) - 1) - # print(f"args:\t{args}") - # print(f"kwargs:\t{kwargs}") - new_p = transform_func(p, path_normal, loop_normal, *args, **kwargs) + loop_normal = i/(len(points) -1) + new_p = transform_func(p, path_normal, loop_normal) result.append(new_p) return result + From 9e39c3db19dd8f7f6f0a131353c914dc50d82c6b Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 12:04:55 -0700 Subject: [PATCH 12/14] cleanup gathering points in open vs closed loop --- solid/splines.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/solid/splines.py b/solid/splines.py index 3464509e..941126ec 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -532,17 +532,16 @@ def hobby_points( # Finally gather up and return the spline points (knots & control points) # as a List of Point2s - res_controls = [(points[i], c0[i], c1[i]) for i in range(2, n - 1)] + res_controls = [(points[i], c0[i], c1[i]) for i in range(1, n - 1)] if close_loop: # Insert the curve from the last input point # to the first input point, - res_controls.insert(0, (points[1], c0[1], c1[1])) res_controls.append( (points[n - 1], c0[n - 1], c1[0]), ) res_controls.append((points[1],)) else: - res_controls.insert(0, (points[1], c0[1], c1[1])) + # Append the last input point only res_controls.append((points[n - 1],)) return [Point3(*p.xyz) for p in flatten(res_controls)] From ac2d58386cf66ce8600c52d3663131ac838c692e Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 12:07:21 -0700 Subject: [PATCH 13/14] fix wording in comment --- solid/splines.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/solid/splines.py b/solid/splines.py index 941126ec..05844069 100644 --- a/solid/splines.py +++ b/solid/splines.py @@ -534,7 +534,7 @@ def hobby_points( # as a List of Point2s res_controls = [(points[i], c0[i], c1[i]) for i in range(1, n - 1)] if close_loop: - # Insert the curve from the last input point + # Append the curve from the last input point # to the first input point, res_controls.append( (points[n - 1], c0[n - 1], c1[0]), From 1d07fa3df2ff8279325a88b7406832ec13504d6a Mon Sep 17 00:00:00 2001 From: cuttlefisch Date: Fri, 10 May 2024 13:14:46 -0700 Subject: [PATCH 14/14] bump version --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 8b31bb73..5fe235dc 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [tool.poetry] name = "solidpython" -version = "1.1.3" +version = "1.1.4" description = "Python interface to the OpenSCAD declarative geometry language" authors = ["Evan Jones "] homepage = "https://github.com/SolidCode/SolidPython"