Coverage for local_installation_linux/mumott/output_handling/projection_viewer.py: 71%
136 statements
« prev ^ index » next coverage.py v7.3.2, created at 2024-08-11 23:08 +0000
1""" Container for the class ProjectionViewer. """
3import logging
4from dataclasses import dataclass
6import numpy as np
7import matplotlib.pyplot as plt
8import colorcet # noqa
10from colorspacious import cspace_converter
11from matplotlib import cm
12from matplotlib.colors import Colormap
13from matplotlib.gridspec import GridSpec
14from numpy.typing import NDArray
15from typing import Union, Tuple
17from .orientation_image_mapper import OrientationImageMapper
18from mumott.data_handling import DataContainer
21@dataclass
22class ProjectionImage:
23 """
24 A dataclass for storing images for instances of
25 class:`ProjectionViewer <mumott.output_handling.ProjectionViewer>`.
27 Attributes
28 ----------
29 mean
30 The mean value of the data.
31 std
32 The standard deviation of the data.
33 angle_color_quad
34 A quadruplet of RGBA colors mapping phase angles
35 that indicate the orientation of the image.
36 """
37 mean: NDArray[float]
38 std: NDArray[float]
39 angle_color_quad: NDArray[float]
42class ProjectionViewer:
43 """
44 Class for viewing data and synthetic projections in a scrollable plot.
45 The viewed projection can be changed with either the arrow keys or
46 with the mouse scrolling wheel.
48 Parameters
49 ----------
50 data_container : DataContainer
51 A :class:`DataContainer
52 <mumott.data_handling.data_container.DataContainer>` object, containing
53 the necessary dimension and rotation information to show orientations correctly.
54 data : numpy.ndarray(dtype=numpy.float64), optional
55 A flat array of data from e.g. reconstruction output. If not given, data
56 will be loaded from the provided ``data_container``. Provide if you e.g.
57 want to look at synthetic projections from a reconstruction.
58 orientation_symmetry : ``'longitudinal'``, ``'transversal'``, optional
59 Specifies the assumed three-dimensional symmetry of the reciprocal
60 space map, which determines whether the minimum or the maximum of the
61 ``cos(phi) ** 2`` component of the of the intensity is interpreted as
62 the orientation angle. The default is ``'longitudinal'``, which interprets the
63 angle of maximum scattering as the orientation, corresponding to
64 an assumed three-dimensional symmetry of polar caps of intensity.
65 The other option is ``transversal``, which corresponds to a symmetry
66 of a great circle of intensity (also known as fiber symmetry), corresponding
67 to a phase shift of ``90`` degrees.
68 mean_colormap : str or matplotlib.colors.Colormap
69 Specifies the colormap of the mean intensity across all angles. The default
70 is ``'cet_linear_bmy_10_95_c71'``, a linear high-contrast with blue, magenta
71 and yellow.
72 std_colormap : str or matplotlib.colors.Colormap
73 Specifies the colormap of the standar deviation of the data across all angles.
74 The default is ``'cet_gouldian'``, a linear high-contrast map with blue,
75 green and yellow.
76 phase_colormap : str or matplotlib.colors.Colormap
77 Specifies the colormap of the phase of the ``cos(detector_angles) ** 2`` fit
78 of the data. The default is ``'cet_CET_C10'``, an isoluminant cyclic map.
79 mean_range : tuple(float, float)
80 Specifies the range of the colormap of the mean.
81 std_range : tuple(float, float)
82 Specifies the range of the colormap of the standard deviation.
83 """
85 def __init__(self,
86 data_container: DataContainer,
87 data: NDArray[np.float64] = None,
88 orientation_symmetry: str = 'longitudinal',
89 mean_colormap: Union[str, Colormap] = 'cet_linear_bmy_10_95_c71',
90 std_colormap: Union[str, Colormap] = 'cet_gouldian',
91 phase_colormap: Union[str, Colormap] = 'cet_CET_C10',
92 mean_range: Tuple[float, float] = None,
93 std_range: Tuple[float, float] = None):
94 detector_angles = data_container.geometry.detector_angles
95 det_vector = \
96 (np.cos(detector_angles[np.newaxis, :]) *
97 data_container.geometry.detector_direction_origin[:, np.newaxis] +
98 np.sin(detector_angles[np.newaxis, :]) *
99 data_container.geometry.detector_direction_positive_90[:, np.newaxis])
100 distance_j = np.dot(data_container.geometry.j_direction_0[np.newaxis, :],
101 det_vector).squeeze()
102 distance_k = np.dot(data_container.geometry.k_direction_0[np.newaxis, :],
103 det_vector).squeeze()
104 if not np.allclose(np.diff(detector_angles).reshape(1, -1), 104 ↛ 107line 104 didn't jump to line 107, because the condition on line 104 was never true
105 np.diff(detector_angles).reshape(-1, 1),
106 rtol=2e-01, atol=8e-02):
107 logging.warning('The provided `detector_angles` do not appear to be '
108 'sorted and equally spaced to within an absolute tolerance '
109 'of 5 degrees and a relative tolerance of twenty percent of the '
110 'segment width.'
111 '\nThe orientation angle may not behave as expected. Please be cautious'
112 'in interpreting the figure.')
113 detector_angles = np.arctan2(distance_k, distance_j)
114 if data is None: 114 ↛ 116line 114 didn't jump to line 116, because the condition on line 114 was never false
115 data = data_container.data
116 if orientation_symmetry == 'transversal': 116 ↛ 117line 116 didn't jump to line 117, because the condition on line 116 was never true
117 detector_angles = np.arctan2(np.sin(detector_angles), np.cos(detector_angles))
118 elif orientation_symmetry == 'longitudinal': 118 ↛ 121line 118 didn't jump to line 121, because the condition on line 118 was never false
119 detector_angles = np.arctan2(-np.cos(detector_angles), np.sin(detector_angles))
120 else:
121 logging.warning(
122 'Unknown `orientation_symmetry` option: "' + orientation_symmetry +
123 '"\nUsing default value: "longitudinal".\n'
124 'Valid options are "transversal" and "longitudinal".')
125 detector_angles = np.arctan2(-np.cos(detector_angles), np.sin(detector_angles))
126 direction = np.arctan2(np.sin(detector_angles[1] - detector_angles[0]),
127 np.cos(detector_angles[1] - detector_angles[0]))
128 a = float(np.round(direction >= 0))
129 b = np.exp((-2j) * detector_angles[0])
130 self._mean_colormap = mean_colormap
131 self._std_colormap = std_colormap
132 self._phase_colormap = phase_colormap
133 self._cielab_converter = cspace_converter('sRGB255', 'JCh')
134 self._srgb_converter = cspace_converter('JCh', 'sRGB255')
135 self._orientation_image_mapper = OrientationImageMapper(phase_colormap)
136 self._wheel_properties = self._orientation_image_mapper.wheel_properties
137 self._phase_colormap_lut = cm.get_cmap(self._phase_colormap, 256)
138 self._projection_number = len(data)
139 self._clims_mean = np.percentile(
140 a=data.mean(axis=2),
141 q=(1., 99.9))
142 if mean_range is not None: 142 ↛ 143line 142 didn't jump to line 143, because the condition on line 142 was never true
143 if mean_range[0] is not None:
144 self._clims_mean[0] = mean_range[0]
145 if mean_range[1] is not None:
146 self._clims_mean[1] = mean_range[1]
147 self._clims_std = np.percentile(
148 a=data.std(axis=2),
149 q=(1., 99.9))
150 if std_range is not None: 150 ↛ 151line 150 didn't jump to line 151, because the condition on line 150 was never true
151 if std_range[0] is not None:
152 self._clims_std[0] = std_range[0]
153 if std_range[1] is not None:
154 self._clims_std[1] = std_range[1]
155 self._projection_images = []
156 for i, frame in enumerate(data):
157 temp_mean = frame[:, :, :detector_angles.size].mean(axis=2).T
158 temp_std = frame[:, :, :detector_angles.size].std(axis=2).T
159 temp_fft = np.fft.rfft(frame[:, :, :detector_angles.size] ** 2, axis=2, norm='ortho')
160 if a: 160 ↛ 163line 160 didn't jump to line 163, because the condition on line 160 was never false
161 temp_angle = (np.sqrt(temp_fft[:, :, 1] * b)).T
162 else:
163 temp_angle = (np.sqrt(np.conj(temp_fft[:, :, 1]) * b)).T
164 temp_angle = np.angle(np.sign(temp_angle.imag) * temp_angle)
165 temp_mask = np.round(
166 np.sqrt(
167 abs(temp_fft[:, :, 0])) > 0.5 * np.sqrt(abs(temp_fft[:, :, 0])).mean()).astype(float).T
168 temp_aniso_factor = np.sqrt(
169 abs(temp_fft[:, :, 1]) / (np.finfo(float).tiny + abs(temp_fft[:, :, 0]))).T
170 temp_shader = temp_aniso_factor.clip(0., 1)
171 phase_rgba = self._phase_colormap_lut(temp_angle / np.pi)
172 phase_cielab = self._cielab_converter(phase_rgba[:, :, :-1])
173 phase_cielab[:, :, :2] *= 1.1 * ((1 - np.exp(-temp_shader)) /
174 (1 - np.exp(-1))).reshape(temp_shader.shape + (1,))
175 phase_rgba[:, :, :-1] = self._srgb_converter(phase_cielab).clip(0., 1.)
176 phase_rgba[:, :, -1] = temp_mask
177 temp_image = ProjectionImage(
178 mean=temp_mean,
179 std=temp_std,
180 angle_color_quad=phase_rgba)
181 self._projection_images.append(temp_image)
182 self._fig = plt.figure()
183 self._grid_spec = GridSpec(6, 9, figure=self._fig)
184 self._axes = []
185 self._axes.append(self._fig.add_subplot(self._grid_spec[0:4, 0:3]))
186 self._axes.append(self._fig.add_subplot(self._grid_spec[0:4, 3:6]))
187 self._axes.append(self._fig.add_subplot(self._grid_spec[0:4, 6:9]))
188 self._axes.append(self._fig.add_subplot(self._grid_spec[4:6, 6:9], projection='polar'))
189 self._colorbar_axes = []
190 self._colorbar_axes.append(self._fig.add_subplot(self._grid_spec[5, 0:3]))
191 self._colorbar_axes.append(self._fig.add_subplot(self._grid_spec[5, 3:6]))
192 for a in self._axes:
193 a.grid(False)
194 a.set_yticks([])
195 a.set_xticks([])
196 a.set_yticklabels([])
197 a.set_xticklabels([])
198 self._im1 = self._axes[0].imshow(
199 self._projection_images[0].mean, origin='lower', cmap=self._mean_colormap)
200 self._axes[0].set_title('Mean')
201 self._fig.colorbar(
202 self._im1, cax=self._colorbar_axes[0], orientation='horizontal')
203 self._im2 = self._axes[1].imshow(
204 self._projection_images[0].std, origin='lower', cmap=self._std_colormap)
205 self._axes[1].set_title('Standard deviation')
206 self._fig.colorbar(self._im2, cax=self._colorbar_axes[1], orientation='horizontal')
208 self._im3 = self._axes[2].imshow(self._projection_images[0].angle_color_quad,
209 origin='lower',)
210 self._axes[2].set_title('Orientation phase')
211 self._axes[2].set_facecolor((254 / 255, 254 / 255, 254 / 255))
212 self._wheel_rgba = self._wheel_properties[2]
213 self._axes[3].set_facecolor((254 / 255, 254 / 255, 254 / 255))
214 wheel_shader = self._wheel_properties[1][:-1, :-1] ** 2
215 wheel_cielab = self._cielab_converter(self._wheel_rgba[:, :, :-1])
216 wheel_cielab[:, :, :2] *= 1.1 * ((1 - np.exp(-wheel_shader)) /
217 (1 - np.exp(-1))).reshape(wheel_shader.shape + (1,))
218 wheel_rgb = self._srgb_converter(wheel_cielab).clip(0., 1.)
219 self._wheel_rgba[:, :, :-1] = wheel_rgb
220 self._axes[3].grid(False)
221 self._color_wheel = self._axes[3].pcolormesh(self._wheel_properties[0],
222 1 - np.arccos(self._wheel_properties[1]) * 2 / np.pi,
223 wheel_shader[:, :],
224 facecolors=self._wheel_rgba.reshape(-1, 4),
225 edgecolors=self._wheel_rgba.reshape(-1, 4))
226 self._axes[3].set_rmax(1.)
227 self._axes[3].set_rmin(0.)
228 self._projection_ind = 0
230 def change_projection(self, ind: int) -> None:
231 """
232 Method for updating the projection shown
233 when the projection slider changes.
235 Parameters
236 ----------
237 ind
238 The index of the projection to be shown.
239 """
240 if ind != self._projection_ind: 240 ↛ 241line 240 didn't jump to line 241, because the condition on line 240 was never true
241 self._projection_ind = ind
242 ind = int(ind)
243 self._im1.set_data(self._projection_images[ind].mean)
244 self._im1.set_clim(*self._clims_mean)
245 self._im2.set_data(self._projection_images[ind].std)
246 self._im2.set_clim(*self._clims_std)
247 self._im3.set_data(self._projection_images[ind].angle_color_quad)
248 self._im3.autoscale()
249 self._axes[0].relim()
250 self._axes[1].relim()
251 self._axes[2].relim()
252 self._axes[0].autoscale_view()
253 self._axes[1].autoscale_view()
254 self._axes[2].autoscale_view()
255 self._axes[0].set_aspect('equal')
256 self._axes[1].set_aspect('equal')
257 self._axes[2].set_aspect('equal')
258 self._fig.canvas.draw_idle()
259 self._fig.canvas.flush_events()