sgt.core module

Core routines for scattering geometry calculations.

class sgt.core.Axis(*values)[source]

Bases: Enum

Spatial axis identifiers

X = 0
Y = 1
Z = 2
sgt.core.make_coords_in_detector_system(px_numbers: tuple[int, int], px_sizes: tuple[float, float], center_coord: tuple[float, float]) tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Makes the matrices of coordinates of each pixel on the 2D detector coordinate system.

The detector coordinate system is a 2D Cartesian coordinate system whose origin is at the image center (= where the direct beam hits the detector plane). Two axes, denoted as u and v, are defined to be parallel to the horizontal and vertical edges of the detector, respectively.

Parameters:
  • px_numbers – Numbers of pixels along the horizontal and vertical axes (hor, ver).

  • px_sizes – Size of a single pixel along the horizontal and vertical axes (hor, ver).

  • center_coord – (x, y) coordinate of the image center, measured from the center of a pixel at index [0,0].

Returns:

2D numpy arrays of the horizontal and vertical coordinates.

Example

>>> u, v = make_coords_in_detector_system( (1475, 1679), (0.172, 0.172), (132.35, 134.39) )
sgt.core.make_rotation_matrix_around_lab_axis(axis: Axis, angle_deg: float) ndarray[tuple[int, ...], dtype[float64]][source]

Makes a rotation matrix representing a rotation around an axis in the lab coordinate system.

Parameters:
  • axis – the axis around which the rotation is performed. Rotation around x means to rotate y axis toward z axis. Rotation around y means to rotate z axis toward x axis. Rotation around z means to rotate x axis toward y axis.

  • angle_deg – Rotation angle in degrees.

Returns:

A 3x3 rotation matrix.

Example

>>> R = make_rotation_matrix_around_lab_axis(sgt.Axis.Y, 19.0)
sgt.core.make_rotation_matrix_euler_xyx(alpha_deg: float, beta_deg: float, gamma_deg: float) ndarray[tuple[int, ...], dtype[float64]][source]

Makes a rotation matrix that defines orientation of the detector.

The rotation matrix is calculated based on the given Euler angles. The Euler angles here is the x-y-x type, that is,

  1. The initial axes \((x, y, z)\) are rotated by alpha_deg around the \(x\) axis.

  2. The resultant axes \((x', y', z')\) are rotated by beta_deg around the \(y'\) axis.

  3. The resultant axes \((x'', y'', z'')\) are rotated by gamma_deg around the \(x''\) axis.

Parameters:
  • alpha_deg – Rotation around the \(x\) axis, in degrees.

  • beta_deg – Rotation around the \(y'\) axis, in degrees.

  • gamma_deg – Rotation around the \(x''\) axis, in degrees.

Returns:

A 3x3 numpy array.

Example

>>> R = make_rotation_matrix(10.0, 10.0, 10.0)
>>> vec = np.array([1.0, 2.0, 3.0])
>>> rotated = np.matmul(R, vec)
sgt.core.make_default_mask(hor_px_num: int, ver_px_num: int) ndarray[tuple[int, ...], dtype[int8]][source]

Makes a default mask array.

A mask array is a 2D array of the same shape as the scattering image but of numpy.uint8 type. Pixels to be masked are assigned with 1 and unmasked pixels are assigned with zero.

Parameters:
  • hor_px_num – Number of pixels along the horizontal edge.

  • ver_px_num – Number of pixels along the horizontal edge.

Returns:

A 2D numpy array of the dtype numpy.uint8.

sgt.core.make_default_array(hor_px_num: int, ver_px_num: int) ndarray[tuple[int, ...], dtype[float64]][source]

Makes a default float array.

Parameters:
  • hor_px_num – Number of pixels along the horizontal edge.

  • ver_px_num – Number of pixels along the horizontal edge.

Returns:

A 2D numpy array of the float type.

sgt.core.make_detector_basis_vectors_in_lab_system(rotation_matrix: ndarray[tuple[int, ...], dtype[float64]]) tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Makes the basis vectors of the detector coordinate system, expressed in the lab coordinate system.

For the definition of the detector coordinate system, refer to make_coords_in_detector_system(). The basis vectors of the detector coordinate system is basically the basis vectors of the lab coordinate system being rotated by the rotation matrix that defines the detector orientation. The input rotation matrix can be created by make_rotation_matrix(). (but actually can be any SO(3) matrix)

Parameters:

rotation_matrix – A 3x3 numpy array representing the detector orientation.

Returns:

Three numpy arrays eu, ev, and en representing the basis vectors. eu and ev are the basis vector along the horizontal and vertical edge of the detector, respectively, and en is the one perpendicular to the detector plane.

sgt.core.transform_detector_to_lab(u: ndarray[tuple[int, ...], dtype[float64]], v: ndarray[tuple[int, ...], dtype[float64]], eu: ndarray[tuple[int, ...], dtype[float64]], ev: ndarray[tuple[int, ...], dtype[float64]], sdd: float) tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Make the lab system coordinate from the detector system coordinate

Parameters:
  • u – horizontal coordinates in the detector system. Must be of the same shape as v.

  • v – vertical coordinates in the detector system. Must be of the same shape as u.

  • eu – basis vector of u axis on the detector, expressed in the lab system.

  • ev – basis vector of v axis on the detector, expressed in the lab system.

  • sdd – sample to detector distance.

Returns:

three arrays, each with the same shape as u and v, representing coordinates in the lab system.

sgt.core.calc_shortest_distance_to_detector(eu: ndarray[tuple[int, ...], dtype[float64]], ev: ndarray[tuple[int, ...], dtype[float64]], en: ndarray[tuple[int, ...], dtype[float64]], sdd: float) float[source]

Computes the shortest distance from the sample to the detector plane.

Let P be the point on the detector such that OP is the shortest distance between the origin and the detector plane. The vector OP must be perpendicular to the detector plane, so the vector \(\vec{OP}\) is proportional to the detector plane normal vector \(\vec{e}_n\). That is,

\[\vec{OP} = (OP)\vec{e}_n\]

The vector OP can also be expressed as

\[\vec{OP} = u_P\vec{e}_u + v_P\vec{e}_v + L\vec{e}_z\]

where \((u_P, v_P)\) is the coordinate of point P on the detector coordinate system and vector \(\vec{e}_z\) is the z basis vector. \(L\) is the sample-to-detector distance.

Equating the two expressions,

\[k\vec{e}_n = u_P\vec{e}_u + v_P\vec{e}_v + L\vec{e}_z\]

which reads

\[\begin{split}u_P (\vec{e}_u)_x + v_P (\vec{e}_v)_x - (OP) (\vec{e}_n)_x &= 0 \\ u_P (\vec{e}_u)_y + v_P (\vec{e}_v)_y - (OP) (\vec{e}_n)_y &= 0 \\ u_P (\vec{e}_u)_z + v_P (\vec{e}_v)_z - (OP) (\vec{e}_n)_z &= -L\end{split}\]

By defining a matrix

\[\begin{split}\mathbf{M} = \begin{pmatrix} (\vec{e}_u)_x & (\vec{e}_v)_x & -(\vec{e}_n)_x \\ (\vec{e}_u)_y & (\vec{e}_v)_y & -(\vec{e}_n)_y \\ (\vec{e}_u)_z & (\vec{e}_v)_z & -(\vec{e}_n)_z \\ \end{pmatrix}\end{split}\]

the equations are simplified to

\[\begin{split}\mathbf{M} \vec{s} &= -L \vec{e}_z \\ \vec{s} &= -L \mathbf{M}^{-1}\vec{e}_z\end{split}\]

where \(\vec{s} = (u_P, v_P, OP)\). This method computes \(\vec{s}\) using the above equation and returns its third component, \(OP\).

Parameters:
  • eu – basis vector of the detector coordinate system along the horizontal edge of the detector.

  • ev – basis vector of the detector coordinate system along the vertical edge of the detector.

  • en – basis vector of the detector coordinate system along the plane normal of the detector.

  • sdd – sample-to-detector distance.

Returns:

The distance in the float value.

Note

The input vectors should be expressed in the lab coordinate system. Use make_basis_vectors_on_detector_in_lab_system() to generate the basis vectors.

sgt.core.make_solid_angle_coverage_correction_factors(x: ndarray[tuple[int, ...], dtype[float64]], y: ndarray[tuple[int, ...], dtype[float64]], z: ndarray[tuple[int, ...], dtype[float64]], shortest_dist_to_detector: float) ndarray[tuple[int, ...], dtype[float64]][source]

Makes an array with the correction factor for solid angle coverage of each pixel.

Based on Equation 28 in Pauw J Phys.: Condens. Matter 25, 383201. DOI: 10.1088/0953-8984/25/38/383201.

Parameters:
  • x – x coordinates of pixels in the lab system.

  • y – y coordinates of pixels in the lab system.

  • z – z coordinates of pixels in the lab system.

  • shortest_dist_to_detector – shortest distance from the sample to the detector plane. see calc_shortest_dist_to_detector().

Returns:

A 2D array of correction factors for each pixel. The correction factor is normalized at the beam center. The correction can be done by multiplying this array to the intensity array.

sgt.core.make_q_vector(x: ndarray[tuple[int, ...], dtype[float64]], y: ndarray[tuple[int, ...], dtype[float64]], z: ndarray[tuple[int, ...], dtype[float64]], wavelength: float) tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Computes q vector.

By definition,

\[\vec{q} = \dfrac{2 \pi}{\lambda}(\vec{e}_\mathrm{s} - \vec{e}_\mathrm{i})\]

where \(\lambda\) is the wavelength, \(\vec{e}_\mathrm{s}\) is the basis vector along the scattered ray, and \(\vec{e}_\mathrm{i}\) is the basis vector along the incident ray.

Here, \(\vec{e}_\mathrm{i}\) is fixed to (0, 0, 1). Since the sample is placed at the origin,

\[\vec{e}_\mathrm{s} = \dfrac{\vec{r}}{|\vec{r}|}\]

where \(\vec{r}\) is the coordinate of the pixel in the lab system.

Parameters:
  • x – x coordinates of pixels in the lab system.

  • y – y coordinates of pixels in the lab system.

  • z – z coordinates of pixels in the lab system.

  • wavelength – wavelength of the incident beam.

Returns:

Three arrays representing x, y, and z components of the q vectors.

sgt.core.make_qabs(qx: ndarray[tuple[int, ...], dtype[float64]], qy: ndarray[tuple[int, ...], dtype[float64]], qz: ndarray[tuple[int, ...], dtype[float64]]) ndarray[tuple[int, ...], dtype[float64]][source]

Computes the magnitude of the q vectors.

Parameters:
  • qx – x components of the q vectors.

  • qy – y components of the q vectors.

  • qz – z components of the q vectors.

Returns:

An array of the same shape as qx, qy, and qz,

sgt.core.make_azimuthal_angle_rad(qx: ndarray[tuple[int, ...], dtype[float64]], qy: ndarray[tuple[int, ...], dtype[float64]]) ndarray[tuple[int, ...], dtype[float64]][source]

Make azimuthal angle of the q vectors.

Parameters:
  • qx – x components of the q vectors.

  • qy – y components of the q vectors.

Returns:

An array of the same shape as qx and qy.

sgt.core.make_bin_indices(bin_edges: ndarray[tuple[int, ...], dtype[float64]], binning_coords: ndarray[tuple[int, ...], dtype[float64]], masked_pixels: ndarray[tuple[int, ...], dtype[any]]) tuple[ndarray[tuple[int, ...], dtype[int64]], ndarray[tuple[int, ...], dtype[int64]]][source]

Make index mapping for binning.

Parameters:
  • bin_edges – Edges of each bin. Must be sorted in an ascending order.

  • binning_coords – Coordinates (positions) used for binning.

  • masked_pixels – Truth array or index array that specifies the pixels to mask.

Returns:

Two int arrays, indices_in_bin and bin_counts. indices_in_bin is an array with the same size as binning_coords. Each element is an index (starting from 1) in the bin to which the pixel belongs. The index 0 means that the pixel is masked or out of the bin range. bin_counts is an array with the same size as bin_edges. bin_counts[i+1] is the number of pixels in i`th bin. `bin_counts[0] is the number of pixels that was masked or out of the bin range.

sgt.core.binning(indices_in_bin: ndarray[tuple[int, ...], dtype[int64]], bin_counts: ndarray[tuple[int, ...], dtype[int64]], intensities: ndarray[tuple[int, ...], dtype[float64]], errors: ndarray[tuple[int, ...], dtype[float64]]) tuple[ndarray[tuple[int, ...], dtype[float64]], ndarray[tuple[int, ...], dtype[float64]]][source]

Bin intensity & error using binning indices

Parameters:
  • indices_in_bin – Binning index array, that can be generated by make_bin_indices().

  • bin_counts – Pixel counts in each bin, that can be generated by make_bin_indices().

  • intensities – Intensity array.

  • errors – Errors of intensity. Can be zero-filled array in case the error output is not necessary.

Returns:

Binned arrays of intensity and error.