python/dpt/volumes.py

1"""Sampled scalar-field metadata; a volume is not a bag of tissue labels.23The sample buffer is contiguous, x-fastest (nz, ny, nx), in inverse mm for4attenuation or dimensionless for an explicitly declared material fraction.5No CT-to-attenuation conversion is inferred here.6"""78from __future__ import annotations910from dataclasses import dataclass11from typing import cast1213from dpt.contracts import ContractError14from dpt.geometry import IDENTITY, Matrix3, Vector3, matvec, rotation_matrix, transpose, vector3151617# region book:volume-sampling-contract18@dataclass(frozen=True, slots=True)19class GridSpec:20    """Oriented anisotropic grid, with origin at the first sample centre.2122    Support extends half a sample spacing outside the outer centres. Within23    that support interpolation clamps to the outer sample; beyond it the field24    is zero. Nonzero boundary samples therefore produce a discontinuous field.25    """2627    shape: tuple[int, int, int]28    spacing_mm: Vector329    origin_mm: Vector3 = (0.0, 0.0, 0.0)30    orientation: Matrix3 = IDENTITY3132    def __post_init__(self) -> None:33        if len(self.shape) != 3 or any(type(x) is not int or x < 1 for x in self.shape):34            raise ContractError("grid shape is (nz, ny, nx), each a positive integer")35        object.__setattr__(self, "shape", tuple(self.shape))36        object.__setattr__(self, "spacing_mm", vector3(self.spacing_mm, "spacing_mm"))37        object.__setattr__(self, "origin_mm", vector3(self.origin_mm, "origin_mm"))38        object.__setattr__(self, "orientation", rotation_matrix(self.orientation, "orientation"))39        if min(self.spacing_mm) <= 0:40            raise ContractError("grid spacing must be positive")41        if self.voxels > 2**31 - 1:42            raise ContractError("grid exceeds signed 32-bit indexing")4344    @property45    def voxels(self) -> int:46        return self.shape[0] * self.shape[1] * self.shape[2]4748    def object_to_grid(self, point_mm: Vector3) -> Vector3:49        displacement: Vector3 = cast(50            Vector3, tuple(point_mm[i] - self.origin_mm[i] for i in range(3))51        )52        aligned = matvec(transpose(self.orientation), displacement)53        return cast(Vector3, tuple(aligned[i] / self.spacing_mm[i] for i in range(3)))5455    def grid_to_object(self, index: Vector3) -> Vector3:56        scaled: Vector3 = cast(Vector3, tuple(index[i] * self.spacing_mm[i] for i in range(3)))57        aligned = matvec(self.orientation, scaled)58        return cast(Vector3, tuple(aligned[i] + self.origin_mm[i] for i in range(3)))5960    @property61    def support(self) -> tuple[Vector3, Vector3]:62        return (-0.5, -0.5, -0.5), (self.shape[2] - 0.5, self.shape[1] - 0.5, self.shape[0] - 0.5)636465# endregion book:volume-sampling-contract66