diff --git a/splipy/io/g2.py b/splipy/io/g2.py
index 536909f..9657877 100644
--- a/splipy/io/g2.py
+++ b/splipy/io/g2.py
@@ -1,11 +1,18 @@
+from __future__ import annotations
+
+from typing import ClassVar, Callable, TextIO, Union, Optional, Type, Literal, Sequence
+from types import TracebackType
+from pathlib import Path
 
 import numpy as np
 from numpy import pi, savetxt
+from typing_extensions import Self
 
 from ..curve import Curve
 from ..surface import Surface
 from ..volume import Volume
 from ..splineobject import SplineObject
+from ..splinemodel import SplineModel
 from ..basis import BSplineBasis
 from ..trimmedsurface import TrimmedSurface
 from ..utils import flip_and_move_plane_geometry, rotate_local_x_axis
@@ -16,13 +23,106 @@
 
 class G2(MasterIO):
 
-    def read_next_non_whitespace(self):
+    fstream: TextIO
+    filename: str
+    mode: Literal['w', 'r']
+    trimming_curves: list[TrimmedSurface]
+
+    g2_type: ClassVar[list[int]] = [100, 200, 700] # curve, surface, volume identifiers
+    g2_generators: ClassVar[dict[int, str]] = {
+        120: 'line',
+        130: 'circle',
+        140: 'ellipse',
+        260: 'cylinder',
+        292: 'disc',
+        270: 'sphere',
+        290: 'torus',
+        250: 'plane',
+        210: 'bounded_surface',
+        261: 'surface_of_linear_extrusion',
+    } #, 280:cone
+
+    def __init__(self, filename: Union[Path, str], mode: Literal["w", "r"] = 'r') -> None:
+        self.filename = str(filename)
+        self.trimming_curves = []
+        self.mode = mode
+
+    def __enter__(self) -> Self:
+        self.fstream = open(self.filename, self.mode).__enter__()
+        return self
+
+    def __exit__(
+        self,
+        exc_type: Optional[Type[BaseException]],
+        exc_val: Optional[BaseException],
+        exc_tb: Optional[TracebackType]
+    ) -> None:
+        self.fstream.__exit__(exc_type, exc_val, exc_tb)
+
+    def _write_obj(self, obj: SplineObject) -> None:
+        for i in range(obj.pardim):
+            if obj.periodic(i):
+                obj = obj.split_periodic(obj.start(i), i)
+
+        self.fstream.write('{} 1 0 0\n'.format(G2.g2_type[obj.pardim-1]))
+        self.fstream.write('{} {}\n'.format(obj.dimension, int(obj.rational)))
+        for b in obj.bases:
+            self.fstream.write('%i %i\n' % (len(b.knots) - b.order, b.order))
+            self.fstream.write(' '.join('%.16g' % k for k in b.knots))
+            self.fstream.write('\n')
+
+        savetxt(self.fstream, obj.controlpoints.reshape(-1, obj.dimension + obj.rational, order='F'),
+                fmt='%.16g', delimiter=' ', newline='\n')
+
+    def write(self, obj: Union[Sequence[SplineObject], SplineObject, SplineModel]) -> None:
+        """Write the object in GoTools format."""
+        if isinstance(obj, (Sequence, SplineModel)): # input SplineModel or list
+            for o in obj:
+                self._write_obj(o)
+            return
+
+        assert isinstance(obj, SplineObject)
+        self._write_obj(obj)
+
+    def read(self) -> list[SplineObject]:
+        result = []
+
+        for line in self.fstream:
+            line = line.strip()
+            if not line:
+                continue
+
+            # read object type
+            objtype, major, minor, patch = map(int, line.split())
+            if (major, minor, patch) != (1, 0, 0):
+                raise IOError('Unknown G2 format')
+
+            # if obj type is in factory methods (cicle, torus etc), create it now
+            if objtype in G2.g2_generators:
+                constructor = getattr(self, G2.g2_generators[objtype])
+                result.append(constructor())
+                continue
+
+            # for "normal" splines (Curves, Surfaces, Volumes) create it now
+            pardim = [i for i in range(len(G2.g2_type)) if G2.g2_type[i] == objtype]
+            if not pardim:
+                raise IOError('Unknown G2 object type {}'.format(objtype))
+            result.append(self.splines(pardim[0] + 1))
+
+        return result
+
+    def read_basis(self) -> BSplineBasis:
+        ncps, order = map(int, next(self.fstream).split())
+        kts = list(map(float, next(self.fstream).split()))
+        return BSplineBasis(order, kts, -1)
+
+    def read_next_non_whitespace(self) -> str:
         line = next(self.fstream).strip()
         while not line:
             line = next(self.fstream).strip()
         return line
 
-    def circle(self):
+    def circle(self) -> Curve:
         self.read_next_non_whitespace()
         # dim   = int(     self.read_next_non_whitespace().strip())
         r     = float(   next(self.fstream).strip())
@@ -38,7 +138,7 @@ def circle(self):
             result.reverse()
         return result
 
-    def ellipse(self):
+    def ellipse(self) -> Curve:
         self.read_next_non_whitespace()
         # dim   = int(     self.read_next_non_whitespace().strip())
         r1    = float(   next(self.fstream).strip())
@@ -55,7 +155,7 @@ def ellipse(self):
             result.reverse()
         return result
 
-    def line(self):
+    def line(self) -> Curve:
         self.read_next_non_whitespace()
         # dim   = int(     self.read_next_non_whitespace().strip())
         start    = np.array(next(self.fstream).split(), dtype=float)
@@ -67,14 +167,13 @@ def line(self):
         s = np.array(start)
         # d /= np.linalg.norm(d)
         if not finite:
-            param = [-state.unlimited, +state.unlimited]
+            param = np.array([-state.unlimited, +state.unlimited])
 
         result = curve_factory.line(s+d*param[0], s+d*param[1])
         if reverse:
             result.reverse()
         return result
 
-
 #   def cone(self):
 #       dim      = int(     self.read_next_non_whitespace().strip())
 #       r        = float(   next(self.fstream).strip())
@@ -87,8 +186,7 @@ def line(self):
 #       if finite:
 #           param_v=np.array(next(self.fstream).split(' '), dtype=float)
 
-
-    def cylinder(self):
+    def cylinder(self) -> Surface:
         self.read_next_non_whitespace()
         # dim   = int(     self.read_next_non_whitespace().strip())
         r        = float(   next(self.fstream).strip())
@@ -98,9 +196,9 @@ def cylinder(self):
         finite   =          next(self.fstream).strip() != '0'
         param_u  = np.array(next(self.fstream).split(), dtype=float)
         if finite:
-            param_v=np.array(next(self.fstream).split(), dtype=float)
+            param_v = np.array(next(self.fstream).split(), dtype=float)
         else:
-            param_v=[-state.unlimited, state.unlimited]
+            param_v = np.array([-state.unlimited, state.unlimited])
         swap     =          next(self.fstream).strip() != '0'
 
         center = center + z_axis*param_v[0]
@@ -111,7 +209,7 @@ def cylinder(self):
             result.swap()
         return result
 
-    def disc(self):
+    def disc(self) -> Surface:
         self.read_next_non_whitespace()
         # dim      = int(     self.read_next_non_whitespace().strip())
         center   = np.array(next(self.fstream).split(), dtype=float)
@@ -135,7 +233,7 @@ def disc(self):
             result.swap()
         return result
 
-    def plane(self):
+    def plane(self) -> Surface:
         self.read_next_non_whitespace()
         # dim        = int(     self.read_next_non_whitespace().strip())
         center     = np.array(next(self.fstream).split(), dtype=float)
@@ -143,22 +241,22 @@ def plane(self):
         x_axis     = np.array(next(self.fstream).split(), dtype=float)
         finite     =          next(self.fstream).strip() != '0'
         if finite:
-            param_u= np.array(next(self.fstream).split(), dtype=float)
-            param_v= np.array(next(self.fstream).split(), dtype=float)
+            param_u = np.array(next(self.fstream).split(), dtype=float)
+            param_v = np.array(next(self.fstream).split(), dtype=float)
         else:
-            param_u= [-state.unlimited, +state.unlimited]
-            param_v= [-state.unlimited, +state.unlimited]
+            param_u = np.array([-state.unlimited, +state.unlimited])
+            param_v = np.array([-state.unlimited, +state.unlimited])
         swap       =          next(self.fstream).strip() != '0'
 
         result = Surface() * [param_u[1]-param_u[0], param_v[1]-param_v[0]] + [param_u[0],param_v[0]]
         result.rotate(rotate_local_x_axis(x_axis, normal))
         result = flip_and_move_plane_geometry(result,center,normal)
         result.reparam(param_u, param_v)
-        if(swap):
+        if swap:
             result.swap()
         return result
 
-    def torus(self):
+    def torus(self) -> Surface:
         self.read_next_non_whitespace()
         # dim      = int(     self.read_next_non_whitespace().strip())
         r2       = float(   next(self.fstream).strip())
@@ -178,7 +276,7 @@ def torus(self):
             result.swap()
         return result
 
-    def sphere(self):
+    def sphere(self) -> Surface:
         self.read_next_non_whitespace()
         # dim      = int(     self.read_next_non_whitespace().strip())
         r        = float(   next(self.fstream).strip())
@@ -195,11 +293,9 @@ def sphere(self):
         result.reparam(param_u, param_v)
         return result
 
-    def splines(self, pardim):
-        cls    = G2.classes[pardim-1]
-
+    def splines(self, pardim: int) -> SplineObject:
         _, rational = self.read_next_non_whitespace().strip().split()
-        rational = bool(int(rational))
+        rational_bool = bool(int(rational))
 
         bases = [self.read_basis() for _ in range(pardim)]
         ncps = 1
@@ -209,10 +305,9 @@ def splines(self, pardim):
         cps = [tuple(map(float, next(self.fstream).split()))
                for _ in range(ncps)]
 
-        args = bases + [cps, rational]
-        return cls(*args)
+        return SplineObject.constructor(pardim)(bases, cps, rational=rational_bool)
 
-    def surface_of_linear_extrusion(self):
+    def surface_of_linear_extrusion(self) -> Surface:
         self.read_next_non_whitespace()
         # dim      = int(      self.read_next_non_whitespace().strip())
         crv      = self.splines(1)
@@ -220,9 +315,9 @@ def surface_of_linear_extrusion(self):
         finite   =          next(self.fstream).strip() != '0'
         param_u  = np.array(next(self.fstream).split(), dtype=float)
         if finite:
-            param_v=np.array(next(self.fstream).split(), dtype=float)
+            param_v = np.array(next(self.fstream).split(), dtype=float)
         else:
-            param_v=[-state.unlimited, +state.unlimited]
+            param_v = np.array([-state.unlimited, +state.unlimited])
         swap     =          next(self.fstream).strip() != '0'
 
         result = surface_factory.extrude(crv + normal*param_v[0], normal*(param_v[1]-param_v[0]))
@@ -232,13 +327,12 @@ def surface_of_linear_extrusion(self):
             result.swap()
         return result
 
-
-    def bounded_surface(self):
-        objtype = int( next(self.fstream).strip() )
+    def bounded_surface(self) -> TrimmedSurface:
+        objtype = int(next(self.fstream).strip())
 
         # create the underlying surface which all trimming curves are to be applied
         if objtype in G2.g2_generators:
-            constructor = getattr(self, G2.g2_generators[objtype].__name__)
+            constructor = getattr(self, G2.g2_generators[objtype])
             surface = constructor()
         elif objtype == 200:
             surface = self.splines(2)
@@ -261,7 +355,7 @@ def bounded_surface(self):
                 two_curves = []
                 for crv_type in [parameter_curve_type, space_curve_type]:
                     if crv_type in G2.g2_generators:
-                        constructor = getattr(self, G2.g2_generators[crv_type].__name__)
+                        constructor = getattr(self, G2.g2_generators[crv_type])
                         crv = constructor()
                     elif crv_type == 100:
                         crv = self.splines(1)
@@ -275,89 +369,3 @@ def bounded_surface(self):
             all_loops.append(one_loop)
 
         return TrimmedSurface(surface.bases[0], surface.bases[1], surface.controlpoints, surface.rational, all_loops, raw=True)
-
-    g2_type = [100, 200, 700] # curve, surface, volume identifiers
-    classes = [Curve, Surface, Volume]
-    g2_generators = {120:line, 130:circle, 140:ellipse,
-                     260:cylinder, 292:disc, 270:sphere, 290:torus, 250:plane,
-                     210:bounded_surface, 261:surface_of_linear_extrusion} #, 280:cone
-
-    def __init__(self, filename):
-        if filename[-3:] != '.g2':
-            filename += '.g2'
-        self.filename = filename
-        self.trimming_curves = []
-
-    def __enter__(self):
-        return self
-
-    def write(self, obj):
-        if not hasattr(self, 'fstream'):
-            self.onlywrite = True
-            self.fstream   = open(self.filename, 'w')
-        if not self.onlywrite:
-            raise IOError('Could not write to file %s' % (self.filename))
-
-        """Write the object in GoTools format. """
-        if isinstance(obj[0], SplineObject): # input SplineModel or list
-            for o in obj:
-                self.write(o)
-            return
-
-        for i in range(obj.pardim):
-            if obj.periodic(i):
-                obj = obj.split(obj.start(i), i)
-
-        self.fstream.write('{} 1 0 0\n'.format(G2.g2_type[obj.pardim-1]))
-        self.fstream.write('{} {}\n'.format(obj.dimension, int(obj.rational)))
-        for b in obj.bases:
-            self.fstream.write('%i %i\n' % (len(b.knots) - b.order, b.order))
-            self.fstream.write(' '.join('%.16g' % k for k in b.knots))
-            self.fstream.write('\n')
-
-        savetxt(self.fstream, obj.controlpoints.reshape(-1, obj.dimension + obj.rational, order='F'),
-                fmt='%.16g', delimiter=' ', newline='\n')
-
-    def read(self):
-        if not hasattr(self, 'fstream'):
-            self.onlywrite = False
-            self.fstream   = open(self.filename, 'r')
-
-        if self.onlywrite:
-            raise IOError('Could not read from file %s' % (self.filename))
-
-        result = []
-
-        for line in self.fstream:
-            line = line.strip()
-            if not line:
-                continue
-
-            # read object type
-            objtype, major, minor, patch = map(int, line.split())
-            if (major, minor, patch) != (1, 0, 0):
-                raise IOError('Unknown G2 format')
-
-            # if obj type is in factory methods (cicle, torus etc), create it now
-            if objtype in G2.g2_generators:
-                constructor = getattr(self, G2.g2_generators[objtype].__name__)
-                result.append( constructor() )
-                continue
-
-            # for "normal" splines (Curves, Surfaces, Volumes) create it now
-            pardim = [i for i in range(len(G2.g2_type)) if G2.g2_type[i] == objtype]
-            if not pardim:
-                raise IOError('Unknown G2 object type {}'.format(objtype))
-            pardim = pardim[0] + 1
-            result.append(self.splines(pardim))
-
-        return result
-
-    def read_basis(self):
-        ncps, order = map(int, next(self.fstream).split())
-        kts = list(map(float, next(self.fstream).split()))
-        return BSplineBasis(order, kts, -1)
-
-    def __exit__(self, exc_type, exc_value, traceback):
-        if hasattr(self, 'fstream'):
-            self.fstream.close()
diff --git a/splipy/io/master.py b/splipy/io/master.py
index 58e16c8..21c2654 100644
--- a/splipy/io/master.py
+++ b/splipy/io/master.py
@@ -1,27 +1,47 @@
-class MasterIO(object):
+from abc import ABC, abstractmethod
+from typing import Union, Sequence, Type, Optional
+from types import TracebackType
+from pathlib import Path
 
-    def __init__(self, filename):
+from typing_extensions import Self
+
+from ..splineobject import SplineObject
+
+
+class MasterIO:
+
+    def __init__(self, filename: Union[str, Path]) -> None:
         """Create an IO object attached to a file.
 
         :param str filename: The file to read from or write to
         """
         raise NotImplementedError()
 
-    def __enter__(self):
-        raise NotImplementedError()
+    @abstractmethod
+    def __enter__(self) -> Self:
+        ...
+
+    @abstractmethod
+    def __exit__(
+        self,
+        exc_type: Optional[Type[BaseException]],
+        exc_val: Optional[BaseException],
+        exc_tb: Optional[TracebackType]
+    ) -> None:
+        ...
 
-    def write(self, obj):
+    @abstractmethod
+    def write(self, obj: Union[SplineObject, Sequence[SplineObject]]) -> None:
         """Write one or more objects to the file.
 
         :param obj: The object(s) to write
         :type obj: [:class:`splipy.SplineObject`] or :class:`splipy.SplineObject`
         """
-        raise NotImplementedError()
 
-    def read(self):
+    @abstractmethod
+    def read(self) -> list[SplineObject]:
         """Reads all the objects from the file.
 
         :return: Objects
         :rtype: [:class:`splipy.SplineObject`]
         """
-        raise NotImplementedError()
diff --git a/splipy/io/svg.py b/splipy/io/svg.py
index ef5d361..5b498f4 100644
--- a/splipy/io/svg.py
+++ b/splipy/io/svg.py
@@ -1,38 +1,43 @@
 import xml.etree.ElementTree as etree
 from xml.dom import minidom
 import re
+from typing import ClassVar, Optional, Type, Union, Sequence, cast
+from types import TracebackType
+from pathlib import Path
 
 import numpy as np
+from typing_extensions import Self
 
 from ..curve import Curve
 from ..surface import Surface
 from ..splineobject import SplineObject
 from ..basis import BSplineBasis
+from ..types import FArray
 from .. import curve_factory, state
 
 from .master import MasterIO
 
 
-def read_number_and_unit(mystring):
+def read_number_and_unit(mystring: str) -> tuple[float, str]:
     unit = ''
     try:
         for i in range(1, len(mystring)+1):
-            number=float(mystring[:i])
+            number = float(mystring[:i])
     except ValueError:
-        unit=mystring[i-1:]
+        unit = mystring[i-1:]
     return (number, unit)
 
 
-def bezier_representation(curve):
-    """ Compute a Bezier representation of a given spline curve. The input
-        curve must be of order less than or equal to 4. The bezier
-        representation is of order 4, and maximal knot multiplicity, i.e.
-        consist of a C0-discretization with knot multiplicity 3.
+def bezier_representation(curve: Curve) -> Curve:
+    """Compute a Bezier representation of a given spline curve. The input
+    curve must be of order less than or equal to 4. The bezier
+    representation is of order 4, and maximal knot multiplicity, i.e.
+    consist of a C0-discretization with knot multiplicity 3.
 
-        :param curve  : Spline curve
-        :type  curve  : Curve
-        :returns      : Bezier curve
-        :rtype        : Curve
+    :param curve  : Spline curve
+    :type  curve  : Curve
+    :returns      : Bezier curve
+    :rtype        : Curve
     """
     # error test input. Another way would be to approximate higher-order curves. Consider looking at Curve.rebuild()
     if curve.order(0) > 4 or curve.rational:
@@ -43,7 +48,7 @@ def bezier_representation(curve):
 
     # make it non-periodic
     if bezier.periodic():
-        bezier = bezier.split(bezier.start(0))
+        bezier = bezier.split_periodic(bezier.start(0))
 
     # make sure it is C0 everywhere
     for k in bezier.knots(0):
@@ -53,23 +58,32 @@ def bezier_representation(curve):
 
 
 class SVG(MasterIO):
-
-    namespace = '{http://www.w3.org/2000/svg}'
-
-    def __init__(self, filename, width=1000, height=1000, margin=0.05):
-        """  Constructor
-          :param filename: Filename to write results to
-          :type  filename: String
-          :param width   : Maximum image width in pixels
-          :type  width   : Int
-          :param height  : Maximum image height in pixels
-          :type  height  : Int
-          :param margin  : White-space around all edges of image, given in percentage of total size (default 5%)
-          :type  margin  : Float
+    namespace: ClassVar[str] = '{http://www.w3.org/2000/svg}'
+
+    filename: str
+    width: float
+    height: float
+    margin: float
+    all_objects: list[SplineObject]
+
+    scale: float
+    center: tuple[float, float]
+    offset: tuple[float, float]
+    xmlRoot: etree.Element
+
+    def __init__(self, filename: Union[Path, str], width: int = 1000, height: int = 1000, margin: float = 0.05) -> None:
+        """Constructor
+
+        :param filename: Filename to write results to
+        :type  filename: String
+        :param width   : Maximum image width in pixels
+        :type  width   : Int
+        :param height  : Maximum image height in pixels
+        :type  height  : Int
+        :param margin  : White-space around all edges of image, given in percentage of total size (default 5%)
+        :type  margin  : Float
         """
-        if filename[-4:] != '.svg':
-            filename += '.svg'
-        self.filename = filename
+        self.filename = str(filename)
 
         self.width  = width
         self.height = height
@@ -77,14 +91,18 @@ def __init__(self, filename, width=1000, height=1000, margin=0.05):
 
         self.all_objects = []
 
-    def __enter__(self):
+    def __enter__(self) -> Self:
         return self
 
-    def __exit__(self, exc_type, exc_value, traceback):
+    def __exit__(
+        self,
+        exc_type: Optional[Type[BaseException]],
+        exc_val: Optional[BaseException],
+        exc_tb: Optional[TracebackType]
+    ) -> None:
         # in case something goes wrong, print error and abandon writing process
         if exc_type is not None:
-            print(exc_type, exc_value, traceback)
-            return False
+            print(exc_type, exc_val, exc_tb)
 
         # compute the bounding box for all geometries
         boundingbox = [np.inf, np.inf, -np.inf, -np.inf]
@@ -97,7 +115,7 @@ def __exit__(self, exc_type, exc_value, traceback):
 
         # compute scaling factors by keeping aspect ratio, and never exceed width or height size (including margins)
         geometryRatio = float(boundingbox[3]-boundingbox[1])/(boundingbox[2]-boundingbox[0])
-        imageRatio    = 1.0*self.height / self.width
+        imageRatio    = 1.0 * self.height / self.width
         if geometryRatio > imageRatio: # scale by y-coordinate
             marginPixels = self.height*self.margin
             self.scale   = self.height*(1-2*self.margin) / (boundingbox[3]-boundingbox[1])
@@ -106,8 +124,8 @@ def __exit__(self, exc_type, exc_value, traceback):
             marginPixels = self.width*self.margin
             self.scale   = self.width*(1-2*self.margin) / (boundingbox[2]-boundingbox[0])
             self.height  = self.width*geometryRatio + 2*marginPixels
-        self.center      = [boundingbox[0], boundingbox[1]]
-        self.offset      = [marginPixels, marginPixels]
+        self.center      = (boundingbox[0], boundingbox[1])
+        self.offset      = (marginPixels, marginPixels)
 
         # create xml root tag
         self.xmlRoot = etree.Element('svg',  {'xmlns':'http://www.w3.org/2000/svg',
@@ -130,21 +148,21 @@ def __exit__(self, exc_type, exc_value, traceback):
             f = open(self.filename, 'w')
             f.write(result)
 
-    def write_curve(self, xmlNode, curve, fill='none', stroke='#000000', width=2):
-        """ Writes a Curve to the xml tree. This will draw a single curve
-
-            :param xmlNode: Node in xml tree
-            :type  xmlNode: Etree.element
-            :param curve  : The spline curve to write
-            :type  curve  : Curve
-            :param fill   : Fill color in hex or none
-            :type  fill   : String
-            :param stroke : Line color written in hex, i.e. '#000000'
-            :type  stroke : String
-            :param width  : Line width, measured in pixels
-            :type  width  : Int
-            :returns: None
-            :rtype  : NoneType
+    def write_curve(self, xmlNode: etree.Element, curve: Curve, fill: str = 'none', stroke: str = '#000000', width: int = 2) -> None:
+        """Write a Curve to the xml tree. This will draw a single curve.
+
+        :param xmlNode: Node in xml tree
+        :type  xmlNode: Etree.element
+        :param curve  : The spline curve to write
+        :type  curve  : Curve
+        :param fill   : Fill color in hex or none
+        :type  fill   : String
+        :param stroke : Line color written in hex, i.e. '#000000'
+        :type  stroke : String
+        :param width  : Line width, measured in pixels
+        :type  width  : Int
+        :returns: None
+        :rtype  : NoneType
         """
         curveNode = etree.SubElement(xmlNode, 'path')
         curveNode.attrib['style'] = 'fill:%s;stroke:%s;stroke-width:%dpx;stroke-linecap:butt;stroke-linejoin:miter;stroke-opacity:1' %(fill,stroke,width)
@@ -159,15 +177,15 @@ def write_curve(self, xmlNode, curve, fill='none', stroke='#000000', width=2):
 
         curveNode.attrib['d'] = pathString
 
-    def write_surface(self, surface, fill='#ffcc99'):
-        """ Writes a Surface to the xml tree. This will draw the surface along with all knot lines
+    def write_surface(self, surface: Surface, fill: str = '#ffcc99') -> None:
+        """Write a Surface to the xml tree. This will draw the surface along with all knot lines.
 
-            :param surface: The spline surface to write
-            :type  surface: Surface
-            :param fill   : Surface color written in hex, i.e. '#ffcc99'
-            :type  fill   : String
-            :returns: None
-            :rtype  : NoneType
+        :param surface: The spline surface to write
+        :type  surface: Surface
+        :param fill   : Surface color written in hex, i.e. '#ffcc99'
+        :type  fill   : String
+        :returns: None
+        :rtype  : NoneType
         """
         # fetch boundary curves and create a connected, oriented bezier loop from it
         bndry_curves = surface.edges()
@@ -196,32 +214,32 @@ def write_surface(self, surface, fill='#ffcc99'):
         for meshline in knotlines:
             self.write_curve(groupNode, meshline, width=1)
 
+    def write(self, obj: Union[SplineObject, Sequence[SplineObject]]) -> None:
+        """Writes a list of planar curves and surfaces to vector graphics SVG file.
+        The image will never be stretched, and the geometry will keep width/height ratio
+        of original geometry, regardless of provided width/height ratio from arguments.
 
-    def write(self, obj):
-        """ Writes a list of planar curves and surfaces to vector graphics SVG file.
-            The image will never be stretched, and the geometry will keep width/height ratio
-            of original geometry, regardless of provided width/height ratio from arguments.
-
-            :param obj: Primitives to write
-            :type  obj: List of Curves and/or Surfaces
-            :returns: None
-            :rtype  : NoneType
+        :param obj: Primitives to write
+        :type  obj: List of Curves and/or Surfaces
+        :returns: None
+        :rtype  : NoneType
         """
         # actually this is a dummy method. It will collect all geometries provided
         # and ton't actually write them to file until __exit__ is called
 
-        if isinstance(obj[0], SplineObject): # input SplineModel or list
+        if isinstance(obj, Sequence): # input SplineModel or list
             for o in obj:
                 self.write(o)
             return
 
+        assert isinstance(obj, SplineObject)
         if obj.dimension != 2:
             raise RuntimeError('SVG files only applicable for 2D geometries')
 
         # have to clone stuff we put here, in case they change on the outside
-        self.all_objects.append( obj.clone() )
+        self.all_objects.append(obj.clone())
 
-    def read(self):
+    def read(self) -> list[SplineObject]:
         tree = etree.parse(self.filename)
         root = tree.getroot()
         parent_map = dict((c, p) for p in tree.iter() for c in p)
@@ -245,12 +263,13 @@ def read(self):
                 result.append(crv)
         return result
 
-    def transform(self, curve, operation):
+    def transform(self, curve: SplineObject, operation: str) -> None:
         # intended input operation string: 'translate(-10,-20) scale(2) rotate(45) translate(5,10)'
         all_operations = re.findall(r'[^\)]*\)', operation.lower())
         all_operations.reverse()
         for one_operation in all_operations:
             parts = re.search(r'([a-z]*)\w*\((.*)\)', one_operation.strip())
+            assert parts is not None
             func = parts.group(1)
             args = [float(d) for d in parts.group(2).split(',')]
             if func == 'translate':
@@ -278,7 +297,7 @@ def transform(self, curve, operation):
                     cp = cp @ M.T
                     curve.controlpoints = np.reshape(np.array(cp), curve.controlpoints.shape)
 
-    def curves_from_path(self, path):
+    def curves_from_path(self, path: str) -> list[SplineObject]:
         # see https://www.w3schools.com/graphics/svg_path.asp for documentation
         # and also https://www.w3.org/TR/SVG/paths.html
 
@@ -289,8 +308,8 @@ def curves_from_path(self, path):
         #     order = 3
         # else:
         #     order = 2
-        last_curve = None
-        result = []
+        last_curve: Optional[SplineObject] = None
+        result: list[SplineObject] = []
 
         # each 'piece' is an operator (M,C,Q,L etc) and accomponying list of argument points
         for piece in re.findall('[a-zA-Z][^a-zA-Z]*', path):
@@ -344,6 +363,7 @@ def curves_from_path(self, path):
                 knot = list(range(int(len(np_pts)/2)+1)) * 3
                 knot += [knot[0], knot[-1]]
                 knot.sort()
+                assert last_curve is not None
                 x0  = np.array(last_curve[-1])
                 xn1 = np.array(last_curve[-2])
                 controlpoints.append(2*x0 -xn1)
@@ -360,6 +380,7 @@ def curves_from_path(self, path):
                 knot = list(range(int(len(np_pts)/2)+1)) * 3
                 knot += [knot[0], knot[-1]]
                 knot.sort()
+                assert last_curve is not None
                 x0  = np.array(last_curve[-1])
                 xn1 = np.array(last_curve[-2])
                 controlpoints.append(2*x0 -xn1)
@@ -381,7 +402,7 @@ def curves_from_path(self, path):
             elif piece[0] == 'Q':
                 # quadratic spline, absolute position
                 controlpoints = [startpoint]
-                knot = list(range(len(np_pts)/2+1)) * 2
+                knot = list(range(len(np_pts)//2+1)) * 2
                 knot += [knot[0], knot[-1]]
                 knot.sort()
                 for cp in np_pts:
@@ -471,9 +492,16 @@ def curves_from_path(self, path):
                                        (rx**2*xp[1]**2 + ry**2*xp[0]**2)) *
                                        np.array([rx*xp[1]/ry, -ry*xp[0]/rx]))
                 center = np.linalg.solve(R.T, cprime) + (startpoint+xend)/2
-                def arccos(vec1, vec2):
-                    return (np.sign(vec1[0]*vec2[1] - vec1[1]*vec2[0]) *
-                            np.arccos(vec1.dot(vec2)/np.linalg.norm(vec1)/np.linalg.norm(vec2)))
+
+                def arccos(vec1: FArray, vec2: FArray) -> float:
+                    return cast(
+                        float,
+                        (
+                            np.sign(vec1[0]*vec2[1] - vec1[1]*vec2[0]) *
+                            np.arccos(vec1.dot(vec2)/np.linalg.norm(vec1)/np.linalg.norm(vec2))
+                        ),
+                    )
+
                 tmp1 = np.divide( xp - cprime, [rx,ry])
                 tmp2 = np.divide(-xp - cprime, [rx,ry])
                 theta1 = arccos(np.array([1,0]), tmp1)
@@ -490,6 +518,7 @@ def arccos(vec1, vec2):
                 # curve_piece = Curve(BSplineBasis(2), [startpoint, last_curve[0]])
                 # curve_piece.reparam([0, curve_piece.length()])
                 # last_curve.append(curve_piece).make_periodic(0)
+                assert last_curve is not None
                 last_curve.make_periodic(0)
                 result.append(last_curve)
                 last_curve = None
@@ -497,12 +526,13 @@ def arccos(vec1, vec2):
             else:
                 raise RuntimeError('Unknown path parameter:' + piece)
 
-            if(curve_piece.length()>state.controlpoint_absolute_tolerance):
-                curve_piece.reparam([0, curve_piece.length()])
+            if curve_piece.length() > state.controlpoint_absolute_tolerance:
+                curve_piece.reparam((0, curve_piece.length()))
                 if last_curve is None:
                     last_curve = curve_piece
                 else:
-                    last_curve.append(curve_piece)
+                    cast(Curve, last_curve).append(curve_piece)
+            assert last_curve is not None
             startpoint = last_curve[-1,:2] # disregard rational weight (if any)
 
         if last_curve is not None:
diff --git a/splipy/splinemodel.py b/splipy/splinemodel.py
index f1ffd79..46b8a25 100644
--- a/splipy/splinemodel.py
+++ b/splipy/splinemodel.py
@@ -912,7 +912,6 @@ def nodes(self, pardim: int) -> list[TopologicalNode]:
 
 # TODO: This class is unfinished, and right now it doesn't do much other than wrap ObjectCatalogue
 class SplineModel:
-
     pardim: int
     dimension: int
     force_right_hand: bool
@@ -967,6 +966,10 @@ def add(
     def __getitem__(self, obj: SplineObject) -> NodeView:
         return self.catalogue[obj]
 
+    def __iter__(self) -> Iterator[SplineObject]:
+        for node in self.catalogue.top_nodes():
+            yield node.obj
+
     def boundary(self, name: Optional[str] = None) -> Iterator[TopologicalNode]:
         for node in self.catalogue.nodes(self.pardim-1):
             if node.nhigher == 1 and (name is None or name == node.name):
@@ -1169,5 +1172,5 @@ def write(self, filename: str) -> None:
 
         # Import here to avoid circular dependencies
         from .io import G2
-        with G2(filename + '.g2') as f:
+        with G2(filename + '.g2', 'w') as f:
             f.write([n.obj for n in self.nodes])
diff --git a/splipy/splineobject.py b/splipy/splineobject.py
index 8c97235..c678c99 100644
--- a/splipy/splineobject.py
+++ b/splipy/splineobject.py
@@ -864,11 +864,11 @@ def refine(self, *ns, **kwargs):  # type: ignore[no-untyped-def]
         return self
 
     @overload
-    def reparam(self, *args: tuple[Scalar, Scalar]) -> Self:
+    def reparam(self, *args: Union[FArray, tuple[Scalar, Scalar]]) -> Self:
         ...
 
     @overload
-    def reparam(self, arg: tuple[Scalar, Scalar], /, direction: Direction) -> Self:
+    def reparam(self, arg: Union[FArray, tuple[Scalar, Scalar]], /, direction: Direction) -> Self:
         ...
 
     @overload
diff --git a/test/io/g2_test.py b/test/io/g2_test.py
index 9be7290..afcc6a9 100644
--- a/test/io/g2_test.py
+++ b/test/io/g2_test.py
@@ -1,5 +1,7 @@
 # -*- coding: utf-8 -*-
 
+from pathlib import Path
+
 from splipy.io import G2
 import splipy.surface_factory as SurfaceFactory
 import splipy.volume_factory as VolumeFactory
@@ -7,12 +9,12 @@
 import unittest
 import os
 
-THIS_DIR = os.path.dirname(os.path.abspath(__file__))
+THIS_DIR = Path(__file__).parent
 
 class TestG2(unittest.TestCase):
 
     def test_read_rational_surf(self):
-        with G2(THIS_DIR + '/geometries/torus.g2') as myfile:
+        with G2(THIS_DIR / 'geometries' / 'torus.g2') as myfile:
             surf = myfile.read()
         self.assertEqual(len(surf), 1)
         surf = surf[0]
@@ -22,7 +24,7 @@ def test_read_rational_surf(self):
     def test_write_and_read_multipatch_surface(self):
         # write teapot to file and test if its there
         teapot = SurfaceFactory.teapot()
-        with G2('teapot.g2') as myfile:
+        with G2('teapot.g2', 'w') as myfile:
             myfile.write(teapot)
         self.assertTrue(os.path.isfile('teapot.g2'))
 
@@ -37,7 +39,7 @@ def test_write_and_read_multipatch_surface(self):
         os.remove('teapot.g2')
 
     def test_read_doublespaced(self):
-        with G2(THIS_DIR + '/geometries/lshape.g2') as myfile: # controlpoints are separated by two spaces
+        with G2(THIS_DIR / 'geometries' / 'lshape.g2') as myfile: # controlpoints are separated by two spaces
             one_surf = myfile.read()
         self.assertEqual(len(one_surf), 1)
         self.assertEqual(one_surf[0].shape[0], 3)
@@ -46,7 +48,7 @@ def test_read_doublespaced(self):
     def test_write_and_read_surface(self):
         # write disc to file and test if its there
         disc = SurfaceFactory.disc(type='square')
-        with G2('disc.g2') as myfile:
+        with G2('disc.g2', 'w') as myfile:
             myfile.write(disc)
         self.assertTrue(os.path.isfile('disc.g2'))
 
@@ -64,7 +66,7 @@ def test_write_and_read_surface(self):
     def test_write_and_read_volume(self):
         # write sphere to file and test if its there
         sphere = VolumeFactory.sphere(type='square')
-        with G2('sphere.g2') as myfile:
+        with G2('sphere.g2', 'w') as myfile:
             myfile.write(sphere)
         self.assertTrue(os.path.isfile('sphere.g2'))
 
@@ -80,7 +82,7 @@ def test_write_and_read_volume(self):
         os.remove('sphere.g2')
 
     def test_read_elementary_curves(self):
-        with G2(THIS_DIR + '/geometries/elementary_curves.g2') as myfile:
+        with G2(THIS_DIR / 'geometries' / 'elementary_curves.g2') as myfile:
             my_curves = myfile.read()
 
         self.assertEqual(len(my_curves), 3)
@@ -106,7 +108,7 @@ def test_read_elementary_curves(self):
         self.assertTrue(np.allclose(line[0], [1,0,0]))
 
     def test_read_elementary_surfaces(self):
-        with G2(THIS_DIR + '/geometries/elementary_surfaces.g2') as myfile:
+        with G2(THIS_DIR / 'geometries' / 'elementary_surfaces.g2') as myfile:
             my_surfaces = myfile.read()
 
         self.assertEqual(len(my_surfaces), 4)
@@ -134,7 +136,7 @@ def test_read_elementary_surfaces(self):
     def test_from_step(self):
         # quite large nasty g2 file which contains cylinders, planes, trimming etc
 
-        with G2(THIS_DIR + '/geometries/winglet_from_step.g2') as myfile:
+        with G2(THIS_DIR / 'geometries' / 'winglet_from_step.g2') as myfile:
             my_surfaces = myfile.read()
             trim_curves = myfile.trimming_curves