Merge branch 'texture-components-in-rotation-class2' into 'development'
Texture components in rotation class2 See merge request damask/DAMASK!225
This commit is contained in:
Franz Roters
commit
190f8a82d4
|
|
@ -438,6 +438,7 @@ class Rotation:
|
|||
if np.any(ax[...,3] < 0.0) or np.any(ax[...,3] > np.pi):
|
||||
raise ValueError('Axis angle rotation angle outside of [0..π].')
|
||||
if not np.all(np.isclose(np.linalg.norm(ax[...,0:3],axis=-1), 1.0)):
|
||||
print(np.linalg.norm(ax[...,0:3],axis=-1))
|
||||
raise ValueError('Axis angle rotation axis is not of unit length.')
|
||||
|
||||
return Rotation(Rotation._ax2qu(ax))
|
||||
|
|
@ -652,6 +653,84 @@ class Rotation:
|
|||
asAxisAngle = as_axis_angle
|
||||
__mul__ = __matmul__
|
||||
|
||||
|
||||
@staticmethod
|
||||
def from_spherical_component(center,sigma,N=500,degrees=True,seed=None):
|
||||
"""
|
||||
Calculate set of rotations with Gaussian distribution around center.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
center : Rotation
|
||||
Central Rotation.
|
||||
sigma : float
|
||||
Standard deviation of (Gaussian) misorientation distribution.
|
||||
N : int, optional
|
||||
Number of samples, defaults to 500.
|
||||
degrees : boolean, optional
|
||||
sigma is given in degrees.
|
||||
seed : {None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
|
||||
A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy
|
||||
will be pulled from the OS.
|
||||
|
||||
"""
|
||||
rng = np.random.default_rng(seed)
|
||||
sigma = np.radians(sigma) if degrees else sigma
|
||||
u,Theta = (rng.random((N,2)) * 2.0 * np.array([1,np.pi]) - np.array([1.0, 0])).T
|
||||
omega = abs(rng.normal(scale=sigma,size=N))
|
||||
p = np.column_stack([np.sqrt(1-u**2)*np.cos(Theta),
|
||||
np.sqrt(1-u**2)*np.sin(Theta),
|
||||
u, omega])
|
||||
|
||||
return Rotation.from_axis_angle(p) @ center
|
||||
|
||||
|
||||
@staticmethod
|
||||
def from_fiber_component(alpha,beta,sigma=0.0,N=500,degrees=True,seed=None):
|
||||
"""
|
||||
Calculate set of rotations with Gaussian distribution around direction.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
alpha : numpy.ndarray of size 2
|
||||
Polar coordinates (phi from x,theta from z) of fiber direction in crystal frame.
|
||||
beta : numpy.ndarray of size 2
|
||||
Polar coordinates (phi from x,theta from z) of fiber direction in sample frame.
|
||||
sigma : float, optional
|
||||
Standard deviation of (Gaussian) misorientation distribution.
|
||||
Defaults to 0.
|
||||
N : int, optional
|
||||
Number of samples, defaults to 500.
|
||||
degrees : boolean, optional
|
||||
sigma, alpha, and beta are given in degrees.
|
||||
seed : {None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
|
||||
A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy
|
||||
will be pulled from the OS.
|
||||
|
||||
"""
|
||||
rng = np.random.default_rng(seed)
|
||||
sigma_,alpha_,beta_ = map(np.radians,(sigma,alpha,beta)) if degrees else (sigma,alpha,beta)
|
||||
|
||||
d_cr = np.array([np.sin(alpha_[0])*np.cos(alpha_[1]), np.sin(alpha_[0])*np.sin(alpha_[1]), np.cos(alpha_[0])])
|
||||
d_lab = np.array([np.sin( beta_[0])*np.cos( beta_[1]), np.sin( beta_[0])*np.sin( beta_[1]), np.cos( beta_[0])])
|
||||
ax_align = np.append(np.cross(d_lab,d_cr), np.arccos(np.dot(d_lab,d_cr)))
|
||||
if np.isclose(ax_align[3],0.0): ax_align[:3] = np.array([1,0,0])
|
||||
R_align = Rotation.from_axis_angle(ax_align if ax_align[3] > 0.0 else -ax_align,normalize=True) # rotate fiber axis from sample to crystal frame
|
||||
|
||||
u,Theta = (rng.random((N,2)) * 2.0 * np.array([1,np.pi]) - np.array([1.0, 0])).T
|
||||
omega = abs(rng.normal(scale=sigma_,size=N))
|
||||
p = np.column_stack([np.sqrt(1-u**2)*np.cos(Theta),
|
||||
np.sqrt(1-u**2)*np.sin(Theta),
|
||||
u, omega])
|
||||
p[:,:3] = np.einsum('ij,...j',np.eye(3)-np.outer(d_lab,d_lab),p[:,:3]) # remove component along fiber axis
|
||||
f = np.column_stack((np.broadcast_to(d_lab,(N,3)),rng.random(N)*np.pi))
|
||||
f[::2,:3] *= -1 # flip half the rotation axes to negative sense
|
||||
|
||||
return R_align.broadcast_to(N) \
|
||||
@ Rotation.from_axis_angle(p,normalize=True) \
|
||||
@ Rotation.from_axis_angle(f)
|
||||
|
||||
|
||||
####################################################################################################
|
||||
# Code below available according to the following conditions on https://github.com/MarDiehl/3Drotations
|
||||
####################################################################################################
|
||||
|
|
|
|||
|
|
@ -2,6 +2,7 @@ import os
|
|||
|
||||
import pytest
|
||||
import numpy as np
|
||||
from scipy import stats
|
||||
|
||||
from damask import Rotation
|
||||
from damask import _rotation
|
||||
|
|
@ -920,3 +921,39 @@ class TestRotation:
|
|||
R_2 = Rotation.from_axis_angle([0,0,1,angle],degrees=True)
|
||||
avg_angle = R_1.average(R_2).as_axis_angle(degrees=True,pair=True)[1]
|
||||
assert np.isclose(avg_angle,10+(angle-10)/2.)
|
||||
|
||||
|
||||
@pytest.mark.parametrize('sigma',[5,10,15,20])
|
||||
@pytest.mark.parametrize('N',[1000,10000,100000])
|
||||
def test_spherical_component(self,N,sigma):
|
||||
c = Rotation.from_random()
|
||||
o = Rotation.from_spherical_component(c,sigma,N)
|
||||
_, angles = c.misorientation(o).as_axis_angle(pair=True,degrees=True)
|
||||
angles[::2] *= -1 # flip angle for every second to symmetrize distribution
|
||||
|
||||
p = stats.normaltest(angles)[1]
|
||||
sigma_out = np.std(angles)
|
||||
print(f'\np: {p}, sigma ratio {sigma/sigma_out}')
|
||||
assert (.9 < sigma/sigma_out < 1.1) and p > 0.001
|
||||
|
||||
|
||||
@pytest.mark.parametrize('sigma',[5,10,15,20])
|
||||
@pytest.mark.parametrize('N',[1000,10000,100000])
|
||||
def test_from_fiber_component(self,N,sigma):
|
||||
"""https://en.wikipedia.org/wiki/Full_width_at_half_maximum."""
|
||||
alpha = np.random.random(2)*np.pi
|
||||
beta = np.random.random(2)*np.pi
|
||||
|
||||
f_in_C = np.array([np.sin(alpha[0])*np.cos(alpha[1]), np.sin(alpha[0])*np.sin(alpha[1]), np.cos(alpha[0])])
|
||||
f_in_S = np.array([np.sin(beta[0] )*np.cos(beta[1] ), np.sin(beta[0] )*np.sin(beta[1] ), np.cos(beta[0] )])
|
||||
ax = np.append(np.cross(f_in_C,f_in_S), - np.arccos(np.dot(f_in_C,f_in_S)))
|
||||
n = Rotation.from_axis_angle(ax if ax[3] > 0.0 else ax*-1.0 ,normalize=True) # rotation to align fiber axis in crystal and sample system
|
||||
|
||||
o = Rotation.from_fiber_component(alpha,beta,np.radians(sigma),N,False)
|
||||
angles = np.arccos(np.clip(np.dot(o@np.broadcast_to(f_in_S,(N,3)),n@f_in_S),-1,1))
|
||||
dist = np.array(angles) * (np.random.randint(0,2,N)*2-1)
|
||||
|
||||
p = stats.normaltest(dist)[1]
|
||||
sigma_out = np.degrees(np.std(dist))
|
||||
print(f'\np: {p}, sigma ratio {sigma/sigma_out}')
|
||||
assert (.9 < sigma/sigma_out < 1.1) and p > 0.001
|
||||
|
|
|
|||
Loading…
Reference in New Issue