Patro/Common/ArgparseAction.py

+####################################################################################################
+#
+# Patro - A Python library to make patterns for fashion design
+# Copyright (C) 2017 Fabrice Salvaire
+#
+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program.  If not, see <http://www.gnu.org/licenses/>.
+#
+####################################################################################################
+
+"""Module to implement argparse actions.
+
+"""
+
+####################################################################################################
+
+__all__ = [
+    'PathAction',
+]
+
+####################################################################################################
+
+import argparse
+from pathlib import Path
+
+####################################################################################################
+
+class PathAction(argparse.Action):
+
+    """Class to implement argparse action for path."""
+
+    ##############################################
+
+    def __call__(self, parser, namespace, values, option_string=None):
+
+        if values is not None:
+            if isinstance(values, list):
+                path = [Path(x) for x in values]
+            else:
+                path = Path(values)
+        else:
+            path = None
+        setattr(namespace, self.dest, path)

Patro/Common/Math/Root.py

@@ -78,27 +78,56 @@ def cubic_root(a, b, c, d):
 
 ####################################################################################################
 
-def cubic_root_sympy(a, b, c, d):
+def x_symbol():
+    return sympy.Symbol('x', real=True)
+
+def real_roots(expression, x):
+    return [i.n() for i in sympy.real_roots(expression, x)]
 
-    x = sympy.Symbol('x', real=True)
+####################################################################################################
+
+def cubic_root_sympy(a, b, c, d):
+    x = x_symbol()
     E = a*x**3 + b*x**2 + c*x + d
+    return real_roots(E, x)
 
-    return [i.n() for i in sympy.real_roots(E, x)]
+####################################################################################################
+
+def cubic_root_normalised(a, b, c):
+    x = x_symbol()
+    E = x**3 + a*x**2 + b*x + c
+    return real_roots(E, x)
 
 ####################################################################################################
 
-def fifth_root_normalised(a, b, c, d, e):
+def fourth_root_normalised(a, b, c, d):
+    x = x_symbol()
+    E = x**4 + a*x**3 + b*x**2 + c*x + d
+    return real_roots(E, x)
 
-    x = sympy.Symbol('x', real=True)
-    E = x**5 + a*x**4 + b*x**3 + c*x**2 + d*x + e
+####################################################################################################
 
-    return [i.n() for i in sympy.real_roots(E, x)]
+def fifth_root_sympy(a, b, c, d, e, f):
+    x = x_symbol()
+    E = a*x**5 + b*x**4 + c*x**3 + d*x**2 + e*x + f
+    return real_roots(E, x)
+
+####################################################################################################
+
+def fifth_root_normalised(a, b, c, d, e):
+    x = x_symbol()
+    E = x**5 + a*x**4 + b*x**3 + c*x**2 + d*x + e
+    return real_roots(E, x)
 
 ####################################################################################################
 
 def fifth_root(*args):
+    # Fixme: RuntimeWarning: divide by zero encountered in double_scalars
     a = args[0]
-    return fifth_root_normalised(*[x/a for x in args[1:]])
+    if a == 0:
+        return fifth_root_sympy(*args)
+    else:
+        return fifth_root_normalised(*[x/a for x in args[1:]])
 
 ####################################################################################################
 

Patro/Common/Xml/XmlFile.py

@@ -39,8 +39,13 @@ class XmlFileMixin:
 
     ##############################################
 
-    def __init__(self, path):
-        self._path = Path(path)
+    def __init__(self, path, data=None):
+
+        if path is not None:
+            self._path = Path(path)
+        else:
+            self._path = None
+        self._data = data
 
     ##############################################
 
@@ -52,8 +57,17 @@ class XmlFileMixin:
 
     def parse(self):
         """Parse a XML file and return the etree"""
-        with open(str(self._path), 'rb') as f:
-            source = f.read()
+
+        data = self._data
+        if data is None:
+            with open(str(self._path), 'rb') as f:
+                source = f.read()
+        else:
+            if isinstance(data, bytes):
+               source = data
+            else:
+                source = bytes(str(self._data).strip(), 'utf-8')
+
         return etree.fromstring(source)
 
     ##############################################

Patro/FileFormat/Dxf/Importer.py

@@ -128,15 +128,15 @@ class DxfImporter:
         major_axis = self._to_vector(item_dxf.major_axis)
         minor_axis = major_axis * item_dxf.ratio
         domain = AngularDomain(item_dxf.start_param, item_dxf.end_param, degrees=False)
-        x_radius, y_radius = major_axis.magnitude, minor_axis.magnitude
+        radius_x, radius_y = major_axis.magnitude, minor_axis.magnitude
         angle = major_axis.orientation
         if angle == 90:
-            x_radius, y_radius = y_radius, x_radius
+            radius_x, radius_y = radius_y, radius_x
             angle = 0
         # Fixme: ...
         ellipse = Ellipse2D(
             center,
-            x_radius, y_radius,
+            radius_x, radius_y,
             angle,
             domain=domain,
         )

Patro/FileFormat/Dxf/Polyline.py

@@ -221,7 +221,7 @@ class Polyline:
                     stop_angle += 360
                 if vertex1.bulge < 0:
                     start_angle, stop_angle = stop_angle, start_angle
-                print('bulb', vertex1, vertex2, vertex1.bulge, start_angle, stop_angle)
+                # print('bulb', vertex1, vertex2, vertex1.bulge, start_angle, stop_angle)
                 arc.domain = AngularDomain(start_angle, stop_angle)
                 # arc = Circle2D(center, vertex1.bulge_radius, domain=AngularDomain(start_angle, stop_angle))
                 items.append(arc)

Patro/FileFormat/Svg/SvgFile.py

@@ -27,6 +27,7 @@ Import Algorithm:
 * line with a small polygon at extremities is a grainline
 * expect pieces are delimited by a path
 * check for paths sharing vertexes and stroke style
+
 """
 
 ####################################################################################################
@@ -35,7 +36,10 @@ import logging
 
 from lxml import etree
 
+from IntervalArithmetic import Interval2D
+
 from Patro.Common.Xml.XmlFile import XmlFileMixin
+from Patro.GeometryEngine.Transformation import AffineTransformation2D
 from . import SvgFormat
 
 ####################################################################################################
@@ -46,19 +50,126 @@ _module_logger = logging.getLogger(__name__)
 
 class RenderState:
 
+    # Fixme: convert type !!!
+
+    STATES = [name for name in SvgFormat.PresentationAttributes.__dict__.keys()
+              if not name.startswith('_')]
+
+    ##############################################
+
+    @classmethod
+    def to_python(cls, value):
+
+        # Fixme: move ???
+
+        if isinstance(value, str):
+            if value == 'none':
+                return None
+            else:
+                try:
+                    float_value = float(value)
+                    if '.' in value:
+                        return float_value
+                    else:
+                        return int(float_value)
+                except ValueError:
+                    pass
+
+        return value
+
+    ##############################################
+
+    def __init__(self, item=None):
+
+        # Init from item else use default value
+        for state in self.STATES:
+            if item is not None and hasattr(item, state):
+                value = self.to_python(getattr(item, state))
+            else:
+                value = getattr(SvgFormat.PresentationAttributes, state)
+            setattr(self, state, value)
+
+    ##############################################
+
+    def clone(self):
+        return self.__class__(self)
+
+    ##############################################
+
+    def to_dict(self, all=False):
+
+        if all:
+            return {state:getattr(self, state) for state in self.STATES}
+        else:
+            d = {}
+            for state in self.STATES:
+                value = getattr(self, state)
+                if value is not None:
+                    d[state] = value
+            return d
+
+    ##############################################
+
+    def merge(self, item):
+
+        for state in self.STATES:
+            if hasattr(item, state):
+                value = getattr(item, state)
+                if state == 'transform':
+                    if value is not None:
+                        # Transform matrix is composed from top to item
+                        # thus left to right
+                        self.transform = self.transform * value
+                elif state == 'style':
+                    pass
+                else:
+                    setattr(self, state, self.to_python(value))
+
+        # Merge style
+        style = getattr(item, 'style', None)
+        if style is not None:
+            for pair in style.split(';'):
+                state, value = [x.strip() for x in pair.split(':')]
+                state = state.replace('-', '_')
+                if state == 'transform':
+                    self.transform = self.transform * value
+                else:
+                    setattr(self, state, self.to_python(value))
+
+        return self
+
+    ##############################################
+
+    def __str__(self):
+        return str(self.to_dict())
+
+####################################################################################################
+
+class RenderStateStack:
+
     ##############################################
 
     def __init__(self):
 
-        self._transformations = []
+        self._stack = [RenderState()]
+
+    ##############################################
+
+    @property
+    def state(self):
+        return self._stack[-1]
 
     ##############################################
 
-    def push_transformation(self, transformation):
-        self._transformations.append(transformation)
+    def push(self, kwargs):
+        new_state = self.state.clone()
+        new_state.merge(kwargs)
+        self._stack.append(new_state)
+
+    ##############################################
 
-    def pop_transformation(self):
-        self._transformations.pop()
+    def pop(self):
+        self._stack.pop()
 
 ####################################################################################################
 
@@ -107,18 +218,30 @@ class SvgDispatcher:
         # 'use',
     ]
 
+    _logger = _module_logger.getChild('SvgDispatcher')
+
     ##############################################
 
-    def __init__(self, root):
+    def __init__(self, reader):
+
+        self._reader =reader
+        self.reset()
+        # self.on_root(root)
 
-        self._state = RenderState()
+    ##############################################
 
-        self.on_root(root)
+    def reset(self):
+        self._state_stack = RenderStateStack()
 
     ##############################################
 
-    def element_tag(self, element):
+    @property
+    def state(self):
+        return self._state_stack.state
+
+    ##############################################
 
+    def element_tag(self, element):
         tag = element.tag
         if '{' in tag:
             tag = tag[tag.find('}')+1:]
@@ -131,7 +254,7 @@ class SvgDispatcher:
         tag = self.element_tag(element)
         tag_class = self.__TAGS__[tag]
         if tag_class is not None:
-            print(element, tag_class)
+            # self._logger.info('\n{}  /  {}'.format(element, tag_class))
             return tag_class(element)
         else:
             raise NotImplementedError
@@ -149,44 +272,54 @@ class SvgDispatcher:
 
     ##############################################
 
-    def on_group(self, group):
+    def on_group(self, element):
 
-        self.on_root(group)
+        group = self.from_xml(element)
+        # self._logger.info('Group: {}\n{}'.format(group.id, group))
+        self._reader.on_group(group)
+
+        self._state_stack.push(group)
+        # self._logger.info('State:\n' + str(self.state))
+
+        self.on_root(element)
 
     ##############################################
 
     def on_graphic_item(self, element):
 
         item = self.from_xml(element)
-        print(item)
+        # self._logger.info('Item: {}\n{}'.format(item.id, item))
+        self._reader.on_graphic_item(item)
 
 ####################################################################################################
 
-class SvgFile(XmlFileMixin):
-
-    """Class to read/write SVG file."""
-
-    _logger = _module_logger.getChild('SvgFile')
+class SvgFileMixin:
 
     SVG_DOCTYPE = '<!DOCTYPE svg PUBLIC "-//W3C//DTD SVG 1.1//EN" "http://www.w3.org/Graphics/SVG/1.1/DTD/svg11.dtd">'
     SVG_xmlns = 'http://www.w3.org/2000/svg'
     SVG_xmlns_xlink = 'http://www.w3.org/1999/xlink'
     SVG_version = '1.1'
 
-    ##############################################
+####################################################################################################
+
+class SvgFileInternal(XmlFileMixin, SvgFileMixin):
 
-    def __init__(self, path=None):
+    """Class to read/write SVG file."""
+
+    _logger = _module_logger.getChild('SvgFile')
+
+    __dispatcher_cls__ = SvgDispatcher
+
+    ##############################################
 
-        # Fixme: path
-        if path is None:
-            path = ''
+    def __init__(self, path, data=None):
 
-        XmlFileMixin.__init__(self, path)
+        super().__init__(path, data)
 
-        # Fixme:
-        # if path is not None:
-        if path != '':
-            self._read()
+        # Fixme: API
+        #  purpose of dispatcher, where must be state ???
+        self._dispatcher = self.__dispatcher_cls__(self)
+        self._read()
 
     ##############################################
 
@@ -200,14 +333,70 @@ class SvgFile(XmlFileMixin):
         #      width="1000.0pt" height="1000.0" viewBox="0 0 1000.0 1000.0"
         # ></svg>
 
-        tree = self._parse()
-        dispatch = SvgDispatcher(tree)
-        # Fixme: ...
+        tree = self.parse()
+
+        svg_root = self._dispatcher.from_xml(tree)
+        self.on_svg_root(svg_root)
+
+        self._dispatcher.on_root(tree)
+
+    ##############################################
+
+    @property
+    def view_box(self):
+        return self._view_box
+
+    @property
+    def width(self):
+        return self._width
+
+    @property
+    def height(self):
+        return self._height
+
+    ##############################################
+
+    def on_svg_root(self, svg_root):
+        x_inf, y_inf, x_sup, y_sup = svg_root.view_box
+        self._view_box = Interval2D((x_inf, x_sup), (y_inf, y_sup))
+        self._width = svg_root.width
+        self._height = svg_root.height
+
+    ##############################################
+
+    def on_group(self, group):
+        self._logger.info('Group: {}\n{}'.format(group.id, group))
+
+    ##############################################
+
+    def on_graphic_item(self, item):
+
+        self._logger.info('Item: {}\n{}'.format(item.id, item))
+        state = self._dispatcher.state.clone().merge(item)
+        self._logger.info('Item State:\n' + str(state))
+
+####################################################################################################
+
+class SvgFileWriter(SvgFileMixin):
+
+    """Class to write a SVF file."""
+
+    _logger = _module_logger.getChild('SvgFileWriter')
+
+    COMMENT = 'Pattern created with Patro (https://github.com/FabriceSalvaire/Patro)'
+
+    ##############################################
+
+    def __init__(self, path, paper, root_tree, transformation=None):
+
+        self._path = str(path)
+
+        self._write(paper, root_tree, transformation)
 
     ##############################################
 
     @classmethod
-    def new_root(cls, paper):
+    def _new_root(cls, paper):
 
         nsmap = {
             None: cls.SVG_xmlns,
@@ -223,15 +412,15 @@ class SvgFile(XmlFileMixin):
         attrib['viewBox'] = '0 0 {:.3f} {:.3f}'.format(paper.width, paper.height)
 
         # Fixme: from conf
-        root.append(etree.Comment('Pattern created with Patro (https://github.com/FabriceSalvaire/Patro)'))
+        root.append(etree.Comment(cls.COMMENT))
 
         return root
 
     ##############################################
 
-    def write(self, paper, root_tree, transformation=None, path=None):
+    def _write(self, paper, root_tree, transformation=None):
 
-        root = self.new_root(paper)
+        root = self._new_root(paper)
 
         # Fixme: implement tree, look at lxml
         if transformation:
@@ -244,14 +433,20 @@ class SvgFile(XmlFileMixin):
         for element in root_tree:
             group.append(element.to_xml())
 
-        if path is None:
-            path = self.path
-
         tree = etree.ElementTree(root)
-        tree.write(str(path),
+        tree.write(self._path,
                    pretty_print=True,
                    xml_declaration=True,
                    encoding='utf-8',
                    standalone=False,
                    doctype=self.SVG_DOCTYPE,
         )
+
+####################################################################################################
+
+class SvgFile:
+
+    ##############################################
+
+    def __init__(self, path):
+        self._interval = SvgFileInternal(path)

Patro/FileFormat/Svg/SvgFormat.py

@@ -66,7 +66,8 @@ __all__ = [
 
 ####################################################################################################
 
-# import logging
+import logging
+
 from Patro.Common.Xml.Objectivity import (
     # BoolAttribute,
     IntAttribute, FloatAttribute,
@@ -76,8 +77,14 @@ from Patro.Common.Xml.Objectivity import (
     TextXmlObjectAdaptator,
 )
 
-# from Patro.GeometryEngine.Vector import Vector2D
+# Fixme: should we mix SVG format and ... ???
+from Patro.GeometryEngine.Path import Path2D
 from Patro.GeometryEngine.Transformation import AffineTransformation2D
+from Patro.GeometryEngine.Vector import Vector2D
+
+####################################################################################################
+
+_module_logger = logging.getLogger(__name__)
 
 ####################################################################################################
 
@@ -127,71 +134,78 @@ class IdMixin:
 #
 ####################################################################################################
 
-####################################################################################################
-#
-# Presentation Attributes
-#  alignment-baseline
-#  baseline-shift
-#  clip
-#  clip-path
-#  clip-rule
-#  color
-#  color-interpolation
-#  color-interpolation-filters
-#  color-profile
-#  color-rendering
-#  cursor
-#  direction
-#  display
-#  dominant-baseline
-#  enable-background
-#  fill
-#  fill-opacity
-#  fill-rule
-#  filter
-#  flood-color
-#  flood-opacity
-#  font-family
-#  font-size
-#  font-size-adjust
-#  font-stretch
-#  font-style
-#  font-variant
-#  font-weight
-#  glyph-orientation-horizontal
-#  glyph-orientation-vertical
-#  image-rendering
-#  kerning
-#  letter-spacing
-#  lighting-color
-#  marker-end
-#  marker-mid
-#  marker-start
-#  mask
-#  opacity
-#  overflow
-#  pointer-events
-#  shape-rendering
-#  stop-color
-#  stop-opacity
-#  stroke
-#  stroke-dasharray
-#  stroke-dashoffset
-#  stroke-linecap
-#  stroke-linejoin
-#  stroke-miterlimit
-#  stroke-opacity
-#  stroke-width
-#  text-anchor
-#  text-decoration
-#  text-rendering
-#  unicode-bidi
-#  visibility
-#  word-spacing
-#  writing-mode
-#
 ####################################################################################################
 
+class PresentationAttributes:
+
+    # https://developer.mozilla.org/en-US/docs/Web/SVG/Attribute
+
+    # Fixme: type !!!
+
+    alignment_baseline = None # check
+    baseline_shift = None # check
+    clip = None # check
+    clip_path = None # check
+    clip_rule = None # check
+    color = None # check
+    color_interpolation = None # check
+    color_interpolation_filters = None # check
+    color_profile = None # check
+    color_rendering = None # check
+    cursor = None # check
+    direction = None # check
+    display = None # check
+    dominant_baseline = None # check
+    enable_background = None # check
+    fill = None
+    fill_opacity = 1
+    fill_rule = 'nonzero'
+    #! filter = None # check
+    flood_color = None # check
+    flood_opacity = None # check
+    font_family = None # check
+    font_size = None # check
+    font_size_adjust = None # check
+    font_stretch = None # check
+    font_style = None # check
+    font_variant = None # check
+    font_weight = None # check
+    glyph_orientation_horizontal = None # check
+    glyph_orientation_vertical = None # check
+    image_rendering = None # check
+    kerning = None # check
+    letter_spacing = None # check
+    lighting_color = None # check
+    marker_end = None # check
+    marker_mid = None # check
+    marker_start = None # check
+    mask = None # check
+    opacity = 1
+    overflow = None # check
+    paint_order = None
+    pointer_events = None # check
+    shape_rendering = None # check
+    stop_color = None # check
+    stop_opacity = None # check
+    stroke = None
+    stroke_dasharray = None # check
+    stroke_dashoffset = None # check
+    stroke_linecap = 'butt'
+    stroke_linejoin = 'miter'
+    stroke_miterlimit = 4
+    stroke_opacity = 1
+    stroke_width = 1 # px
+    style = None
+    text_anchor = None # check
+    text_decoration = None # check
+    text_rendering = None # check
+    transform = AffineTransformation2D.Identity()
+    unicode_bidi = None # check
+    vector_effect = None
+    visibility = None # check
+    word_spacing = None # check
+    writing_mode = None # check
+
 ####################################################################################################
 
 class InheritAttribute(StringAttribute):
@@ -439,13 +453,59 @@ class TransformAttribute(StringAttribute):
                 values = [float(x) for x in value[pos0+1:-1].split(',')]
                 transforms.append((transform_type, values))
                 # Fixme:
-            return transforms
+
+            # return transforms
+            return cls.to_python(transforms, concat=True)
 
     ##############################################
 
     @classmethod
     def to_xml(cls, value):
-        return 'matrix({})'.format(' '.join([str(x) for x in value.to_list()])) # Fixme: to func
+        # Fixme: wrong if value is AffineTransformation2D !!!
+        # Fixme: to func
+        return 'matrix({})'.format(' '.join([str(x) for x in value.to_list()]))
+
+    ##############################################
+
+    @classmethod
+    def to_python(cls, transforms, concat=True):
+
+        def complete(values, size):
+            return values + [0]*(size - len(values))
+
+        global_transform = AffineTransformation2D.Identity()
+        py_transforms = []
+        for name, values in transforms:
+            transform = None
+            if name == 'matrix':
+                array = [values[i] for i in (0, 2, 4, 1, 3, 5)] + [0, 0, 1]
+                transform = AffineTransformation2D(array)
+            elif name == 'translate':
+                vector = Vector2D(complete(values, 2))
+                transform = AffineTransformation2D.Translation(vector)
+            elif name == 'scale':
+                transform = AffineTransformation2D.Scale(*values)
+            elif name == 'rotate':
+                angle, *vector = complete(values, 3)
+                vector = Vector2D(vector)
+                transform = AffineTransformation2D.RotationAt(vector, angle)
+            elif name == 'skewX':
+                angle = values[0]
+                raise NotImplementedError
+            elif name == 'skewY':
+                angle = values[0]
+                raise NotImplementedError
+            else:
+                raise NotImplementedError
+            if concat:
+                global_transform = transform * global_transform
+            else:
+                py_transforms.append(transform)
+
+        if concat:
+            return global_transform
+        else:
+            return py_transforms
 
 ####################################################################################################
 
@@ -891,23 +951,27 @@ class PathDataAttribute(StringAttribute):
 
     COMMANDS = ''.join(NUMBER_OF_ARGS.keys())
 
+    _logger = _module_logger.getChild('PathDataAttribute')
+
     ##############################################
 
     @classmethod
-    def from_xml(cls, value):
+    def from_xml(cls, svg_path):
+
+        # cls._logger.info('SVG path:\n'+ svg_path)
 
         # Replace comma separator by space
-        value = value.replace(',', ' ')
+        cleaned_svg_path = svg_path.replace(',', ' ')
         # Add space after letter
         data_path = ''
-        for c in value:
+        for c in cleaned_svg_path:
             data_path += c
             if c.isalpha:
                 data_path += ' '
-        # Convert float
+        # Convert float values
         parts = []
-        for part in split_space_list(value):
-            if not(len(part) == 1 and part.isalpha):
+        for part in split_space_list(cleaned_svg_path):
+            if not(len(part) == 1 and part.isalpha()):
                 part = float(part)
             parts.append(part)
 
@@ -921,17 +985,23 @@ class PathDataAttribute(StringAttribute):
                 command = part
                 command_lower = command.lower()
                 if command_lower not in cls.COMMANDS:
-                    raise ValueError('Invalid path instruction')
+                    raise ValueError("Invalid path instruction: '{}' in\n{}".format(command, svg_path))
                 number_of_args = cls.NUMBER_OF_ARGS[command_lower]
+                i += 1 # move to first arg
             # else repeated instruction
-            next_i = i + number_of_args + 1
-            values = parts[i+1:next_i]
-            #! points = [Vector2D(values[2*i], values[2*i+1]) for i in range(number_of_args / 2)]
-            points = values
-            commands.append((command, points))
+            next_i = i + number_of_args
+            args = parts[i:next_i]
+            commands.append((command, args))
             i = next_i
+            # for implicit line to
+            if command == 'm':
+                command = 'l'
+            elif command == 'M':
+                command = 'L'
 
-        return commands
+        # return commands
+        # Fixme: do later ???
+        return cls.to_geometry(commands)
 
     ##############################################
 
@@ -945,6 +1015,62 @@ class PathDataAttribute(StringAttribute):
             path_data += ' '.join(list(command[0]) + [str(x) for x in command[1]])
         return path_data
 
+    ##############################################
+
+    @classmethod
+    def as_vector(cls, args):
+
+        number_of_args = len(args)
+        number_of_vectors = number_of_args // 2
+        if number_of_args != number_of_vectors * 2:
+            raise ValueError('len(args) is not // 2: {}'.format(number_of_args))
+        return [Vector2D(args[i:i+2]) for i in range(0, number_of_args, 2)]
+
+    ##############################################
+
+    @classmethod
+    def to_geometry(cls, commands):
+
+        # cls._logger.info('Path:\n' + str(commands).replace('), ', '),\n '))
+        path = None
+        for command, args in commands:
+            command_lower = command.lower()
+            absolute = command_lower != command # Upper case means absolute
+            # if is_lower:
+            #     cls._logger.warning('incremental command')
+            #     raise NotImplementedError
+            if path is None:
+                if command_lower != 'm':
+                    raise NameError('Path must start with m')
+                path = Path2D(args) # Vector2D()
+            else:
+                if command_lower == 'l':
+                    path.line_to(args, absolute=absolute)
+                elif command == 'h':
+                    path.horizontal_to(*args, absolute=False)
+                elif command == 'H':
+                    path.absolute_horizontal_to(*args)
+                elif command_lower == 'v':
+                    path.vertical_to(*args, absolute=absolute)
+                elif command == 'V':
+                    path.absolute_vertical_to(*args)
+                elif command_lower == 'c':
+                    path.cubic_to(*cls.as_vector(args), absolute=absolute)
+                elif command_lower == 's':
+                    path.stringed_quadratic_to(*cls.as_vector(args), absolute=absolute)
+                elif command_lower == 'q':
+                    path.quadratic_to(*cls.as_vector(args), absolute=absolute)
+                elif command_lower == 't':
+                    path.stringed_cubic_to(*cls.as_vector(args), absolute=absolute)
+                elif command_lower == 'a':
+                    radius_x, radius_y, angle, large_arc, sweep, x, y = args
+                    point = Vector2D(x, y)
+                    path.arc_to(point, radius_x, radius_y, angle, bool(large_arc), bool(sweep), absolute=absolute)
+                elif command_lower == 'z':
+                    path.close()
+
+        return path
+
 ####################################################################################################
 
 class Path(PathMixin, SvgElementMixin, XmlObjectAdaptator):
@@ -1028,6 +1154,32 @@ class Rect(PositionMixin, RadiusMixin, SizeMixin, PathMixin, SvgElementMixin, Xm
 
     __tag__ = 'rect'
 
+    ##############################################
+
+    @property
+    def geometry(self):
+
+        # Fixme: width is str
+        width = float(self.width)
+        height = float(self.height)
+
+        # Fixme: which one ???
+        radius_x = self.rx
+        radius_y = self.ry
+        if radius_y == 0:
+            radius = None
+        else:
+            radius = radius_y
+
+        point = Vector2D(self.x, self.y)
+        path = Path2D(point)
+        path.horizontal_to(width)
+        path.vertical_to(height, radius=radius)
+        path.horizontal_to(-width, radius=radius)
+        path.close(radius=radius, close_radius=radius)
+
+        return path
+
 ####################################################################################################
 
 class Stop(XmlObjectAdaptator):

Patro/FileFormat/Valentina/Measurement.py

@@ -75,7 +75,7 @@ class VitFileInternal(XmlFileMixin):
 
     def __init__(self, path):
 
-        XmlFileMixin.__init__(self, path)
+        super().__init__(path)
         self.measurements = ValentinaMeasurements()
         self.read()
 

Patro/GeometryEngine/Bezier.py

@@ -20,185 +20,20 @@
 
 r"""Module to implement Bézier curve.
 
-Definitions
------------
-
-A Bézier curve is defined by a set of control points :math:`\mathbf{P}_0` through
-:math:`\mathbf{P}_n`, where :math:`n` is called its order (:math:`n = 1` for linear, 2 for
-quadratic, 3 for cubic etc.). The first and last control points are always the end points of the
-curve;
-
-In the following :math:`0 \le t \le 1`.
-
-Linear Bézier Curves
----------------------
-
-Given distinct points :math:`\mathbf{P}_0` and :math:`\mathbf{P}_1`, a linear Bézier curve is simply
-a straight line between those two points. The curve is given by
-
-.. math::
-    \begin{align}
-    \mathbf{B}(t) &= \mathbf{P}_0 + t (\mathbf{P}_1 - \mathbf{P}_0) \\
-                  &= (1-t) \mathbf{P}_0 + t \mathbf{P}_1
-    \end{align}
-
-and is equivalent to linear interpolation.
-
-Quadratic Bézier Curves
------------------------
-
-A quadratic Bézier curve is the path traced by the function :math:`\mathbf{B}(t)`, given points
-:math:`\mathbf{P}_0`, :math:`\mathbf{P}_1`, and :math:`\mathbf{P}_2`,
-
-.. math::
-    \mathbf{B}(t) = (1 - t)[(1 - t) \mathbf{P}_0 + t \mathbf{P}_1] + t [(1 - t) \mathbf{P}_1 + t \mathbf{P}_2]
-
-which can be interpreted as the linear interpolant of corresponding points on the linear Bézier
-curves from :math:`\mathbf{P}_0` to :math:`\mathbf{P}_1` and from :math:`\mathbf{P}_1` to
-:math:`\mathbf{P}_2` respectively.
-
-Rearranging the preceding equation yields:
-
-.. math::
-    \mathbf{B}(t) = (1 - t)^{2} \mathbf{P}_0 + 2(1 - t)t \mathbf{P}_1 + t^{2} \mathbf{P}_2
-
-This can be written in a way that highlights the symmetry with respect to :math:`\mathbf{P}_1`:
-
-.. math::
-    \mathbf{B}(t) = \mathbf{P}_1 + (1 - t)^{2} ( \mathbf{P}_0 - \mathbf{P}_1) + t^{2} (\mathbf{P}_2 - \mathbf{P}_1)
-
-Which immediately gives the derivative of the Bézier curve with respect to `t`:
-
-.. math::
-    \mathbf{B}'(t) = 2(1 - t) (\mathbf{P}_1 - \mathbf{P}_0) + 2t (\mathbf{P}_2 - \mathbf{P}_1)
-
-from which it can be concluded that the tangents to the curve at :math:`\mathbf{P}_0` and
-:math:`\mathbf{P}_2` intersect at :math:`\mathbf{P}_1`. As :math:`t` increases from 0 to 1, the
-curve departs from :math:`\mathbf{P}_0` in the direction of :math:`\mathbf{P}_1`, then bends to
-arrive at :math:`\mathbf{P}_2` from the direction of :math:`\mathbf{P}_1`.
-
-The second derivative of the Bézier curve with respect to :math:`t` is
-
-.. math::
-    \mathbf{B}''(t) = 2 (\mathbf{P}_2 - 2 \mathbf{P}_1 + \mathbf{P}_0)
-
-Cubic Bézier Curves
--------------------
-
-Four points :math:`\mathbf{P}_0`, :math:`\mathbf{P}_1`, :math:`\mathbf{P}_2` and
-:math:`\mathbf{P}_3` in the plane or in higher-dimensional space define a cubic Bézier curve.  The
-curve starts at :math:`\mathbf{P}_0` going toward :math:`\mathbf{P}_1` and arrives at
-:math:`\mathbf{P}_3` coming from the direction of :math:`\mathbf{P}_2`. Usually, it will not pass
-through :math:`\mathbf{P}_1` or :math:`\mathbf{P}_2`; these points are only there to provide
-directional information. The distance between :math:`\mathbf{P}_1` and :math:`\mathbf{P}_2`
-determines "how far" and "how fast" the curve moves towards :math:`\mathbf{P}_1` before turning
-towards :math:`\mathbf{P}_2`.
-
-Writing :math:`\mathbf{B}_{\mathbf P_i,\mathbf P_j,\mathbf P_k}(t)` for the quadratic Bézier curve
-defined by points :math:`\mathbf{P}_i`, :math:`\mathbf{P}_j`, and :math:`\mathbf{P}_k`, the cubic
-Bézier curve can be defined as an affine combination of two quadratic Bézier curves:
-
-.. math::
-    \mathbf{B}(t) = (1-t) \mathbf{B}_{\mathbf P_0,\mathbf P_1,\mathbf P_2}(t) +
-                        t \mathbf{B}_{\mathbf P_1,\mathbf P_2,\mathbf P_3}(t)
-
-The explicit form of the curve is:
-
-.. math::
-    \mathbf{B}(t) = (1-t)^3 \mathbf{P}_0 + 3(1-t)^2t \mathbf{P}_1 + 3(1-t)t^2 \mathbf{P}_2 + t^3\mathbf{P}_3
-
-For some choices of :math:`\mathbf{P}_1` and :math:`\mathbf{P}_2` the curve may intersect itself, or
-contain a cusp.
-
-The derivative of the cubic Bézier curve with respect to :math:`t` is
-
-.. math::
-    \mathbf{B}'(t) = 3(1-t)^2 (\mathbf{P}_1 - \mathbf{P}_0) + 6(1-t)t (\mathbf{P}_2 - \mathbf{P}_1) + 3t^2 (\mathbf{P}_3 - \mathbf{P}_2)
-
-The second derivative of the Bézier curve with respect to :math:`t` is
-
-.. math::
-    \mathbf{B}''(t) = 6(1-t) (\mathbf{P}_2 - 2 \mathbf{P}_1 + \mathbf{P}_0) +  6t (\mathbf{P}_3 - 2 \mathbf{P}_2 + \mathbf{P}_1)
-
-Recursive definition
---------------------
-
-A recursive definition for the Bézier curve of degree :math:`n` expresses it as a point-to-point
-linear combination of a pair of corresponding points in two Bézier curves of degree :math:`n-1`.
-
-Let :math:`\mathbf{B}_{\mathbf{P}_0\mathbf{P}_1\ldots\mathbf{P}_n}` denote the Bézier curve
-determined by any selection of points :math:`\mathbf{P}_0`, :math:`\mathbf{P}_1`, :math:`\ldots`,
-:math:`\mathbf{P}_{n-1}`.
-
-The recursive definition is
-
-.. math::
-    \begin{align}
-    \mathbf{B}_{\mathbf{P}_0}(t) &= \mathbf{P}_0 \\[1em]
-    \mathbf{B}(t) &= \mathbf{B}_{\mathbf{P}_0\mathbf{P}_1\ldots\mathbf{P}_n}(t) \\
-                  &= (1-t) \mathbf{B}_{\mathbf{P}_0\mathbf{P}_1\ldots\mathbf{P}_{n-1}}(t) +
-                         t \mathbf{B}_{\mathbf{P}_1\mathbf{P}_2\ldots\mathbf{P}_n}(t)
-    \end{align}
-
-The formula can be expressed explicitly as follows:
-
-.. math::
-    \begin{align}
-    \mathbf{B}(t) &= \sum_{i=0}^n b_{i,n}(t) \mathbf{P}_i \\
-                  &= \sum_{i=0}^n {n\choose i}(1 - t)^{n - i}t^i \mathbf{P}_i \\
-                  &= (1 - t)^n \mathbf{P}_0 +
-                     {n\choose 1}(1 - t)^{n - 1}t \mathbf{P}_1 +
-                     \cdots +
-                     {n\choose n - 1}(1 - t)t^{n - 1} \mathbf{P}_{n - 1} +
-                     t^n \mathbf{P}_n
-    \end{align}
-
-where :math:`b_{i,n}(t)` are the Bernstein basis polynomials of degree :math:`n` and :math:`n
-\choose i` are the binomial coefficients.
-
-Degree elevation
-----------------
-
-A Bézier curve of degree :math:`n` can be converted into a Bézier curve of degree :math:`n + 1` with
-the same shape.
-
-To do degree elevation, we use the equality
-
-.. math::
-       \mathbf{B}(t) = (1-t) \mathbf{B}(t) + t \mathbf{B}(t)`
-
-Each component :math:`\mathbf{b}_{i,n}(t) \mathbf{P}_i` is multiplied by :math:`(1-t)` and
-:math:`t`, thus increasing a degree by one, without changing the value.
-
-For arbitrary :math:`n`, we have
-
-.. math::
-    \begin{align}
-    \mathbf{B}(t) &= (1 - t) \sum_{i=0}^n \mathbf{b}_{i,n}(t) \mathbf{P}_i +
-                           t \sum_{i=0}^n \mathbf{b}_{i,n}(t) \mathbf{P}_i \\
-                  &= \sum_{i=0}^n \frac{n + 1 - i}{n + 1} \mathbf{b}_{i, n + 1}(t) \mathbf{P}_i +
-                     \sum_{i=0}^n \frac{i + 1}{n + 1} \mathbf{b}_{i + 1, n + 1}(t) \mathbf{P}_i \\
-                  &= \sum_{i=0}^{n + 1} \mathbf{b}_{i, n + 1}(t)
-                                        \left(\frac{i}{n + 1}         \mathbf{P}_{i - 1} +
-                                              \frac{n + 1 - i}{n + 1} \mathbf{P}_i\right) \\
-                  &= \sum_{i=0}^{n + 1} \mathbf{b}_{i, n + 1}(t) \mathbf{P'}_i
-    \end{align}
-
-Therefore the new control points are
-
-.. math::
-      \mathbf{P'}_i = \frac{i}{n + 1} \mathbf{P}_{i - 1} + \frac{n + 1 - i}{n + 1} \mathbf{P}_i
-
-It introduces two arbitrary points :math:`\mathbf{P}_{-1}` and :math:`\mathbf{P}_{n+1}` which are
-cancelled in :math:`\mathbf{P'}_i`.
+For resources on Bézier curve see :ref:`this section <bezier-geometry-ressources-page>`.
 
 """
 
+####################################################################################################
+#
+# Notes: algorithm details are on bezier.rst
+#
+####################################################################################################
+
 # Fixme:
 #   max distance to the chord for linear approximation
 #   fitting
 
-
 # C0 = continuous
 # G1 = geometric continuity
 #    Tangents point to the same direction
@@ -216,6 +51,8 @@ __all__ = [
 
 ####################################################################################################
 
+import logging
+
 from math import log, sqrt
 
 import numpy as np
@@ -229,12 +66,18 @@ from .Vector import Vector2D
 
 ####################################################################################################
 
+_module_logger = logging.getLogger(__name__)
+
+####################################################################################################
+
 class BezierMixin2D(Primitive2DMixin):
 
     """Mixin to implements 2D Bezier Curve."""
 
     LineInterpolationPrecision = 0.05
 
+    _logger = _module_logger.getChild('BezierMixin2D')
+
     ##############################################
 
     def interpolated_length(self, dt=None):
@@ -368,16 +211,6 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
 
     """Class to implements 2D Quadratic Bezier Curve."""
 
-    # Q(t) = Transformation * Control * Basis * T(t)
-    #
-    #           / P1x P2x P3x \  / 1 -2  1 \ / 1    \
-    # Q(t) = Tr | P1y P2x P3x |  | 0  2 -2 | | t    |
-    #           \  1   1   1  /  \ 0  0  1 / \ t**2 /
-    #
-    # Q(t) =    P0 * (1 - 2*t + t**2) +
-    #           P1 * (    2*t - t**2) +
-    #           P2 *            t**2
-
     BASIS = np.array((
         (1, -2,  1),
         (0,  2, -2),
@@ -390,6 +223,8 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
         (-1, -1, -2),
     ))
 
+    _logger = _module_logger.getChild('QuadraticBezier2D')
+
     ##############################################
 
     def __init__(self, p0, p1, p2):
@@ -405,28 +240,11 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
     @property
     def length(self):
 
-        # Algorithm:
-        #
-        # http://www.gamedev.net/topic/551455-length-of-a-generalized-quadratic-bezier-curve-in-3d
-        # Dave Eberly Posted October 25, 2009
-        #
-        # The quadratic Bezier is
-        #   (x(t),y(t)) = (1-t)^2*(x0,y0) + 2*t*(1-t)*(x1,y1) + t^2*(x2,y2)
-        #
-        # The derivative is
-        #   (x'(t),y'(t)) = -2*(1-t)*(x0,y0) + (2-4*t)*(x1,y1) + 2*t*(x2,y2)
-        #
-        # The length of the curve for 0 <= t <= 1 is
-        #   Integral[0,1] sqrt((x'(t))^2 + (y'(t))^2) dt
-        # The integrand is of the form sqrt(c*t^2 + b*t + a)
-        #
-        # You have three separate cases: c = 0, c > 0, or c < 0.
-        # * The case c = 0 is easy.
-        # * For the case c > 0, an antiderivative is
-        #     (2*c*t+b)*sqrt(c*t^2+b*t+a)/(4*c) + (0.5*k)*log(2*sqrt(c*(c*t^2+b*t+a)) + 2*c*t + b)/sqrt(c)
-        #   where k = 4*c/q with q = 4*a*c - b*b.
-        # * For the case c < 0, an antiderivative is
-        #    (2*c*t+b)*sqrt(c*t^2+b*t+a)/(4*c) - (0.5*k)*arcsin((2*c*t+b)/sqrt(-q))/sqrt(-c)
+        r"""Compute the length of the curve.
+
+        For more details see :ref:`this section <bezier-curve-length-section>`.
+
+        """
 
         A0 = self._p1 - self._p0
         A1 = self._p0 - self._p1 * 2 + self._p2
@@ -505,19 +323,10 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
 
         """Find the intersections of the curve with a line.
 
-        Algorithm
-
-        * Apply a transformation to the curve that maps the line onto the X-axis.
-        * Then we only need to test the Y-values for a zero.
+        For more details see :ref:`this section <bezier-curve-line-intersection-section>`.
 
         """
 
-        # t, p0, p1, p2, p3 = symbols('t p0 p1 p2 p3')
-        # u = 1 - t
-        # B = p0 * u**2  +  p1 * 2*t*u  +  p2 * t**2
-        # collect(expand(B), t)
-        # solveset(B, t)
-
         curve = self._map_to_line(line)
 
         p0 = curve.p0.y
@@ -566,48 +375,10 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
 
         """Return the closest point on the curve to the given *point*.
 
-        Reference
-
-        * https://hal.archives-ouvertes.fr/inria-00518379/document
-          Improved Algebraic Algorithm On Point Projection For Bézier Curves
-          Xiao-Diao Chen, Yin Zhou, Zhenyu Shu, Hua Su, Jean-Claude Paul
+        For more details see :ref:`this section <bezier-curve-closest-point-section>`.
 
         """
 
-        # Condition:
-        #   (P - B(t)) . B'(t) = 0   where t in [0,1]
-        #
-        #   P. B'(t) - B(t). B'(t) = 0
-
-        # A = P1 - P0
-        # B = P2 - P1 - A
-        # M = P0 - P
-
-        # Q(t)  = P0*(1-t)**2 + P1*2*t*(1-t) + P2*t**2
-        # Q'(t) = -2*P0*(1 - t) + 2*P1*(1 - 2*t) + 2*P2*t
-        #       = 2*(A + B*t)
-
-        # P0, P1, P2, P, t = symbols('P0 P1 P2 P t')
-        # Q = P0 * (1-t)**2  +  P1 * 2*t*(1-t)  +  P2 * t**2
-        # Qp = simplify(Q.diff(t))
-        # collect(expand((P*Qp - Q*Qp)/-2), t)
-
-        # (P0**2 - 4*P0*P1 + 2*P0*P2 + 4*P1**2 - 4*P1*P2 + P2**2) * t**3
-        # (-3*P0**2 + 9*P0*P1 - 3*P0*P2 - 6*P1**2 + 3*P1*P2) * t**2
-        # (-P*P0 + 2*P*P1 - P*P2 + 3*P0**2 - 6*P0*P1 + P0*P2 + 2*P1**2) * t
-        # P*P0 - P*P1 - P0**2 + P0*P1
-
-        # factorisation
-        # (P0 - 2*P1 + P2)**2 * t**3
-        # 3*(P1 - P0)*(P0 - 2*P1 + P2) * t**2
-        # ...
-        # (P0 - P)*(P1 - P0)
-
-        # B**2 * t**3
-        # 3*A*B * t**2
-        # (2*A**2 + M*B) * t
-        # M*A
-
         A = self._p1 - self._p0
         B = self._p2 - self._p1 - A
         M = self._p0 - point
@@ -622,7 +393,8 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
         if not t:
             return None
         elif len(t) > 1:
-            raise NameError("Found more than on root")
+            self._logger.warning("Found more than one root {} for {} and point {}".format(t, self, point))
+            return None
         else:
             return self.point_at_t(t)
 
@@ -632,22 +404,14 @@ class QuadraticBezier2D(BezierMixin2D, Primitive3P):
 
         r"""Elevate the quadratic Bézier curve to a cubic Bézier cubic with the same shape.
 
-        The new control points are
-
-        .. math::
-            \begin{align}
-            \mathbf{P'}_0 &= \mathbf{P}_0 \\
-            \mathbf{P'}_1 &= \mathbf{P}_0 + \frac{2}{3} (\mathbf{P}_1 - \mathbf{P}_0) \\
-            \mathbf{P'}_1 &= \mathbf{P}_2 + \frac{2}{3} (\mathbf{P}_1 - \mathbf{P}_2) \\
-            \mathbf{P'}_2 &= \mathbf{P}_2
-            \end{align}
+        For more details see :ref:`this section <bezier-curve-degree-elevation-section>`.
 
         """
 
         p1 = (self._p0 + self._p1 * 2) / 3
         p2 = (self._p2 + self._p1 * 2) / 3
 
-        return CubicBezier2D(self._p0, p1, p2, self._p3)
+        return CubicBezier2D(self._p0, p1, p2, self._p2)
 
 ####################################################################################################
 
@@ -662,13 +426,6 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
 
     InterpolationPrecision = 0.001
 
-    # Q(t) = Transformation * Control * Basis * T(t)
-    #
-    #           / P1x P2x P3x P4x \  / 1 -3  3 -1 \ / 1    \
-    # Q(t) = Tr | P1y P2x P3x P4x |  | 0  3 -6  3 | | t    |
-    #           |  0   0   0   0  |  | 0  0  3 -3 | | t**2 |
-    #           \  1   1   1   1  /  \ 0  0  0  1 / \ t**3 /
-
     BASIS = np.array((
         (1, -3,  3, -1),
         (0,  3, -6,  3),
@@ -683,6 +440,8 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
         (0,    0,   0, 1),
     ))
 
+    _logger = _module_logger.getChild('CubicMixin2D')
+
     #######################################
 
     def __init__(self, p0, p1, p2, p3):
@@ -698,7 +457,7 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
     def to_spline(self):
         from .Spline import CubicUniformSpline2D
         basis = np.dot(self.BASIS, CubicUniformSpline2D.INVERSE_BASIS)
-        points = np.dot(self.geometry_matrix, basis).transpose()
+        points = np.dot(self.point_array, basis).transpose()
         return CubicUniformSpline2D(*points)
 
     ##############################################
@@ -804,7 +563,6 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
 
         # Algorithm: same as for quadratic
 
-        # t, p0, p1, p2, p3, p4 = symbols('t p0 p1 p2 p3 p4')
         # u = 1 - t
         # B = p0 *     u**3        +
         #     p1 * 3 * u**2 * t    +
@@ -977,10 +735,15 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
             locations=[],
     ) :
 
-        # Code inspired from
-        #   https://github.com/paperjs/paper.js/blob/master/src/path/Curve.js
-        #   http://nbviewer.jupyter.org/gist/hkrish/0a128f21a5b9e5a7a914 The Bezier Clipping Algorithm
-        #   https://gist.github.com/hkrish/5ef0f2da7f9882341ee5 hkrish/bezclip_manual.py
+        """Compute the intersection of two Bézier curves.
+
+        Code inspired from
+
+          *  https://github.com/paperjs/paper.js/blob/master/src/path/Curve.js
+          *  http://nbviewer.jupyter.org/gist/hkrish/0a128f21a5b9e5a7a914 The Bezier Clipping Algorithm
+          *  https://gist.github.com/hkrish/5ef0f2da7f9882341ee5 hkrish/bezclip_manual.py
+
+        """
 
         # Note:
         #   see https://github.com/paperjs/paper.js/issues/565
@@ -1071,58 +834,25 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
 
     def is_flat_enough(self, flatness):
 
-         """Determines if a curve is sufficiently flat, meaning it appears as a straight line and has
+         r"""Determines if a curve is sufficiently flat, meaning it appears as a straight line and has
          curve-time that is enough linear, as specified by the given *flatness* parameter.
 
-         *flatness* is the maximum error allowed for the straight line to deviate from the curve.
+         For more details see :ref:`this section <bezier-curve-flatness-section>`.
+
+         """
 
-         Reference
+         u = 3*self._p1 - 2*self._p0 - self._p3
+         v = 3*self._p2 - 2*self._p3 - self._p0
 
-         * Kaspar Fischer and Roger Willcocks  http://hcklbrrfnn.files.wordpress.com/2012/08/bez.pdf
-         * PostScript Language Reference. Addison- Wesley, third edition, 1999
+         criterion = max(u.x**2, v.x**2) + max(u.y**2, v.y**2)
+         threshold = 16 * flatness**2
 
-         """
+         self._logger.warning("is flat {} <= {} with flatness {}".format(criterion, threshold, flatness))
 
-         # We define the flatness of the curve as the argmax of the distance from the curve to the
-         # line passing by the start and stop point.
-         #
-         # flatness = argmax(d(t)) for t in [0, 1] where d(t) = | B(t) - L(t) |
-         #
-         # L = (1-t)*P0 + t*P1
-         #
-         # Let
-         #   u = 3*P1 - 2*P0 - P3
-         #   v = 3*P2 - P0 - 2*P3
-         #
-         # d(t) = (1-t)**2 * t * (3*P1 - 2*P0 - P3)  +  (1-t) * t**2 * (3*P2 - P0 - 2*P3)
-         #      = (1-t)**2 * t * u  +  (1-t) * t**2 * v
-         #
-         # d(t)**2 = (1 - t)**2 * t**2 * (((1 - t)*ux + t*vx)**2 + ((1 - t)*uy + t*vy)**2
-         #
-         # argmax((1 - t)**2 * t**2) = 1/16
-         # argmax((1 - t)*a + t*b)   = argmax(a, b)
-         #
-         # flatness**2 = argmax(d(t)**2) <= 1/16 * (argmax(ux**2, vx**2) + argmax(uy**2, vy**2))
-         #
-         # argmax(ux**2, vx**2) + argmax(uy**2, vy**2) is thus an upper bound of 16 * flatness**2
-
-         # x0, y0 = list(self._p0)
-         # x1, y1 = list(self._p1)
-         # x2, y2 = list(self._p2)
-         # x3, y3 = list(self._p3)
-
-         # ux = 3*x1 - 2*x0 - x3
-         # uy = 3*y1 - 2*y0 - y3
-         # vx = 3*x2 - 2*x3 - x0
-         # vy = 3*y2 - 2*y3 - y0
-
-         u = 3*P1 - 2*P0 - P3
-         v = 3*P2 - 2*P3 - P0
-
-         return max(u.x**2, v.x**2) + max(u.y**2, v.y**2) <= 16 * flatness**2
+         return criterion <= threshold
 
     ##############################################
-
+    
     @property
     def area(self):
 
@@ -1145,31 +875,11 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
 
     def closest_point(self, point):
 
-        # Q(t) = (P3 - 3*P2 + 3*P1 - P0) * t**3 +
-        #        3*(P2 - 2*P1 + P0) * t**2 +
-        #        3*(P1 - P0) * t  +
-        #        P0
-
-        # n = P3 - 3*P2 + 3*P1 - P0
-        # r = 3*(P2 - 2*P1 + P0
-        # s = 3*(P1 - P0)
-        # v = P0
-
-        # Q(t)  = n*t**3 + r*t**2 + s*t + v
-        # Q'(t) = 3*n*t**2 + 2*r*t + s
+        """Return the closest point on the curve to the given *point*.
 
-        # P0, P1, P2, P3, P, t = symbols('P0 P1 P2 P3 P t')
-        # n, r, s, v = symbols('n r s v')
-        # Q = n*t**3 + r*t**2 + s*t + v
-        # Qp = simplify(Q.diff(t))
-        # collect(expand((P*Qp - Q*Qp)), t)
+        For more details see :ref:`this section <bezier-curve-closest-point-section>`.
 
-        # -3*n**2 * t**5
-        # -5*n*r * t**4
-        # -2*(2*n*s + r**2) * t**3
-        # 3*(P*n - n*v - r*s) * t**2
-        # (2*P*r - 2*r*v - s**2) * t
-        # P*s - s*v
+        """
 
         n = self._p3 - self._p2*3 + self._p1*3 - self._p0
         r = (self._p2 - self._p1*2 + self._p0)*3
@@ -1189,6 +899,18 @@ class CubicBezier2D(BezierMixin2D, Primitive4P):
         if not t:
             return None
         elif len(t) > 1:
-            raise NameError("Found more than on root")
+            # Fixme:
+            # Found more than one root [0, 0.516373783749732]
+            # for CubicBezier2D(
+            #   Vector2D[1394.4334 1672.0004], Vector2D[1394.4334 1672.0004],
+            #   Vector2D[1585.0004 1624.9634], Vector2D[1585.0004 1622.0004])
+            # and point Vector2D[1495.11502887 1649.7386517 ]
+            # raise NameError("Found more than one root: {}".format(t))
+            self._logger.warning("Found more than one root {} for {} and point {}".format(t, self, point))
+            # self._logger.warning("is flat {}".format(self.is_flat_enough(.1)))
+            if len(t) == 2 and t[0] == 0:
+                return self.point_at_t(t[1])
+            else:
+                return None
         else:
             return self.point_at_t(t[0])

Patro/GeometryEngine/Conic.py

@@ -20,16 +20,16 @@
 
 """Module to implement conic geometry like circle and ellipse.
 
-Valentina Requirements
+*Valentina Requirements*
 
-* circle with angular domain
-* circle with start angle and arc length
-* curvilinear distance on circle
-* line-circle intersection
-* circle-circle intersection
-* point constructed from a virtual circle and a point on a tangent : right triangle
-* point from tangent circle and segment ???
-* ellipse with angular domain and rotation
+  * circle with angular domain
+  * circle with start angle and arc length
+  * curvilinear distance on circle
+  * line-circle intersection
+  * circle-circle intersection
+  * point constructed from a virtual circle and a point on a tangent : right triangle
+  * point from tangent circle and segment ???
+  * ellipse with angular domain and rotation
 
 """
 
@@ -48,6 +48,8 @@ __all__ = [
 
 ####################################################################################################
 
+import logging
+
 import math
 from math import fabs, sqrt, radians, pi, cos, sin # , degrees
 
@@ -59,6 +61,11 @@ from .Line import Line2D
 from .Mixin import AngularDomainMixin, CenterMixin, AngularDomain
 from .Primitive import Primitive, Primitive2DMixin
 from .Segment import Segment2D
+from .Transformation import Transformation2D
+
+####################################################################################################
+
+_module_logger = logging.getLogger(__name__)
 
 ####################################################################################################
 
@@ -436,35 +443,144 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
     \ge e_1 > 0`. The ellipse points are
 
     .. math::
-         P = C + x_0 U_0 + x_1 U_1
+        \begin{equation}
+        P = C + x_0 U_0 + x_1 U_1
+        \end{equation}
 
     where
 
     .. math::
-
+        \begin{equation}
         \left(\frac{x_0}{e_0}\right)^2 + \left(\frac{x_1}{e_1}\right)^2 = 1
+        \end{equation}
+
+    If :math:`e_0 = e_1`, then the ellipse is a circle with center `C` and radius :math:`e_0`.
 
-    If :math:`e_0 = e_1`, then the ellipse is a circle with center `C` and radius :math:`e_0`. The
-    orthonormality of the axis directions and Equation (1) imply :math:`x_i = U_i \dot (P −
+    The orthonormality of the axis directions and Equation (1) imply :math:`x_i = U_i \dot (P −
     C)`. Substituting this into Equation (2) we obtain
 
     .. math::
-
         (P − C)^T M (P − C) = 1
 
     where :math:`M = R D R^T`, `R` is an orthogonal matrix whose columns are :math:`U_0` and
     :math:`U_1` , and `D` is a diagonal matrix whose diagonal entries are :math:`1/e_0^2` and
     :math:`1/e_1^2`.
 
+    An ellipse can also be parameterised by an angle :math:`\theta`
+
+    .. math::
+        \begin{pmatrix} x \\ y \end{pmatrix} =
+        \begin{bmatrix}
+        \cos\phi & \sin\phi \\
+        -\sin\phi & \cos\phi
+        \end{bmatrix}
+        \begin{pmatrix} r_x \cos\theta \\ r_y \sin\theta \end{pmatrix}
+        + \begin{pmatrix} C_x \\ C_y \end{pmatrix}
+
+    where :math:`\phi` is the angle from the x-axis, :math:`r_x` is the semi-major and :math:`r_y`
+    semi-minor axes.
+
     """
 
+    _logger = _module_logger.getChild('Ellipse2D')
+
+    ##############################################
+
+    @classmethod
+    def svg_arc(cls, point1, point2, radius_x, radius_y, angle, large_arc, sweep):
+
+        """Implement SVG Arc.
+
+        Parameters
+
+          * *point1* is the start point and *point2* is the end point.
+          * *radius_x* and *radius_y* are the radii of the ellipse, also known as its semi-major and
+            semi-minor axes.
+          * *angle* is the angle from the x-axis of the current coordinate system to the x-axis of the ellipse.
+          * if the *large arc* flag is unset then arc spanning less than or equal to 180 degrees is
+            chosen, else an arc spanning greater than 180 degrees is chosen.
+          * if the *sweep* flag is unset then the line joining centre to arc sweeps through decreasing
+            angles, else if it sweeps through increasing angles.
+
+        References
+
+          * https://www.w3.org/TR/SVG/implnote.html#ArcConversionEndpointToCenter
+          * https://www.w3.org/TR/SVG/implnote.html#ArcCorrectionOutOfRangeRadii
+
+        """
+
+        # Ensure radii are non-zero
+        if radius_x == 0 or radius_y == 0:
+            return Segment2D(point1, point2)
+
+        # Ensure radii are positive
+        radius_x = abs(radius_x)
+        radius_y = abs(radius_y)
+
+        # step 1
+
+        radius_x2 = radius_x**2
+        radius_y2 = radius_y**2
+
+        # We define a new referential with the origin is set to the middle of P1 — P2
+        origin_prime = (point1 + point2)/2
+
+        # P1 is exprimed in this referential where the ellipse major axis line up with the x axis
+        point1_prime = Transformation2D.Rotation(-angle) * (point1 - point2)/2
+
+        # Ensure radii are large enough
+        radii_scale = point1_prime.x**2/radius_x2 + point1_prime.y**2/radius_y2
+        if radii_scale > 1:
+            self._logger.warning('SVG Arc: radii must be scale')
+            radii_scale = math.sqrt(radii_scale)
+            radius_x = radii_scale * radius_x
+            radius_y = radii_scale * radius_y
+            radius_x2 = radius_x**2
+            radius_y2 = radius_y**2
+
+        # step 2
+
+        den = radius_x2 * point1_prime.y**2 + radius_y2 * point1_prime.x**2
+        num = radius_x2*radius_y2 - den
+
+        ratio = radius_x/radius_y
+
+        sign = 1 if large_arc != sweep else -1
+        # print(point1_prime)
+        # print(point1_prime.anti_normal)
+        # print(ratio)
+        # print(point1_prime.anti_normal.scale(ratio, 1/ratio))
+        sign *= -1 # Fixme: solve mirroring artefacts for y-axis pointing to the top
+        center_prime = sign * math.sqrt(num / den) * point1_prime.anti_normal.scale(ratio, 1/ratio)
+
+        center = Transformation2D.Rotation(angle) * center_prime + origin_prime
+
+        vector1 =   (point1_prime - center_prime).divide(radius_x, radius_y)
+        vector2 = - (point1_prime + center_prime).divide(radius_x, radius_y)
+        theta = cls.__vector_cls__(1, 0).angle_with(vector1)
+        delta_theta = vector1.angle_with(vector2)
+        # if theta < 0:
+        #     theta = 180 + theta
+        # if delta_theta < 0:
+        #     delta_theta = 180 + delta_theta
+        delta_theta = delta_theta % 360
+        # print('theta', theta, delta_theta)
+        if not sweep and delta_theta > 0:
+            delta_theta -= 360
+        elif sweep and delta_theta < 0:
+            delta_theta += 360
+        # print('theta', theta, delta_theta, theta + delta_theta)
+        domain = domain = AngularDomain(theta, theta + delta_theta)
+
+        return cls(center, radius_x, radius_y, angle, domain)
+
     #######################################
 
-    def __init__(self, center, x_radius, y_radius, angle, domain=None):
+    def __init__(self, center, radius_x, radius_y, angle, domain=None):
 
         self.center = center
-        self.x_radius = x_radius
-        self.y_radius = y_radius
+        self.radius_x = radius_x
+        self.radius_y = radius_y
         self.angle = angle
         self.domain = domain
 
@@ -475,7 +591,7 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
     def clone(self):
         return self.__class__(
             self._center,
-            self._x_radius, self._y_radius,
+            self._radius_x, self._radius_y,
             self._angle,
             self._domain,
         )
@@ -484,32 +600,32 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
 
     def apply_transformation(self, transformation):
         self._center = transformation * self._center
-        self._x_radius = transformation * self._x_radius
-        self._y_radius = transformation * self._y_radius
+        self._radius_x = transformation * self._radius_x
+        self._radius_y = transformation * self._radius_y
         self._bounding_box = None
 
     ##############################################
 
     def __repr__(self):
-        return '{0}({1._center}, {1._x_radius}, {1._x_radius}, {1._angle})'.format(self.__class__.__name__, self)
+        return '{0}({1._center}, {1._radius_x}, {1._radius_x}, {1._angle})'.format(self.__class__.__name__, self)
 
     ##############################################
 
     @property
-    def x_radius(self):
-        return self._x_radius
+    def radius_x(self):
+        return self._radius_x
 
-    @x_radius.setter
-    def x_radius(self, value):
-        self._x_radius = float(value)
+    @radius_x.setter
+    def radius_x(self, value):
+        self._radius_x = float(value)
 
     @property
-    def y_radius(self):
-        return self._y_radius
+    def radius_y(self):
+        return self._radius_y
 
-    @y_radius.setter
-    def y_radius(self, value):
-        self._y_radius = float(value)
+    @radius_y.setter
+    def radius_y(self, value):
+        self._radius_y = float(value)
 
     @property
     def angle(self):
@@ -522,21 +638,21 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
     @property
     def major_vector(self):
         # Fixme: x < y
-        return self.__vector_cls__.from_polar(self._angle, self._x_radius)
+        return self.__vector_cls__.from_polar(self._angle, self._radius_x)
 
     @property
     def minor_vector(self):
         # Fixme: x < y
-        return self.__vector_cls__.from_polar(self._angle + 90, self._y_radius)
+        return self.__vector_cls__.from_polar(self._angle + 90, self._radius_y)
 
     ##############################################
 
     @property
     def eccentricity(self):
         # focal distance
-        # c = sqrt(self._x_radius**2 - self._y_radius**2)
+        # c = sqrt(self._radius_x**2 - self._radius_y**2)
         # e = c / a
-        return sqrt(1 - (self._y_radius/self._x_radius)**2)
+        return sqrt(1 - (self._radius_y/self._radius_x)**2)
 
     ##############################################
 
@@ -550,8 +666,8 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
         c2 = c**2
         s2 = s**2
 
-        a = self._x_radius
-        b = self._y_radius
+        a = self._radius_x
+        b = self._radius_y
         a2 = a**2
         b2 = b**2
 
@@ -584,9 +700,9 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
     ##############################################
 
     def point_at_angle(self, angle):
-        # point = self.__vector_cls__.from_ellipse(self._x_radius, self._y_radius, angle)
+        # point = self.__vector_cls__.from_ellipse(self._radius_x, self._radius_y, angle)
         # return self.point_from_ellipse_frame(point)
-        point = self.__vector_cls__.from_ellipse(self._x_radius, self._y_radius, self._angle + angle)
+        point = self.__vector_cls__.from_ellipse(self._radius_x, self._radius_y, self._angle + angle)
         return self._center + point
 
     ##############################################
@@ -595,21 +711,21 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
     def bounding_box(self):
 
         if self._bounding_box is None:
-            x_radius, y_radius = self._x_radius, self._y_radius
+            radius_x, radius_y = self._radius_x, self._radius_y
             if self._angle == 0:
                 bounding_box = self._center.bounding_box
-                bounding_box.x.enlarge(x_radius)
-                bounding_box.y.enlarge(y_radius)
+                bounding_box.x.enlarge(radius_x)
+                bounding_box.y.enlarge(radius_y)
                 self._bounding_box = bounding_box
             else:
                 angle_x = self._angle
                 angle_y = angle_x + 90
                 Vector2D = self.__vector_cls__
                 points = [self._center + offset for offset in (
-                    Vector2D.from_polar(angle_x,  x_radius),
-                    Vector2D.from_polar(angle_x, -x_radius),
-                    Vector2D.from_polar(angle_y,  y_radius),
-                    Vector2D.from_polar(angle_y, -y_radius),
+                    Vector2D.from_polar(angle_x,  radius_x),
+                    Vector2D.from_polar(angle_x, -radius_x),
+                    Vector2D.from_polar(angle_y,  radius_y),
+                    Vector2D.from_polar(angle_y, -radius_y),
                 )]
                 self._bounding_box = bounding_box_from_points(points)
 
@@ -672,7 +788,7 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
         # Fixme: make a 3D plot to check the algorithm on a 2D grid and rotated ellipse
 
         y0, y1 = point
-        e0, e1 = self._x_radius, self._y_radius
+        e0, e1 = self._radius_x, self._radius_y
 
         if y1 > 0:
             if  y0 > 0:
@@ -749,11 +865,11 @@ class Ellipse2D(Primitive2DMixin, CenterMixin, AngularDomainMixin, Primitive):
         # Fixme: to be checked
 
         # Map segment in ellipse frame and scale y axis so as to transform the ellipse to a circle
-        y_scale = self._x_radius / self._y_radius
+        y_scale = self._radius_x / self._radius_y
         points = [self.point_in_ellipse_frame(point) for point in segment.points]
         points = [self.__vector_cls__(point.x, point.y * y_scale) for point in points]
         segment_in_frame = Segment2D(*points)
-        circle = Circle2D(self.__vector_cls__(0, 0), self._x_radius)
+        circle = Circle2D(self.__vector_cls__(0, 0), self._radius_x)
 
         points = circle.intersect_segment(segment_in_frame)
         points = [self.__vector_cls__(point.x, point.y / y_scale) for point in points]

Patro/GeometryEngine/Line.py

@@ -111,21 +111,21 @@ class Line2D(Primitive2DMixin, Primitive):
 
     def get_y_from_x(self, x):
         """Return y corresponding to x"""
-        return self.v.tan() * (x - self.p.x) + self.p.y
+        return self.v.tan * (x - self.p.x) + self.p.y
 
     ##############################################
 
     def get_x_from_y(self, y):
         """Return x corresponding to y"""
-        return self.v.inverse_tan() * (y - self.p.y) + self.p.x
+        return self.v.inverse_tan * (y - self.p.y) + self.p.x
 
     ##############################################
 
     # Fixme: is_parallel_to
 
-    def is_parallel(self, other, cross=False):
+    def is_parallel(self, other, return_cross=False):
         """Self is parallel to other"""
-        return self.v.is_parallel(other.v, cross)
+        return self.v.is_parallel(other.v, return_cross)
 
     ##############################################
 
@@ -139,7 +139,7 @@ class Line2D(Primitive2DMixin, Primitive):
 
         """Return the shifted parallel line"""
 
-        n = self.v.normal()
+        n = self.v.normal
         n.normalise()
         point = self.p + n*shift
 
@@ -152,7 +152,7 @@ class Line2D(Primitive2DMixin, Primitive):
         """Return the orthogonal line at abscissa s"""
 
         point = self.interpolate(s)
-        vector = self.v.normal()
+        vector = self.v.normal
 
         return self.__class__(point, vector)
 
@@ -170,7 +170,7 @@ class Line2D(Primitive2DMixin, Primitive):
         # delta x v1 = - s2 * v2 x v1 = s2 * v1 x v2
         # delta x v2                  = s1 * v1 x v2
 
-        test, cross = l1.is_parallel(l2, cross=True)
+        test, cross = l1.is_parallel(l2, return_cross=True)
         if test:
             return (None, None)
         else:
@@ -255,7 +255,7 @@ class Line2D(Primitive2DMixin, Primitive):
 
         left, bottom, right, top = interval.bounding_box()
         vb = Vector2D(interval.size())
-        if abs(self.v.tan()) > vb.tan():
+        if abs(self.v.tan) > vb.tan:
             x_min, y_min = self.get_x_from_y(bottom), bottom
             x_max, y_max = self.get_x_from_y(top), top
         else:

Patro/GeometryEngine/Path.py

@@ -20,6 +20,8 @@
 
 """Module to implement path.
 
+For resources on path see :ref:`this section <path-geometry-ressources-page>`.
+
 """
 
 ####################################################################################################
@@ -33,24 +35,35 @@ __all__ = [
 
 ####################################################################################################
 
+import logging
 import math
 
+from Patro.Common.Math.Functions import sign
 from .Primitive import Primitive1P, Primitive2DMixin
 from .Bezier import QuadraticBezier2D, CubicBezier2D
-from .Conic import Circle2D, AngularDomain
+from .Conic import AngularDomain, Circle2D, Ellipse2D
 from .Segment import Segment2D
 from .Vector import Vector2D
 
 ####################################################################################################
 
+_module_logger = logging.getLogger(__name__)
+
+####################################################################################################
+
 class PathPart:
 
     ##############################################
 
-    def __init__(self, path, position):
+    def __init__(self, path, index):
 
         self._path = path
-        self._position = position
+        self._index = index
+
+    ##############################################
+
+    def _init_absolute(self, absolute):
+        self._absolute = bool(absolute)
 
     ##############################################
 
@@ -60,7 +73,7 @@ class PathPart:
     ##############################################
 
     def __repr__(self):
-        return self.__class__.__name__
+        return '{0}(@{1._index})'.format(self.__class__.__name__, self)
 
     ##############################################
 
@@ -69,22 +82,22 @@ class PathPart:
         return self._path
 
     @property
-    def position(self):
-        return self._position
+    def index(self):
+        return self._index
 
-    @position.setter
-    def position(self, value):
-        self._position = int(value)
+    @index.setter
+    def index(self, value):
+        self._index = int(value)
 
     ##############################################
 
     @property
     def prev_part(self):
-        return self._path[self._position -1]
+        return self._path[self._index -1]
 
     @property
     def next_part(self):
-        return self._path[self._position +1]
+        return self._path[self._index +1]
 
     ##############################################
 
@@ -103,6 +116,17 @@ class PathPart:
     def stop_point(self):
         raise NotImplementedError
 
+    ##############################################
+
+    def to_absolute_point(self, point):
+        # Fixme: cache ???
+        if self._absolute:
+            return point
+        else:
+            return point + self.start_point
+
+    ##############################################
+
     @property
     def geometry(self):
         raise NotImplementedError
@@ -115,36 +139,87 @@ class PathPart:
 
 ####################################################################################################
 
-class LinearSegment(PathPart):
+class OnePointMixin:
+
+    ##############################################
+
+    @property
+    def point(self):
+        return self._point
+
+    @point.setter
+    def point(self, value):
+        self._point = Vector2D(value) # self._path.__vector_cls__
+
+    ##############################################
+
+    @property
+    def stop_point(self):
+        return self.to_absolute_point(self._point)
+
+    ##############################################
+
+    def apply_transformation(self, transformation):
+        # Fixme: right for relative ???
+        self._point = transformation * self._point
+
+####################################################################################################
+
+class TwoPointMixin:
+
+    ##############################################
+
+    @property
+    def point1(self):
+        return self.to_absolute_point(self._point1)
 
-    r"""
+    @point1.setter
+    def point1(self, value):
+        self._point1 = Vector2D(value) # self._path.__vector_cls__
 
-    Bulge
+    ##############################################
 
-    Let `P0`, `P1`, `P2` the vertices and `R` the bulge radius.
+    @property
+    def point2(self):
+        return self.to_absolute_point(self._point2)
 
-    The deflection :math:`\theta = 2 \alpha` at the corner is
+    @point2.setter
+    def point2(self, value):
+        self._point2 = Vector2D(value)
 
-    .. math::
+    ##############################################
 
-       D_1 \cdot D_0 = (P_2 - P_1) \cdot (P_1 - P_0) = \cos \theta
+    def apply_transformation(self, transformation):
+        # Fixme: right for relative ???
+        self._point1 = transformation * self._point1
+        self._point2 = transformation * self._point2
 
-    The bisector direction is
+####################################################################################################
 
-    .. math::
+class ThreePointMixin(TwoPointMixin):
 
-       Bis = D_1 - D_0 = (P_2 - P_1) - (P_1 - P_0) = P_2 -2 P_1 + P_0
+    ##############################################
 
-    Bulge Center is
+    @property
+    def point3(self):
+        return self.to_absolute_point(self._point3)
 
-    .. math::
+    @point3.setter
+    def point3(self, value):
+        self._point3 = Vector2D(value) # self._path.__vector_cls__
 
-        C = P_1 + Bis \times \frac{R}{\sin \alpha}
+    ##############################################
+
+    def apply_transformation(self, transformation):
+        # Fixme: right for relative ???
+        TwoPointMixin.apply_transformation(self, transformation)
+        self._point3 = transformation * self._point3
 
-    Extremities are
+####################################################################################################
 
-        \prime P_1 = P_1 - d_0 \times \frac{R}{\tan \alpha}
-        \prime P_1 = P_1 + d_1 \times \frac{R}{\tan \alpha}
+class LinearSegment(PathPart):
+
+    """Class to implement a linear segment.
 
     """
 
@@ -153,15 +228,19 @@ class LinearSegment(PathPart):
     #    If two successive vertices share the same circle, then it should be merged to one.
     #
 
+    _logger = _module_logger.getChild('LinearSegment')
+
     ##############################################
 
-    def __init__(self, path, position, radius):
+    def __init__(self, path, index, radius, closing=False):
 
-        super().__init__(path, position)
+        super().__init__(path, index)
 
         self._bissector = None
         self._direction = None
 
+        self._start_bulge = False
+        self._closing = bool(closing)
         self.radius = radius
         if self._radius is not None:
             if not isinstance(self.prev_part, LinearSegment):
@@ -179,6 +258,31 @@ class LinearSegment(PathPart):
 
     ##############################################
 
+    def close(self, radius):
+        """Set the bulge radius at the closure"""
+        self.radius = radius
+        self._reset_cache()
+        self._start_bulge = True
+
+    ##############################################
+
+    @property
+    def prev_part(self):
+        if self._start_bulge:
+            return self._path.stop_segment # or [-1] don't work
+        else:
+            # Fixme: super
+            return self._path[self._index -1]
+
+    @property
+    def next_part(self):
+        if self._closing:
+            return self._path.start_segment # or [0]
+        else:
+            return self._path[self._index +1]
+
+    ##############################################
+
     @property
     def points(self):
 
@@ -204,9 +308,10 @@ class LinearSegment(PathPart):
     @radius.setter
     def radius(self, value):
         if value is not None:
-            self._radius = float(value)
-        else:
-            self._radius = None
+            value = abs(float(value))
+            if value == 0:
+                radius = None
+        self._radius = value
 
     ##############################################
 
@@ -231,7 +336,12 @@ class LinearSegment(PathPart):
     @property
     def bulge_angle_rad(self):
         if self._bulge_angle is None:
-            angle = self.direction.orientation_with(self.prev_part.direction)
+            # Fixme: rad vs degree
+            angle = self.direction.angle_with(self.prev_part.direction)
+            if angle >= 0:
+                angle = 180 - angle
+            else:
+                angle = -(180 + angle)
             self._bulge_angle = math.radians(angle)
         return self._bulge_angle
 
@@ -241,15 +351,16 @@ class LinearSegment(PathPart):
 
     @property
     def half_bulge_angle(self):
-        return abs(self.bulge_angle_rad / 2)
+        return self.bulge_angle_rad / 2
 
     ##############################################
 
     @property
     def bulge_center(self):
         if self._bulge_center is None:
-            offset = self.bissector * self._radius / math.sin(self.half_bulge_angle)
+            offset = self.bissector * self._radius / math.sin(abs(self.half_bulge_angle))
             self._bulge_center = self.start_point + offset
+            # Note: -offset create external loop
         return self._bulge_center
 
     ##############################################
@@ -257,15 +368,18 @@ class LinearSegment(PathPart):
     @property
     def bulge_start_point(self):
         if self._start_bulge_point is None:
-            offset = self.prev_part.direction * self._radius / math.tan(self.half_bulge_angle)
-            self._start_bulge_point = self.start_point - offset
+            angle = self.half_bulge_angle
+            offset = self.prev_part.direction * self._radius / math.tan(angle)
+            self._start_bulge_point = self.start_point - sign(angle) *offset
+            # Note: -offset create internal loop
         return self._start_bulge_point
 
     @property
     def bulge_stop_point(self):
         if self._stop_bulge_point is None:
+            angle = self.half_bulge_angle
             offset = self.direction * self._radius / math.tan(self.half_bulge_angle)
-            self._stop_bulge_point = self.start_point + offset
+            self._stop_bulge_point = self.start_point + sign(angle) * offset
         return self._stop_bulge_point
 
     ##############################################
@@ -279,80 +393,201 @@ class LinearSegment(PathPart):
         if self.bulge_angle < 0:
             start_angle, stop_angle = stop_angle, start_angle
         arc.domain = AngularDomain(start_angle, stop_angle)
+        # self._dump_bulge(arc)
         return arc
 
+    ##############################################
+
+    def _dump_bulge(self, arc):
+        self._logger.info(
+            'Bulge @{}\n'.format(self._index) +
+            str(arc)
+        )
+
 ####################################################################################################
 
-class PathSegment(LinearSegment):
+class PathSegment(OnePointMixin, LinearSegment):
 
     ##############################################
 
-    def __init__(self, path, position, point, radius=None, absolute=False):
-        super().__init__(path, position, radius)
+    def __init__(self, path, index, point, radius=None, absolute=False, closing=False):
+        super().__init__(path, index, radius, closing)
         self.point = point
-        self._absolute = bool(absolute)
+        self._init_absolute(absolute)
 
     ##############################################
 
     def clone(self, path):
-        return self.__class__(path, self._position, self._point, self._radius, self._absolute)
+
+        # Fixme: check
+        if obj._start_bulge:
+            radius = None
+        else:
+            radius = self._radius
+        obj = self.__class__(path, self._index, self._point, radius, self._absolute, self._closing)
+        if obj._start_bulge:
+            self.close(self._radius)
+
+        return obj
+
+    ##############################################
+
+    def to_absolute(self):
+        self._point = self.stop_point
+        self._absolute = True
 
     ##############################################
 
     def apply_transformation(self, transformation):
-        self._point = transformation * self._point
+        OnePointMixin.apply_transformation(self, transformation)
         if self._radius is not None:
             self._radius = transformation * self._radius
 
     ##############################################
 
     @property
-    def point(self):
-        return self._point
+    def geometry(self):
+        # Fixme: cache ???
+        return Segment2D(*self.points)
 
-    @point.setter
-    def point(self, value):
-        self._point = Vector2D(value) # self._path.__vector_cls__
+####################################################################################################
+
+class DirectionalSegmentMixin(LinearSegment):
 
     ##############################################
 
-    @property
-    def stop_point(self):
-        if self._absolute:
-            return self._point
-        else:
-            return self._point + self.start_point
+    def apply_transformation(self, transformation):
+        # Since a rotation will change the direction
+        # DirectionalSegment must be casted to PathSegment
+        raise NotImplementedError
 
     ##############################################
 
     @property
     def geometry(self):
         # Fixme: cache ???
-        return Segment2D(*self.points)
+        return Segment2D(self.start_point, self.stop_point)
+
+####################################################################################################
+
+class AbsoluteHVSegment(DirectionalSegmentMixin):
+
+    ##############################################
+
+    def to_path_segment(self):
+
+        # Fixme: duplicted
+        if self._index == 0 and self._radius is not None:
+            radius = None
+            close = True
+        else:
+            radius = self._radius
+            close = False
+
+        path = PathSegment(self._path, self._index, self.stop_point, radius, absolute=True)
+
+        if close:
+            path.close(self._radius)
+
+        return path
 
 ####################################################################################################
 
-class DirectionalSegment(LinearSegment):
+class AbsoluteHorizontalSegment(AbsoluteHVSegment):
+
+    ##############################################
+
+    def __init__(self, path, index, x, radius=None):
+
+        super().__init__(path, index, radius)
+        self.x = x
+        self._init_absolute(False) # Fixme: mix
+
+    ##############################################
+
+    def __repr__(self):
+        return '{0}(@{1._index}, {1._x})'.format(self.__class__.__name__, self)
+
+    ##############################################
+
+    def clone(self, path):
+        return self.__class__(path, self._index, self._x, self._radius)
+
+    ##############################################
+
+    @property
+    def x(self):
+        return self._x
+
+    @x.setter
+    def x(self, value):
+        self._x = float(value)
+
+    ##############################################
+
+    @property
+    def stop_point(self):
+        return Vector2D(self._x, self.start_point.y)
+
+####################################################################################################
+
+class AbsoluteVerticalSegment(AbsoluteHVSegment):
+
+    ##############################################
+
+    def __init__(self, path, index, y, radius=None):
+
+        super().__init__(path, index, radius)
+        self.y = y
+        self._init_absolute(False) # Fixme: mix
+
+    ##############################################
+
+    def __repr__(self):
+        return '{0}(@{1._index}, {1._y})'.format(self.__class__.__name__, self)
+
+    ##############################################
+
+    def clone(self, path):
+        return self.__class__(path, self._index, self._y, self._radius)
+
+    ##############################################
+
+    @property
+    def y(self):
+        return self._y
+
+    @y.setter
+    def y(self, value):
+        self._y = float(value)
+
+    ##############################################
+
+    @property
+    def stop_point(self):
+        return Vector2D(self.start_point.x, self._y)
+
+####################################################################################################
+
+class DirectionalSegment(DirectionalSegmentMixin):
 
     __angle__ = None
 
     ##############################################
 
-    def __init__(self, path, position, length, radius=None):
-        super().__init__(path, position, radius)
+    def __init__(self, path, index, length, radius=None):
+        super().__init__(path, index, radius)
         self.length = length
 
     ##############################################
 
-    def apply_transformation(self, transformation):
-        # Since a rotation will change the direction
-        # DirectionalSegment must be casted to PathSegment
-        raise NotImplementedError
+    def __repr__(self):
+        return '{0}(@{1._index}, {1.offset})'.format(self.__class__.__name__, self)
 
     ##############################################
 
     def clone(self, path):
-        return self.__class__(path, self._position, self._length, self._radius)
+        return self.__class__(path, self._index, self._length, self._radius)
 
     ##############################################
 
@@ -378,15 +613,21 @@ class DirectionalSegment(LinearSegment):
 
     ##############################################
 
-    @property
-    def geometry(self):
-        # Fixme: cache ???
-        return Segment2D(self.start_point, self.stop_point)
+    def to_path_segment(self):
 
-    ##############################################
+        if self._index == 0 and self._radius is not None:
+            radius = None
+            close = True
+        else:
+            radius = self._radius
+            close = False
 
-    def to_path_segment(self):
-        return PathSegment(self._path, self._position, self.offset, self._radius, absolute=False)
+        path = PathSegment(self._path, self._index, self.offset, radius, absolute=False)
+
+        if close:
+            path.close(self._radius)
+
+        return path
 
 ####################################################################################################
 
@@ -422,118 +663,208 @@ class SouthWestSegment(DirectionalSegment):
 
 ####################################################################################################
 
-class TwoPointsMixin:
+class QuadraticBezierSegment(PathPart, TwoPointMixin):
 
-    @property
-    def point1(self):
-        return self._point1
+    # Fixme: abs / inc
 
-    @point1.setter
-    def point1(self, value):
-        self._point1 = Vector2D(value) # self._path.__vector_cls__
+    ##############################################
+
+    def __init__(self, path, index, point1, point2, absolute=False):
+
+        PathPart.__init__(self, path, index)
+        self._init_absolute(absolute)
+
+        self.point1 = point1
+        self.point2 = point2
+
+    ##############################################
+
+    def clone(self, path):
+        return self.__class__(path, self._index, self._point1, self._point2, self._absolute)
+
+    ##############################################
+
+    def to_absolute(self):
+        self._point1 = self.point1
+        self._point2 = self.point2
+        self._absolute = True
+
+    ##############################################
 
     @property
-    def point2(self):
-        return self._point2
+    def stop_point(self):
+        return self.point2
 
-    @point2.setter
-    def point2(self, value):
-        self._point2 = Vector2D(value)
+    @property
+    def points(self):
+        return (self.start_point, self.point1, self.point2)
 
     ##############################################
 
-    def apply_transformation(self, transformation):
-        self._point1 = transformation * self._point1
-        self._point2 = transformation * self._point2
+    @property
+    def geometry(self):
+        # Fixme: cache ???
+        return QuadraticBezier2D(*self.points)
 
 ####################################################################################################
 
-class QuadraticBezierSegment(PathPart, TwoPointsMixin):
-
-    # Fixme: abs / inc
+class CubicBezierSegment(PathPart, ThreePointMixin):
 
     ##############################################
 
-    def __init__(self, path, position, point1, point2):
+    def __init__(self, path, index, point1, point2, point3, absolute=False):
 
-        PathPart.__init__(self, path, position)
+        PathPart.__init__(self, path, index)
+        self._init_absolute(absolute)
 
         self.point1 = point1
         self.point2 = point2
+        self.point3 = point3
 
     ##############################################
 
     def clone(self, path):
-        return self.__class__(path, self._position, self._point1, self._point2)
+        return self.__class__(path, self._index, self._point1, self._point2, self._point3, absolute)
+
+    ##############################################
+
+    def to_absolute(self):
+        self._point1 = self.point1
+        self._point2 = self.point2
+        self._point3 = self.point3
+        self._absolute = True
 
     ##############################################
 
     @property
     def stop_point(self):
-        return self._point2
+        return self.point3
 
     @property
     def points(self):
-        return (self.start_point, self._point1, self._point2)
+        return (self.start_point, self.point1, self.point2, self.point3)
 
     ##############################################
 
     @property
     def geometry(self):
         # Fixme: cache ???
-        return QuadraticBezier2D(self.start_point, self._point1, self._point2)
+        return CubicBezier2D(*self.points)
 
 ####################################################################################################
 
-class CubicBezierSegment(PathPart, TwoPointsMixin):
+class StringedQuadtraticBezierSegment(PathPart, TwoPointMixin):
 
     ##############################################
 
-    def __init__(self, path, position, point1, point2, point3):
+    def __init__(self, path, index, point1, absolute=False):
 
-        PathPart.__init__(self, path, position)
+        PathPart.__init__(self, path, index)
+        self._init_absolute(absolute)
 
         self.point1 = point1
-        self.point2 = point2
-        self.point3 = point3
 
     ##############################################
 
     def clone(self, path):
-        return self.__class__(path, self._position, self._point1, self._point2, self._point3)
+        return self.__class__(path, self._index, self._point1, absolute)
 
     ##############################################
 
-    def apply_transformation(self, transformation):
-        TwoPointsMixin.transform(self, transformation)
-        self._point3 = transformation * self._point3
+    @property
+    def geometry(self):
+        # Fixme: cache ???
+        # Fixme: !!!
+        return Segment2D(self.start_point, self._point2)
+        # return CubicBezier2D(self.start_point, self._point1, self._point2, self._point3)
+
+####################################################################################################
+
+class StringedCubicBezierSegment(PathPart, TwoPointMixin):
 
     ##############################################
 
-    @property
-    def point3(self):
-        return self._point3
+    def __init__(self, path, index, point1, point2, absolute=False):
 
-    @point3.setter
-    def point3(self, value):
-        self._point3 = Vector2D(value) # self._path.__vector_cls__
+        PathPart.__init__(self, path, index)
+        self._init_absolute(absolute)
+
+        # self.point1 = point1
+
+    ##############################################
+
+    def clone(self, path):
+        return self.__class__(path, self._index, self._point1, self._point2, absolute)
 
     ##############################################
 
     @property
-    def stop_point(self):
-        return self._point3
+    def geometry(self):
+        # Fixme: cache ???
+        # Fixme: !!!
+        return Segment2D(self.start_point, self._point2)
+        # return CubicBezier2D(self.start_point, self._point1, self._point2, self._point3)
+
+####################################################################################################
+
+class ArcSegment(OnePointMixin, PathPart):
+
+    ##############################################
+
+    def __init__(self, path, index, point, radius_x, radius_y, angle, large_arc, sweep, absolute=False):
+
+        PathPart.__init__(self, path, index)
+        self._init_absolute(absolute)
+
+        self.point = point
+
+        self._large_arc = bool(large_arc)
+        self._sweep = bool(sweep)
+        self._radius_x = radius_x
+        self._radius_y = radius_y
+        self._angle = angle
+
+    ##############################################
+
+    def clone(self, path):
+        return self.__class__(
+            path,
+            self._index,
+            self._point,
+            self._radius_x, self._radius_y,
+            self._angle,
+            self._large_arc, self._sweep,
+            self._absolute,
+        )
+
+    ##############################################
+
+    def __repr__(self):
+        template = '{0}(@{1._index} {1._point} rx={1._radius_x} ry={1._radius_y} a={1._angle} la={1._large_arc} s={1._sweep})'
+        return template.format(self.__class__.__name__, self)
+
+    ##############################################
+
+    def to_absolute(self):
+        self._point = self.stop_point
+        self._absolute = True
+
+    ##############################################
 
     @property
     def points(self):
-        return (self.start_point, self._point1, self._point2, self._point3)
+        return self.start_point, self.stop_point
 
     ##############################################
 
     @property
     def geometry(self):
-        # Fixme: cache ???
-        return CubicBezier2D(self.start_point, self._point1, self._point2, self._point3)
+        return Ellipse2D.svg_arc(
+            self.start_point, self.stop_point,
+            self._radius_x, self._radius_y,
+            self._angle,
+            self._large_arc, self._sweep,
+        )
 
 ####################################################################################################
 
@@ -541,13 +872,16 @@ class Path2D(Primitive2DMixin, Primitive1P):
 
     """Class to implements 2D Path."""
 
+    _logger = _module_logger.getChild('Path2D')
+
     ##############################################
 
     def __init__(self, start_point):
 
         Primitive1P.__init__(self, start_point)
 
-        self._parts = []
+        self._parts = [] # Fixme: segment ???
+        self._is_closed = False
 
     ##############################################
 
@@ -570,39 +904,89 @@ class Path2D(Primitive2DMixin, Primitive1P):
     def __iter__(self):
         return iter(self._parts)
 
-    def __getitem__(self, position):
+    def __getitem__(self, index):
         # try:
         #     return self._parts[slice_]
         # except IndexError:
         #     return None
-        position = int(position)
-        if 0 <= position < len(self._parts):
-            return self._parts[position]
-        else:
-            return None
+        index = int(index)
+        number_of_parts = len(self._parts)
+        if 0 <= index < number_of_parts:
+            return self._parts[index]
+        # elif self._is_closed:
+        #     if index == -1:
+        #         return self.start_segment
+        #     elif index == number_of_parts:
+        #         return self.stop_segment
+        return None
+
+    ##############################################
+
+    @property
+    def start_segment(self):
+        # Fixme: start_part ???
+        return self._parts[0]
+
+    @property
+    def stop_segment(self):
+        return self._parts[-1]
+
+    ##############################################
+
+    @property
+    def is_closed(self):
+        return self._is_closed
 
     ##############################################
 
     def _add_part(self, part_cls, *args, **kwargs):
-        obj = part_cls(self, len(self._parts), *args, **kwargs)
-        self._parts.append(obj)
-        return obj
+        if not self._is_closed:
+            obj = part_cls(self, len(self._parts), *args, **kwargs)
+            self._parts.append(obj)
+            return obj
 
     ##############################################
 
     def apply_transformation(self, transformation):
 
-        self._p0 = transformation * self._p0
+        # self._logger.info(str(self) + '\n' + str(transformation.type) + '\n' + str(transformation) )
 
-        for i, part in enumerate(self._parts):
+        for part in self._parts:
+            # print(part)
             if isinstance(part, PathSegment):
                 part._reset_cache()
-            if isinstance(part, DirectionalSegment):
+            if isinstance(part, DirectionalSegmentMixin):
                 # Since a rotation will change the direction
                 # DirectionalSegment must be casted to PathSegment
                 part = part.to_path_segment()
-                self._parts[i] = part
+                self._parts[part.index] = part
+            if part._absolute is False:
+                # print('incremental', part, part.points)
+                part.to_absolute()
+                # print('->', part.points)
+
+        # print()
+        self._p0 = transformation * self._p0
+        # print('p0', self._p0)
+        for part in self._parts:
+            # print(part)
             part.apply_transformation(transformation)
+            # print('->', part.points)
+
+    ##############################################
+
+    @property
+    def bounding_box(self):
+
+        # Fixme: cache
+        bounding_box = None
+        for item in self._parts:
+            interval = item.geometry.bounding_box
+            if bounding_box is None:
+                bounding_box = interval
+            else:
+                bounding_box |= interval
+        return bounding_box
 
     ##############################################
 
@@ -611,52 +995,124 @@ class Path2D(Primitive2DMixin, Primitive1P):
 
     ##############################################
 
-    def horizontal_to(self, distance, radius=None):
-        return self._add_part(HorizontalSegment, distance, radius)
+    def horizontal_to(self, length, radius=None, absolute=False):
+        if absolute:
+            return self._add_part(PathSegment, self.__vector_cls__(length, 0), radius,
+                                  absolute=True)
+        else:
+            return self._add_part(HorizontalSegment, length, radius)
+
+    ##############################################
+
+    def vertical_to(self, length, radius=None, absolute=False):
+        if absolute:
+            return self._add_part(PathSegment, self.__vector_cls__(0, length), radius,
+                                  absolute=True)
+        else:
+            return self._add_part(VerticalSegment, length, radius)
+
+    ##############################################
+
+    def absolute_horizontal_to(self, x, radius=None):
+        return self._add_part(AbsoluteHorizontalSegment, x, radius)
+
+    def absolute_vertical_to(self, y, radius=None):
+        return self._add_part(AbsoluteVerticalSegment, y, radius)
+
+    ##############################################
+
+    def north_to(self, length, radius=None):
+        return self._add_part(NorthSegment, length, radius)
+
+    def south_to(self, length, radius=None):
+        return self._add_part(SouthSegment, length, radius)
+
+    def west_to(self, length, radius=None):
+        return self._add_part(WestSegment, length, radius)
+
+    def east_to(self, length, radius=None):
+        return self._add_part(EastSegment, length, radius)
+
+    def north_east_to(self, length, radius=None):
+        return self._add_part(NorthEastSegment, length, radius)
+
+    def south_east_to(self, length, radius=None):
+        return self._add_part(SouthEastSegment, length, radius)
+
+    def north_west_to(self, length, radius=None):
+        return self._add_part(NorthWestSegment, length, radius)
+
+    def south_west_to(self, length, radius=None):
+        return self._add_part(SouthWestSegment, length, radius)
+
+    ##############################################
+
+    def line_to(self, point, radius=None, absolute=False):
+        return self._add_part(PathSegment, point, radius, absolute=absolute)
 
-    def vertical_to(self, distance, radius=None):
-        return self._add_part(VerticalSegment, distance, radius)
+    ##############################################
 
-    def north_to(self, distance, radius=None):
-        return self._add_part(NorthSegment, distance, radius)
+    def close(self, radius=None, close_radius=None):
 
-    def south_to(self, distance, radius=None):
-        return self._add_part(SouthSegment, distance, radius)
+        # Fixme: identify as close for SVG export <-- meaning ???
 
-    def west_to(self, distance, radius=None):
-        return self._add_part(WestSegment, distance, radius)
+        closing = close_radius is not None
+        segment = self._add_part(PathSegment, self._p0, radius, absolute=True, closing=closing)
+        if closing:
+            self.start_segment.close(close_radius)
+        self._is_closed = True
 
-    def east_to(self, distance, radius=None):
-        return self._add_part(EastSegment, distance, radius)
+        return segment
 
-    def north_east_to(self, distance, radius=None):
-        return self._add_part(NorthEastSegment, distance, radius)
+    ##############################################
 
-    def south_east_to(self, distance, radius=None):
-        return self._add_part(SouthEastSegment, distance, radius)
+    def quadratic_to(self, point1, point2, absolute=False):
+        return self._add_part(QuadraticBezierSegment, point1, point2, absolute=absolute)
 
-    def north_west_to(self, distance, radius=None):
-        return self._add_part(NorthWestSegment, distance, radius)
+    ##############################################
 
-    def south_west_to(self, distance, radius=None):
-        return self._add_part(SouthWestSegment, distance, radius)
+    def cubic_to(self, point1, point2, point3, absolute=False):
+        return self._add_part(CubicBezierSegment, point1, point2, point3, absolute=absolute)
 
     ##############################################
 
-    def line_to(self, point, radius=None):
-        return self._add_part(PathSegment, point, radius)
+    def stringed_quadratic_to(self, point, absolute=False):
+        return self._add_part(StringedQuadraticBezierSegment, point, absolute=absolute)
+
+    ##############################################
+
+    def stringed_cubic_to(self, point1, point2, absolute=False):
+        return self._add_part(StringedCubicBezierSegment, point1, point2, absolute=absolute)
+
+    ##############################################
 
-    def close(self, radius=None):
-        # Fixme: identify as close for SVG export
-        # Fixme: radius must apply to start and stop
-        return self._add_part(PathSegment, self._p0, radius, absolute=True)
+    def arc_to(self, point, radius_x, radius_y, angle, large_arc, sweep, absolute=False):
+        return self._add_part(ArcSegment, point, radius_x, radius_y, angle, large_arc, sweep,
+                              absolute=absolute)
 
     ##############################################
 
-    def quadratic_to(self, point1, point2):
-        return self._add_part(QuadraticBezierSegment, point1, point2)
+    @classmethod
+    def rounded_rectangle(cls, point, width, height, radius=None):
+
+        path = cls(point)
+        path.horizontal_to(width)
+        path.vertical_to(height, radius=radius)
+        path.horizontal_to(-width, radius=radius)
+        path.close(radius=radius, close_radius=radius)
+
+        return path
 
     ##############################################
 
-    def cubic_to(self, point1, point2, point3):
-        return self._add_part(CubicBezierSegment, point1, point2, point3)
+    @classmethod
+    def circle(cls, point, radius):
+
+        diameter = 2*float(radius)
+        path = cls(point)
+        path.horizontal_to(diameter)
+        path.vertical_to(diameter, radius=radius)
+        path.horizontal_to(-diameter, radius=radius)
+        path.close(radius=radius, close_radius=radius)
+
+        return path

Patro/GeometryEngine/Polygon.py

@@ -19,6 +19,7 @@
 ####################################################################################################
 
 """Module to implement polygon.
+
 """
 
 ####################################################################################################
@@ -29,14 +30,18 @@ __all__ = ['Polygon2D']
 
 import math
 
-from Patro.Common.Math.Functions import sign
-from .Primitive import PrimitiveNP, Primitive2DMixin
+import numpy as np
+
+from .Primitive import PrimitiveNP, ClosedPrimitiveMixin, PathMixin, Primitive2DMixin
 from .Segment import Segment2D
 from .Triangle import Triangle2D
+from Patro.Common.Math.Functions import sign
 
 ####################################################################################################
 
-class Polygon2D(Primitive2DMixin, PrimitiveNP):
+# Fixme: PrimitiveNP last ???
+
+class Polygon2D(Primitive2DMixin, ClosedPrimitiveMixin, PathMixin, PrimitiveNP):
 
     """Class to implements 2D Polygon."""
 
@@ -86,11 +91,14 @@ class Polygon2D(Primitive2DMixin, PrimitiveNP):
         self._is_simple = None
         self._is_convex = None
 
-    ##############################################
+        self._area = None
+        # self._cross = None
+        # self._barycenter = None
 
-    @property
-    def is_closed(self):
-        return True
+        # self._major_axis_angle = None
+        self._major_axis = None
+        # self._minor_axis = None
+        # self._axis_ratio = None
 
     ##############################################
 
@@ -106,20 +114,17 @@ class Polygon2D(Primitive2DMixin, PrimitiveNP):
 
     ##############################################
 
-    # barycenter
-    # momentum
-
-    ##############################################
-
     @property
     def edges(self):
 
         if self._edges is None:
+            edges = []
             N = self.number_of_points
             for i in range(N):
                 j = (i+1) % N
                 edge = Segment2D(self._points[i], self._points[j])
-                self._edges.append(edge)
+                edges.append(edge)
+            self._edges = edges
 
         return iter(self._edges)
 
@@ -147,6 +152,7 @@ class Polygon2D(Primitive2DMixin, PrimitiveNP):
                             # two edge intersect
                             # intersections.append(intersection)
                             return False
+        return True
 
     ##############################################
 
@@ -192,29 +198,258 @@ class Polygon2D(Primitive2DMixin, PrimitiveNP):
 
     @property
     def perimeter(self):
-        return sum([edge.magnitude for edge in self.edges])
+        return sum([edge.length for edge in self.edges])
 
     ##############################################
 
     @property
-    def area(self):
+    def point_barycenter(self):
+        center = self.start_point
+        for point in self.iter_from_second_point():
+            center += point
+        return center / self.number_of_points
+
+    ##############################################
+
+    def _compute_area_barycenter(self):
+
+        r"""Compute polygon area and barycenter.
+
+        Polygon area is determined by
+
+        .. math::
+            \begin{align}
+             \mathbf{A} &= \frac{1}{2} \sum_{i=0}^{n-1} P_i \otimes P_{i+1} \\
+                        &= \frac{1}{2} \sum_{i=0}^{n-1}
+                        \begin{vmatrix}
+                         x_i & x_{i+1} \\
+                         y_i & y_{i+1}
+                        \end{vmatrix} \\
+                        &= \frac{1}{2} \sum_{i=0}^{n-1} x_i y_{i+1} - x_{i+1} y_i
+            \end{align}
+
+        where :math:`x_n = x_0`
+
+        Polygon barycenter is determined by
+
+        .. math::
+            \begin{align}
+            \mathbf{C} &= \frac{1}{6\mathbf{A}} \sum_{i=0}^{n-1}
+                          (P_i + P_{i+1}) \times (P_i \otimes P_{i+1}) \\
+                       &= \frac{1}{6\mathbf{A}} \sum_{i=0}^{n-1}
+                          \begin{pmatrix}
+                              (x_i + x_{i+1}) (x_i y_{i+1} - x_{i+1} y_i) \\
+                              (y_i + y_{i+1}) (x_i y_{i+1} - x_{i+1} y_i)
+                          \end{pmatrix}
+            \end{align}
+
+        References
+
+          * On the Calculation of Arbitrary Moments of Polygons,
+            Carsten Steger,
+            Technical Report FGBV–96–05,
+            October 1996
+          * http://mathworld.wolfram.com/PolygonArea.html
+          * https://en.wikipedia.org/wiki/Polygon#Area_and_centroid
+
+        """
 
         if not self.is_simple:
             return None
 
-        # http://mathworld.wolfram.com/PolygonArea.html
+        # area = self._points[-1].cross(self._points[0])
+        # for i in range(self.number_of_points):
+        #     area *= self._points[i].cross(self._points[i+1])
+
+        # P0, P1, Pn-1, P0
+        points = self.closed_point_array
+
+        # from 0 to n-1 : P0, ..., Pn-1
+        xi  = points[0,:-1]
+        yi  = points[1,:-1]
+        # from 1 to n : P1, ..., Pn-1, P0
+        xi1 = points[0,1:]
+        yi1 = points[1,1:]
+
+        # Fixme: np.cross ???
+        cross = xi * yi1 - xi1 * yi
+        self._cross = cross
+
+        area = .5 * np.sum(cross)
+
+        if area == 0:
+            # print('Null area')
+            self._area = 0
+            self._barycenter = self.start_point
+        else:
+            factor = 1 / (6*area)
+            x = factor * np.sum((xi + xi1) * cross)
+            y = factor * np.sum((yi + yi1) * cross)
+
+            # area of a convex polygon is defined to be positive if the points are arranged in a
+            # counterclockwise order, and negative if they are in clockwise order (Beyer 1987).
+            self._area = abs(area)
+            self._barycenter = self.__vector_cls__(x, y)
+
+    ##############################################
+
+    def _compute_inertia_moment(self):
+
+        r"""Compute inertia moment on vertices.
+
+        .. warning:: untrusted formulae
+
+        .. math::
+             \begin{align}
+             I_x    &= \frac{1}{12} \sum (y_i^2 + y_i y_{i+1} + y_{i+1}^2) (x_i y_{i+1} - x_{i+1} y_i) \\
+             I_y    &= \frac{1}{12} \sum (x_i^2 + x_i x_{i+1} + x_{i+1}^2) (x_i y_{i+1} - x_{i+1} y_i) \\
+             I_{xy} &= \frac{1}{24} \sum (x_i y_{i+1} + 2 x_i y_i + 2 x_{i+1} y_{i+1} + x_{i+1} y_i) (x_i y_{i+1} - x_{i+1} y_i)
+             \end{align}
+
+        Reference
+
+          * https://en.wikipedia.org/wiki/Second_moment_of_area#Any_cross_section_defined_as_polygon
+
+        """
+
+        # self.recenter()
+
+        # Fixme: duplicated code
+        # P0, P1, Pn-1, P0
+        points = self.closed_point_array
 
-        # A = 1/2 (x1*y2 - x2*y1 + x2*y3 - x3*y2 + ... + x(n-1)*yn - xn*y(n-1) + xn*y1 - x1*yn)
-        #           determinant
+        # from 0 to n-1 : P0, ..., Pn-1
+        xi  = points[0,:-1]
+        yi  = points[1,:-1]
+        # from 1 to n : P1, ..., Pn-1, P0
+        xi1 = points[0,1:]
+        yi1 = points[1,1:]
 
-        area = self._points[-1].cross(self._points[0])
-        for i in range(self.number_of_points):
-            area *= self._points[i].cross(self._points[i+1])
+        # computation on vertices
+        number_of_points = self.number_of_points
+        Ix =    np.sum(yi**2) / number_of_points
+        Iy =    np.sum(xi**2) / number_of_points
+        Ixy = - np.sum(xi*yi) / number_of_points
 
-        # area of a convex polygon is defined to be positive if the points are arranged in a
-        # counterclockwise order, and negative if they are in clockwise order (Beyer 1987).
+        # cross = xi * yi1 - xi1 * yi
+        # cross = self._cross
+        # Ix  = 1/(12*self._area) * np.sum((yi**2 + yi*yi1 + yi1**2) * cross)
+        # Iy  = 1/(12*self._area) * np.sum((xi**2 + xi*xi1 + xi1**2) * cross)
+        # Ixy = 1/(24*self._area) * np.sum((xi*yi1 + 2*(xi*yi + xi1*yi1) + xi1*yi) * cross)
+        # cx, cy = self._barycenter
+        # Ix -= cy**2
+        # Iy -= cx**2
+        # Ixy -= cx*cy
+        # Ix = -Ix
+        # Iy = -Iy
+        # print(Ix, Iy, Ixy)
 
-        return abs(area) / 2
+        if Ixy == 0:
+            if Iy >= Ix:
+                self._major_axis_angle = 0
+
+                lambda1 = Iy
+                lambda2 = Ix
+
+                vx  = 0
+                v1y = 1
+                v2y = 0
+            else:
+                self._major_axis_angle = 90
+
+                lambda1 = Ix
+                lambda2 = Iy
+
+                vx  = 1
+                v1y = 0
+                v2y = 1
+        else:
+            Is = Iy + Ix
+            Id = Ix - Iy
+
+            sqrt0 = math.sqrt(Id*Id + 4*Ixy*Ixy)
+
+            lambda1 = (Is + sqrt0) / 2
+            lambda2 = (Is - sqrt0) / 2
+
+            vx  = Ixy
+            v1y = (Id + sqrt0) / 2
+            v2y = (Id - sqrt0) / 2
+
+            if lambda1 < lambda2:
+                v1y, v2y = v2y, v1y
+                lambda1, lambda2 = lambda2, lambda1
+
+            self._major_axis_angle = - math.degrees(math.atan(v1y/vx))
+
+        self._major_axis = 4 * math.sqrt(math.fabs(lambda1))
+        self._minor_axis = 4 * math.sqrt(math.fabs(lambda2))
+
+        if self._minor_axis != 0:
+            self._axis_ratio = self._major_axis / self._minor_axis
+        else:
+            self._axis_ratio = 0
+
+    ##############################################
+
+    def _check_area(self):
+        if self.is_simple and self._area is None:
+            self._compute_area_barycenter()
+
+    ##############################################
+
+    @property
+    def area(self):
+        """Return polygon area."""
+        self._check_area()
+        return self._area
+
+    ##############################################
+
+    @property
+    def barycenter(self):
+        """Return polygon barycenter."""
+        self._check_area()
+        return self._barycenter
+
+    ##############################################
+
+    def recenter(self):
+        """Recenter the polygon to the barycenter."""
+        # if self._centred:
+        #     return
+        barycenter = self._barycenter
+        for point in self._points:
+            point -= barycenter
+        # self._centred = True
+
+    ##############################################
+
+    def _check_moment(self):
+        if self.is_simple and self._major_axis is None:
+            self._compute_inertia_moment()
+
+    ##############################################
+
+    @property
+    def major_axis_angle(self):
+        self._check_moment()
+        return self._major_axis_angle
+
+    @property
+    def major_axis(self):
+        self._check_moment()
+        return self._major_axis
+
+    @property
+    def minor_axis(self):
+        self._check_moment()
+        return self._minor_axis
+
+    @property
+    def axis_ratio(self):
+        self._check_moment()
+        return self._axis_ratio
 
     ##############################################
 

Patro/GeometryEngine/Polyline.py

@@ -28,6 +28,7 @@ __all__ = ['Polyline2D']
 
 ####################################################################################################
 
+from .Path import Path2D
 from .Primitive import PrimitiveNP, Primitive2DMixin
 from .Segment import Segment2D
 
@@ -74,3 +75,13 @@ class Polyline2D(Primitive2DMixin, PrimitiveNP):
             if distance is None or edge_distance < distance:
                 distance = edge_distance
         return distance
+
+    ##############################################
+
+    def to_path(self):
+
+        path = Path2D(self.start_point)
+        for point in self.iter_from_second_point():
+            path.line_to(point)
+
+        return path

Patro/GeometryEngine/Primitive.py

@@ -22,6 +22,10 @@
 
 """
 
+# Fixme:
+#  length, interpolate path
+#  area
+
 ####################################################################################################
 
 __all__ = [
@@ -35,7 +39,7 @@ __all__ = [
 
 ####################################################################################################
 
-import collections
+from collections import abc as collections
 
 import numpy as np
 
@@ -45,18 +49,10 @@ from .BoundingBox import bounding_box_from_points
 
 ####################################################################################################
 
-# Fixme:
-#  length, interpolate path
-#  area
-
-####################################################################################################
-
 class Primitive:
 
     """Base class for geometric primitive"""
 
-    # __dimension__ = None # in [2, 3] for 2D / 3D primitive
-
     __vector_cls__ = None
 
     ##############################################
@@ -134,9 +130,7 @@ class Primitive:
 
         Return None if primitive is infinite.
         """
-
         # Fixme: cache
-
         if self.is_infinite:
             return None
         else:
@@ -156,9 +150,8 @@ class Primitive:
 
         """
         obj = self.clone() if clone else self
-        # for point in obj.points:
-        #     point *= transformation # don't work
-        obj.apply_transformation(transformation)
+        if not transformation.is_identity:
+            obj.apply_transformation(transformation)
         return obj
 
     ##############################################
@@ -198,29 +191,83 @@ class Primitive:
     ##############################################
 
     @property
-    def geometry_matrix(self):
+    def point_array(self):
+        r"""Return the geometry matrix as a Numpy array.
+
+        .. math::
+           \mathrm{Geometry\ Matrix} =
+           \begin{bmatrix}
+           x_0  & x_1 & \ldots & x_{n-1} \\
+           y_0  & y_1 & \ldots & y_{n-1}
+           \end{bmatrix}
+
+        """
+        # Fixme: geometry_matrix vs point_array
+        # Fixme: cache ??? but point set and init
+        # if self._point_array is None:
+        #     self._point_array = np.array(list(self.points)).transpose()
+        # return self._point_array
         return np.array(list(self.points)).transpose()
 
     ##############################################
 
-    def is_close(self, other):
+    def is_point_close(self, other):
         # Fixme: verus is_closed
-        return np.allclose(self.geometry_matrix, other.geometry_matrix)
+        #  is_similar
+        return np.allclose(self.point_array, other.point_array)
 
 ####################################################################################################
 
 class Primitive2DMixin:
 
-    # __dimension__ = 2
-
     __vector_cls__ = None # Fixme: due to import, done in module's __init__.py
 
+    # __dimension__ = 2
     @property
     def dimension(self):
         return 2
 
 ####################################################################################################
 
+class ClosedPrimitiveMixin:
+
+    # Fixme: should be reversible
+
+    ##############################################
+
+    @property
+    def is_closed(self):
+        """True if the primitive is a closed path."""
+        return True
+
+    ##############################################
+
+    @property
+    def closed_points(self):
+        points = list(self.points)
+        points.append(self.start_point)
+        return points
+
+    ##############################################
+
+    @property
+    def closed_point_array(self):
+        r"""Return the geometry matrix as a Numpy array for a closed primitive.
+
+        .. math::
+           \mathrm{Geometry\ Matrix} =
+           \begin{bmatrix}
+           x_0  & x_1 & \ldots & x_{n-1} & x_0 \\
+           y_0  & y_1 & \ldots & y_{n-1} & y_0
+           \end{bmatrix}
+
+        """
+        # Fixme: place, func for closed_point
+        # Fixme: cache ???
+        return np.array(self.closed_points).transpose()
+
+####################################################################################################
+
 class ReversiblePrimitiveMixin:
 
     ##############################################
@@ -258,7 +305,6 @@ class Primitive1P(Primitive):
     ##############################################
 
     def __init__(self, p0):
-
         self.p0 = p0
 
     ##############################################
@@ -297,7 +343,6 @@ class Primitive2P(Primitive1P, ReversiblePrimitiveMixin):
     ##############################################
 
     def __init__(self, p0, p1):
-
         # We don't call super().__init__(p0) for speed
         self.p0 = p0
         self.p1 = p1
@@ -356,7 +401,6 @@ class Primitive3P(Primitive2P):
     ##############################################
 
     def __init__(self, p0, p1, p2):
-
         self.p0 = p0
         self.p1 = p1
         self.p2 = p2
@@ -389,8 +433,13 @@ class Primitive3P(Primitive2P):
 
     @property
     def reversed_points(self):
+        # Fixme: share code ???
         return iter((self._p2, self._p1, self._p0))
 
+    def iter_from_second_point(self):
+        # Fixme: share code ???
+        return iter(self._p1, self._p2)
+
     ##############################################
 
     def _set_points(self, points):
@@ -403,7 +452,6 @@ class Primitive4P(Primitive3P):
     ##############################################
 
     def __init__(self, p0, p1, p2, p3):
-
         self.p0 = p0
         self.p1 = p1
         self.p2 = p2
@@ -439,6 +487,9 @@ class Primitive4P(Primitive3P):
     def reversed_points(self):
         return iter((self._p3, self._p2, self._p1, self._p0))
 
+    def iter_from_second_point(self):
+        return iter(self._p1, self._p2, self._p3)
+
     ##############################################
 
     def _set_points(self, points):
@@ -490,13 +541,8 @@ class PrimitiveNP(Primitive, ReversiblePrimitiveMixin):
     def reversed_points(self):
         return reversed(self._points)
 
-    ##############################################
-
-    @property
-    def point_array(self):
-        if self._point_array is None:
-            self._point_array = np.array([point for point in self._points])
-        return self._point_array
+    def iter_from_second_point(self):
+        return iter(self._points[1:])
 
     ##############################################
 
@@ -518,3 +564,33 @@ class PrimitiveNP(Primitive, ReversiblePrimitiveMixin):
 
         for i in range(self.number_of_points - size +1):
             yield self._points[i:i+size]
+
+####################################################################################################
+
+class PolygonMixin:
+
+    ##############################################
+
+    def to_polygon(self):
+        from .Polygon import Polygon2D
+        return Polygon2D(self.points)
+
+####################################################################################################
+
+class PathMixin:
+
+    ##############################################
+
+    def to_path(self):
+
+        from .Path import Path2D
+
+        path = Path2D(self.start_point)
+        for point in self.iter_from_second_point():
+            path.line_to(point)
+        if self.is_closed:
+            path.close()
+
+        return path
+
+

Patro/GeometryEngine/Rectangle.py

@@ -30,12 +30,13 @@ __all__ = ['Rectangle2D']
 
 import math
 
-from .Primitive import Primitive2P, Primitive2DMixin
+from .Path import Path2D
+from .Primitive import Primitive2P, ClosedPrimitiveMixin, PathMixin, PolygonMixin, Primitive2DMixin
 from .Segment import Segment2D
 
 ####################################################################################################
 
-class Rectangle2D(Primitive2DMixin, Primitive2P):
+class Rectangle2D(Primitive2DMixin, ClosedPrimitiveMixin, PathMixin, PolygonMixin, Primitive2P):
 
     """Class to implements 2D Rectangle."""
 

Patro/GeometryEngine/Segment.py

@@ -76,7 +76,7 @@ class Segment2D(Primitive2DMixin, Primitive2P):
 
     def to_line(self):
         # Fixme: cache
-        return Line2D.from_two_points(self._p1, self._p0)
+        return Line2D.from_two_points(self._p0, self._p1)
 
     ##############################################
 
@@ -116,7 +116,7 @@ class Segment2D(Primitive2DMixin, Primitive2P):
 
         s1, s2 = line1.intersection_abscissae(line2)
         if s1 is None:
-            return None
+            return None, None
         else:
             intersect = (0 <= s1 <= 1) and (0 <= s2 <= 1)
             return self.interpolate(s1), intersect

Patro/GeometryEngine/Spline.py

@@ -20,256 +20,15 @@
 
 r"""Module to implement Spline curve.
 
-B-spline Basis
---------------
-
-A nonuniform, nonrational B-spline of order `k` is a piecewise polynomial function of degree
-:math:`k - 1` in a variable `t`.
-
-.. check: k+1 knots ???
-
-.. It is defined over :math:`k + 1` locations :math:`t_i`, called knots, which must be in
-   non-descending order :math:`t_i \leq t_{i+1}`. This series defines a knot vector :math:`T = (t_0,
-   \ldots, t_{k})`.
-
-A set of non-descending breaking points, called knot, :math:`t_0 \le t_1 \le \ldots \le t_m` defines
-a knot vector :math:`T = (t_0, \ldots, t_{m})`.
-
-If each knot is separated by the same distance `h` (where :math:`h = t_{i+1} - t_i`) from its
-predecessor, the knot vector and the corresponding B-splines are called "uniform".
-
-Given a knot vector `T`, the associated B-spline basis functions, :math:`B_i^k(t)` are defined as:
-
-.. t \in [t_i, t_{i+1}[
-
-.. math::
-    B_i^1(t) =
-    \left\lbrace
-    \begin{array}{l}
-       1 \;\textrm{if}\; t_i \le t < t_{i+1} \\
-       0 \;\textrm{otherwise}
-    \end{array}
-    \right.
-
-.. math::
-   \begin{split}
-    B_i^k(t) &= \frac{t - t_i}{t_{i+k-1} - t_i} B_i^{k-1}(t)
-              + \frac{t_{i+k} - t}{t_{i+k} - t_{i+1}} B_{i+1}^{k-1}(t) \\
-             &= w_i^{k-1}(t) B_i^{k-1}(t) + [1 - w_{i+1}^{k-1}(t)] B_{i+1}^{k-1}(t)
-   \end{split}
-
-where
-
-.. math::
-    w_i^k(t) =
-    \left\lbrace
-    \begin{array}{l}
-       \frac{t - t_i}{t_{i+k} - t_i} \;\textrm{if}\; t_i < t_{i+k} \\
-       0 \;\textrm{otherwise}
-    \end{array}
-    \right.
-
-These equations have the following properties, for :math:`k > 1` and :math:`i = 0, 1, \ldots, n` :
-
-* Positivity: :math:`B_i^k(t) > 0`, for :math:`t_i < t < t_{i+k}`
-* Local Support: :math:`B_i^k(t) = 0`, for :math:`t_0 \le t \le t_i` and :math:`t_{i+k} \le t \le t_{n+k}`
-* Partition of unity: :math:`\sum_{i=0}^n B_i^k(t)= 1`, for :math:`t \in [t_0, t_m]`
-* Continuity: :math:`B_i^k(t)` as :math:`C^{k-2}` continuity at each simple knot
-
-.. The B-spline contributes only in the range between the first and last of these knots and is zero
-   elsewhere.
-
-B-spline Curve
---------------
-
-A B-spline curve of order `k` is defined as a linear combination of control points :math:`p_i` and
-B-spline basis functions :math:`B_i^k(t)` given by
-
-.. math::
-    S^k(t) = \sum_{i=0}^{n} p_i\; B_i^k(t) ,\quad n \ge k - 1,\; t \in [t_{k-1}, t_{n+1}]
-
-In this context the control points are called De Boor points. The basis functions :math:`B_i^k(t)`
-is defined on a knot vector
-
-.. math::
-    T = (t_0, t_1, \ldots, t_{k-1}, t_k, t_{k+1}, \ldots, t_{n-1}, t_n, t_{n+1}, \ldots, t_{n+k})
-
-where there are :math:`n+k+1` elements, i.e. the number of control points :math:`n+1` plus the order
-of the curve `k`.  Each knot span :math:`t_i \le t \le t_{i+1}` is mapped onto a polynomial curve
-between two successive joints :math:`S(t_i)` and :math:`S(t_{i+1})`.
-
-Unlike Bézier curves, B-spline curves do not in general pass through the two end control points.
-Increasing the multiplicity of a knot reduces the continuity of the curve at that knot.
-Specifically, the curve is :math:`(k-p-1)` times continuously differentiable at a knot with
-multiplicity :math:`p (\le k)`, and thus has :math:`C^{k-p-1}` continuity.  Therefore, the control
-polygon will coincide with the curve at a knot of multiplicity :math:`k-1`, and a knot with
-multiplicity `k` indicates :math:`C^{-1}` continuity, or a discontinuous curve.  Repeating the knots
-at the end `k` times will force the endpoints to coincide with the control polygon.  Thus the first
-and the last control points of a curve with a knot vector described by
-
-.. math::
-    \begin{eqnarray}
-    T = (
-    \underbrace{t_0, t_1, \ldots, t_{k-1},}_{\mbox{$k$ equal knots}}
-    \quad
-    \underbrace{t_k, t_{k+1}, \ldots, t_{n-1}, t_n,}_{\mbox{$n$-$k$+1 internal knots}}
-    \quad
-    \underbrace{t_{n+1}, \ldots, t_{n+k}}_{\mbox{$k$ equal knots}})
-    \end{eqnarray}
-
-coincide with the endpoints of the curve.  Such knot vectors and curves are known as *clamped*. In
-other words, *clamped/unclamped* refers to whether both ends of the knot vector have multiplicity
-equal to `k` or not.
-
-**Local support property**: A single span of a B-spline curve is controlled only by `k` control
-points, and any control point affects `k` spans.  Specifically, changing :math:`p_i` affects the
-curve in the parameter range :math:`t_i < t < t_{i+k}` and the curve at a point where :math:`t_r < t
-< t_{r+1}` is determined completely by the control points :math:`p_{r-(k-1)}, \ldots, p_r`.
-
-**B-spline to Bézier property**: From the discussion of end points geometric property, it can be
-seen that a Bézier curve of order `k` (degree :math:`k-1`) is a B-spline curve with no internal
-knots and the end knots repeated `k` times.  The knot vector is thus
-
-.. math::
- \begin{eqnarray}
-    T = (
-    \underbrace{t_0, t_1, \ldots, t_{k-1}}_{\mbox{$k$ equal knots}}
-    ,\quad
-    \underbrace{t_{n+1}, \ldots, t_{n+k}}_{\mbox{$k$ equal knots}}
-    )
-    \end{eqnarray}
-
-where :math:`n+k+1 = 2k` or :math:`n = k-1`.
-
-Algorithms for B-spline curves
-------------------------------
-
-Evaluation and subdivision algorithm
-~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
-A B-spline curve can be evaluated at a specific parameter value `t` using the de Boor algorithm,
-which is a generalization of the de Casteljau algorithm.  The repeated substitution of the recursive
-definition of the B-spline basis function into the previous definition and re-indexing leads to the
-following de Boor algorithm:
-
-..  math::
-   S(t) = \sum_{i=0}^{n+j} p_i^j B_i^{k-j}(t) ,\quad j = 0, 1, \ldots, k-1
-
-where
-
-.. math::
-   p_i^j = \Big[1 - w_i^j\Big] p_{i-1}^{j-1} + w_i^j p_i^{j-1}, \; j > 0
-
-with
-
-.. math::
-
-  w_i^j = \frac{t - t_i}{t_{i+k-j} - t_i} \quad \textrm{and} \; p_j^0 = p_j
-
-For :math:`j = k-1`, the B-spline basis function reduces to :math:`B_l^1` for :math:`t \in [t_l,
-t_{l+1}]`, and :math:`p_l^{k-1}` coincides with the curve :math:`S(t) = p_l^{k-1}`.
-
-The de Boor algorithm is a generalization of the de Casteljau algorithm. The de Boor algorithm also
-permits the subdivision of the B-spline curve into two segments of the same order.
-
-De Boor Algorithm
-~~~~~~~~~~~~~~~~~
-
-Let the index `l` define the knot interval that contains the position, :math:`t \in [t_l ,
-t_{l+1}]`.  We can see in the recursion formula that only B-splines with :math:`i = l-K, \dots, l`
-are non-zero for this knot interval, where :math:`K = k - 1` is the degree.  Thus, the sum is
-reduced to:
-
-.. math::
-    S^k(t) = \sum _{i=l-K}^{l} p_{i} B_i^k(t)
-
-The algorithm does not compute the B-spline functions :math:`B_i^k(t)` directly.  Instead it
-evaluates :math:`S(t)` through an equivalent recursion formula.
-
-Let :math:`d _i^r` be new control points with :math:`d_i^1 = p_i` for :math:`i = l-K, \dots, l`.
-
-For :math:`r = 2, \dots, k` the following recursion is applied:
-
-.. math::
-    d_i^r = (1 - w_i^r) d_{i-1}^{r-1} + w_i^r d_i^{r-1} \quad i = l-K+r, \dots, l
-
-    w_i^r = \frac{t - t_i}{t_{i+1+l-r} - t_{i}}
-
-Once the iterations are complete, we have :math:`S^k(t) = d_l^k`.
-
-.. , meaning that :math:`d_l^k` is the desired result.
-
-De Boor's algorithm is more efficient than an explicit calculation of B-splines :math:`B_i^k(t)`
-with the Cox-de Boor recursion formula, because it does not compute terms which are guaranteed to be
-multiplied by zero.
-
-..
-    :math:`S(t) = p_j^k` for :math:`t \in [t_j , t_{j+1}[` for :math:`k \le j \le n` with the following relation:
-    .. math::
-        \begin{split}
-        p_i^{r+1} &= \frac{t - t_i}{t_{i+k-r} - t} p_i^r + \frac{t_{i+k-r} - t_i}{t_{i+k-r} - t_i} p_{i-1}^r \\
-                  &= w_i^{k-r}(t) p_i^r +  (1 - w_i^{k-r}(t)) p_{i-1}^r
-        \end{split}
-
-Knot insertion
-~~~~~~~~~~~~~~
-
-A knot can be inserted into a B-spline curve without changing the
-geometry of the curve.  The new curve is identical to
-
-.. math::
-    \begin{array}{lcl}
-    \sum_{i=0}^n p_i B_i^k(t) & \textrm{becomes} & \sum_{i=0}^{n+1} \bar{p}_i \bar B_i^k(t) \\
-    \mbox{over}\; T = (t_0, t_1, \ldots, t_l, t_{l+1}, \ldots) & &
-    \mbox{over}\; T = (t_0, t_1, \ldots, t_l, \bar t, t_{l+1}, \ldots) & &
-    \end{array}
-
-when a new knot :math:`\bar t` is inserted between knots :math:`t_l` and :math:`t_{l+1}`.  The new
-de Boor points are given by
-
-.. math::
-    \bar{p}_i = (1 - w_i) p_{i-1} + w_i p_i
-
-where
-
-.. math::
-    w_i =
-    \left\{ \begin{array}{ll}
-    1 & i \le l-k+1 \\
-    0 & i \ge l+1 \\
-    \frac{\bar{t} - t_i}{t_{l+k-1} - t_i} & l-k+2 \le i \leq l
-    \end{array}
-    \right.
-
-The above algorithm is also known as **Boehm's algorithm**.  A more general (but also more complex)
-insertion algorithm permitting insertion of several (possibly multiple) knots into a B-spline knot
-vector, known as the Oslo algorithm, was developed by Cohen et al.
-
-A B-spline curve is :math:`C^{\infty}` continuous in the interior of a span.  Within exact
-arithmetic, inserting a knot does not change the curve, so it does not change the continuity.
-However, if any of the control points are moved after knot insertion, the continuity at the knot
-will become :math:`C^{k-p-1}`, where `p` is the multiplicity of the knot.
-
-The B-spline curve can be subdivided into Bézier segments by knot insertion at each internal knot
-until the multiplicity of each internal knot is equal to `k`.
-
-References
-----------
-
-* Computer Graphics, Principle and Practice, Foley et al., Adison Wesley
-* http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node15.html
+For resources on Spline curve see :ref:`this section <spline-geometry-ressources-page>`.
 
 """
 
-# The DeBoor-Cox algorithm permits to evaluate recursively a B-Spline in a similar way to the De
-# Casteljaud algorithm for Bézier curves.
+####################################################################################################
 #
-# Given `k` the degree of the B-spline, `n + 1` control points :math:`p_0, \ldots, p_n`, and an
-# increasing series of scalars :math:`t_0 \le t_1 \le \ldots \le t_m` with :math:`m = n + k + 1`,
-# called knots.
+# Notes: algorithm details are on spline.rst
 #
-# The number of points must respect the condition :math:`n + 1 \le k`, e.g. a B-spline of degree 3
-# must have 4 control points.
+####################################################################################################
 
 ####################################################################################################
 
@@ -316,7 +75,7 @@ class QuadraticUniformSpline2D(Primitive2DMixin, Primitive3P):
 
     def to_bezier(self):
         basis = np.dot(self.BASIS, QuadraticBezier2D.INVERSE_BASIS)
-        points = np.dot(self.geometry_matrix, basis).transpose()
+        points = np.dot(self.point_array, basis).transpose()
         return QuadraticBezier2D(*points)
 
     ##############################################
@@ -376,7 +135,7 @@ class CubicUniformSpline2D(Primitive2DMixin, Primitive4P):
 
     def to_bezier(self):
         basis = np.dot(self.BASIS, CubicBezier2D.INVERSE_BASIS)
-        points = np.dot(self.geometry_matrix, basis).transpose()
+        points = np.dot(self.point_array, basis).transpose()
         if self._start:
             # list(self.points)[:2]
             points[:2] = self._p0, self._p1

Patro/GeometryEngine/Transformation.py

@@ -18,7 +18,9 @@
 #
 ####################################################################################################
 
-"""Module to implement transformations like scale, rotation and translation.
+r"""Module to implement transformations like scale, rotation and translation.
+
+For resources on transformations see :ref:`this section <transformation-geometry-ressources-page>`.
 
 """
 
@@ -56,10 +58,17 @@ class TransformationType(Enum):
 
     Rotation = auto()
 
+    Translation = auto()
+
     Generic = auto()
 
 ####################################################################################################
 
+class IncompatibleArrayDimension(ValueError):
+    pass
+
+####################################################################################################
+
 class Transformation:
 
     __dimension__ = None
@@ -79,17 +88,26 @@ class Transformation:
             if self.same_dimension(obj):
                 array = obj.array # *._m
             else:
-                raise ValueError
+                raise IncompatibleArrayDimension
         elif isinstance(obj, np.ndarray):
             if obj.shape == (self.__size__, self.__size__):
                 array = obj
             else:
-                raise ValueError
+                raise IncompatibleArrayDimension
+        elif isinstance(obj, (list, tuple)):
+            if len(obj) == self.__size__ **2:
+                array = np.array(obj)
+                array.shape = (self.__size__, self.__size__)
+            else:
+                raise IncompatibleArrayDimension
         else:
             array = np.array((self.__size__, self.__size__))
             array[...] = obj
 
         self._m = np.array(array)
+
+        if transformation_type == TransformationType.Generic:
+            transformation_type = self._check_type()
         self._type = transformation_type
 
     ##############################################
@@ -110,6 +128,10 @@ class Transformation:
     def type(self):
         return self._type
 
+    @property
+    def is_identity(self):
+        return self._type == TransformationType.Identity
+
     ##############################################
 
     def __repr__(self):
@@ -125,13 +147,45 @@ class Transformation:
     def same_dimension(self, other):
         return self.__size__ == other.dimension
 
+    ##############################################
+
+    def _check_type(self):
+        raise NotImplementedError
+
+    ##############################################
+
+    def _mul_type(self, obj):
+
+        # Fixme: check matrix value ???
+        #   usage identity/rotation, scale/parity test
+        #   metric test ?
+        # if t in (parity, xparity, yparity) t*t = Id
+        # if t in (rotation, scale) t*t = t
+
+        if self._type == obj._type:
+            if self._type in (
+                    TransformationType.Parity,
+                    TransformationType.XParity,
+                    TransformationType.YParity
+            ):
+                return TransformationType.Identity
+            elif self._type not in (TransformationType.Rotation, TransformationType.Scale):
+                return TransformationType.Generic
+            else:
+                return self._type
+        else: # shear, generic
+            return TransformationType.Generic
+
     #######################################
 
     def __mul__(self, obj):
 
+        """Return self * obj composition."""
+
         if isinstance(obj, Transformation):
+            # T = T1 * T2
             array = np.matmul(self._m, obj.array)
-            return self.__class__(array)
+            return self.__class__(array, self._mul_type(obj))
         elif isinstance(obj, Vector2D):
             array = np.matmul(self._m, np.transpose(obj.v))
             return Vector2D(array)
@@ -150,23 +204,19 @@ class Transformation:
 
     def __imul__(self, obj):
 
+        """Set transformation to obj * self composition."""
+
+        # Fixme: order ???
+
         if isinstance(obj, Transformation):
-            if obj._type != TransformationType.Identity:
-                self._m = np.matmul(self._m, obj.array)
-                # Fixme: check matrix value ???
-                #   usage identity/rotation, scale/parity test
-                #   metric test ?
-                # if t in (parity, xparity, yparity) t*t = Id
-                # if t in (rotation, scale) t*t = t
-                if self._type == obj._type:
-                    if self._type in (TransformationType.Parity,
-                                      TransformationType.XParity,
-                                      TransformationType.YParity):
-                        self._type = TransformationType.Identity
-                    elif self._type not in (TransformationType.Rotation, TransformationType.Scale):
-                        self._type = TransformationType.Generic
-                else: # shear, generic
-                    self._type = TransformationType.Generic
+            if obj.type != TransformationType.Identity:
+                # (T = T1) *= T2
+                # T = T2 * T1
+                # order is inverted !
+                self._m = np.matmul(obj.array, self._m)
+                self._type = self._mul_type(obj)
+                if self._type == TransformationType.Generic:
+                    self._type = self._check_type()
         else:
             raise ValueError
 
@@ -181,17 +231,6 @@ class Transformation2D(Transformation):
 
     ##############################################
 
-    @classmethod
-    def Rotation(cls, angle):
-
-        angle = radians(angle)
-        c = cos(angle)
-        s = sin(angle)
-
-        return cls(np.array(((c, -s), (s,  c))), TransformationType.Rotation)
-
-    ##############################################
-
     @classmethod
     def type_for_scale(cls, x_scale, y_scale):
 
@@ -214,6 +253,51 @@ class Transformation2D(Transformation):
 
     ##############################################
 
+    @classmethod
+    def check_matrix_type(self, matrix):
+
+        # Fixme: check
+        m00, m01, m10, m11 = matrix
+        if m01 == 0 and m10 == 0:
+            if m00 == 1:
+                if m11 == 1:
+                    return TransformationType.Identity
+                elif m11 == -1:
+                    return TransformationType.YParity
+            elif m00 == -1:
+                if m11 == 1:
+                    return TransformationType.XParity
+                elif m11 == -1:
+                    return TransformationType.Parity
+            elif m00 == m11:
+                return TransformationType.Scale
+
+        # Fixme: check for rotation
+        return TransformationType.Generic
+
+    ##############################################
+
+    def _check_type(self):
+        return self._check_matrix_type(self.to_list())
+
+    ##############################################
+
+    @classmethod
+    def Rotation(cls, angle):
+
+        angle = radians(angle)
+        c = cos(angle)
+        s = sin(angle)
+
+        return cls(
+            np.array((
+                (c, -s),
+                (s,  c))),
+            TransformationType.Rotation,
+        )
+
+    ##############################################
+
     @classmethod
     def Scale(cls, x_scale, y_scale=None):
         if y_scale is None:
@@ -250,6 +334,7 @@ class AffineTransformation(Transformation):
 
         transformation = cls.Identity()
         transformation.translation_part[...] = vector.v[...]
+        transformation._type = TransformationType.Translation
         return transformation
 
     ##############################################
@@ -257,9 +342,11 @@ class AffineTransformation(Transformation):
     @classmethod
     def RotationAt(cls, center, angle):
 
-        transformation = cls.Translation(center)
+        # return cls.Translation(center) * cls.Rotation(angle) * cls.Translation(-center)
+
+        transformation = cls.Translation(-center)
         transformation *= cls.Rotation(angle)
-        transformation *= cls.Translation(-center)
+        transformation *= cls.Translation(center)
         return transformation
 
     ##############################################
@@ -281,6 +368,12 @@ class AffineTransformation2D(AffineTransformation):
 
     ##############################################
 
+    def _check_type(self):
+        matrix_type = Transformation2D.check_matrix_type(self.matrix_part.flat)
+        # Fixme: translation etc. !!!
+
+    ##############################################
+
     @classmethod
     def Rotation(cls, angle):
 
@@ -294,13 +387,25 @@ class AffineTransformation2D(AffineTransformation):
     @classmethod
     def Scale(cls, x_scale, y_scale):
 
-        # Fixme: others, use *= ?
+        # Fixme: others, use *= ? (comment means ???)
 
         transformation = cls.Identity()
         transformation.matrix_part[...] = Transformation2D.Scale(x_scale, y_scale).array
         transformation._type = cls.type_for_scale(x_scale, y_scale)
         return transformation
 
+    ##############################################
+
+    @classmethod
+    def Screen(cls, y_height):
+
+        transformation = cls.Identity()
+        # Fixme: better ?
+        transformation.matrix_part[...] = Transformation2D.YReflection().array
+        transformation.translation_part[...] = Vector2D(0, y_height).v[...]
+        transformation._type = TransformationType.Generic
+        return transformation
+
     #######################################
 
     def __mul__(self, obj):
@@ -310,7 +415,8 @@ class AffineTransformation2D(AffineTransformation):
             return obj.__class__(array)
         elif isinstance(obj, Vector2D):
             array = np.matmul(self._m, HomogeneousVector2D(obj).v)
-            return HomogeneousVector2D(array)
+            # return HomogeneousVector2D(array).to_vector()
+            return Vector2D(array[:2])
         else:
             return super(AffineTransformation, self).__mul__(obj)
 

Patro/GeometryEngine/Triangle.py

@@ -30,7 +30,7 @@ __all__ = ['Triangle2D']
 
 import math
 
-from .Primitive import Primitive3P, Primitive2DMixin
+from .Primitive import Primitive3P, ClosedPrimitiveMixin, PathMixin, PolygonMixin, Primitive2DMixin
 from .Line import Line2D
 
 ####################################################################################################
@@ -78,7 +78,7 @@ def same_side(p1, p2, a, b):
 
 ####################################################################################################
 
-class Triangle2D(Primitive2DMixin, Primitive3P):
+class Triangle2D(Primitive2DMixin, ClosedPrimitiveMixin, PathMixin, PolygonMixin, Primitive3P):
 
     """Class to implements 2D Triangle."""
 
@@ -215,7 +215,6 @@ class Triangle2D(Primitive2DMixin, Primitive3P):
 
     @property
     def is_isosceles(self):
-
         self._cache_length()
         # two sides of equal length
         return not(self.is_equilateral) and not(self.is_scalene)
@@ -224,7 +223,6 @@ class Triangle2D(Primitive2DMixin, Primitive3P):
 
     @property
     def is_right(self):
-
         self._cache_angle()
         # one angle = 90
         raise NotImplementedError
@@ -233,7 +231,6 @@ class Triangle2D(Primitive2DMixin, Primitive3P):
 
     @property
     def is_obtuse(self):
-
         self._cache_angle()
         # one angle > 90
         return max(self._a10, self._a21, self._a02) > 90
@@ -242,7 +239,6 @@ class Triangle2D(Primitive2DMixin, Primitive3P):
 
     @property
     def is_acute(self):
-
         self._cache_angle()
         # all angle < 90
         return max(self._a10, self._a21, self._a02) < 90
@@ -251,7 +247,6 @@ class Triangle2D(Primitive2DMixin, Primitive3P):
 
     @property
     def is_oblique(self):
-
         return not self.is_equilateral
 
     ##############################################