Coverage for local_installation_linux/mumott/methods/basis_sets/gaussian_kernels.py: 93%
201 statements
« prev ^ index » next coverage.py v7.3.2, created at 2024-08-11 23:08 +0000
« prev ^ index » next coverage.py v7.3.2, created at 2024-08-11 23:08 +0000
1import logging
2from typing import Any, Dict, Iterator, Tuple
4import numpy as np
5from numpy.typing import NDArray
7from mumott import ProbedCoordinates, SphericalHarmonicMapper
8from mumott.core.hashing import list_to_hash
9from mumott.methods.utilities.tensor_operations import (framewise_contraction,
10 framewise_contraction_transpose)
11from .base_basis_set import BasisSet
12from .spherical_harmonics import SphericalHarmonics
15logger = logging.getLogger(__name__)
18class GaussianKernels(BasisSet):
19 r""" Basis set class for gaussian kernels, a simple local representation on the sphere.
20 The kernels follow a pseudo-even distribution similar to that described by
21 Y. Kurihara `in 1965 <https://doi.org/10.1175/1520-0493%281965%29093%3C0399%3ANIOTPE%3E2.3.CO%3B2>`_,
22 except with offsets added at the poles.
24 Notes
25 -----
26 The Gaussian kernel at location :math:`\rho_i` is given by
28 .. math ::
29 N_i \exp\left[ -\frac{1}{2} \left(\frac{d(\rho_i, r)}{\sigma}\right)^2 \right]
31 .. math ::
32 \sigma = \frac{\nu \pi}{2 (g + 1)}
34 where :math:`\nu` is the kernel scale parameter and :math:`g` is the grid scale, and
36 .. math ::
37 d(\rho, r) = \arctan_2(\Vert \rho \times r \Vert, \rho \cdot r),
39 that is, the great circle distance from the kernel location :math:`\rho` to the
40 probed location :math:`r`. If Friedel symmetry is assumed, the expression is instead
42 .. math ::
43 d(\rho, r) = \arctan_2(\Vert \rho \times r \Vert, \vert \rho \cdot r \vert)
45 The normalization factor :math:`\rho_i` is given by
47 .. math ::
48 N_i = \sum_j \exp\left[ -\frac{1}{2} \left( \frac{d(\rho_i, \rho_j)}{\sigma} \right)^2 \right]
50 where the sum goes over the coordinates of all grid points. This leads to an
51 approximately even spherical function, such that a set of coefficients which are all equal
52 is approximately isotropic, to the extent possible with respect to restrictions
53 imposed by grid resolution and scale parameter.
55 Parameters
56 ----------
57 probed_coordinates : ProbedCoordinates
58 Optional. A container with the coordinates on the sphere probed at each detector segment by the
59 experimental method. Its construction from the system geometry is method-dependent.
60 By default, an empty instance of :class:`mumott.ProbedCoordinates` is created.
61 grid_scale : int
62 The size of the coordinate grid on the sphere. Denotes the number of azimuthal rings between the
63 pole and the equator, where each ring has between ``2`` and ``2 * grid_scale`` points
64 along the azimuth.
65 kernel_scale_parameter : float
66 The scale parameter of the kernel in units of :math:`\frac{\pi}{2 (g + 1)}`, where
67 :math:`g` is ``grid_scale``.
68 enforce_friedel_symmetry : bool
69 If set to ``True``, Friedel symmetry will be enforced, using the assumption that points
70 on opposite sides of the sphere are equivalent.
71 kwargs
72 Miscellaneous arguments which relate to segment integrations can be
73 passed as keyword arguments:
75 integration_mode
76 Mode to integrate line segments on the reciprocal space sphere. Possible options are
77 ``'simpson'``, ``'midpoint'``, ``'romberg'``, ``'trapezoid'``.
78 ``'simpson'``, ``'trapezoid'``, and ``'romberg'`` use adaptive
79 integration with the respective quadrature rule from ``scipy.integrate``.
80 ``'midpoint'`` uses a single mid-point approximation of the integral.
81 Default value is ``'simpson'``.
82 n_integration_starting_points
83 Number of points used in the first iteration of the adaptive integration.
84 The number increases by the rule ``N`` ← ``2 * N - 1`` for each iteration.
85 Default value is 3.
86 integration_tolerance
87 Tolerance for the maximum relative error between iterations before the integral
88 is considered converged. Default is ``1e-5``.
89 integration_maxiter
90 Maximum number of iterations. Default is ``10``.
91 enforce_sparsity
92 If ``True``, makes matrix sparse by limiting the number of basis set elements
93 that can map to each segment. Default is ``False``.
94 sparsity_count
95 Number of basis set elements that can map to each segment,
96 if ``enforce_sparsity`` is set to ``True``. Default is ``3``.
97 """
98 def __init__(self,
99 probed_coordinates: ProbedCoordinates = None,
100 grid_scale: int = 4,
101 kernel_scale_parameter: float = 1.,
102 enforce_friedel_symmetry: bool = True,
103 **kwargs):
104 super().__init__(probed_coordinates, **kwargs)
105 self._probed_coordinates_hash = hash(self.probed_coordinates)
106 self._grid_scale = grid_scale
107 self._kernel_scale_parameter = kernel_scale_parameter
108 self._enforce_friedel_symmetry = enforce_friedel_symmetry
109 self._projection_matrix = self._get_integrated_projection_matrix()
111 def _get_kurihara_mesh(self, N) -> Tuple[NDArray, NDArray]:
112 phi = []
113 theta = []
114 for i in np.arange(N, -1, -1):
115 for j in np.arange(0, (2 * (i + 0.5))):
116 phi.append((j + 0.5) / (2 * (i + 0.5)) * np.pi)
117 phi.append((j + 0.5) / (2 * (i + 0.5)) * -np.pi)
118 theta.append((i + 0.5) / (N + 1) * np.pi / 2)
119 theta.append((i + 0.5) / (N + 1) * np.pi / 2)
120 theta = np.array(theta)
121 phi = np.mod(phi, 2 * np.pi)
123 if not self._enforce_friedel_symmetry: 123 ↛ 124line 123 didn't jump to line 124, because the condition on line 123 was never true
124 theta = np.concatenate((theta, np.arccos(-np.cos(theta))))
125 phi = np.concatenate((phi, phi))
126 return theta, phi
128 def get_inner_product(self,
129 u: NDArray,
130 v: NDArray) -> NDArray:
131 r""" Retrieves the inner product of two coefficient arrays, that is to say,
132 the sum-product over the last axis.
134 Parameters
135 ----------
136 u
137 The first coefficient array, of arbitrary shape and dimension, so long as
138 the number of coefficients equals the length of this :class:`GaussianKernels` instance.
139 v
140 The second coefficient array, of the same shape as :attr:`u`.
141 """
142 assert u.shape[-1] == len(self)
143 assert u.shape == v.shape
144 return np.einsum('...i, ...i -> ...', u, v, optimize='greedy')
146 def _get_spherical_distances(self,
147 theta_1: NDArray[float], theta_2: NDArray[float],
148 phi_1: NDArray[float], phi_2: NDArray[float],
149 radius: float = 1.,
150 enforce_friedel_symmetry: bool = False) -> NDArray[float]:
151 r""" Function for obtaining the distances between two point sets
152 on a sphere, possibly with Friedel symmetry enforced.
153 Arrays can have any shape, but they must all be broadcastable, and the polar angles of
154 each set must have the same shape as the azimuthal angles.
155 If the first and second set of points have the same shape, then the distances will be
156 computed pointwise. Otherwise, the distances will be computed by broadcasting.
158 Parameters
159 ----------
160 theta_1
161 The polar angle of the first set of points, defined as :math:`\arccos(z_1)`.
162 theta_2
163 The polar angle fo the second set of points, defined as :math:`\arccos(z_2)`.
164 phi_2
165 The azimuthal angle of the first set of points, defined as :math:`\arctan_2(y_1, x_1)`.
166 phi_2
167 The azimuthal angle of the second set of points, defined as :math:`\arctan_2(y_2, x_2)`.
168 radius
169 The radius of the sphere. Default is `1`, i.e., describing a unit sphere.
170 enforce_friedel_symmetry
171 If ``True`` (default), the point :math:`(x, y, z)` will be considered as
172 equivalent to the point :math:`(-x, -y, -z)` and the maximum possible distance on the sphere will
173 be :math:`\frac{\sqrt{x^2 + y^2 + z^2} \pi}{2}`, i.e., at a 90-degree angle.
174 Otherwise, the two points will be considered distinct and the maximum
175 distance will be :math:`\sqrt{x^2 + y^2 + z^2} \pi`.
176 """
177 phi_diff = abs(phi_1 - phi_2)
178 sine_factor = np.sin(theta_1) * np.sin(theta_2) * np.cos(phi_diff)
179 cosine_factor = np.cos(theta_1) * np.cos(theta_2)
180 if enforce_friedel_symmetry: 180 ↛ 183line 180 didn't jump to line 183, because the condition on line 180 was never false
181 return radius * np.arccos(abs(sine_factor + cosine_factor).clip(0., 1.))
182 else:
183 return radius * np.arccos((sine_factor + cosine_factor).clip(-1., 1.))
185 def _get_projection_matrix(self, probed_coordinates: ProbedCoordinates = None) -> NDArray[float]:
186 """ Computes the matrix necessary for forward and gradient calculations.
187 Called when the coordinate system has been updated, or one of
188 ``kernel_scale_parameter`` or ``grid_scale`` has been changed."""
189 if probed_coordinates is None: 189 ↛ 190line 189 didn't jump to line 190, because the condition on line 189 was never true
190 probed_coordinates = self.probed_coordinates
191 theta, phi = self._get_kurihara_mesh(self._grid_scale)
192 phi = phi.reshape(1, 1, 1, -1)
193 theta = theta.reshape(1, 1, 1, -1)
194 _, probed_polar_angles, probed_azim_angles = probed_coordinates.to_spherical
195 probed_polar_angles = probed_polar_angles[..., np.newaxis]
196 probed_azim_angles = probed_azim_angles[..., np.newaxis]
197 # Find distances to all probed detector points on sphere
198 distances = self._get_spherical_distances(theta, probed_polar_angles,
199 phi, probed_azim_angles,
200 enforce_friedel_symmetry=self._enforce_friedel_symmetry)
202 # Probe at location of each kernel function to normalize over sphere
203 mesh_distances = self._get_spherical_distances(
204 theta.reshape(-1, 1), theta.reshape(1, -1),
205 phi.reshape(-1, 1), phi.reshape(1, -1),
206 enforce_friedel_symmetry=self._enforce_friedel_symmetry)
207 std = (self._kernel_scale_parameter * np.sqrt(2 * np.pi)) / (2 * (self._grid_scale + 1))
208 matrix = np.exp(-(1 / 2) * (distances / std) ** 2)
209 norm_matrix = np.exp(-(1 / 2) * (mesh_distances / std) ** 2)
210 # The normalization factor is the inverse of the unnormalized function value at each grid point
211 norm_factors = np.reciprocal(norm_matrix.sum(-1).reshape(1, 1, 1, -1))
212 return matrix * norm_factors
214 def forward(self,
215 coefficients: NDArray,
216 indices: NDArray = None) -> NDArray:
217 """ Carries out a forward computation of projections from Gaussian kernel space
218 into detector space, for one or several tomographic projections.
220 Parameters
221 ----------
222 coefficients
223 An array of coefficients, of arbitrary shape so long as the last
224 axis has the same size as :attr:`~.GaussianKernels.kernel_scale_parameter`, and if
225 :attr:`indices` is ``None`` or greater than one, the first axis should have the
226 same length as :attr:`indices`
227 indices
228 Optional. Indices of the tomographic projections for which the forward
229 computation is to be performed. If ``None``, the forward computation will
230 be performed for all projections.
232 Returns
233 -------
234 An array of values on the detector corresponding to the :attr:`coefficients` given.
235 If :attr:`indices` contains exactly one index, the shape is ``(coefficients.shape[:-1], J)``
236 where ``J`` is the number of detector segments. If :attr:`indices` is ``None`` or contains
237 several indices, the shape is ``(N, coefficients.shape[1:-1], J)`` where ``N``
238 is the number of tomographic projections for which the computation is performed.
240 Notes
241 -----
242 The assumption is made in this implementation that computations over several
243 indices act on sets of images from different projections. For special usage
244 where multiple projections of entire fields are desired, it may be better
245 to use :attr:`projection_matrix` directly. This also applies to
246 :meth:`gradient`.
247 """
248 assert coefficients.shape[-1] == len(self)
249 self._update()
250 output = np.zeros(coefficients.shape[:-1] + (self._projection_matrix.shape[1],),
251 coefficients.dtype)
252 if indices is None:
253 framewise_contraction_transpose(self._projection_matrix,
254 coefficients,
255 output)
256 elif indices.size == 1: 256 ↛ 264line 256 didn't jump to line 264, because the condition on line 256 was never false
257 np.einsum('ijk, ...k -> ...j',
258 self._projection_matrix[indices],
259 coefficients,
260 out=output,
261 optimize='greedy',
262 casting='unsafe')
263 else:
264 framewise_contraction_transpose(self._projection_matrix[indices],
265 coefficients,
266 output)
267 return output
269 def gradient(self,
270 coefficients: NDArray,
271 indices: NDArray = None) -> NDArray:
272 """ Carries out a gradient computation of projections from Gaussian kernel space
273 into detector space for one or several tomographic projections.
275 Parameters
276 ----------
277 coefficients
278 An array of coefficients (or residuals) of arbitrary shape so long as the last
279 axis has the same size as the number of detector segments.
280 indices
281 Optional. Indices of the tomographic projections for which the gradient
282 computation is to be performed. If ``None``, the gradient computation will
283 be performed for all projections.
285 Returns
286 -------
287 An array of gradient values based on the :attr:`coefficients` given.
288 If :attr:`indices` contains exactly one index, the shape is ``(coefficients.shape[:-1], J)``
289 where ``J`` is the number of detector segments. If indices is ``None`` or contains
290 several indices, the shape is ``(N, coefficients.shape[1:-1], J)`` where ``N``
291 is the number of tomographic projections for which the computation is performed.
293 Notes
294 -----
295 When solving an inverse problem, one should not attempt to optimize the
296 coefficients directly using the gradient one obtains by applying this method to the data.
297 Instead, one must either take the gradient of the residual between the
298 :meth:`~.GaussianKernels.forward` computation of the coefficients and the data.
299 Alternatively one can apply both the forward and the gradient computation to the
300 coefficients to be optimized, and the gradient computation to the data, and treat
301 the residual of the two as the gradient of the optimization coefficients. The approaches
302 are algebraically equivalent, but one may be more efficient than the other in some
303 circumstances. However, normally, the projection between detector and
304 ``GaussianKernel`` space is only a small part of the overall computation,
305 so there is typically not much to be gained from optimizing it.
306 """
307 self._update()
308 output = np.zeros(coefficients.shape[:-1] + (self._projection_matrix.shape[2],),
309 coefficients.dtype)
310 if indices is None:
311 framewise_contraction(self._projection_matrix,
312 coefficients,
313 output)
314 elif indices.size == 1: 314 ↛ 322line 314 didn't jump to line 322, because the condition on line 314 was never false
315 np.einsum('ikj, ...k -> ...j',
316 self._projection_matrix[indices],
317 coefficients,
318 out=output,
319 optimize='greedy',
320 casting='unsafe')
321 else:
322 framewise_contraction(self._projection_matrix[indices],
323 coefficients,
324 output)
325 return output
327 def get_amplitudes(self, coefficients: NDArray[float],
328 probed_coordinates: ProbedCoordinates = None) -> NDArray[float]:
329 """ Computes the amplitudes of the spherical function represented by the provided
330 :attr:`coefficients` at the :attr:`probed_coordinates`.
332 Parameters
333 ----------
334 coefficients
335 An array of coefficients of arbitrary shape, provided that the
336 last dimension contains the coefficients for one spherical function.
337 probed_coordinates
338 An instance of :class:`mumott.core.ProbedCoordinates` with its :attr:`vector`
339 attribute indicating the points of the sphere for which to evaluate the amplitudes.
340 """
341 if probed_coordinates is None: 341 ↛ 342line 341 didn't jump to line 342, because the condition on line 341 was never true
342 probed_coordinates = self._probed_coordinates
343 matrix = self._get_projection_matrix(probed_coordinates)
344 self._make_projection_matrix_sparse(matrix)
345 return np.einsum('ij, ...j', matrix.squeeze(), coefficients, optimize='greedy')
347 def get_spherical_harmonic_coefficients(self, coefficients: NDArray, ell_max: int = None):
348 """ Computes the spherical harmonic coefficients of the spherical function
349 represented by the provided :attr:`coefficients` using a Driscoll-Healy grid.
351 For details on the Driscoll-Healy grid, see
352 `the SHTools page <https://shtools.github.io/SHTOOLS/grid-formats.html>`_ for a
353 comprehensive overview.
355 Parameters
356 ----------
357 coefficients
358 An array of coefficients of arbitrary shape, provided that the
359 last dimension contains the coefficients for one function.
360 ell_max
361 The bandlimit of the spherical harmonic expansion. By default, it is ``2 * grid_scale``.
363 """
364 if ell_max is None: 364 ↛ 366line 364 didn't jump to line 366, because the condition on line 364 was never false
365 ell_max = 2 * self._grid_scale
366 dh_grid_size = 4 * (self._grid_scale + 1)
367 mapper = SphericalHarmonicMapper(ell_max=ell_max, polar_resolution=dh_grid_size,
368 azimuthal_resolution=dh_grid_size,
369 enforce_friedel_symmetry=self._enforce_friedel_symmetry)
370 coordinates = ProbedCoordinates(mapper.unit_vectors.reshape((1, -1, 1, 3)))
371 amplitudes = self.get_amplitudes(coefficients, coordinates)\
372 .reshape(coefficients.shape[:-1] + (dh_grid_size, dh_grid_size))
373 spherical_harmonics_coefficients = mapper.get_harmonic_coefficients(amplitudes)
374 return spherical_harmonics_coefficients
376 def get_output(self,
377 coefficients: NDArray) -> Dict[str, Any]:
378 r""" Returns a dictionary of output data for a given array of basis set coefficients.
380 Parameters
381 ----------
382 coefficients
383 An array of coefficients of arbitrary shape and dimensions, except
384 its last dimension must be the same length as the :attr:`len` of this instance.
385 Computations only operate over the last axis of :attr:`coefficients`, so derived
386 properties in the output will have the shape ``(*coefficients.shape[:-1], ...)``.
388 Returns
389 -------
390 A dictionary containing two sub-dictionaries, ``basis_set`` and ``spherical_harmonic_analysis``.
391 ``basis_set`` contains information particular to :class:`GaussianKernels`, whereas
392 ``spherical_harmonic_analysis`` contains an analysis of the spherical function using
393 a spherical harmonic transform.
395 Notes
396 -----
397 In detail, the two sub-dictionaries ``basis_set`` and ``spherical_harmonic_analysis``
398 have the following members:
400 basis_set
401 name
402 The name of the basis set, i.e., ``'GaussianKernels'``
403 coefficients
404 A copy of :attr:`coefficients`.
405 grid_scale
406 A copy of :attr:`~.GaussianKernels.grid_scale`.
407 kernel_Scale_paramter
408 A copy of :attr:`~.GaussianKernels.kernel_scale_parameter`.
409 enforce_friedel_symmetry
410 A copy of :attr:`~.GaussianKernels.enforce_friedel_symmetry`.
411 projection_matrix
412 A copy of :attr:`~.GaussianKernels.projection_matrix`.
413 spherical_harmonic_analysis
414 An analysis of the spherical function in terms of
415 spherical harmonics. See
416 :meth:`SphericalHarmonics.get_output <SphericalHarmonics.get_output()>` for details.
417 """
418 assert coefficients.shape[-1] == len(self)
419 # Update to ensure non-dirty output state.
420 self._update()
421 sph_coefficients = self.get_spherical_harmonic_coefficients(coefficients)
422 spherical_harmonics_set = SphericalHarmonics(ell_max=2 * self._grid_scale,
423 enforce_friedel_symmetry=self._enforce_friedel_symmetry)
424 sph_output = spherical_harmonics_set.get_output(sph_coefficients)
425 output_dictionary = {}
427 # basis set-specific information
428 basis_set = {}
429 output_dictionary['basis_set'] = basis_set
430 basis_set['name'] = type(self).__name__
431 basis_set['coefficients'] = coefficients.copy()
432 basis_set['grid_scale'] = self._grid_scale
433 basis_set['grid'] = self.grid
434 basis_set['kernel_scale_parameter'] = self._kernel_scale_parameter
435 basis_set['enforce_friedel_symmetry'] = self._enforce_friedel_symmetry
436 basis_set['projection_matrix'] = self._projection_matrix.copy()
437 basis_set['hash'] = hex(hash(self))
439 # Analysis is more easily done in spherical harmonic space.
440 output_dictionary['spherical_harmonic_analysis'] = sph_output
442 return output_dictionary
444 def __iter__(self) -> Iterator[Tuple[str, Any]]:
445 """ Allows class to be iterated over and in particular be cast as a dictionary.
446 """
447 yield 'name', type(self).__name__
448 yield 'grid_scale', self._grid_scale
449 yield 'kernel_scale_parameter', self._kernel_scale_parameter
450 yield 'enforce_friedel_symmetry', self._enforce_friedel_symmetry
451 yield 'projection_matrix', self._projection_matrix
452 yield 'hash', hex(hash(self))[2:]
454 def __len__(self) -> int:
455 return self._projection_matrix.shape[-1]
457 def __hash__(self) -> int:
458 """Returns a hash reflecting the internal state of the instance.
460 Returns
461 -------
462 A hash of the internal state of the instance,
463 cast as an ``int``.
464 """
465 to_hash = [self._grid_scale,
466 self.grid,
467 self._kernel_scale_parameter,
468 self._enforce_friedel_symmetry,
469 self._projection_matrix,
470 self._probed_coordinates_hash]
471 return int(list_to_hash(to_hash), 16)
473 def _update(self) -> None:
474 # We only run updates if the hashes do not match.
475 if self.is_dirty:
476 self._projection_matrix = self._get_integrated_projection_matrix()
477 self._probed_coordinates_hash = hash(self._probed_coordinates)
479 @property
480 def is_dirty(self) -> bool:
481 return hash(self._probed_coordinates) != self._probed_coordinates_hash
483 @property
484 def projection_matrix(self) -> NDArray:
485 """ The matrix used to project spherical functions from the unit sphere onto the detector.
486 If ``v`` is a vector of gaussian kernel coefficients, and ``M`` is the ``projection_matrix``,
487 then ``M @ v`` gives the corresponding values on the detector segments associated with
488 each projection. ``M[i] @ v`` gives the values on the detector segments associated with
489 projection ``i``.
491 If ``r`` is a residual between a projection from Gaussian kernel to detector space and data from
492 projection ``i``, then ``M[i].T @ r`` gives the associated gradient in Gaussian kernel space.
493 """
494 self._update()
495 return self._projection_matrix
497 @property
498 def grid_scale(self) -> int:
499 """ The number of azimuthal rings from each pole to the equator in the
500 spherical grid.
501 """
502 return self._grid_scale
504 @grid_scale.setter
505 def grid_scale(self, val: int) -> None:
506 if val < 0 or val != round(val):
507 raise ValueError('grid_scale must be a non-negative integer,'
508 f' but a value of {val} was given!')
509 self._grid_scale = val
510 self._projection_matrix = self._get_integrated_projection_matrix()
512 @property
513 def kernel_scale_parameter(self) -> float:
514 """ The scale parameter for each kernel.
515 """
516 return self._kernel_scale_parameter
518 @kernel_scale_parameter.setter
519 def kernel_scale_parameter(self, val: float) -> float:
520 self._kernel_scale_parameter = val
521 self._projection_matrix = self._get_integrated_projection_matrix()
523 @property
524 def enforce_friedel_symmetry(self) -> bool:
525 """ If ``True``, Friedel symmetry is enforced, i.e., the point
526 :math:`-r` is treated as equivalent to :math:`r`. """
527 return self._enforce_friedel_symmetry
529 @property
530 def grid(self) -> Tuple[NDArray['float'], NDArray['float']]:
531 r""" Returns the polar and azimuthal angles of the grid used by the basis.
533 Returns
534 -------
535 A ``Tuple`` with contents ``(polar_angle, azimuthal_angle)``, where the
536 polar angle is defined as :math:`\arccos(z)`.
537 """
538 return self._get_kurihara_mesh(self._grid_scale)
540 @property
541 def grid_hash(self) -> str:
542 """ Returns a hash of :attr:`grid`.
543 """
544 return list_to_hash([self.grid])
546 @property
547 def projection_matrix_hash(self) -> str:
548 """ Returns a hash of :attr:`projection_matrix`.
549 """
550 return list_to_hash([self.projection_matrix])
552 def __str__(self) -> str:
553 wdt = 74
554 s = []
555 s += ['-' * wdt]
556 s += ['GaussianKernels'.center(wdt)]
557 s += ['-' * wdt]
558 with np.printoptions(threshold=4, edgeitems=2, precision=5, linewidth=60):
559 s += ['{:18} : {}'.format('grid_scale', self.grid_scale)]
560 s += ['{:18} : {}'.format('grid_hash', self.grid_hash)]
561 s += ['{:18} : {}'.format('enforce_friedel_symmetry', self.enforce_friedel_symmetry)]
562 s += ['{:18} : {}'.format('kernel_scale_parameter', self.kernel_scale_parameter)]
563 s += ['{:18} : {}'.format('projection_matrix_hash', self.projection_matrix_hash)]
564 s += ['{:18} : {}'.format('hash', hex(hash(self))[2:8])]
565 s += ['-' * wdt]
566 return '\n'.join(s)
568 def _repr_html_(self) -> str:
569 s = []
570 s += [f'<h3>{self.__class__.__name__}</h3>']
571 s += ['<table border="1" class="dataframe">']
572 s += ['<thead><tr><th style="text-align: left;">Field</th><th>Size</th><th>Data</th></tr></thead>']
573 s += ['<tbody>']
574 with np.printoptions(threshold=4, edgeitems=2, precision=2, linewidth=40):
575 s += ['<tr><td style="text-align: left;">grid_scale</td>']
576 s += [f'<td>1</td><td>{self.grid_scale}</td></tr>']
577 s += ['<tr><td style="text-align: left;">grid_hash</td>']
578 s += [f'<td>{len(self.grid_hash)}</td><td>{self.grid_hash[:6]}</td></tr>']
579 s += ['<tr><td style="text-align: left;">kernel_scale_parameter</td>']
580 s += [f'<td>1</td><td>{self.kernel_scale_parameter}</td></tr>']
581 s += ['<tr><td style="text-align: left;">enforce_friedel_symmetry</td>']
582 s += [f'<td>1</td>'
583 f'<td>{self.enforce_friedel_symmetry}</td></tr>']
584 s += ['<tr><td style="text-align: left;">projection_matrix</td>']
585 s += [f'<td>{len(self.projection_matrix_hash)}</td>'
586 f'<td>{self.projection_matrix_hash[:6]}</td></tr>']
587 s += ['<tr><td style="text-align: left;">hash</td>']
588 s += [f'<td>{len(hex(hash(self)))}</td><td>{hex(hash(self))[2:8]}</td></tr>']
589 s += ['</tbody>']
590 s += ['</table>']
591 return '\n'.join(s)