Unverified Commit 95aeea76 authored by Gigg, Martyn Anthony's avatar Gigg, Martyn Anthony Committed by GitHub
Browse files

Merge pull request #30445 from...

Merge pull request #30445 from mantidproject/30433_correct_transform_ellipsoid_when_axes_swap_sliceviewer

Apply correct transform to ellipsoid eigenvectors when axes swapped in sliceviewer
parents 4819915b 68981f74
......@@ -73,6 +73,7 @@ Bugfixes
- Fixed bug in plotting elliptical shell of integrated peaks in sliceviewer - the inner background radius is now correct.
- Fixed a bug in error bars tab in plot settings where the Error Every property was not being shown correctly
- Fixed a bug where the fit action (Fit > Fit) in the fit browser wasn't disabled if all the functions were individually removed.
- Fixed a bug in sliveviewer that wasn't transforming ellipsoid axes of integrated peaks correctly when axes swapped.
:ref:`Release 6.0.0 <v6.0.0>`
......@@ -33,11 +33,11 @@ class EllipsoidalIntergratedPeakRepresentation():
"""
shape_info = json_loads(peak_shape.toJSON())
# use signal ellipse to see if it is valid at this slice
axes, signal_radii = _signal_ellipsoid_info(shape_info, slice_info.transform)
axes, signal_radii = _signal_ellipsoid_info(shape_info)
peak_origin = slice_info.transform(peak_origin)
slice_origin, major_radius, minor_radius, angle = slice_ellipsoid(
peak_origin, *axes, *signal_radii, slice_info.z_value)
peak_origin, *axes, *signal_radii, slice_info.z_value, slice_info.transform)
if not np.any(np.isfinite((major_radius, minor_radius))):
# slice not possible
return None
......@@ -69,11 +69,11 @@ class EllipsoidalIntergratedPeakRepresentation():
]
# add background shell
a, b, c, shell_thick = _bkgd_ellipsoid_info(shape_info, slice_info.transform)
a, b, c, shell_thick = _bkgd_ellipsoid_info(shape_info)
# only add background shell if it is greater than 0
if shell_thick > 0:
_, major_radius, minor_radius, angle = slice_ellipsoid(peak_origin, *axes, a, b, c,
slice_info.z_value)
slice_info.z_value, slice_info.transform)
bkgd_width, bkgd_height = 2 * major_radius, 2 * minor_radius
artists.append(
......@@ -93,7 +93,7 @@ class EllipsoidalIntergratedPeakRepresentation():
return Painted(painter, artists)
def _signal_ellipsoid_info(shape_info, transform):
def _signal_ellipsoid_info(shape_info):
"""
Retrieve axes and radii from the PeakShape in the slice frame
:param shape_info: A dictionary of ellipsoid properties
......@@ -108,27 +108,28 @@ def _signal_ellipsoid_info(shape_info, transform):
"direction2")
a, b, c = float(shape_info["radius0"]), float(shape_info["radius1"]), float(
shape_info["radius2"])
return transform((axis_a, axis_b, axis_c)), transform((a, b, c))
return (axis_a, axis_b, axis_c), (a, b, c)
def _bkgd_ellipsoid_info(shape_info, transform):
def _bkgd_ellipsoid_info(shape_info):
"""
Retrieve background radii and width from the PeakShape in the slice frame
:param shape_info: A dictionary of ellipsoid properties
:param transform: Callable to move to the slice frame
"""
a, b, c = float(shape_info["background_outer_radius0"]), float(
shape_info["background_outer_radius1"]), float(shape_info["background_outer_radius2"])
inner_a, inner_b, inner_c = float(shape_info["background_inner_radius0"]), float(
shape_info["background_inner_radius1"]), float(shape_info["background_inner_radius2"])
a, b, c = transform((a, b, c))
inner_a, inner_b, inner_c = transform((inner_a, inner_b, inner_c))
width = (max((a, b, c)) - max((inner_a, inner_b, inner_c))) / max((a, b, c)) # fractional width of shell
if min([a, b, c]) > 1e-15:
inner_a, inner_b, inner_c = float(shape_info["background_inner_radius0"]), float(
shape_info["background_inner_radius1"]), float(shape_info["background_inner_radius2"])
width = (max((a, b, c)) - max((inner_a, inner_b, inner_c))) / max((a, b, c)) # fractional width of shell
else:
# no background specified (inner bg defaults to peak radius) assign unphysical width
width = -1
return (a, b, c, width)
return a, b, c, width
def slice_ellipsoid(origin, axis_a, axis_b, axis_c, a, b, c, zp):
def slice_ellipsoid(origin, axis_a, axis_b, axis_c, a, b, c, zp, transform):
"""Return semi-axes lengths and angle of ellipse formed
by cutting the ellipsoid with a plane at z=zp
:param origin: Origin of ellipsoid
......@@ -139,6 +140,7 @@ def slice_ellipsoid(origin, axis_a, axis_b, axis_c, a, b, c, zp):
:param b: Axis 1 radius
:param c: Axis 2 radius
:param zp: Slice point in slice dimension
:param transform: Callable to move to the slice frame
:return (slice_origin, major_radius, minor_radius, angle):
"""
# From https://en.wikipedia.org/wiki/Ellipsoid#As_quadric:
......@@ -170,10 +172,10 @@ def slice_ellipsoid(origin, axis_a, axis_b, axis_c, a, b, c, zp):
#
# c = m22*zk^2
return slice_ellipsoid_matrix(origin, zp,
calculate_ellipsoid_matrix(axis_a, axis_b, axis_c, a, b, c))
calculate_ellipsoid_matrix(axis_a, axis_b, axis_c, a, b, c, transform))
def calculate_ellipsoid_matrix(axis_a, axis_b, axis_c, a, b, c):
def calculate_ellipsoid_matrix(axis_a, axis_b, axis_c, a, b, c, transform):
"""
Compute matrix A defining ellipsoid
https://en.wikipedia.org/wiki/Ellipsoid#As_quadric
......@@ -183,17 +185,16 @@ def calculate_ellipsoid_matrix(axis_a, axis_b, axis_c, a, b, c):
:param a: Axis 0 radius
:param b: Axis 1 radius
:param c: Axis 2 radius
:param transform: Callable to move to the slice frame
"""
# Create matrix whose eigenvalues are squares of semi-axes lengths
# and eigen vectors define principle axis vectors
axes_lengths = np.diag((1 / a ** 2, 1 / b ** 2, 1 / c ** 2))
# yapf: disable
axes_dir = np.array((
(axis_a[0], axis_b[0], axis_c[0]),
(axis_a[1], axis_b[1], axis_c[1]),
(axis_a[2], axis_b[2], axis_c[2])
))
# yapf: enable
# transformation matrix (inverse gives transform from MD to view basis)
R = np.vstack((transform((1, 0, 0)), transform((0, 1, 0)), transform((0, 0, 1))))
# matrix with eigenvectors (in view basis) in columns
axes_dir = linalg.inv(R) @ np.hstack((axis_a[:, None], axis_b[:, None], axis_c[:, None]))
return axes_dir @ axes_lengths @ np.transpose(axes_dir)
......@@ -201,7 +202,7 @@ def slice_ellipsoid_matrix(origin, zp, ellipMatrix):
"""
:param origin: Origin of the ellipsoid
:param zp: Slice point in the slicing dimension
:param M: Ellipsoid Matrix. See `calculate_ellipsoid_matrix` for definition
:param ellipMatrix: Ellipsoid Matrix. See `calculate_ellipsoid_matrix` for definition
"""
# slice point in frame with ellipsoid at origin
z = zp - origin[2]
......@@ -235,11 +236,14 @@ def slice_ellipsoid_matrix(origin, zp, ellipMatrix):
MM = A / (0.25 * np.transpose(B) @ linalg.inv(A) @ B - (c - 1))
try:
eigvalues, eigvectors = linalg.eig(MM)
isort = np.argsort(eigvalues) # smallest eigenvalue corresponds to largest axis length
eigvalues = eigvalues[isort]
eigvectors = eigvectors[:, isort]
except np.linalg.LinAlgError:
# Sometimes the integration volume may not intersect the slices of data
return origin, np.nan, np.nan, 0
minor_radius, major_radius = 1 / np.sqrt(eigvalues[0]), 1 / np.sqrt(eigvalues[1])
major_axis = eigvectors[:, 1]
minor_radius, major_radius = 1 / np.sqrt(eigvalues[1]), 1 / np.sqrt(eigvalues[0])
major_axis = eigvectors[:, 0]
angle = np.rad2deg(np.arctan2(major_axis[1], major_axis[0]))
slice_origin_xy = -0.5 * linalg.inv(A) @ B
slice_origin = np.array(
......
......@@ -24,9 +24,10 @@ class EllipticalShell(Patch):
"""
def __str__(self):
return f"EllipticalShell(center={self.center}, width={self.width}, height={self.height}, thick={self.thick}, angle={self.angle})"
return f"EllipticalShell(center={self.center}, width={self.width}, height={self.height}, " \
f"frac_thick={self.frac_thick}, angle={self.angle})"
def __init__(self, center, width, height, thick, angle=0.0, **kwargs):
def __init__(self, center, width, height, frac_thick, angle=0.0, **kwargs):
"""
Draw an elliptical ring centered at *x*, *y* center with outer width (horizontal diameter)
*width* and outer height (vertical diameter) *height* with a fractional ring thickness of *thick*
......@@ -37,7 +38,7 @@ class EllipticalShell(Patch):
super().__init__(**kwargs)
self.center = center
self.height, self.width = height, width
self.thick = thick
self.frac_thick = frac_thick
self.angle = angle
self._recompute_path()
# Note: This cannot be calculated until this is added to an Axes
......@@ -48,7 +49,7 @@ class EllipticalShell(Patch):
arc = Path.arc(theta1=0.0, theta2=360.0)
# Draw the outer unit circle followed by a reversed and scaled inner circle
v1 = arc.vertices
v2 = arc.vertices[::-1] * float(1.0 - self.thick) # self.thick is fractional thickness
v2 = arc.vertices[::-1] * float(1.0 - self.frac_thick)
v = np.vstack([v1, v2, v1[0, :], (0, 0)])
c = np.hstack([arc.codes, arc.codes, Path.MOVETO, Path.CLOSEPOLY])
c[len(arc.codes)] = Path.MOVETO
......@@ -140,7 +141,7 @@ class MplPainter():
"""
return self.axes.add_patch(Ellipse((x, y), width, height, angle, **kwargs))
def elliptical_shell(self, x, y, outer_width, outer_height, thick, angle=0.0, **kwargs):
def elliptical_shell(self, x, y, outer_width, outer_height, frac_thick, angle=0.0, **kwargs):
"""Draw an ellipse at the given location
:param x: X coordinate of the center
:param y: Y coordinate of the center
......@@ -151,7 +152,7 @@ class MplPainter():
:param kwargs: Additional matplotlib properties to pass to the call
"""
return self.axes.add_patch(
EllipticalShell((x, y), outer_width, outer_height, thick, angle, **kwargs))
EllipticalShell((x, y), outer_width, outer_height, frac_thick, angle, **kwargs))
def shell(self, x, y, outer_radius, thick, **kwargs):
"""Draw a wedge on the Axes
......
......@@ -19,11 +19,19 @@ from mantidqt.widgets.sliceviewer.peaksviewer.representation.test.shapetesthelpe
import FuzzyMatch, create_ellipsoid_info, draw_representation
def create_test_ellipsoid():
def create_test_ellipsoid(bg_shell=True):
radii = [1.5, 1.3, 1.4]
if bg_shell:
frac_thick = 0.1 # fractional thickness of major axis
bg_outer = [2 * r for r in radii]
bg_inner = [(1 - frac_thick) * r for r in bg_outer]
else:
bg_outer = 3 * [0]
bg_inner = radii
return create_ellipsoid_info(
radii=(1.5, 1.3, 1.4),
radii=radii,
axes=("1 0 0", "0 1 0", "0 0 1"),
bkgd_radii=(((2.5, 2.3, 2.4)), ((2.6, 2.4, 2.5))))
bkgd_radii=((bg_inner), (bg_outer)))
@patch("mantidqt.widgets.sliceviewer.peaksviewer.representation.ellipsoid.compute_alpha")
......@@ -42,7 +50,7 @@ class EllipsoidalIntergratedPeakRepresentationTest(unittest.TestCase):
painter.ellipse.assert_not_called()
painter.elliptical_shell.assert_not_called()
def test_draw_creates_ellipse_with_expected_properties_with_nonzero_alpha(
def test_draw_creates_ellipse_with_expected_properties_with_nonzero_alpha_and_background(
self, compute_alpha_mock):
peak_center = (1, 2, 3)
ellipsoid = create_test_ellipsoid()
......@@ -58,17 +66,39 @@ class EllipsoidalIntergratedPeakRepresentationTest(unittest.TestCase):
self._assert_painter_calls(
painter,
peak_center[:2],
cross_width=0.26,
signal_width=2.6,
signal_height=3.0,
angle=90,
cross_width=0.3,
signal_width=3.0,
signal_height=2.6,
angle=0,
alpha=fake_alpha,
fg_color=fg_color,
bkgd_width=4.8,
bkgd_width=6.0,
bkgd_height=5.2,
thickness=0.1 / 2.6,
thickness=0.1,
bg_color=bg_color)
def test_draw_creates_ellipse_with_expected_properties_with_nonzero_alpha_no_background(self, compute_alpha_mock):
peak_center = (1, 2, 3)
ellipsoid = create_test_ellipsoid(bg_shell=False)
painter = MagicMock()
fg_color, bg_color = 'r', 'g'
fake_alpha = 0.5
compute_alpha_mock.return_value = fake_alpha
painted = draw_representation(EllipsoidalIntergratedPeakRepresentation, peak_center,
ellipsoid, painter, fg_color, bg_color)
self.assertTrue(painted is not None)
self._assert_painter_calls(
painter,
peak_center[:2],
cross_width=0.3,
signal_width=3.0,
signal_height=2.6,
angle=0,
alpha=fake_alpha,
fg_color=fg_color)
def test_draw_respects_transform(self, compute_alpha_mock):
def slice_transform(x):
# set slice(x)=data(z)
......@@ -87,15 +117,15 @@ class EllipsoidalIntergratedPeakRepresentationTest(unittest.TestCase):
self.assertTrue(painted is not None)
self._assert_painter_calls(
painter, (peak_center[2], peak_center[0]),
cross_width=0.182,
signal_width=1.82,
signal_height=2.1,
cross_width=0.192,
signal_width=1.92,
signal_height=1.79,
angle=90,
alpha=fake_alpha,
fg_color=fg_color,
bkgd_width=4.4,
bkgd_height=4.77,
thickness=0.1 / 2.6,
bkgd_width=5.54,
bkgd_height=5.17,
thickness=0.1,
bg_color=bg_color)
# private
......@@ -147,58 +177,63 @@ class EllipsoidalIntergratedPeakRepresentationTest(unittest.TestCase):
class EllipsoidalIntergratedPeakRepresentationSliceEllipsoidTest(unittest.TestCase):
def test_slice_ellipsoid_zp(self):
origin = (0, 0, 0)
axis_a = (1, 0, 0)
axis_b = (0, 1, 0)
axis_c = (0, 0, 1)
axis_a = np.array([1, 0, 0])
axis_b = np.array([0, 1, 0])
axis_c = np.array([0, 0, 1])
a = 1
b = 1
b = 2
c = 1
zp = 0
def slice_transform(x):
return (x[1], x[0], x[2])
expected_slice_origin = (0, 0, 0)
expected_major_radius = 1
expected_minar_radius = 1
expected_angle = 90
expected_major_radius = 2
expected_minor_radius = 1
expected_angle = 0
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp),
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp, slice_transform),
(*expected_slice_origin,
expected_major_radius,
expected_minar_radius,
expected_minor_radius,
expected_angle))
zp = np.sin(np.pi / 3)
expected_slice_origin = (0, 0, zp)
expected_major_radius = 0.5 # cos(pi/3)
expected_minar_radius = 0.5 # cos(pi/3)
expected_major_radius = 1.0 # 2*cos(pi/3)
expected_minor_radius = 0.5 # cos(pi/3)
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp),
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp, slice_transform),
(*expected_slice_origin,
expected_major_radius,
expected_minar_radius,
expected_minor_radius,
expected_angle))
# This causes negative eignevalues there np.sqrt(eignevalues) gives NaN radius
zp = 2
expected_slice_origin = (0, 0, zp)
expected_major_radius = np.nan
expected_minar_radius = np.nan
expected_minor_radius = np.nan
expected_angle = 90 # flips to 90 as x-axis (Y on view) has negative eigenvalue and is first in sorted list
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp),
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp, slice_transform),
(*expected_slice_origin,
expected_major_radius,
expected_minar_radius,
expected_minor_radius,
expected_angle))
# This causes `eigvalues, eigvectors = linalg.eig(MM)` to throw np.linalg.LinAlgError
zp = 1
expected_slice_origin = (0, 0, 0)
expected_major_radius = np.nan
expected_minar_radius = np.nan
expected_angle = 0
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp),
expected_minor_radius = np.nan
expected_angle = 0 # returns 0 is could not diag ellipMatrix
self._run_slice_ellipsoid_and_compare((origin, axis_a, axis_b, axis_c, a, b, c, zp, slice_transform),
(*expected_slice_origin,
expected_major_radius,
expected_minar_radius,
expected_minor_radius,
expected_angle))
def _run_slice_ellipsoid_and_compare(self, input_values, expectecd):
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment