/[escript]/trunk/pycad/py_src/transformations.py
ViewVC logotype

Contents of /trunk/pycad/py_src/transformations.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3259 - (show annotations)
Mon Oct 11 01:48:14 2010 UTC (8 years, 11 months ago) by jfenwick
File MIME type: text/x-python
File size: 4827 byte(s)
Merging dudley and scons updates from branches

1
2 ########################################################
3 #
4 # Copyright (c) 2003-2010 by University of Queensland
5 # Earth Systems Science Computational Center (ESSCC)
6 # http://www.uq.edu.au/esscc
7 #
8 # Primary Business: Queensland, Australia
9 # Licensed under the Open Software License version 3.0
10 # http://www.opensource.org/licenses/osl-3.0.php
11 #
12 ########################################################
13
14 __copyright__="""Copyright (c) 2003-2010 by University of Queensland
15 Earth Systems Science Computational Center (ESSCC)
16 http://www.uq.edu.au/esscc
17 Primary Business: Queensland, Australia"""
18 __license__="""Licensed under the Open Software License version 3.0
19 http://www.opensource.org/licenses/osl-3.0.php"""
20 __url__="https://launchpad.net/escript-finley"
21
22 """
23 transformations
24
25 :var __author__: name of author
26 :var __copyright__: copyrights
27 :var __license__: licence agreement
28 :var __url__: url entry point on documentation
29 :var __version__: version
30 :var __date__: date of the version
31 :var DEG: unit of degree
32 :var RAD: unit of radiant
33 """
34
35 __author__="Lutz Gross, l.gross@uq.edu.au"
36
37 import numpy
38 import math
39
40 _TYPE=numpy.float64
41 DEG=math.pi/180.
42 RAD=1.
43 class Transformation(object):
44 """
45 General class to define an affine transformation *x->Ax+b*.
46 """
47 def __init__(self):
48 """
49 Creates a linear transformation.
50 """
51 pass
52
53 def __call__(self,x=numpy.zeros((3,))):
54 """
55 Applies transformation to ``x``.
56 """
57 raise NotImplementeError()
58
59 class Translation(Transformation):
60 """
61 Defines a translation *x->x+b*.
62 """
63 def __init__(self,b=numpy.zeros((3,),dtype=_TYPE)):
64 """
65 Creates the linear transformation *x->x+b*.
66 """
67 super(Translation, self).__init__()
68 self.__b=numpy.array(b,_TYPE)
69
70 def __call__(self,x=numpy.zeros((3,))):
71 """
72 Applies translation to ``x``.
73 """
74 return numpy.array(x,_TYPE)+self.__b
75
76 class Rotatation(Transformation):
77 """
78 Defines a rotation.
79 """
80 def __init__(self,axis=numpy.ones((3,),dtype=_TYPE),point=numpy.zeros((3,),dtype=_TYPE),angle=0.*RAD):
81 """
82 Creates a rotation using an axis and a point on the axis.
83 """
84 self.__axis=numpy.array(axis,dtype=_TYPE)
85 self.__point=numpy.array(point,dtype=_TYPE)
86 lax=numpy.dot(self.__axis,self.__axis)
87 if not lax>0:
88 raise ValueError("points must be distinct.")
89 self.__axis/=math.sqrt(lax)
90 self.__angle=float(angle)
91
92 def __call__(self,x=numpy.zeros((3,))):
93 """
94 Applies the rotation to ``x``.
95 """
96 x=numpy.array(x,_TYPE)
97 z=x-self.__point
98 z0=numpy.dot(z,self.__axis)
99 z_per=z-z0*self.__axis
100 lz_per=numpy.dot(z_per,z_per)
101 if lz_per>0:
102 axis1=z_per/math.sqrt(lz_per)
103 axis2=_cross(axis1,self.__axis)
104 lax2=numpy.dot(axis2,axis2)
105 if lax2>0:
106 axis2/=math.sqrt(lax2)
107 return z0*self.__axis+math.sqrt(lz_per)*(math.cos(self.__angle)*axis1-math.sin(self.__angle)*axis2)+self.__point
108 else:
109 return x
110 else:
111 return x
112
113 def _cross(x, y):
114 """
115 Returns the cross product of ``x`` and ``y``.
116 """
117 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)
118
119 class Dilation(Transformation):
120 """
121 Defines a dilation.
122 """
123 def __init__(self,factor=1.,center=numpy.zeros((3,),dtype=_TYPE)):
124 """
125 Creates a dilation with a center and a given expansion/contraction
126 factor.
127 """
128 if not abs(factor)>0:
129 raise ValueError("factor must be non-zero.")
130 self.__factor=factor
131 self.__center=numpy.array(center,dtype=_TYPE)
132
133 def __call__(self,x=numpy.zeros((3,))):
134 """
135 Applies dilation to ``x``.
136 """
137 x=numpy.array(x,_TYPE)
138 return self.__factor*(x-self.__center)+self.__center
139
140 class Reflection(Transformation):
141 """
142 Defines a reflection on a plane.
143 """
144 def __init__(self,normal=numpy.ones((3,),dtype=_TYPE),offset=0.):
145 """
146 Defines a reflection on a plane defined in normal form.
147 """
148 self.__normal=numpy.array(normal,dtype=_TYPE)
149 ln=math.sqrt(numpy.dot(self.__normal,self.__normal))
150 if not ln>0.:
151 raise ValueError("normal must have positive length.")
152 self.__normal/=ln
153 if isinstance(offset,float) or isinstance(offset,int):
154 self.__offset=offset/ln
155 else:
156 self.__offset=numpy.dot(numpy.array(offset,dtype=_TYPE),self.__normal)
157
158 def __call__(self,x=numpy.zeros((3,))):
159 """
160 Applies reflection to ``x``.
161 """
162 x=numpy.array(x,_TYPE)
163 return x - 2*(numpy.dot(x,self.__normal)-self.__offset)*self.__normal
164

  ViewVC Help
Powered by ViewVC 1.1.26