Merge branch 'vectorize_rotation' into 'development'

Vectorize rotation

See merge request damask/DAMASK!171
This commit is contained in:
Karo 2020-05-28 15:25:40 +02:00
commit f63e412bcf
11 changed files with 1288 additions and 925 deletions

View File

@ -2,6 +2,7 @@
import os
import sys
from io import StringIO
from optparse import OptionParser
import numpy as np
@ -24,14 +25,8 @@ Additional (globally fixed) rotations of the lab frame and/or crystal frame can
""", version = scriptID)
representations = {
'quaternion': ['qu',4],
'rodrigues': ['ro',4],
'Rodrigues': ['Ro',3],
'eulers': ['eu',3],
'matrix': ['om',9],
'angleaxis': ['ax',4],
}
representations = ['quaternion', 'rodrigues', 'eulers', 'matrix', 'axisangle']
parser.add_option('-o',
'--output',
@ -93,10 +88,10 @@ parser.set_defaults(output = [],
)
(options, filenames) = parser.parse_args()
if filenames == []: filenames = [None]
#options.output = list(map(lambda x: x.lower(), options.output))
if options.output == [] or (not set(options.output).issubset(set(representations))):
parser.error('output must be chosen from {}.'.format(', '.join(representations)))
parser.error('output must be chosen from {}.'.format(', '.join(representations)))
input = [options.eulers is not None,
options.rodrigues is not None,
@ -109,95 +104,47 @@ input = [options.eulers is not None,
if np.sum(input) != 1: parser.error('needs exactly one input format.')
(label,dim,inputtype) = [(options.eulers,representations['eulers'][1],'eulers'),
(options.rodrigues,representations['rodrigues'][1],'rodrigues'),
([options.x,options.y,options.z],[3,3,3],'frame'),
(options.matrix,representations['matrix'][1],'matrix'),
(options.quaternion,representations['quaternion'][1],'quaternion'),
][np.where(input)[0][0]] # select input label that was requested
r = damask.Rotation.fromAxisAngle(np.array(options.crystalrotation),options.degrees,normalise=True)
R = damask.Rotation.fromAxisAngle(np.array(options.labrotation),options.degrees,normalise=True)
# --- loop over input files ------------------------------------------------------------------------
if filenames == []: filenames = [None]
r = damask.Rotation.from_axis_angle(np.array(options.crystalrotation),options.degrees,normalise=True)
R = damask.Rotation.from_axis_angle(np.array(options.labrotation),options.degrees,normalise=True)
for name in filenames:
try:
table = damask.ASCIItable(name = name)
except IOError:
continue
damask.util.report(scriptName,name)
damask.util.report(scriptName,name)
# ------------------------------------------ read header ------------------------------------------
table = damask.Table.from_ASCII(StringIO(''.join(sys.stdin.read())) if name is None else name)
table.head_read()
if options.eulers is not None:
label = options.eulers
print(np.max(table.get(options.eulers),axis=0))
o = damask.Rotation.from_Eulers(table.get(options.eulers), options.degrees)
elif options.rodrigues is not None:
label = options.rodrigues
o = damask.Rotation.from_Rodrigues(table.get(options.rodrigues))
elif options.matrix is not None:
label = options.matrix
o = damask.Rotation.from_matrix(table.get(options.matrix).reshape(-1,3,3))
elif options.x is not None:
label = '<{},{},{}>'.format(options.x,options.y,options.z)
M = np.block([table.get(options.x),table.get(options.y),table.get(options.z)]).reshape(-1,3,3)
o = damask.Rotation.from_matrix(M/np.linalg.norm(M,axis=0))
elif options.quaternion is not None:
label = options.quaternion
o = damask.Rotation.from_quaternion(table.get(options.quaternion))
# ------------------------------------------ sanity checks -----------------------------------------
o = r.broadcast_to(o.shape) @ o @ R.broadcast_to(o.shape)
errors = []
remarks = []
if not np.all(table.label_dimension(label) == dim): errors.append('input {} does not have dimension {}.'.format(label,dim))
else: column = table.label_index(label)
#if options.lattice is not None:
# o = damask.Orientation(rotation = o,lattice = options.lattice).reduced().rotation
if remarks != []: damask.util.croak(remarks)
if errors != []:
damask.util.croak(errors)
table.close(dismiss = True)
continue
# ------------------------------------------ assemble header ---------------------------------------
if 'rodrigues' in options.output:
table.add('ro({})'.format(label),o.as_Rodrigues(), scriptID+' '+' '.join(sys.argv[1:]))
if 'eulers' in options.output:
table.add('eu({})'.format(label),o.as_Eulers(options.degrees), scriptID+' '+' '.join(sys.argv[1:]))
if 'quaternion' in options.output:
table.add('qu({})'.format(label),o.as_quaternion(), scriptID+' '+' '.join(sys.argv[1:]))
if 'matrix' in options.output:
table.add('om({})'.format(label),o.as_matrix(), scriptID+' '+' '.join(sys.argv[1:]))
if 'axisangle' in options.output:
table.add('om({})'.format(label),o.as_axisangle(options.degrees), scriptID+' '+' '.join(sys.argv[1:]))
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
for output in options.output:
if output in representations:
table.labels_append(['{}_{}({})'.format(i+1,representations[output][0],label) \
for i in range(representations[output][1])])
table.head_write()
# ------------------------------------------ process data ------------------------------------------
outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table
if inputtype == 'eulers':
d = representations['eulers'][1]
o = damask.Rotation.fromEulers(list(map(float,table.data[column:column+d])),options.degrees)
elif inputtype == 'rodrigues':
d = representations['rodrigues'][1]
o = damask.Rotation.fromRodrigues(list(map(float,table.data[column:column+d])))
elif inputtype == 'matrix':
d = representations['matrix'][1]
o = damask.Rotation.fromMatrix(np.array(list(map(float,table.data[column:column+d]))).reshape(3,3))
elif inputtype == 'frame':
M = np.array(list(map(float,table.data[column[0]:column[0]+3] + \
table.data[column[1]:column[1]+3] + \
table.data[column[2]:column[2]+3]))).reshape(3,3).T
o = damask.Rotation.fromMatrix(M/np.linalg.norm(M,axis=0))
elif inputtype == 'quaternion':
d = representations['quaternion'][1]
o = damask.Rotation.fromQuaternion(list(map(float,table.data[column:column+d])))
o = r*o*R # apply additional lab and crystal frame rotations
if options.lattice is not None:
o = damask.Orientation(rotation = o,lattice = options.lattice).reduced().rotation
for output in options.output:
if output == 'quaternion': table.data_append(o.asQuaternion())
elif output == 'rodrigues': table.data_append(o.asRodrigues())
elif output == 'Rodrigues': table.data_append(o.asRodrigues(vector=True))
elif output == 'eulers': table.data_append(o.asEulers(degrees=options.degrees))
elif output == 'matrix': table.data_append(o.asMatrix())
elif output == 'angleaxis': table.data_append(o.asAxisAngle(degrees=options.degrees))
outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization -----------------------------------
table.close() # close ASCII tables
table.to_ASCII(sys.stdout if name is None else name)

View File

@ -37,27 +37,26 @@ parser.add_option('--degrees',
parser.set_defaults(rotation = (1.,1.,1.,0), # no rotation about (1,1,1)
degrees = False,
)
(options,filenames) = parser.parse_args()
if filenames == []: filenames = [None]
if options.data is None:
parser.error('no data column specified.')
parser.error('no data column specified.')
r = damask.Rotation.fromAxisAngle(options.rotation,options.degrees,normalise=True)
r = damask.Rotation.from_axis_angle(options.rotation,options.degrees,normalise=True)
for name in filenames:
damask.util.report(scriptName,name)
table = damask.Table.from_ASCII(StringIO(''.join(sys.stdin.read())) if name is None else name)
for data in options.data:
d = table.get(data)
if table.shapes[data] == (9,): d=d.reshape(-1,3,3)
for i,l in enumerate(d):
d[i] = r*l
d = r.broadcast_to(d.shape[0:1]) @ d
if table.shapes[data] == (9,): d=d.reshape(-1,9)
table.set(data,d,scriptID+' '+' '.join(sys.argv[1:]))
table.to_ASCII(sys.stdout if name is None else name)

View File

@ -84,9 +84,9 @@ if [options.angleaxis,options.quaternion].count(None) == 0:
parser.error('more than one rotation specified.')
if options.angleaxis is not None:
rotation = damask.Rotation.fromAxisAngle(np.array(options.angleaxis),options.degrees,normalise=True)
rotation = damask.Rotation.from_axis_angle(np.array(options.angleaxis),options.degrees,normalise=True)
elif options.quaternion is not None:
rotation = damask.Rotation.fromQuaternion(options.quaternion)
rotation = damask.Rotation.from_quaternion(options.quaternion)
else:
rotation = damask.Rotation()

View File

@ -97,7 +97,7 @@ for name in filenames:
dataset = os.path.join(group_pointwise,options.quaternion)
try:
quats = np.reshape(inFile[dataset][...],(np.product(grid),4))
rot = [damask.Rotation.fromQuaternion(q,True,P=+1) for q in quats]
rot = [damask.Rotation.from_quaternion(q,True,P=+1) for q in quats]
except KeyError:
errors.append('Pointwise orientation (quaternion) data ({}) not readable'.format(dataset))
@ -123,7 +123,7 @@ for name in filenames:
dataset = os.path.join(group_average,options.quaternion)
try:
rot = [damask.Rotation.fromQuaternion(q,True,P=+1) for q in inFile[dataset][...][1:]] # skip first entry (unindexed)
rot = [damask.Rotation.from_quaternion(q,True,P=+1) for q in inFile[dataset][...][1:]] # skip first entry (unindexed)
except KeyError:
errors.append('Average orientation data ({}) not readable'.format(dataset))
@ -140,7 +140,7 @@ for name in filenames:
config_header = ['<texture>']
for i in range(np.nanmax(microstructure)):
config_header += ['[{}{}]'.format(label,i+1),
'(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 {:.2f}'.format(*rot[i].asEulers(degrees = True)),
'(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 {:.2f}'.format(*rot[i].as_Eulers(degrees = True)),
]
config_header += ['<microstructure>']
for i in range(np.nanmax(microstructure)):

View File

@ -89,7 +89,7 @@ for name in filenames:
for i,data in enumerate(unique):
ori = damask.Rotation(data[0:4])
config_header += ['[Grain{}]'.format(i+1),
'(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 {:.2f}'.format(*ori.asEulers(degrees = True)),
'(gauss)\tphi1 {:.2f}\tPhi {:.2f}\tphi2 {:.2f}'.format(*ori.as_Eulers(degrees = True)),
]
if options.axes is not None: config_header += ['axes\t{} {} {}'.format(*options.axes)]

View File

@ -59,15 +59,15 @@ if [options.rotation,options.eulers,options.matrix,options.quaternion].count(Non
parser.error('no rotation specified.')
if options.quaternion is not None:
rot = damask.Rotation.fromQuaternion(np.array(options.quaternion)) # we might need P=+1 here, too...
rot = damask.Rotation.from_quaternion(np.array(options.quaternion)) # we might need P=+1 here, too...
if options.rotation is not None:
rot = damask.Rotation.fromAxisAngle(np.array(options.rotation),degrees=options.degrees,normalise=True,P=+1)
rot = damask.Rotation.from_axis_angle(np.array(options.rotation),degrees=options.degrees,normalise=True,P=+1)
if options.matrix is not None:
rot = damask.Rotation.fromMatrix(np.array(options.Matrix))
rot = damask.Rotation.from_matrix(np.array(options.Matrix))
if options.eulers is not None:
rot = damask.Rotation.fromEulers(np.array(options.eulers),degrees=options.degrees)
rot = damask.Rotation.from_Eulers(np.array(options.eulers),degrees=options.degrees)
eulers = rot.asEulers(degrees=True)
eulers = rot.as_Eulers(degrees=True)
if filenames == []: filenames = [None]

View File

@ -635,6 +635,6 @@ class Lattice:
otherDir = miller[otherDir_id]/ np.linalg.norm(miller[otherDir_id])
otherMatrix = np.array([otherDir,np.cross(otherPlane,otherDir),otherPlane])
r['rotations'].append(Rotation.fromMatrix(np.dot(otherMatrix.T,myMatrix)))
r['rotations'].append(Rotation.from_matrix(np.dot(otherMatrix.T,myMatrix)))
return r

View File

@ -37,7 +37,7 @@ class Orientation:
if isinstance(rotation, Rotation):
self.rotation = rotation
else:
self.rotation = Rotation.fromQuaternion(rotation) # assume quaternion
self.rotation = Rotation.from_quaternion(rotation) # assume quaternion
if self.rotation.quaternion.shape != (4,):
raise NotImplementedError('Support for multiple rotations missing')
@ -68,8 +68,8 @@ class Orientation:
r = b*aInv
for k in range(2):
r.inverse()
breaker = self.lattice.symmetry.inFZ(r.asRodrigues(vector=True)) \
and (not SST or other.lattice.symmetry.inDisorientationSST(r.asRodrigues(vector=True)))
breaker = self.lattice.symmetry.inFZ(r.as_Rodrigues(vector=True)) \
and (not SST or other.lattice.symmetry.inDisorientationSST(r.as_Rodrigues(vector=True)))
if breaker: break
if breaker: break
if breaker: break
@ -78,7 +78,7 @@ class Orientation:
# ... own sym, other sym,
# self-->other: True, self<--other: False
def inFZ(self):
return self.lattice.symmetry.inFZ(self.rotation.asRodrigues(vector=True))
return self.lattice.symmetry.inFZ(self.rotation.as_Rodrigues(vector=True))
def equivalentOrientations(self,members=[]):
@ -100,7 +100,7 @@ class Orientation:
def reduced(self):
"""Transform orientation to fall into fundamental zone according to symmetry."""
for me in self.equivalentOrientations():
if self.lattice.symmetry.inFZ(me.rotation.asRodrigues(vector=True)): break
if self.lattice.symmetry.inFZ(me.rotation.as_Rodrigues(vector=True)): break
return self.__class__(me.rotation,self.lattice)

File diff suppressed because it is too large Load Diff

View File

@ -8,13 +8,8 @@ import damask
from damask import Rotation
from damask import Orientation
from damask import Lattice
n = 1000
@pytest.fixture
def default():
"""A set of n random rotations."""
return [Rotation.fromRandom() for r in range(n)]
n = 1000
@pytest.fixture
def reference_dir(reference_dir_base):
@ -28,15 +23,15 @@ class TestOrientation:
{'label':'green','RGB':[0,1,0],'direction':[0,1,1]},
{'label':'blue', 'RGB':[0,0,1],'direction':[1,1,1]}])
@pytest.mark.parametrize('lattice',['fcc','bcc'])
def test_IPF_cubic(self,default,color,lattice):
def test_IPF_cubic(self,color,lattice):
cube = damask.Orientation(damask.Rotation(),lattice)
for direction in set(permutations(np.array(color['direction']))):
assert np.allclose(cube.IPFcolor(direction),np.array(color['RGB']))
assert np.allclose(cube.IPFcolor(np.array(direction)),np.array(color['RGB']))
@pytest.mark.parametrize('lattice',Lattice.lattices)
def test_IPF(self,lattice):
direction = np.random.random(3)*2.0-1
for rot in [Rotation.fromRandom() for r in range(n//100)]:
for rot in [Rotation.from_random() for r in range(n//100)]:
R = damask.Orientation(rot,lattice)
color = R.IPFcolor(direction)
for equivalent in R.equivalentOrientations():
@ -45,7 +40,7 @@ class TestOrientation:
@pytest.mark.parametrize('model',['Bain','KS','GT','GT_prime','NW','Pitsch'])
@pytest.mark.parametrize('lattice',['fcc','bcc'])
def test_relationship_forward_backward(self,model,lattice):
ori = Orientation(Rotation.fromRandom(),lattice)
ori = Orientation(Rotation.from_random(),lattice)
for i,r in enumerate(ori.relatedOrientations(model)):
ori2 = r.relatedOrientations(model)[i]
misorientation = ori.rotation.misorientation(ori2.rotation)
@ -56,8 +51,8 @@ class TestOrientation:
def test_relationship_reference(self,update,reference_dir,model,lattice):
reference = os.path.join(reference_dir,'{}_{}.txt'.format(lattice,model))
ori = Orientation(Rotation(),lattice)
eu = np.array([o.rotation.asEulers(degrees=True) for o in ori.relatedOrientations(model)])
if update:
eu = np.array([o.rotation.as_Eulers(degrees=True) for o in ori.relatedOrientations(model)])
if update:
coords = np.array([(1,i+1) for i,x in enumerate(eu)])
table = damask.Table(eu,{'Eulers':(3,)})
table.add('pos',coords)

File diff suppressed because it is too large Load Diff