Feat: change Segment to accept analytic functions

src/shapepy/bool2d/boolean.py

@@ -145,7 +145,7 @@ def pursue_path(
             if (index_jordan, index_segment) in matrix:
                 break
             matrix.append((index_jordan, index_segment))
-            last_point = segment.ctrlpoints[-1]
+            last_point = segment(1)
             possibles = []
             for i, jordan in enumerate(jordans):
                 if i == index_jordan:
@@ -157,7 +157,7 @@ def pursue_path(
                 continue
             index_jordan = possibles[0]
             for j, segj in enumerate(all_segments[index_jordan]):
-                if segj.ctrlpoints[0] == last_point:
+                if segj(0) == last_point:
                     index_segment = j
                     break
         return CyclicContainer(matrix)
@@ -221,7 +221,7 @@ def midpoints_one_shape(
 
         """
         for i, jordan in enumerate(shapea.jordans):
-            for j, segment in enumerate(jordan.piecewise):
+            for j, segment in enumerate(jordan.parametrize()):
                 mid_point = segment(Fraction(1, 2))
                 density = shapeb.density(mid_point)
                 mid_point_in = (density > 0 and closed) or density == 1

src/shapepy/bool2d/primitive.py

@@ -6,17 +6,15 @@
 
 """
 
+from __future__ import annotations
+
 import math
-from copy import copy
 from typing import Tuple
 
 import numpy as np
 
 from ..geometry.factory import FactoryJordan
-from ..geometry.jordancurve import JordanCurve
-from ..geometry.point import Point2D, cartesian
-from ..geometry.segment import Segment
-from ..geometry.unparam import USegment
+from ..geometry.point import Point2D
 from ..loggers import debug
 from ..tools import Is, To
 from .base import EmptyShape, WholeShape
@@ -225,25 +223,5 @@ def circle(
             raise ValueError
         if not Is.integer(ndivangle) or ndivangle < 4:
             raise ValueError
-        center = To.point(center)
-
-        angle = math.tau / ndivangle
-        height = np.tan(angle / 2)
-
-        start_point = radius * cartesian(1, 0)
-        middle_point = radius * cartesian(1, height)
-        beziers = []
-        for _ in range(ndivangle - 1):
-            end_point = copy(start_point).rotate(angle)
-            new_bezier = Segment([start_point, middle_point, end_point])
-            beziers.append(new_bezier)
-            start_point = end_point
-            middle_point = copy(middle_point).rotate(angle)
-        end_point = beziers[0].ctrlpoints[0]
-        new_bezier = Segment([start_point, middle_point, end_point])
-        beziers.append(new_bezier)
-
-        jordan_curve = JordanCurve(map(USegment, beziers))
-        jordan_curve.move(center)
-        circle = SimpleShape(jordan_curve)
-        return circle
+        jordan_curve = FactoryJordan.circle(ndivangle)
+        return SimpleShape(jordan_curve).scale(radius).move(center)

src/shapepy/bool2d/shape.py

@@ -45,7 +45,7 @@ def __deepcopy__(self, memo) -> SimpleShape:
 
     def __str__(self) -> str:  # pragma: no cover  # For debug
         area = float(self.area)
-        vertices = tuple(map(tuple, self.jordan.vertices))
+        vertices = tuple(map(tuple, self.jordan.vertices()))
         return f"SimpleShape[{area:.2f}]:[{vertices}]"
 
     def __eq__(self, other: SubSetR2) -> bool:

src/shapepy/geometry/base.py

@@ -67,7 +67,6 @@ def parametrize(self) -> IParametrizedCurve:
     def __or__(self, other: IGeometricCurve) -> IGeometricCurve:
         return Future.concatenate((self, other))
 
-    @abstractmethod
     def move(self, vector: Point2D) -> IGeometricCurve:
         """
         Moves/translate entire shape by an amount

src/shapepy/geometry/concatenate.py

@@ -4,103 +4,93 @@
 
 from typing import Iterable, Union
 
-import pynurbs
-
-from ..tools import Is, NotExpectedError, To
-from .base import IParametrizedCurve
+from ..analytic import Bezier
+from ..tools import Is, NotExpectedError
+from .base import IGeometricCurve
 from .piecewise import PiecewiseCurve
-from .point import cross, inner
+from .point import cross
 from .segment import Segment
+from .unparam import UPiecewiseCurve, USegment
 
 
-def concatenate(curves: Iterable[IParametrizedCurve]) -> IParametrizedCurve:
+def concatenate(curves: Iterable[IGeometricCurve]) -> IGeometricCurve:
     """
     Concatenates the given curves.
 
     Ignores all the curves parametrization
     """
     curves = tuple(curves)
-    if not all(Is.instance(curve, IParametrizedCurve) for curve in curves):
+    if not all(Is.instance(curve, IGeometricCurve) for curve in curves):
         raise ValueError
-    # Check if the curves are connected
-    for i, curvei in enumerate(curves[:-1]):
-        curvej = curves[i + 1]
-        if curvei(curvei.knots[-1]) != curvej(curvej.knots[0]):
-            raise ValueError("The curves to concatenate are not connected")
-    segments = []
-    for curve in curves:
-        if Is.instance(curve, Segment):
-            segments.append(curve)
-        elif Is.instance(curve, PiecewiseCurve):
-            segments += list(curve)
-        else:
-            raise NotExpectedError(f"Unknown type: {type(curve)}")
-    return concatenate_segments(segments)
+    if all(Is.instance(curve, Segment) for curve in curves):
+        return concatenate_segments(curves)
+    if all(Is.instance(curve, USegment) for curve in curves):
+        return concatenate_usegments(curves)
+    raise NotExpectedError(str(tuple(str(type(c)) for c in curves)))
 
 
-def concatenate_segments(segments: Iterable[Segment]) -> IParametrizedCurve:
+def concatenate_usegments(
+    usegments: Iterable[USegment],
+) -> Union[USegment, UPiecewiseCurve]:
+    """
+    Concatenates all the unparametrized segments
+    """
+    usegments = tuple(usegments)
+    assert all(Is.instance(useg, USegment) for useg in usegments)
+    union = concatenate_segments(useg.parametrize() for useg in usegments)
+    return (
+        USegment(union)
+        if Is.instance(union, Segment)
+        else UPiecewiseCurve(map(USegment, union))
+    )
+
+
+def concatenate_segments(
+    segments: Iterable[Segment],
+) -> Union[Segment, PiecewiseCurve]:
     """
     Concatenates all the segments
     """
-    segments = list(segments)
-    assert all(map(Is.segment, segments))
+    segments = tuple(segments)
+    if len(segments) == 0:
+        raise ValueError(f"Number sizes: {len(segments)}")
     filtsegments = []
-    segment0 = segments.pop(0)
-    while len(segments) > 0:
-        segment1 = segments.pop(0)
+    segments = iter(segments)
+    segmenti = next(segments)
+    for segmentj in segments:
         try:
-            segment0 = bezier_and_bezier(segment0, segment1)
+            segmenti = bezier_and_bezier(segmenti, segmentj)
         except ValueError:
-            filtsegments.append(segment0)
-            segment0 = segment1
-    filtsegments.append(segment0)
-    try:
-        union = bezier_and_bezier(filtsegments[-1], filtsegments[0])
-        filtsegments.pop(0)
-        filtsegments.pop()
-        filtsegments.append(union)
-    except ValueError:
-        pass
+            filtsegments.append(segmenti)
+            segmenti = segmentj
+    filtsegments.append(segmenti)
+    print("filtsegments = ", filtsegments)
     return (
         PiecewiseCurve(filtsegments)
         if len(filtsegments) > 1
         else filtsegments[0]
     )
 
 
-def bezier_and_bezier(
-    curvea: Segment, curveb: Segment
-) -> Union[Segment, PiecewiseCurve]:
+def bezier_and_bezier(curvea: Segment, curveb: Segment) -> Segment:
     """Computes the union of two bezier curves"""
-    assert Is.instance(curvea, Segment)
-    assert Is.instance(curveb, Segment)
-    assert curvea.degree == curveb.degree
-    if curvea.ctrlpoints[-1] != curveb.ctrlpoints[0]:
+    if not Is.instance(curvea, Segment):
+        raise TypeError(f"Invalid type: {type(curvea)}")
+    if not Is.instance(curveb, Segment):
+        raise TypeError(f"Invalid type: {type(curveb)}")
+    if abs(cross(curvea(1, 1), curveb(0, 1))) > 1e-6:
         raise ValueError
-    # Last point of first derivative
-    dapt = curvea.ctrlpoints[-1] - curvea.ctrlpoints[-2]
-    # First point of first derivative
-    dbpt = curveb.ctrlpoints[1] - curveb.ctrlpoints[0]
-    if abs(cross(dapt, dbpt)) > 1e-6:
-        node = To.rational(1, 2)
-    else:
-        dsumpt = dapt + dbpt
-        denomin = inner(dsumpt, dsumpt)
-        node = inner(dapt, dsumpt) / denomin
-    knotvectora = pynurbs.GeneratorKnotVector.bezier(
-        curvea.degree, To.rational
-    )
-    knotvectora.scale(node)
-    knotvectorb = pynurbs.GeneratorKnotVector.bezier(
-        curveb.degree, To.rational
-    )
-    knotvectorb.scale(1 - node).shift(node)
-    newknotvector = tuple(knotvectora) + tuple(
-        knotvectorb[curvea.degree + 1 :]
-    )
-    finalcurve = pynurbs.Curve(newknotvector)
-    finalcurve.ctrlpoints = tuple(curvea.ctrlpoints) + tuple(curveb.ctrlpoints)
-    finalcurve.knot_clean((node,))
-    if finalcurve.degree + 1 != finalcurve.npts:
-        raise ValueError("Union is not a bezier curve!")
-    return Segment(finalcurve.ctrlpoints)
+    if curvea.xfunc.degree != curveb.xfunc.degree:
+        raise ValueError
+    if curvea.yfunc.degree != curveb.yfunc.degree:
+        raise ValueError
+    if curvea.xfunc.degree > 1 or curvea.yfunc.degree > 1:
+        raise ValueError
+
+    if curvea(1) != curveb(0):
+        raise ValueError
+    startpoint = curvea(0)
+    endpoint = curveb(1)
+    nxfunc = Bezier((startpoint[0], endpoint[0]))
+    nyfunc = Bezier((startpoint[1], endpoint[1]))
+    return Segment(nxfunc, nyfunc)

src/shapepy/geometry/factory.py

@@ -2,16 +2,45 @@
 Defines the Factory to create Jordan Curves
 """
 
-from typing import Tuple
+from __future__ import annotations
 
+import math
+from copy import copy
+from typing import Iterable, Tuple
+
+import numpy as np
+
+from ..analytic import Bezier
 from ..loggers import debug
 from ..tools import To
 from .jordancurve import JordanCurve
-from .point import Point2D
-from .segment import Segment, clean_segment
+from .point import Point2D, cartesian
+from .segment import Segment
 from .unparam import USegment
 
 
+# pylint: disable=too-few-public-methods
+class FactorySegment:
+    """
+    Define functions to create Segments
+    """
+
+    @staticmethod
+    @debug("shapepy.geometry.factory")
+    def bezier(ctrlpoints: Iterable[Point2D]) -> Segment:
+        """Initialize a bezier segment from a list of control points
+
+        :param ctrlpoints: The list of control points
+        :type ctrlpoints: Iterable[Point2D]
+        :return: The created segment
+        :rtype: Segment
+        """
+        ctrlpoints = tuple(map(To.point, ctrlpoints))
+        xfunc = Bezier((pt[0] for pt in ctrlpoints))
+        yfunc = Bezier((pt[1] for pt in ctrlpoints))
+        return Segment(xfunc, yfunc)
+
+
 class FactoryJordan:
     """
     Define functions to create Jordan Curves
@@ -43,7 +72,7 @@ def polygon(vertices: Tuple[Point2D, ...]) -> JordanCurve:
         beziers = [None] * nverts
         for i in range(nverts):
             ctrlpoints = vertices[i : i + 2]
-            new_bezier = Segment(ctrlpoints)
+            new_bezier = FactorySegment.bezier(ctrlpoints)
             beziers[i] = USegment(new_bezier)
         return JordanCurve(beziers)
 
@@ -76,6 +105,28 @@ def spline_curve(spline_curve) -> JordanCurve:
         """
         beziers = spline_curve.split(spline_curve.knots)
         segments = (
-            clean_segment(Segment(bezier.ctrlpoints)) for bezier in beziers
+            FactorySegment.bezier(bezier.ctrlpoints).clean()
+            for bezier in beziers
         )
         return JordanCurve(map(USegment, segments))
+
+    @staticmethod
+    @debug("shapepy.geometry.factory")
+    def circle(ndivangle: int):
+        """Creates a jordan curve that represents a unit circle"""
+        angle = math.tau / ndivangle
+        height = np.tan(angle / 2)
+
+        start_point = cartesian(1, 0)
+        middle_point = cartesian(1, height)
+        all_ctrlpoints = []
+        for _ in range(ndivangle - 1):
+            end_point = copy(start_point).rotate(angle)
+            all_ctrlpoints.append([start_point, middle_point, end_point])
+            start_point = end_point
+            middle_point = copy(middle_point).rotate(angle)
+        end_point = all_ctrlpoints[0][0]
+        all_ctrlpoints.append([start_point, middle_point, end_point])
+        return JordanCurve(
+            map(USegment, map(FactorySegment.bezier, all_ctrlpoints))
+        )

src/shapepy/geometry/intersection.py

@@ -214,6 +214,11 @@ def param_and_param(
     raise NotExpectedError
 
 
+def segment_is_linear(segment: Segment) -> bool:
+    """Tells if the segment is a polynomial linear"""
+    return segment.xfunc.degree <= 1 and segment.yfunc.degree <= 1
+
+
 def segment_and_segment(
     curvea: Segment, curveb: Segment
 ) -> Tuple[SubSetR1, SubSetR1]:
@@ -224,11 +229,12 @@ def segment_and_segment(
         return Empty(), Empty()
     if curvea == curveb:
         return Interval(0, 1), Interval(0, 1)
-    if curvea.degree == 1 and curveb.degree == 1:
-        subseta, subsetb = IntersectionSegments.lines(curvea, curveb)
-        return subseta, subsetb
-    usample = list(NodeSampleFactory.closed_linspace(curvea.npts + 3))
-    vsample = list(NodeSampleFactory.closed_linspace(curveb.npts + 3))
+    if segment_is_linear(curvea) and segment_is_linear(curveb):
+        return IntersectionSegments.lines(curvea, curveb)
+    nptsa = max(curvea.xfunc.degree, curvea.yfunc.degree) + 4
+    nptsb = max(curveb.xfunc.degree, curveb.yfunc.degree) + 4
+    usample = list(NodeSampleFactory.closed_linspace(nptsa))
+    vsample = list(NodeSampleFactory.closed_linspace(nptsb))
     pairs = []
     for ui in usample:
         pairs += [(ui, vj) for vj in vsample]
@@ -264,11 +270,9 @@ class IntersectionSegments:
     @staticmethod
     def lines(curvea: Segment, curveb: Segment) -> Tuple[SubSetR1, SubSetR1]:
         """Finds the intersection of two line segments"""
-        assert curvea.degree == 1
-        assert curveb.degree == 1
         empty = Empty()
-        A0, A1 = curvea.ctrlpoints
-        B0, B1 = curveb.ctrlpoints
+        A0, A1 = curvea(0), curvea(1)
+        B0, B1 = curveb(0), curveb(1)
         dA = A1 - A0
         dB = B1 - B0
         B0mA0 = B0 - A0