added gitignore and some derived functions

.gitignore

+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+env/
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+*.egg-info/
+.installed.cfg
+*.egg
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*,cover
+.hypothesis/
+
+# pyc files:
+*.pyc 
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+
+# Flask instance folder
+instance/
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+target/
+
+# IPython Notebook
+.ipynb_checkpoints
+
+# pyenv
+.python-version
+
+# celery beat schedule file
+celerybeat-schedule
+
+# dotenv
+.env
+
+# virtualenv
+venv/
+ENV/
+
+# Spyder project settings
+.spyderproject
+
+# vim files
+*.sw*
+
+# folder
+tmp/
+.idea

tfields/planes3D.py

+#!/usr/bin/env
+# encoding: utf-8
+"""
+Author:     Daniel Boeckenhoff
+Mail:       daniel.boeckenhoff@ipp.mpg.de
+
+part of npTools library
+"""
+import npTools
+import sympy
+import mplTools as mpt
+import loggingTools
+logger = loggingTools.Logger(__name__)
+
+
+class Planes3D(npTools.VectorField3D):
+    """
+    Point-NormVector representaion of planes
+    """
+
+    def symbolic(self):
+        """
+        Returns:
+            list: list with sympy.Plane objects
+        """
+        return [sympy.Plane(point, normal_vector=vector)
+                for point, vector in zip(self.points, self.vectors)]
+
+    def plot(self, **kwargs):
+        """
+        forward to Mesh3D plotting
+        """
+        artists = []
+        for i in range(self.points.shape[0]):
+            artists.append(mpt.plotPlane(self.points[0],
+                                         self.vectors[0],
+                                         **kwargs))
+        # symbolic = self.symbolic()
+        # planeMeshes = [npTools.Mesh3D([pl.arbitrary_point(t=(i + 1) * 1. / 2 * np.pi)
+        #                                for i in range(4)],
+        #                               faces=[[0, 1, 2], [1, 2, 3], [2, 3, 0]]) for pl in symbolic]
+        # artists = [m.plot(**kwargs) for m in planeMeshes]
+        return artists
+
+
+if __name__ == '__main__':
+    import doctest
+    doctest.testmod()

tfields/scalarField3D.py

+#!/usr/bin/env
+# encoding: utf-8
+"""
+Author:     Daniel Boeckenhoff
+Mail:       daniel.boeckenhoff@ipp.mpg.de
+
+part of npTools library
+"""
+import numpy as np
+import npTools
+import loggingTools
+logger = loggingTools.Logger(__name__)
+
+
+class ScalarField3D(npTools.Points3D):
+    """Class to Store scalar Fields in 3D.
+    Multiple scalars can be saved. Thus scalars is a multidimensional array.
+
+    Note:
+
+    Args:
+        inputArray: see Points3D
+        scalars (list like object)
+
+    Attributes:
+        x1 (list): coordinate x
+        x2 (list): coordinate y
+        x3 (list): coordinate z
+        scalars (list, can be list of lists.)
+
+    Examples:
+        No need to fille a scalar to the field
+        >>> empty = npTools.ScalarField3D([[1,2,3], [4,5,6], [7,8,9]])
+        >>> empty.scalars
+        array([], shape=(3, 0), dtype=float64)
+
+        Multiple scalars also ok.
+        >>> s = npTools.ScalarField3D.createTest()
+        >>> s
+        points: [[ 3.  0.  6.]
+         [ 3.  0.  9.]]
+        scalars: [[ 0.  2.  0.]
+         [ 1.  2.  3.]]
+
+        Saving and reading works
+        >>> from tempfile import NamedTemporaryFile
+        >>> outFile = NamedTemporaryFile(suffix='.npz')
+        >>> s.savez(outFile.name)
+        >>> _ = outFile.seek(0)
+        >>> s1 = npTools.ScalarField3D.createFromFile(outFile.name)
+        >>> bool(np.all(s == s1))
+        True
+        >>> s1.scalars
+        array([[ 0.,  2.,  0.],
+               [ 1.,  2.,  3.]])
+
+        Adding scalar fields
+        >>> s1[:, 2] = 42
+        >>> s1.scalars[:, 0] = 24
+        >>> s2 = npTools.ScalarField3D([s, s1, s])
+        >>> s2
+        points: [[  3.   0.   6.]
+         [  3.   0.   9.]
+         [  3.   0.  42.]
+         [  3.   0.  42.]
+         [  3.   0.   6.]
+         [  3.   0.   9.]]
+        scalars: [[  0.   2.   0.]
+         [  1.   2.   3.]
+         [ 24.   2.   0.]
+         [ 24.   2.   3.]
+         [  0.   2.   0.]
+         [  1.   2.   3.]]
+
+        >>> s3 = npTools.ScalarField3D([s, s1, s], scalars = [[0], [1], [2], [3], [4], [5]])
+        >>> s3
+        points: [[  3.   0.   6.]
+         [  3.   0.   9.]
+         [  3.   0.  42.]
+         [  3.   0.  42.]
+         [  3.   0.   6.]
+         [  3.   0.   9.]]
+        scalars: [[0]
+         [1]
+         [2]
+         [3]
+         [4]
+         [5]]
+
+    """
+    __slots__ = ['coordSys', 'scalars']
+
+    def __new__(cls, inputArray, scalars=None, **kwargs):
+
+        if any([isinstance(a, ScalarField3D) for a in inputArray]) and scalars is None:
+            if not all([isinstance(a, ScalarField3D) for a in inputArray]):
+                # TODO: could allow if all scalars are none
+                raise TypeError("Copy constructor of ScalarField3D only accepts"
+                                "ScalarField3D instances.")
+            if not len(set([o.scalars.shape for o in inputArray])) == 1:
+                raise TypeError("Shape of scalars that should be concatenated must be the same.")
+            scalars = np.concatenate([o.scalars for o in inputArray], axis=0)
+
+        obj = npTools.Points3D.__new__(cls, inputArray, **kwargs)
+
+        # add scalars
+        if scalars is None:
+            scalars = np.empty((len(obj), 0))
+        scalars = np.array(scalars)
+        if not len(scalars) == len(obj):
+            raise ValueError("Length of scalars ({0}) has to be the same as "
+                             "points length ({1}) but is not.".format(len(scalars), len(obj)))
+        obj.scalars = scalars
+        return obj
+
+    def __repr__(self):
+        return "points: %s\nscalars: %s" % (self, self.scalars)
+
+    def cut(self, cutExpression, coordSys=None):
+        """
+        Examples:
+            >>> import npTools
+            >>> x, y, z = sympy.symbols('x y z')
+            >>> p = npTools.ScalarField3D([[1., 2., 3.], [4., 5., 6.], [1, 2, -6],
+            ...                    [-5, -5, -5], [1,0,-1], [0,1,-1]], scalars = [1,2,3,4,5,6])
+            >>> p.cut(x > 0)
+            points: [[ 1.  2.  3.]
+             [ 4.  5.  6.]
+             [ 1.  2. -6.]
+             [ 1.  0. -1.]]
+            scalars: [1 2 3 5]
+
+            cut combination
+            >>> p.cut((x > 0) & (z < 0))
+            points: [[ 1.  2. -6.]
+             [ 1.  0. -1.]]
+            scalars: [3 5]
+
+        """
+        mask = self.getMask(cutExpression, coordSys=coordSys)
+
+        inst = self[mask].copy()
+        inst.scalars = inst.scalars[mask]
+
+        return inst
+
+    @classmethod
+    def createTest(cls):
+        p = npTools.MeshGrid3D(3, 4, 1, 0, 2, 1, 6, 9, 2, iterOrder=[1, 0, 2], coordSys=cls.CYLINDER)
+        scalars = [[i * 1., 2, 3 * i] for i in range(len(p))]
+        return cls(p, scalars)
+
+    def evalAtPoint(self, point):
+        """
+        Examples:
+            >>> s = npTools.ScalarField3D.createTest()
+            >>> s.evalAtPoint(npTools.Points3D([[3, 0, 6]]))
+            array([[ 0.,  2.,  0.]])
+            >>> s.evalAtPoint([3, 0, 9])
+            array([[ 1.,  2.,  3.]])
+            >>> s.evalAtPoint([[0, 0, 0]])
+            array([], shape=(0, 3), dtype=float64)
+
+        """
+        return self.scalars[np.equal(self, point).all(axis=1)]
+
+    def getScalarDepth(self):
+        return self.scalars.shape[1]
+
+    def getScalarArrays(self):
+        """
+        Examples:
+            >>> s = npTools.ScalarField3D.createTest()
+            >>> s.getScalarArrays()
+            array([[ 0.,  1.],
+                   [ 2.,  2.],
+                   [ 0.,  3.]])
+
+        """
+        return self.scalars.T
+
+    def setScalarArray(self, scalarIndex, scalarArray):
+        """
+        Examples:
+            >>> s = npTools.ScalarField3D.createTest()
+            >>> s.setScalarArray(1, [42, 84])
+            >>> s.getScalarArrays()
+            array([[  0.,   1.],
+                   [ 42.,  84.],
+                   [  0.,   3.]])
+
+            >>> s.setScalarArray(3, [-1, -1])
+            >>> s.getScalarArrays()
+            array([[  0.,   1.],
+                   [ 42.,  84.],
+                   [  0.,   3.],
+                   [ -1.,  -1.]])
+
+        """
+        if scalarIndex == self.getScalarDepth():
+            self.appendScalars(scalarArray)
+        self.scalars[:, scalarIndex] = np.array(scalarArray)
+
+    def appendScalars(self, scalars):
+        """
+        Examples:
+            >>> s = npTools.ScalarField3D.createTest()
+            >>> s.appendScalars([42, 84])
+            >>> s.getScalarArrays()
+            array([[  0.,   1.],
+                   [  2.,   2.],
+                   [  0.,   3.],
+                   [ 42.,  84.]])
+
+        """
+        self.scalars = np.concatenate([self.scalars, np.array([scalars]).T], 1)
+
+    def dropScalars(self, indices=None):
+        """
+        Examples:
+            >>> s = npTools.ScalarField3D.createTest()
+            >>> s.appendScalars([42, 84])
+
+            give index
+            >>> s.dropScalars(-2)
+            >>> s.getScalarArrays()
+            array([[  0.,   1.],
+                   [  2.,   2.],
+                   [ 42.,  84.]])
+
+            or indices
+            >>> s.dropScalars([0, 1])
+            >>> s.getScalarArrays()
+            array([[ 42.,  84.]])
+
+        """
+        if indices:
+            self.scalars = np.delete(self.scalars, indices, 1)
+
+    def triangulate(self, method="areaConserving", **kwargs):
+        """
+        Args:
+            method (str): triangulation method
+
+                valid values:
+                  - "areaConserving"
+                      kwargs:
+                          vicinity (list of list of int): map of point
+                              neigbourship relations.
+
+        Examples:
+            Imagine you have a camera with some pixels.
+            Each pixel has a center point which can be mapped to the area viewed in 3d.
+            >>> points = np.array([[0,0,0],[1,0,0],[2,0,0],[0,1,0],[1,1,0],[2,1,0],[3,1,0],
+            ...                    [0,2,0],[1,2,0],[1.7,1.7,0],[3,2,1],[1,3,1],[2,3,1],[3,2,1]])
+
+            Every pixel and thus every point has a value belonging to it.
+            >>> scalars = points * 42.1
+            >>> sf = npTools.ScalarField3D(points, scalars=scalars)
+
+            It is known which point belongs to which pixel
+            >>> pixels = np.array([[0,0],[1,0],[2,0],[0,1],[1,1],[2,1],[3,1],
+            ...                    [0,2],[1,2],[2,2],[3,2],[1,3],[2,3],[3,3]])
+
+            Triangulation:
+            >>> t = sf.triangulate(method="areaConserving", vicinity=pixels)
+            >>> t.getNTriangles()
+            84
+            >>> t[9:12]  # 3rd triangle
+            Triangles3D([[ 1.        ,  1.        ,  0.        ],
+                         [ 1.        ,  0.5       ,  0.        ],
+                         [ 0.66666667,  0.66666667,  0.        ]])
+            >>> t.triangleScalars[3]
+            array([ 42.1,  42.1,   0. ])
+
+        """
+        if method == "areaConserving":
+            vicinity = kwargs.get('vicinity')
+            if vicinity is None:
+                raise AttributeError("kwarg vicinity is mandatory for areaConserving method.")
+            pointsVicinityMap = ScalarField3D(self, scalars=vicinity)
+            pixelTriangleCornerPoints = []
+            if self.scalars.size == 0:
+                pixelTriangleScalars = None
+            else:
+                pixelTriangleScalars = []
+            relPixels = [[1, 0], [1, -1], [0, -1]]
+            for i, pixel in enumerate(pointsVicinityMap.scalars):
+                locNeighbours = []
+                for relP in relPixels:
+                    mask = (pointsVicinityMap.scalars == pixel + relP).all(axis=1)
+                    locNeighbours.append(mask)
+
+                pointI = pointsVicinityMap[i]
+                # locNeighbours.append(locNeighbours[0])
+                locNeighbours = np.array(locNeighbours)
+                # if not locNeighbours.sum() == 9:
+                #     continue
+                for n in range(2):
+                    if not locNeighbours[n: n + 2].any(axis=1).all():
+                        # not 2 neighbour points for triangle n
+                        continue
+                        
+                    j = np.where(locNeighbours[n])[0][0]
+                    k = np.where(locNeighbours[n + 1])[0][0]
+
+                    pointJ = pointsVicinityMap[j]
+                    pointK = pointsVicinityMap[k]
+            
+                    # points bisecting pointI and closePoints
+                    halfPointIJ = 0.5 * (pointI + pointJ)
+                    halfPointIK = 0.5 * (pointI + pointK)
+                    halfPointJK = 0.5 * (pointJ + pointK)
+            
+                    triangle = npTools.Triangles3D([pointI, pointJ, pointK])
+                    centroidPoint = triangle.getCentroids()
+                    
+                    pixelTriangleCornerPoints.extend([pointI, halfPointIJ, centroidPoint])
+                    pixelTriangleCornerPoints.extend([pointI, halfPointIK, centroidPoint])
+                    pixelTriangleCornerPoints.extend([pointJ, halfPointIJ, centroidPoint])
+                    pixelTriangleCornerPoints.extend([pointJ, halfPointJK, centroidPoint])
+                    pixelTriangleCornerPoints.extend([pointK, halfPointIK, centroidPoint])
+                    pixelTriangleCornerPoints.extend([pointK, halfPointJK, centroidPoint])
+
+                    if pixelTriangleScalars is not None:
+                        pixelTriangleScalars.extend([self.scalars[i]] * 2)
+                        pixelTriangleScalars.extend([self.scalars[j]] * 2)
+                        pixelTriangleScalars.extend([self.scalars[k]] * 2)
+
+            return npTools.Triangles3D(pixelTriangleCornerPoints,
+                                       triangleScalars=pixelTriangleScalars)
+
+        else:
+            raise NotImplementedError("Method {0} not implemented.".format(method))
+
+    def plot(self, **kwargs):
+        scalarIndex = kwargs.pop('scalarIndex', 0)
+        if not (self.scalars.size == 0 or scalarIndex is None):
+            kwargs['c'] = kwargs.pop('c', self.scalars.T[scalarIndex])
+        return super(ScalarField3D, self).plot(**kwargs)
+
+
+if __name__ == '__main__':
+    import doctest
+    import sympy  # NOQA: F401
+    doctest.testmod()

tfields/vectorField3D.py

+#!/usr/bin/env
+# encoding: utf-8
+"""
+Author:     Daniel Boeckenhoff
+Mail:       daniel.boeckenhoff@ipp.mpg.de
+
+part of npTools library
+"""
+import numpy as np
+import npTools
+import mplTools as mpt
+import loggingTools
+logger = loggingTools.Logger(__name__)
+
+
+class VectorField3D(object):
+    """
+    Container class for two points3D instances points and vectors
+    Attributes:
+        points
+        vectors
+    """
+    def __init__(self, points, vectors, **kwargs):
+        """
+        Examples:
+            >>> p = npTools.Points3D([[1, 1],[2, 2],[3, 3]], coordSys=npTools.Points3D.CYLINDER)
+            >>> v = npTools.Points3D([[4, 3],[5, 3],[6, 3]])
+            >>> vf = npTools.VectorField3D(p, v, coordSys=npTools.Points3D.CARTESIAN)
+            >>> vf
+            points: [[-0.41614684  0.90929743  3.        ]
+             [-0.41614684  0.90929743  3.        ]], vectors: [[ 4.  5.  6.]
+             [ 3.  3.  3.]]
+            >>> vf2 = npTools.VectorField3D(p, v, coordSys=npTools.Points3D.CYLINDER)
+            >>> vf2 == vf
+            True
+
+        """
+        coordSys = kwargs.pop('coordSys', None)
+        if not isinstance(points, npTools.Points3D):
+            self.points = npTools.Points3D(points, coordSys=coordSys)
+        else:
+            self.points = points.copy()
+            if coordSys is None:
+                coordSys = points.getCoordSys()
+        if not isinstance(vectors, npTools.Points3D):
+            self.vectors = npTools.Points3D(vectors, coordSys=coordSys)
+        else:
+            self.vectors = vectors.copy()
+            if coordSys is None:
+                coordSys = vectors.getCoordSys()
+        if coordSys is not None:  # meaning a coordinate system has been specified
+            self.coordinateTransform(coordSys)
+        if not len(self.vectors) == len(self.points):
+            raise ValueError("Vectors and points must be of same dimension!")
+
+    @classmethod
+    def createTest(cls):
+        """
+        A test factory method
+        """
+        p = npTools.MeshGrid3D(3, 4, 2, 0, 2, 3, 6, 9, 4, iterOrder=[1, 0, 2],
+                               coordSys=npTools.Points3D.CYLINDER)
+        v = npTools.MeshGrid3D(3, 4, 2, 0, 2, 3, 6, 9, 4, iterOrder=[1, 0, 2],
+                               coordSys=npTools.Points3D.CYLINDER)
+        return cls(p, v)
+
+    def coordinateTransform(self, coordSys):
+        """
+        Args:
+            coordSys (str): see Points3D.
+        """
+        self.points.coordinateTransform(coordSys)
+        self.vectors.coordinateTransform(coordSys)
+
+    def mirrorCoordinate(self, coordinate, condition=None):
+        """
+        Args:
+            coordSys (str): see Points3D.
+        """
+        if condition is not None:
+            raise NotImplementedError("Mirror with condition not defined in vectorfield")
+        self.points.mirrorCoordinate(coordinate)
+        self.vectors.mirrorCoordinate(coordinate)
+
+    def changeIterOrder(self, iterOrder):
+        """
+        Set a new iterorder
+        """
+        if not issubclass(self.points, npTools.MeshGrid3D):
+            raise TypeError("Can only change iterOrder for MeshGrid3D")
+        self.points.changeIterOrder(iterOrder)
+        self.vectors.changeGridIterOrder(iterOrder, *self.points.baseLengths,
+                                         actualIterOrder=self.points.iterOrder)
+
+    def getCoordSys(self):
+        """
+        Returns:
+            str: coordSys
+        """
+        return self.points.getCoordSys()
+
+    def copy(self):
+        """
+        Returns:
+            copy of object
+        """
+        return self.__class__(self.points, self.vectors)
+
+    def __getitem__(self, index):
+        return self.points[index], self.vectors[index]
+
+    def __repr__(self):
+        return "points: %s, vectors: %s" % (self.points, self.vectors)
+
+    def __eq__(self, other):
+        if self.getCoordSys() != other.getCoordSys():
+            other = other.copy()
+            other.coordinateTransform(self.getCoordSys())
+        return bool(np.all(np.isclose(self.points, other.points)) and
+                    np.all(np.isclose(self.vectors, other.vectors)))
+
+    def getVectors(self):
+        """
+        Returns:
+            vectors
+        """
+        return self.vectors.copy()
+
+    def getPoints(self):
+        """
+        Returns:
+            point bases
+        """
+        return self.points.copy()