
##############################################################################
#
# Copyright (c) 2003-2020 by The University of Queensland
# http://www.uq.edu.au
#
# Primary Business: Queensland, Australia
# Licensed under the Apache License, version 2.0
# http://www.apache.org/licenses/LICENSE-2.0
#
# Development until 2012 by Earth Systems Science Computational Center (ESSCC)
# Development 2012-2013 by School of Earth Sciences
# Development from 2014 by Centre for Geoscience Computing (GeoComp)
# Development from 2019 by School of Earth and Environmental Sciences
#
##############################################################################

from __future__ import print_function, division

__copyright__="""Copyright (c) 2003-2020 by The University of Queensland
http://www.uq.edu.au
Primary Business: Queensland, Australia"""
__license__="""Licensed under the Apache License, version 2.0
http://www.apache.org/licenses/LICENSE-2.0"""
__url__="https://launchpad.net/escript-finley"

"""
transformations

:var __author__: name of author
:var __copyright__: copyrights
:var __license__: licence agreement
:var __url__: url entry point on documentation
:var __version__: version
:var __date__: date of the version
:var DEG: unit of degree
:var RAD: unit of radiant
"""

__author__="Lutz Gross, l.gross@uq.edu.au"

import numpy
import math

_TYPE=numpy.float64
DEG=math.pi/180.
RAD=1.
class Transformation(object):
   """
   General class to define an affine transformation *x->Ax+b*.
   """
   def __init__(self):
       """
       Creates a linear transformation.
       """
       pass

   def __call__(self,x=numpy.zeros((3,))):
       """
       Applies transformation to ``x``.
       """
       raise NotImplementeError()

class Translation(Transformation):
    """
    Defines a translation *x->x+b*.
    """
    def __init__(self,b=numpy.zeros((3,),dtype=_TYPE)):
       """
       Creates the linear transformation *x->x+b*.
       """
       super(Translation, self).__init__()
       self.__b=numpy.array(b,_TYPE)

    def __call__(self,x=numpy.zeros((3,))):
       """
       Applies translation to ``x``.
       """
       return numpy.array(x,_TYPE)+self.__b

class Rotatation(Transformation):
    """
    Defines a rotation.
    """
    def __init__(self,axis=numpy.ones((3,),dtype=_TYPE),point=numpy.zeros((3,),dtype=_TYPE),angle=0.*RAD):
       """
       Creates a rotation using an axis and a point on the axis.
       """
       self.__axis=numpy.array(axis,dtype=_TYPE)
       self.__point=numpy.array(point,dtype=_TYPE)
       lax=numpy.dot(self.__axis,self.__axis)
       if not lax>0:
          raise ValueError("points must be distinct.")
       self.__axis/=math.sqrt(lax)
       self.__angle=float(angle)

    def __call__(self,x=numpy.zeros((3,))):
       """
       Applies the rotation to ``x``.
       """
       x=numpy.array(x,_TYPE)
       z=x-self.__point
       z0=numpy.dot(z,self.__axis)
       z_per=z-z0*self.__axis
       lz_per=numpy.dot(z_per,z_per)
       if lz_per>0:
         axis1=z_per/math.sqrt(lz_per)
         axis2=_cross(axis1,self.__axis)
         lax2=numpy.dot(axis2,axis2)
         if lax2>0:
            axis2/=math.sqrt(lax2)
            return z0*self.__axis+math.sqrt(lz_per)*(math.cos(self.__angle)*axis1-math.sin(self.__angle)*axis2)+self.__point
         else:
            return x
       else:
         return x

def _cross(x, y):
    """
    Returns the cross product of ``x`` and ``y``.
    """
    return numpy.array([x[1] * y[2] - x[2] * y[1], x[2] * y[0] - x[0] * y[2], x[0] * y[1] - x[1] * y[0]], _TYPE)

class Dilation(Transformation):
    """
    Defines a dilation.
    """
    def __init__(self,factor=1.,center=numpy.zeros((3,),dtype=_TYPE)):
       """
       Creates a dilation with a center and a given expansion/contraction
       factor.
       """
       if not abs(factor)>0:
          raise ValueError("factor must be non-zero.")
       self.__factor=factor
       self.__center=numpy.array(center,dtype=_TYPE)

    def __call__(self,x=numpy.zeros((3,))):
       """
       Applies dilation to ``x``.
       """
       x=numpy.array(x,_TYPE)
       return self.__factor*(x-self.__center)+self.__center

class Reflection(Transformation):
    """
    Defines a reflection on a plane.
    """
    def __init__(self,normal=numpy.ones((3,),dtype=_TYPE),offset=0.):
       """
       Defines a reflection on a plane defined in normal form.
       """
       self.__normal=numpy.array(normal,dtype=_TYPE)
       ln=math.sqrt(numpy.dot(self.__normal,self.__normal))
       if not ln>0.:
          raise ValueError("normal must have positive length.")
       self.__normal/=ln
       if isinstance(offset,float) or isinstance(offset,int):
          self.__offset=offset/ln
       else:
          self.__offset=numpy.dot(numpy.array(offset,dtype=_TYPE),self.__normal)

    def __call__(self,x=numpy.zeros((3,))):
       """
       Applies reflection to ``x``.
       """
       x=numpy.array(x,_TYPE)
       return x - 2*(numpy.dot(x,self.__normal)-self.__offset)*self.__normal


[(122, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (58, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (68, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (73, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (75, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (85, 'numpy.ones', 'numpy.ones', False, 'import numpy\n'), (85, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (89, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (90, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (91, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (94, 'math.sqrt', 'math.sqrt', False, 'import math\n'), (97, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (101, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (103, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (105, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (128, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (136, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (138, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (142, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (149, 'numpy.ones', 'numpy.ones', False, 'import numpy\n'), (153, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (163, 'numpy.zeros', 'numpy.zeros', False, 'import numpy\n'), (167, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (79, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (109, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (154, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (107, 'math.sqrt', 'math.sqrt', False, 'import math\n'), (111, 'math.sqrt', 'math.sqrt', False, 'import math\n'), (161, 'numpy.array', 'numpy.array', False, 'import numpy\n'), (168, 'numpy.dot', 'numpy.dot', False, 'import numpy\n'), (112, 'math.sqrt', 'math.sqrt', False, 'import math\n'), (112, 'math.cos', 'math.cos', False, 'import math\n'), (112, 'math.sin', 'math.sin', False, 'import math\n')]