[cig-commits] r12487 - short/3D/PyLith/trunk/playpen/euler
willic3 at geodynamics.org
willic3 at geodynamics.org
Tue Jul 29 11:44:32 PDT 2008
Author: willic3
Date: 2008-07-29 11:44:31 -0700 (Tue, 29 Jul 2008)
New Revision: 12487
Modified:
short/3D/PyLith/trunk/playpen/euler/euler.py
short/3D/PyLith/trunk/playpen/euler/transform.py
Log:
Slight change to transform.py to (hopefully) correctly determine the updip
direction using results from Meade & Hager.
Explicitly make all numpy arrays float64.
Modified: short/3D/PyLith/trunk/playpen/euler/euler.py
===================================================================
--- short/3D/PyLith/trunk/playpen/euler/euler.py 2008-07-29 17:32:04 UTC (rev 12486)
+++ short/3D/PyLith/trunk/playpen/euler/euler.py 2008-07-29 18:44:31 UTC (rev 12487)
@@ -139,7 +139,7 @@
self.numPoints = 0
self.pointsUTM = []
self.normals = []
- self.eulerPole = numpy.array([0.0, 0.0, 0.0], dtype=float)
+ self.eulerPole = numpy.array([0.0, 0.0, 0.0], dtype=float64)
# Note that we use a mean radius since we are doing rotations on a
# spherical Earth.
self.earthRad = 6372795.477598
@@ -199,11 +199,11 @@
self.eulerPole[2] = self.earthRad * rot * slat
self.upVec = numpy.array([float(self.upDir[0]), float(self.upDir[1]),
- float(self.upDir[2])], dtype=float)
+ float(self.upDir[2])], dtype=float64)
self.normalVec = numpy.array([float(self.normalDir[0]),
float(self.normalDir[1]),
- float(self.normalDir[2])], dtype=float)
+ float(self.normalDir[2])], dtype=float64)
self.dipCutoffProj = abs(math.sin(self.dipCutoff.value))
@@ -244,12 +244,14 @@
self.destCoordSys.initialize()
from spatialdata.geocoords.Converter import convert
- pointsLL = numpy.array(self.pointsUTM, dtype=float).reshape(self.numPoints,
- self.spaceDim)
+ pointsLL = numpy.array(self.pointsUTM,
+ dtype=float64).reshape(self.numPoints,
+ self.spaceDim)
convert(pointsLL, self.destCoordSys, self.srcCoordSys)
- normalsArr = numpy.array(self.normals, dtype=float).reshape(self.numPoints,
- self.spaceDim)
+ normalsArr = numpy.array(self.normals,
+ dtype=float64).reshape(self.numPoints,
+ self.spaceDim)
iCount = 0
velocity = [0.0, 0.0, 0.0]
@@ -355,13 +357,13 @@
pX = clat * clon
pY = clat * slon
pZ = slat
- pointPole = numpy.array([pX, pY, pZ], dtype=float)
+ pointPole = numpy.array([pX, pY, pZ], dtype=float64)
velGlobal = numpy.cross(self.eulerPole, pointPole)
rotMatrix = self._makeRot(latPoint, lonPoint)
velNED = numpy.dot(rotMatrix, velGlobal)
# NOTE: I should rearrange the rotation matrix so this transformation
# isn't necessary.
- velENU = numpy.array([velNED[1], velNED[0], -velNED[2]], dtype=float)
+ velENU = numpy.array([velNED[1], velNED[0], -velNED[2]], dtype=float64)
return velENU
@@ -376,7 +378,7 @@
vec1 = [ -slat * clon, -slat * slon, clat]
vec2 = [ -slon, clon, 0.0]
vec3 = [ -clat * clon, -clat * slon, -slat]
- rotMatrix = numpy.array([vec1, vec2, vec3], dtype=float)
+ rotMatrix = numpy.array([vec1, vec2, vec3], dtype=float64)
return rotMatrix
Modified: short/3D/PyLith/trunk/playpen/euler/transform.py
===================================================================
--- short/3D/PyLith/trunk/playpen/euler/transform.py 2008-07-29 17:32:04 UTC (rev 12486)
+++ short/3D/PyLith/trunk/playpen/euler/transform.py 2008-07-29 18:44:31 UTC (rev 12487)
@@ -163,11 +163,11 @@
self.defaultValues = self.inventory.defaultValues
self.upVec = numpy.array([float(self.upDir[0]), float(self.upDir[1]),
- float(self.upDir[2])], dtype=float)
+ float(self.upDir[2])], dtype=float64)
self.normalVec = numpy.array([float(self.normalDir[0]),
float(self.normalDir[1]),
- float(self.normalDir[2])], dtype=float)
+ float(self.normalDir[2])], dtype=float64)
self.dipCutoffProj = abs(math.sin(self.dipCutoff.value))
@@ -207,11 +207,12 @@
spatial DB.
"""
- normalsArr = numpy.array(self.normals, dtype=float).reshape(self.numPoints,
+ normalsArr = numpy.array(self.normals,
+ dtype=float64).reshape(self.numPoints,
+ self.spaceDim)
+
+ points = numpy.array(self.pointsUTM, dtype=float64).reshape(self.numPoints,
self.spaceDim)
-
- points = numpy.array(self.pointsUTM, dtype=float).reshape(self.numPoints,
- self.spaceDim)
iCount = 0
velocity = [0.0, 0.0, 0.0]
@@ -311,23 +312,25 @@
ts = self.segInfo[iSeg,7]
dx = sx2 - sx1
dy = sy2 - sy1
- asVec = numpy.array([dx, dy, 0.0], dtype=float)
+ asVec = numpy.array([dx, dy, 0.0], dtype=float64)
mag = math.sqrt(numpy.dot(asVec,asVec))
asVec /= mag
+ horPerp = numpy.cross(self.upVec, asVec)
+ mag = math.sqrt(numpy.dot(horPerp,horPerp))
+ horPerp /= mag
+ if dip > 90.0:
+ horPerp *= -1.0
testVec = numpy.cross(asVec, self.normalVec)
if testVec[2] > 0.0:
asVec *= -1.0
- horPerp = numpy.cross(self.upVec, asVec)
- mag = math.sqrt(numpy.dot(horPerp,horPerp))
- horPerp /= mag
dipCos = math.sin(math.radians(dip))
if math.fabs(dipCos) != 1.0:
r = math.sqrt(1.0/(1.0-dipCos*dipCos))
- udVec = numpy.array([horPerp[0]/r, horPerp[1]/r, dipCos], dtype=float)
+ udVec = numpy.array([horPerp[0]/r, horPerp[1]/r, dipCos], dtype=float64)
else:
- udVec = numpy.array([0.0, 0.0, dipCos], dtype=float)
+ udVec = numpy.array([0.0, 0.0, dipCos], dtype=float64)
normVec = numpy.cross(asVec, udVec)
- slipVec = numpy.array([ss, ds,ts], dtype=float)
+ slipVec = numpy.array([ss, ds,ts], dtype=float64)
rot1 = numpy.vstack((asVec, udVec, normVec))
rot = rot1.transpose()
velocity = numpy.dot(rot,slipVec)
@@ -426,7 +429,7 @@
self.numSegs += 1
iCount += 1
f.close()
- self.segInfo = numpy.array(segtmp, dtype=float).reshape(self.numSegs, 8)
+ self.segInfo = numpy.array(segtmp, dtype=float64).reshape(self.numSegs, 8)
return
# ----------------------------------------------------------------------
More information about the cig-commits
mailing list