Source

geoalg / geoalg / vector3d.py

#!/usr/bin/env python
#coding:utf-8
# Author:   mozman
# Created: 27.03.2010
# Purpose: 3d vector math
# module belongs to package geoalg
# Copyright (C) 2010, Manfred Moitzi
# License: GPLv3

from __future__ import division
from __future__ import print_function
from __future__ import unicode_literals
from __future__ import absolute_import

__all__= ['Vec3']

X, Y, Z = (0, 1, 2)

class _Vec3Factory():
    def __call__(self, vector):
        """ Vector creation, returns the NULLVEC for vector == (0, 0, 0).
        """
        v = _Vec3(vector)
        return v if not v.isNULL else NULLVEC

    @staticmethod
    def from_points(p1, p2):
        x = p2[X] - p1[X]
        y = p2[Y] - p1[Y]
        z = p2[Z] - p1[Z]
        return Vec3((x, y, z))

    @staticmethod
    def from_floats(x, y, z):
        return Vec3((x, y, z))

# Vec3 is a factory used like a class
# like: v = Vec3((1, 2, 3))
Vec3 = _Vec3Factory()

class _Vec3(object):
    __slots__ = ['x', 'y', 'z']

    def __init__(self, vector):
        """Constructor

        :param vector: 2-tuple or 3-tuple
        """
        self.x = float(vector[X])
        self.y = float(vector[Y])
        try: # accept also (x,y) tuples -> z=0
            self.z = float(vector[Z])
        except IndexError:
            self.z = 0.

    def __len__(self):
        return 3

    def __hash__(self):
        return hash((self.x, self.y, self.z))

    def copy(self):
        return Vec3((self.x, self.y, self.z))
    __copy__ = copy

    def __getitem__(self, key):
        """Access vector elements by index.
        """
        # this is important to use a Vec3 object as replacement for a tuple and
        # vice verse.
        return self.__getattribute__(('x', 'y', 'z')[key])

    def __iter__(self):
        return iter((self.x, self.y, self.z))

    def __setitem__(self, key, value):
        """Set vector elements by index.
        """
        self.__setattr__(('x', 'y', 'z')[key], float(value))

    def __abs__(self):
        """Magnitude"""
        return self.length2()**0.5

    length = __abs__

    def length2(self):
        return self.x**2 + self.y**2 + self.z**2

    @property
    def isNULL(self):
        return not any(self)

    def normalized(self):
        m = self.__abs__()
        return Vec3( (self.x/m, self.y/m, self.z/m) )
    normalised = normalized

    def __neg__(self):
        return Vec3( (-self.x, -self.y, -self.z) )

    def __nonzero__(self):
        return bool(self.x or self.y or self.z)

    def __eq__(self, other):
        return self.x == other.x and \
               self.y == other.y and \
               self.z == other.z

    def __lt__(self, other):
        if self.x == other.x:
            if self.y == other.y:
                return self.z < other.z
            else:
                return self.y < other.y
        else:
            return self.x < other.x

    def __add__(self, other):
        if isinstance(other, _Vec3):
            return Vec3( (self.x+other.x, self.y+other.y, self.z+other.z) )
        else:
            scalar = float(other)
            return Vec3( (self.x+scalar, self.y+scalar, self.z+scalar) )

    def __radd__(self, other):
        return self.__add__(other)

    def __sub__(self, other):
        if isinstance(other, _Vec3):
            return Vec3( (self.x-other.x, self.y-other.y, self.z-other.z) )
        else:
            scalar = float(other)
            return Vec3( (self.x-scalar, self.y-scalar, self.z-scalar) )

    def __rsub__(self, other):
        if isinstance(other, _Vec3):
            return Vec3( (other.x-self.x, other.y-self.y, other.z-self.z) )
        else:
            scalar = float(other)
            return Vec3( (scalar-self.x, scalar-self.y, scalar-self.z) )

    def __mul__(self, other):
        scalar = float(other)
        if scalar:
            return Vec3( (self.x*scalar, self.y*scalar, self.z*scalar) )
        else:
            return NULLVEC

    def __rmul__(self, other):
        return self.__mul__(other)

    def __truediv__(self, other):
        scalar = float(other)
        return Vec3( (self.x/scalar, self.y/scalar, self.z/scalar) )
    __div__ = __truediv__

    def __rtruediv__(self, other):
        scalar = float(other)
        return Vec3( (scalar/self.x, scalar/self.y, scalar/self.z) )
    __rdiv__ = __rtruediv__

    def dot(self, other):
        """Scalar product"""
        return self.x*other[X] + self.y*other[Y] + self.z*other[Z]

    def cross(self, other):
        """ Cross product"""
        return Vec3((self.y*other[Z] - self.z*other[Y],
                     self.z*other[X] - self.x*other[Z],
                     self.x*other[Y] - self.y*other[X]))

    def distance(self, other):
        return self.distance2(other)**0.5

    def distance2(self, other):
        if not isinstance(other, _Vec3):
            other = Vec3(other)
        dx = other[X] - self.x
        dy = other[Y] - self.y
        dz = other[Z] - self.z
        return dx*dx + dy*dy + dz*dz

    def mid(self, other):
        if not isinstance(other, _Vec3):
            other = Vec3(other)
        return Vec3.from_floats((self.x + other[X]) / 2.,
                                (self.y + other[Y]) / 2.,
                                (self.z + other[Z]) / 2.)
class _NullVec(_Vec3):
    def angle(self):
        raise ValueError("for NULLVEC not defined")

    def normalized(self):
        return NULLVEC

    def perpendicular(self):
        return NULLVEC

    def __mul__(self, other):
        return NULLVEC

    def __truediv__(self, other):
        return NULLVEC
    __div__ = __truediv__

    def __rtruediv__(self, other):
        raise ZeroDivisionError("division by NULLVEC")
    __rdiv__ = __rtruediv__

NULLVEC = _NullVec((0, 0, 0))


def distance(point1, point2):
    """Calc distance between two 3d points """
    return ((point1[X]-point2[X])**2 +
            (point1[Y]-point2[Y])**2 +
            (point1[Z]-point2[Z])**2)**0.5

def midpoint(point1, point2):
    """Calc midpoint between point1 and point2 """
    return ((point1[X]+point2[X])*.5,
            (point1[Y]+point2[Y])*.5,
            (point1[Z]+point2[Z])*.5)