diff --git a/doc/changelog.rst b/doc/changelog.rst index 4f7c4c6a..1834401b 100644 --- a/doc/changelog.rst +++ b/doc/changelog.rst @@ -17,11 +17,15 @@ Unreleased Contributors ------------ - Håkon Wiik Ånes +- Lars Andreas Hastad Lervik Added ----- - Reader for EMsoft's simulated EBSD patterns returned by their ``EMEBSD.f90`` program. (`#202 `_) +- Modified Lambert mapping, and its inverse, from points on the unit sphere to a + 2D square grid, as implemented in Callahan and De Graef (2013). + (`#214 `_) Changed ------- diff --git a/kikuchipy/projections/__init__.py b/kikuchipy/projections/__init__.py index 33bdc298..e531d11b 100644 --- a/kikuchipy/projections/__init__.py +++ b/kikuchipy/projections/__init__.py @@ -20,12 +20,14 @@ ebsd_projections, gnomonic_projection, hesse_normal_form, + lambert_projection, spherical_projection, ) __all__ = [ "ebsd_projections", - "gnomonic_projection.py", + "gnomonic_projection", "hesse_normal_form", - "spherical_projection.py", + "lambert_projection", + "spherical_projection", ] diff --git a/kikuchipy/projections/gnomonic_projection.py b/kikuchipy/projections/gnomonic_projection.py index 09c66ca0..8c71c8bf 100644 --- a/kikuchipy/projections/gnomonic_projection.py +++ b/kikuchipy/projections/gnomonic_projection.py @@ -27,7 +27,9 @@ class GnomonicProjection(SphericalProjection): """Gnomonic projection of a vector as implemented in MTEX.""" + @classmethod def project(self, v: Union[Vector3d, np.ndarray]) -> np.ndarray: + """Convert from Cartesian to the Gnomonic projection.""" polar = super().project(v) theta, phi = polar[..., 0], polar[..., 1] @@ -53,6 +55,7 @@ def project(self, v: Union[Vector3d, np.ndarray]) -> np.ndarray: @staticmethod def iproject(x: np.ndarray, y: np.ndarray) -> Vector3d: + """Convert from the Gnomonic projection to Cartesian coordinates.""" theta = np.arctan(np.sqrt(x ** 2 + y ** 2)) - phi = np.atan2(y, x) + phi = np.arctan2(y, x) return Vector3d.from_polar(theta=theta, phi=phi) diff --git a/kikuchipy/projections/lambert_projection.py b/kikuchipy/projections/lambert_projection.py new file mode 100644 index 00000000..5c7d0527 --- /dev/null +++ b/kikuchipy/projections/lambert_projection.py @@ -0,0 +1,117 @@ +# -*- coding: utf-8 -*- +# Copyright 2019-2020 The kikuchipy developers +# +# This file is part of kikuchipy. +# +# kikuchipy is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# kikuchipy is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with kikuchipy. If not, see . + +from typing import Union + +import numpy as np +from orix.vector import Vector3d + +from kikuchipy.projections.gnomonic_projection import GnomonicProjection + + +class LambertProjection: + """Lambert projection of a vector [Callahan2013]_.""" + + @classmethod + def project(cls, v: Union[Vector3d, np.ndarray]) -> np.ndarray: + """Convert (n, 3) vector from Cartesian to the Lambert projection.""" + if isinstance(v, Vector3d): + w = v.unit.data + x = w[..., 0] + y = w[..., 1] + z = w[..., 2] + else: + norm = np.sqrt( + np.sum(np.square([v[..., 0], v[..., 1], v[..., 2]]), axis=0) + ) + x = v[..., 0] / norm + y = v[..., 1] / norm + z = v[..., 2] / norm + + # Arrays used in both setting X and Y + sqrt_z = np.sqrt(2 * (1 - z)) + sign_x = np.sign(x) + sign_y = np.sign(y) + abs_yx = abs(y) <= abs(x) + + # Reusable constants + sqrt_pi = np.sqrt(np.pi) + sqrt_pi_half = sqrt_pi / 2 + two_over_sqrt_pi = 2 / sqrt_pi + + # Equations 10a and 10b from Callahan and De Graef (2013) + X = np.where( + abs_yx, + sign_x * sqrt_z * sqrt_pi_half, + sign_y * sqrt_z * (two_over_sqrt_pi * np.arctan2(x, y)), + ) + Y = np.where( + abs_yx, + sign_x * sqrt_z * (two_over_sqrt_pi * np.arctan2(y, x)), + sign_y * sqrt_z * sqrt_pi_half, + ) + + return np.column_stack((X, Y)) + + @staticmethod + def iproject(xy: np.ndarray) -> Vector3d: + """Convert (n, 2) array from Lambert to Cartesian coordinates.""" + X = xy[..., 0] + Y = xy[..., 1] + + # Arrays used in setting x and y + true_term = Y * np.pi / (4 * X) + false_term = X * np.pi / (4 * Y) + abs_yx = abs(Y) <= abs(X) + c_x = _eq_c(X) + c_y = _eq_c(Y) + + # Equations 8a and 8b from Callahan and De Graef (2013) + x = np.where(abs_yx, c_x * np.cos(true_term), c_y * np.sin(false_term)) + y = np.where(abs_yx, c_x * np.sin(true_term), c_y * np.cos(false_term)) + z = np.where( + abs_yx, 1 - (2 * (X ** 2)) / np.pi, 1 - (2 * (Y ** 2)) / np.pi, + ) + + return Vector3d(np.column_stack((x, y, z))) + + @staticmethod + def lambert_to_gnomonic(xy: np.ndarray) -> np.ndarray: + """Convert (n,2) array from Lambert via Cartesian coordinates to + Gnomonic.""" + # These two functions could probably be combined into 1 to decrease + # runtime + vec = LambertProjection.iproject(xy) + xy = GnomonicProjection.project(vec) + return xy + + @staticmethod + def gnomonic_to_lambert(xy: np.ndarray) -> np.ndarray: + """Convert (n,2) array from Gnomonic via Cartesian coordinates to + Lambert.""" + # These two functions could probably be combined into 1 to decrease + # runtime + vec = GnomonicProjection.iproject(xy[..., 0], xy[..., 1]) + xy = LambertProjection.project(vec) + return xy + + +def _eq_c(p: Union[np.ndarray, float, int]) -> Union[np.ndarray, float, int]: + """Private function used inside LambertProjection.iproject to increase + readability.""" + return (2 * p / np.pi) * np.sqrt(np.pi - p ** 2) diff --git a/kikuchipy/projections/spherical_projection.py b/kikuchipy/projections/spherical_projection.py index d8ff3c65..d0de463d 100644 --- a/kikuchipy/projections/spherical_projection.py +++ b/kikuchipy/projections/spherical_projection.py @@ -28,6 +28,7 @@ class SphericalProjection: spherical_region = SphericalRegion([0, 0, 1]) + @classmethod def project(self, v: Union[Vector3d, np.ndarray]) -> np.ndarray: """Convert from cartesian to spherical coordinates.""" polar = get_polar(v) diff --git a/kikuchipy/projections/tests/__init__.py b/kikuchipy/projections/tests/__init__.py new file mode 100644 index 00000000..17cc82a1 --- /dev/null +++ b/kikuchipy/projections/tests/__init__.py @@ -0,0 +1,17 @@ +# -*- coding: utf-8 -*- +# Copyright 2019-2020 The kikuchipy developers +# +# This file is part of kikuchipy. +# +# kikuchipy is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# kikuchipy is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with kikuchipy. If not, see . diff --git a/kikuchipy/projections/tests/test_lambert_projection.py b/kikuchipy/projections/tests/test_lambert_projection.py new file mode 100644 index 00000000..706ddd23 --- /dev/null +++ b/kikuchipy/projections/tests/test_lambert_projection.py @@ -0,0 +1,92 @@ +# -*- coding: utf-8 -*- +# Copyright 2019-2020 The kikuchipy developers +# +# This file is part of kikuchipy. +# +# kikuchipy is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# kikuchipy is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with kikuchipy. If not, see . + +import numpy as np +from orix.vector import Vector3d +import pytest + +from kikuchipy.projections.lambert_projection import LambertProjection, _eq_c + + +class TestLambertProjection: + def test_project_vector3d(self): + """Works for Vector3d objects with single and multiple vectors""" + + vector_one = Vector3d((0.578, 0.578, 0.578)) + output_a = LambertProjection.project(vector_one) + expected_a = np.array((0.81417, 0.81417)) + + assert output_a[..., 0][0] == pytest.approx(expected_a[..., 0], rel=1e4) + assert output_a[..., 1][0] == pytest.approx(expected_a[..., 1], rel=1e4) + + vector_two = Vector3d( + np.array( + [[0.578, 0.578, 0.578], [0, 0.707, 0.707], [0.707, 0, 0.707]] + ) + ) + output_b = LambertProjection.project(vector_two) + + expected_b = np.array( + ( + (0.814172, 0.814172, 0.814172), + (0, 0, 0), + (0.6784123, 0, 0.6784123), + ) + ) + + assert output_a[..., 0][0] == pytest.approx(expected_a[..., 0], rel=1e4) + assert output_a[..., 1][0] == pytest.approx(expected_a[..., 1], rel=1e4) + + def test_project_ndarray(self): + "Works for numpy ndarrays" + ipt = np.array((0.578, 0.578, 0.578)) + output = LambertProjection.project(ipt) + expected = np.array((0.81417, 0.81417)) + assert output[..., 0][0] == pytest.approx(expected[..., 0], rel=1e4) + + def test_iproject(self): + """Conversion from Lambert to Cartesian coordinates works""" + vec = np.array((0.81417, 0.81417)) + expected = Vector3d((0.57702409, 0.577, 0.578)) + output = LambertProjection.iproject(vec) + assert np.allclose(output.x.data[0], expected.x.data[0], rtol=1e-05) + + def test_eq_c(self): + """Helper function works""" + arr_1 = np.array((0.578, 0.578, 0.578)) + value = arr_1[..., 0] + output_arr = np.array((0.61655, 0.61655, 0.61655)) + expected = output_arr[..., 0] + output = _eq_c(value) + assert output == pytest.approx(expected, rel=1e-5) + + def test_lambert_to_gnomonic(self): + """Conversion from Lambert to Gnomonic works""" + vec = np.array( + (0.81417, 0.81417) + ) # Should give x,y,z = 1/sqrt(3) (1, 1, 1) + output = LambertProjection.lambert_to_gnomonic(vec) + expected = np.array((1, 1)) + assert output[..., 0][0] == pytest.approx(expected[..., 0], abs=1e-2) + + def test_gnomonic_to_lambert(self): + """Conversion from Gnomonic to Lambert works""" + vec = np.array((1, 1)) # Similar to the case above + output = LambertProjection.gnomonic_to_lambert(vec) + expected = np.array((0.81417, 0.81417)) + assert output[..., 0][0] == pytest.approx(expected[..., 0], rel=1e-3) diff --git a/kikuchipy/release.py b/kikuchipy/release.py index 963de205..46ddee21 100644 --- a/kikuchipy/release.py +++ b/kikuchipy/release.py @@ -18,7 +18,11 @@ author = "kikuchipy developers" copyright = "Copyright 2019-2020, kikuchipy" -credits = ["Håkon Wiik Ånes", "Tina Bergh"] +credits = [ + "Håkon Wiik Ånes", + "Tina Bergh", + "Lars Andreas Hastad Lervik", +] license = "GPLv3+" maintainer = "Håkon Wiik Ånes" maintainer_email = "hakon.w.anes@ntnu.no"