///|
/// A data-axis slice, with the two remaining axes in ascending order.
pub(all) struct Slice {
  width : Int
  height : Int
  values : Array[Double]
}

///|
pub fn Image::slice(
  self : Image,
  axis : Int,
  index : Int,
  t? : Int = 0,
  scaled? : Bool = true,
) -> Slice raise NiftiError {
  require(axis >= 0 && axis < 3, "slice_axis")
  require(index >= 0 && index < self.shape[axis], "slice_index")
  require(t >= 0 && t < self.shape[3], "time_index")
  let axes = [0, 1, 2].filter(a => a != axis)
  let width = self.shape[axes[0]]
  let height = self.shape[axes[1]]
  let values = []
  for y in 0.. Statistics raise NiftiError {
  let accumulator = SampleAccumulator::new()
  for i in 0.. GridComparison raise NiftiError {
  require(finite(tolerance_mm) && tolerance_mm >= 0, "grid_tolerance")
  for a in 0..<3 {
    if self.shape[a] != other.shape[a] {
      return {
        compatible: false,
        reason: "spatial_shape",
        maximum_corner_error_mm: None,
      }
    }
  }
  fn code(image : Image, space : Space) -> Int {
    match space {
      QForm => image.qcode
      SForm => image.scode
      PreferSForm => if image.scode > 0 { image.scode } else { image.qcode }
    }
  }
  let ac = code(self, space)
  let bc = code(other, space)
  if ac == 0 || bc == 0 {
    return {
      compatible: false,
      reason: "unknown_space",
      maximum_corner_error_mm: None,
    }
  }
  if ac != bc {
    return {
      compatible: false,
      reason: "coordinate_system_code",
      maximum_corner_error_mm: None,
    }
  }
  fn unit(code : Int) -> Double {
    match code {
      1 => 1000.0
      2 => 1.0
      3 => 0.001
      _ => 0.0
    }
  }
  let au = unit(self.spatial_unit)
  let bu = unit(other.spatial_unit)
  if au == 0 || bu == 0 {
    return {
      compatible: false,
      reason: "unknown_units",
      maximum_corner_error_mm: None,
    }
  }
  let mut error = 0.0
  for c in 0..<8 {
    let p = Array::makei(3, a => {
      if ((c >> a) & 1) == 0 {
        0.0
      } else {
        // A singleton dimension still has a voxel basis; compare one step.
        (if self.shape[a] > 1 { self.shape[a] - 1 } else { 1 }).to_double()
      }
    })
    let av = self.world(p[0], p[1], p[2], space~)
    let bv = other.world(p[0], p[1], p[2], space~)
    // Maximum coordinate-component displacement, not Euclidean distance.
    for a in 0..<3 {
      let delta = (av[a] * au - bv[a] * bu).abs()
      require(finite(delta), "coordinate_overflow")
      if delta > error {
        error = delta
      }
    }
  }
  {
    compatible: error <= tolerance_mm,
    reason: if error <= tolerance_mm {
      "same_spatial_grid"
    } else {
      "world_coordinates"
    },
    maximum_corner_error_mm: Some(error),
  }
}