Repository: pyfar/pyfar (python)

Issue:
Arithmetic operations with pf.Coordinates 
It would be nice to be able to use arithmetic operations with pf.Coordinates, similar to audio objects. 
ideas would be:
- [x] addition of 2 Coordinates objects
- [x] substraction of 2 Coordinates objects
- [x] dot product
- [x] cross product

Pull request:
add arithmetics for coordinates
### Which issue(s) are closed by this pull request?

Closes #632

### Changes proposed in this pull request:

- [x] addition of 2 Coordinates objects
- [x] substraction of 2 Coordinates objects
- [x] dot product
- [x] cross product
- [x] multiply and devide with numbers
- [x] add doc at the beginning


Patch:
diff --git a/pyfar/__init__.py b/pyfar/__init__.py
index f205bb0b..08912f95 100644
--- a/pyfar/__init__.py
+++ b/pyfar/__init__.py
@@ -11,7 +11,7 @@ from .classes.audio import Signal, TimeData, FrequencyData
 from .classes.audio import (add, subtract, multiply, divide, power,
                             matrix_multiplication)
 from .classes.coordinates import Coordinates
-from .classes.coordinates import (deg2rad, rad2deg)
+from .classes.coordinates import (deg2rad, rad2deg, dot, cross)
 from .classes.orientations import Orientations
 from .classes.filter import FilterFIR, FilterIIR, FilterSOS
 
@@ -45,4 +45,7 @@ __all__ = [
     'io',
     'dsp',
     'signals',
-    'utils']
+    'utils',
+    'dot',
+    'cross',
+    ]
diff --git a/pyfar/classes/coordinates.py b/pyfar/classes/coordinates.py
index ce2f9889..a93408e9 100644
--- a/pyfar/classes/coordinates.py
+++ b/pyfar/classes/coordinates.py
@@ -10,8 +10,12 @@ Different coordinate systems are frequently used in acoustics research and
 handling sampling points and different systems can be cumbersome. The
 Coordinates class was designed with this in mind. It stores coordinates in
 cartesian coordinates internally and can convert to all coordinate systems
-listed below. Additionally, the the class can  query and plot coordinates
-points. Functions for converting coordinates not stored in a Coordinates object
+listed below. Additionally, the class can query and plot coordinates
+points. Addition and subtraction are supported with numbers and Coordinates
+objects, while multiplication and division are supported with numbers only.
+All arithmetic operations are performed element-wise on Cartesian coordinates
+using the appropriate operator.
+Functions for converting coordinates not stored in a Coordinates object
 are available for convenience. However, it is strongly recommended to use the
 Coordinates class for all conversions.
 
@@ -2650,6 +2654,11 @@ class Coordinates():
         y = np.broadcast_to(y, shapes)
         z = np.broadcast_to(z, shapes)
 
+        # set writeable, to make sure that the class does not become read-only
+        x.setflags(write=True)
+        y.setflags(write=True)
+        z.setflags(write=True)
+
         # set values
         self._x = x
         self._y = y
@@ -2807,12 +2816,195 @@ class Coordinates():
         return (eq_x & eq_y & eq_z).all() & eq_weights & eq_comment \
             & eq_sh_order & eq_system
 
+    def __add__(self, other):
+        """Add two numbers/Coordinates objects."""
+        return _arithmetics(self, other, 'add')
+
+    def __radd__(self, other):
+        """Add two numbers/Coordinates objects."""
+        return _arithmetics(other, self, 'add')
+
+    def __sub__(self, other):
+        """Subtract two numbers/Coordinates objects."""
+        return _arithmetics(self, other, 'sub')
+
+    def __rsub__(self, other):
+        """Subtract two numbers/Coordinates objects."""
+        return _arithmetics(other, self, 'sub')
+
+    def __mul__(self, other):
+        """Multiply Coordinates object with number."""
+        return _arithmetics(self, other, 'mul')
+
+    def __rmul__(self, other):
+        """Multiply number with Coordinates object."""
+        return _arithmetics(other, self, 'mul')
+
+    def __div__(self, other):
+        """Divide Coordinates object with number."""
+        return _arithmetics(self, other, 'div')
+
+    def __truediv__(self, other):
+        """Divide Coordinates object with number."""
+        return _arithmetics(self, other, 'div')
+
+    def __rtruediv__(self, other):
+        """Divide number with Coordinates object."""
+        return _arithmetics(other, self, 'div')
+
+    def __rdiv__(self, other):
+        """Divide number with Coordinates object."""
+        return _arithmetics(other, self, 'div')
+
+
     def _check_empty(self):
         """check if object is empty"""
         if self.cshape == (0,):
             raise ValueError('Object is empty.')
 
 
+def dot(a, b):
+    """Dot product of two Coordinates objects.
+
+    .. math::
+        \\vec{a} \\cdot \\vec{b}
+        = a_x \\cdot b_x + a_y \\cdot b_y + a_z \\cdot b_z
+
+    Parameters
+    ----------
+    a : pf.Coordinates
+        first argument, must be broadcastable with b
+    b : pf.Coordinates
+        second argument, much be broadcastable with a
+
+    Returns
+    -------
+    result : np.ndarray
+        array with the dot product of the two objects
+
+    Examples
+    --------
+    >>> import pyfar as pf
+    >>> a = pf.Coordinates(1, 0, 0)
+    >>> b = pf.Coordinates(1, 0, 0)
+    >>> pf.dot(a, b)
+    array([1.])
+    """
+
+    if not isinstance(a, Coordinates) or not isinstance(b, Coordinates):
+        raise TypeError(
+            "Dot product is only possible with Coordinates objects.")
+
+    return a.x * b.x + a.y * b.y + a.z * b.z
+
+
+def cross(a, b):
+    """Cross product of two Coordinates objects
+
+    .. math::
+        \\vec{a} \\times \\vec{b}
+        = (a_y \\cdot b_z - a_z \\cdot b_y) \\cdot \\hat{x}
+        + (a_z \\cdot b_x - a_x \\cdot b_z) \\cdot \\hat{y}
+        + (a_x \\cdot b_y - a_y \\cdot b_x) \\cdot \\hat{z}
+
+    Parameters
+    ----------
+    a : pf.Coordinates
+        first argument, must be broadcastable with b
+    b : pf.Coordinates
+        second argument, much be broadcastable with a
+
+    Returns
+    -------
+    result : pf.Coordinates
+        new Coordinates object with the cross product of the two objects
+
+    Examples
+    --------
+    >>> import pyfar as pf
+    >>> a = pf.Coordinates(1, 0, 0)
+    >>> b = pf.Coordinates(0, 1, 0)
+    >>> result = pf.cross(a, b)
+    >>> result.cartesian
+    array([0., 0., 1.])
+    """
+
+    if not isinstance(a, Coordinates) or not isinstance(b, Coordinates):
+        raise TypeError(
+            "Dot product is only possible with Coordinates objects.")
+
+    new = Coordinates()
+    new.cartesian = np.zeros(np.broadcast_shapes(
+        a.cartesian.shape, b.cartesian.shape))
+
+    # apply cross product
+    new.x = a.y * b.z - a.z * b.y
+    new.y = a.z * b.x - a.x * b.z
+    new.z = a.x * b.y - a.y * b.x
+
+    return new
+
+def _arithmetics(first, second, operation):
+    """Add or Subtract two Coordinates objects, numbers or arrays.
+
+    Parameters
+    ----------
+    first : Coordinates, number, array
+        first operand
+    second : Coordinates, number, array
+        second operand
+    operation : 'add', 'sub', 'mul', 'div'
+        whether to add or subtract the two objects
+
+    Returns
+    -------
+    new : Coordinates
+        result of the operation
+
+    """
+    # convert data
+    data = []
+    num_objects = 0
+    for obj in [first, second]:
+        if isinstance(obj, Coordinates):
+            data.append(obj.cartesian)
+            num_objects += 1
+        elif isinstance(obj, (int, float)):
+            data.append(np.array(obj))
+        else:
+            if operation == 'add':
+                op = 'Addition'
+            elif operation == 'sub':
+                op = 'Subtraction'
+            elif operation == 'mul':
+                op = 'Multiplication'
+            elif operation == 'div':
+                op = 'Division'
+            raise TypeError(
+                f"{op} is only possible with Coordinates or number.")
+
+    if operation in ['mul', 'div'] and num_objects > 1:
+        raise TypeError(
+            "Multiplication and division are only possible with one "
+            "Coordinates object.")
+
+    # broadcast shapes
+    shape = np.broadcast_shapes(data[0].shape, data[1].shape)
+    new = pf.Coordinates()
+    new.cartesian = np.zeros(shape)
+
+    # perform operation
+    if operation == 'add':
+        new.cartesian = data[0] + data[1]
+    elif operation == 'sub':
+        new.cartesian = data[0] - data[1]
+    elif operation == 'mul':
+        new.cartesian = data[0] * data[1]
+    elif operation == 'div':
+        new.cartesian = data[0] / data[1]
+    return new
+
+
 def cart2sph(x, y, z):
     """
     Transforms from Cartesian to spherical coordinates.


Test patch:
diff --git a/tests/test_coordiantes_arithmetic.py b/tests/test_coordiantes_arithmetic.py
new file mode 100644
index 00000000..8187699f
--- /dev/null
+++ b/tests/test_coordiantes_arithmetic.py
@@ -0,0 +1,144 @@
+import pytest
+import pyfar as pf
+import numpy.testing as npt
+import numpy as np
+
+
+@pytest.mark.parametrize(
+    "other",
+    (
+        # 1d coords
+        pf.Coordinates(1, 1, 1),
+        # 2d coords
+        pf.Coordinates([1, 1], [1, 1], [1, 1]),
+        # numbers
+        1,
+        1.0,
+    ),
+)
+@pytest.mark.parametrize(
+    "operator",
+    (
+        "-",
+        "+",
+    ),
+)
+@pytest.mark.parametrize(
+    "order",
+    (
+        True,
+        False,
+    ),
+)
+def test_arithmetic_coordinates(other, operator, order):
+    coords = pf.Coordinates([0, 1], [0, 1], [0, 1])
+    part1 = coords if order else other
+    part2 = other if order else coords
+    if operator == "+":
+        new = part1 + part2
+        desired = pf.Coordinates([1, 2], [1, 2], [1, 2])
+    elif operator == "-":
+        new = part1 - part2
+        if order:
+            desired = pf.Coordinates([-1, 0], [-1, 0], [-1, 0])
+        else:
+            desired = pf.Coordinates([1, 0], [1, 0], [1, 0])
+    assert isinstance(new, pf.Coordinates)
+    npt.assert_array_equal(new.cartesian, desired.cartesian)
+
+
+
+@pytest.mark.parametrize(
+    "other",
+    (
+        5,
+        5.0,
+    ),
+)
+@pytest.mark.parametrize(
+    "operator",
+    (
+        "*",
+        "/",
+    ),
+)
+@pytest.mark.parametrize(
+    "order",
+    (
+        True,
+        False,
+    ),
+)
+def test_arithmetic_coordinates_mul_div(other, operator, order):
+    coords = pf.Coordinates([0, 1], [0, 1], [0, 1])
+    part1 = coords if order else other
+    part2 = other if order else coords
+    if operator == "*":
+        new = part1 * part2
+        desired = pf.Coordinates([0, 5], [0, 5], [0, 5])
+    elif operator == "/":
+        new = part1 / part2
+        if order:
+            desired = pf.Coordinates([0, .2], [0, .2], [0, .2])
+        else:
+            desired = pf.Coordinates(
+                [np.inf, 5], [np.inf, 5], [np.inf, 5])
+    assert isinstance(new, pf.Coordinates)
+    npt.assert_array_equal(new.cartesian, desired.cartesian)
+
+
+@pytest.mark.parametrize(
+    "other",
+    (
+        "wrong",
+        1 + 1j * 1,
+        np.array([1, 1]),
+    ),
+)
+@pytest.mark.parametrize(
+    "operator",
+    (
+        "-",
+        "+",
+        "dot",
+        "cross",
+    ),
+)
+def test_arithmetic_coordinates_error(other, operator):
+    coords = pf.Coordinates([0, 1], [0, 1], [0, 1])
+    if operator == "+":
+        match = "Addition is only possible with Coordinates or number."
+        with pytest.raises(TypeError, match=match):
+            coords + other
+    elif operator == "-":
+        match = "Subtraction is only possible with Coordinates or number."
+        with pytest.raises(TypeError, match=match):
+            coords - other
+    elif operator == "dot":
+        match = "Dot product is only possible with Coordinates."
+        with pytest.raises(TypeError, match=match):
+            pf.dot(coords, other)
+    elif operator == "cross":
+        match = "Dot product is only possible with Coordinates."
+        with pytest.raises(TypeError, match=match):
+            pf.cross(coords, other)
+
+
+@pytest.mark.parametrize(
+    "other",
+    (
+        pf.Coordinates(1, 1, 1),
+        pf.Coordinates([1, 1], [1, 1], [1, 1]),
+    ),
+)
+def test_dot_product(other):
+    coords = pf.Coordinates(1, 1, 1)
+    dot = pf.dot(coords, other)
+    npt.assert_array_equal(dot, 3 * other.x)
+
+
+def test_cross_product():
+    coords_1 = pf.Coordinates(1, 0, 0)
+    coords_2 = pf.Coordinates(0, 1, 0)
+    cross = pf.cross(coords_1, coords_2)
+    npt.assert_array_equal(cross.cartesian, [[0, 0, 1]])