wim/ndbioimage
1
0 Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

- testing 2d transforms

Browse Source
This commit is contained in:
Wim Pomp committed 2026-10-10 12:44:37 +02:00
1 parent 4815690621
commit 453237cae9
20 files changed
+798 -304
Download Patch File Download Diff File Expand all files Collapse all files

No files matched your search

Cargo.toml
+2 -2
View File
@@ -3,9 +3,9 @@ resolver = "3"
members = ["ndbioimage", "image-tiff", "image-registration"]
[workspace.package]
version = "0.2.0"
version = "0.3.0"
edition = "2024"
rust-version = "1.94.0"
rust-version = "1.94.1"
authors = ["Wim Pomp <w.pomp@nki.nl>"]
license = "MIT OR Apache-2.0"
description = "Read bio image formats using the bio-formats java package."
image-registration
+1 -1
Submodule image-registration updated: 7eb39a2fbf...5c580dc80f.
image-tiff
+1 -1
Submodule image-tiff updated: 9b452e93b7...a6e8ed772d.
ndbioimage/Cargo.toml
+1 -1
View File
@@ -37,7 +37,7 @@ ome-metadata = "0.6"
ordered-float = { version = "5", optional = true }
phf = { version = "0.14", features = ["macros"] }
postcard = { version = "1", features = ["use-std"], optional = true }
pyo3 = { version = "0.29", features = ["abi3-py310", "anyhow", "generate-import-lib"], optional = true }
pyo3 = { version = "0.29", features = ["abi3-py310", "anyhow", "generate-import-lib", "multiple-pymethods"], optional = true }
pyo3-stub-gen = { version = "0.23", optional = true }
rayon = { version = "1", optional = true }
regex = "1"
ndbioimage/pyproject.toml
+1 -1
View File
@@ -4,7 +4,7 @@ build-backend = "maturin"
[project]
name = "ndbioimage"
version = "2027.0.5"
version = "2027.0.6"
requires-python = ">=3.10"
classifiers = [
"License :: OSI Approved :: MIT License",
ndbioimage/src/axes.rs
+3 -3
View File
@@ -198,9 +198,9 @@ impl<'de> DeserializeAs<'de, SliceInfoElem> for SliceInfoElemDef {
#[derive(Clone, Debug)]
pub(crate) struct Slice {
start: isize,
end: isize,
step: isize,
pub(crate) start: isize,
pub(crate) end: isize,
pub(crate) step: isize,
}
impl Slice {
ndbioimage/src/cache.rs
+86 -12
View File
@@ -1,14 +1,22 @@
use crate::axes::{Axis, Operation};
use crate::readers::Frame;
#[cfg(feature = "transforms")]
use image_registration::transform::Transform;
use indexmap::{Equivalent, IndexMap};
#[cfg(feature = "transforms")]
use ndarray::IxDyn;
use ndarray::{ArrayD, SliceInfoElem};
use std::any::Any;
#[cfg(feature = "transforms")]
use std::borrow::Cow;
use std::cmp::PartialEq;
use std::fmt::{Debug, Formatter};
use std::hash::{Hash, Hasher};
use std::path::{Path, PathBuf};
use std::sync::{Arc, Mutex, OnceLock};
use std::time::SystemTime;
use crate::readers::Frame;
#[cfg(feature = "transforms")]
use crate::transforms::Transforms;
/// maximum number of frames held in the cache
const DEFAULT_FRAME_CACHE_SIZE: usize = 128;
@@ -132,7 +140,15 @@ impl ReaderKeyRef<'_> {
}
}
pub(crate) type FrameKey = (ReaderKey, usize, usize, usize);
#[derive(Debug, PartialEq, Eq, Hash)]
pub(crate) struct FrameKey {
pub(crate) reader_key: ReaderKey,
pub(crate) c: usize,
pub(crate) z: usize,
pub(crate) t: usize,
#[cfg(feature = "transforms")]
pub(crate) transform: Option<(Transform<IxDyn>, usize)>,
}
/// borrowed view of [`FrameKey`] for cache lookups without allocation.
#[derive(Debug, PartialEq, Eq)]
@@ -141,6 +157,8 @@ pub(crate) struct FrameKeyRef<'a> {
pub(crate) c: usize,
pub(crate) z: usize,
pub(crate) t: usize,
#[cfg(feature = "transforms")]
pub(crate) transform: Option<(Cow<'a, Transform<IxDyn>>, usize)>,
}
impl Hash for FrameKeyRef<'_> {
@@ -149,19 +167,50 @@ impl Hash for FrameKeyRef<'_> {
self.c.hash(state);
self.z.hash(state);
self.t.hash(state);
#[cfg(feature = "transforms")]
self.transform.hash(state);
}
}
impl Equivalent<FrameKey> for FrameKeyRef<'_> {
fn equivalent(&self, key: &FrameKey) -> bool {
let (rk, c, z, t) = key;
self.reader.equivalent(rk) && self.c == *c && self.z == *z && self.t == *t
let FrameKey {
reader_key,
c,
z,
t,
#[cfg(feature = "transforms")]
transform,
} = key;
#[cfg(feature = "transforms")]
{
self.reader.equivalent(reader_key)
&& self.c == *c
&& self.z == *z
&& self.t == *t
&& self.transform.as_ref().map(|(t, i)| (t.as_ref(), i))
== transform.as_ref().map(|(t, i)| (t, i))
}
#[cfg(not(feature = "transforms"))]
{
self.reader.equivalent(reader_key) && self.c == *c && self.z == *z && self.t == *t
}
}
}
impl FrameKeyRef<'_> {
pub(crate) fn to_owned(&self) -> FrameKey {
(self.reader.to_owned(), self.c, self.z, self.t)
FrameKey {
reader_key: self.reader.to_owned(),
c: self.c,
z: self.z,
t: self.t,
#[cfg(feature = "transforms")]
transform: self
.transform
.as_ref()
.map(|(t, i)| (t.clone().into_owned(), *i)),
}
}
}
@@ -173,6 +222,8 @@ pub(crate) struct ArrayKey {
pub(crate) slice: Vec<SliceInfoElem>,
pub(crate) axes: Vec<Axis>,
pub(crate) operations: Vec<(Axis, Operation)>,
#[cfg(feature = "transforms")]
pub(crate) transforms: Transforms,
}
/// borrowed view of [`ArrayKey`] for cache lookups without allocation.
@@ -183,6 +234,8 @@ pub(crate) struct ArrayKeyRef<'a> {
pub(crate) slice: &'a [SliceInfoElem],
pub(crate) axes: &'a [Axis],
pub(crate) operations: &'a IndexMap<Axis, Operation>,
#[cfg(feature = "transforms")]
pub(crate) transforms: &'a Transforms,
}
impl Hash for ArrayKeyRef<'_> {
@@ -196,21 +249,40 @@ impl Hash for ArrayKeyRef<'_> {
ax.hash(state);
op.hash(state);
}
#[cfg(feature = "transforms")]
self.transforms.hash(state);
}
}
impl Equivalent<ArrayKey> for ArrayKeyRef<'_> {
fn equivalent(&self, key: &ArrayKey) -> bool {
self.reader.equivalent(&key.reader)
&& self.dtype == key.dtype
&& self.slice == key.slice
&& self.axes == key.axes
&& self.operations.len() == key.operations.len()
&& self
#[cfg(feature = "transforms")]
{
self.reader.equivalent(&key.reader)
&& self.dtype == key.dtype
&& self.slice == key.slice
&& self.axes == key.axes
&& self.operations.len() == key.operations.len()
&& self
.operations
.iter()
.zip(&key.operations)
.all(|((ax, op), (key_ax, key_op))| ax == key_ax && op == key_op)
&& self.transforms == key.transforms
}
#[cfg(not(feature = "transforms"))]
{
self.reader.equivalent(&key.reader)
&& self.dtype == key.dtype
&& self.slice == key.slice
&& self.axes == key.axes
&& self.operations.len() == key.operations.len()
&& self
.operations
.iter()
.zip(&key.operations)
.all(|((ax, op), (key_ax, key_op))| ax == key_ax && op == key_op)
}
}
}
@@ -226,6 +298,8 @@ impl ArrayKeyRef<'_> {
.iter()
.map(|(ax, op)| (*ax, op.clone()))
.collect(),
#[cfg(feature = "transforms")]
transforms: self.transforms.to_owned(),
}
}
}
ndbioimage/src/error.rs
+4
View File
@@ -179,6 +179,10 @@ pub enum Error {
it cannot use. read the file before forking, or start the child with spawn/exec instead."
)]
ForkedJvm(u32, u32),
#[error("{0}D transformation is not implemented")]
TransformationNotImplemented(usize),
#[error("no channel {0} present")]
NoChannel(usize),
}
impl Error {
ndbioimage/src/lib.rs
+4 -4
View File
@@ -1,4 +1,4 @@
#![cfg_attr(docsrs, feature(doc_cfg))]
#![cfg_attr(docsrs, feature(doc_cfg, doc_auto_cfg))]
//! The ndbioimage crate exposes (bio) images as a struct that can be sliced like an ndarray Array
//! (Rust), but without loading the whole image into memory, reading from the file only when needed.
//! Some metadata is read
@@ -16,7 +16,7 @@
//! # fn main() -> Result<(), ndbioimage::error::Error> {
//! let path = "/path/to/file";
//! let reader = DynReader::new(&path, 0, 0)?;
//! println!("shape: {}", reader.shape());
//! println!("shape: {}", reader.raw_shape());
//! let frame = reader.get_frame(0, 0, 0)?;
//! if let Ok(arr) = <Frame as TryInto<Array2<i8>>>::try_into(frame) {
//! println!("{:?}", arr);
@@ -70,7 +70,7 @@ pub mod readers;
pub mod tiffwrite;
#[cfg(feature = "transforms")]
/// image registration and warping images using affine transforms
mod transforms;
pub mod transforms;
mod utils;
/// main entry point for the application
@@ -338,7 +338,7 @@ mod tests {
fn read_ser() -> Result<(), Error> {
let file = "czi/Experiment-2029.czi";
let reader = open(file)?;
println!("shape: {}", reader.shape());
println!("shape: {}", reader.raw_shape());
let frame = reader.get_frame(0, 0, 0)?;
if let Ok(arr) = <Frame as TryInto<Array2<i8>>>::try_into(frame) {
println!("{:?}", arr);
ndbioimage/src/py.rs
+1
View File
@@ -67,6 +67,7 @@ fn main() -> PyResult<()> {
#[pymodule]
#[pyo3(name = "ndbioimage_rs")]
mod ndbioimage_rs {
#[cfg(feature = "transforms")]
use pyo3::prelude::*;
#[pymodule_export]
ndbioimage/src/py/imread.rs
+3 -139
View File
@@ -3,8 +3,6 @@ use crate::error::Error;
use crate::metadata::Metadata;
use crate::movie::MovieOptions;
#[cfg(feature = "transforms")]
use crate::py::transforms::{PyTransform, PyTransforms};
use crate::readers::{DynReader, PixelType, Reader};
use crate::view::{Item, View};
use itertools::Itertools;
@@ -27,7 +25,7 @@ pub(crate) use pyo3::{
use pyo3_stub_gen::derive::*;
use pyo3_stub_gen::inventory::submit;
use serde::{Deserialize, Serialize};
use std::collections::{HashMap, HashSet};
use std::collections::HashSet;
use std::path::PathBuf;
/// class to read image files, while taking good care of important metadata,
@@ -61,7 +59,7 @@ use std::path::PathBuf;
#[pyclass(subclass, from_py_object, name = "Imread", module = "ndbioimage")]
#[derive(Clone, Debug, Serialize, Deserialize)]
pub(crate) struct PyView {
view: View<IxDyn, DynReader>,
pub(crate) view: View<IxDyn, DynReader>,
dtype: PixelType,
#[serde(skip)]
ome: Ome,
@@ -296,141 +294,6 @@ impl PyView {
todo!()
}
/// return a new view with transformations applied (channel alignment, drift correction)
/// backward compatibility
#[cfg(feature = "transforms")]
#[pyo3(signature = (channels = true, drift = false, file = None, bead_files = None, main_channel = None, default_transform = None))]
fn with_transform<'py>(
&self,
channels: bool,
drift: bool,
file: Option<Bound<'py, PyAny>>,
bead_files: Option<Bound<'py, PyAny>>,
main_channel: Option<usize>,
default_transform: Option<Vec<f64>>,
) -> PyResult<PyView> {
self.with_transforms(
channels,
drift,
file,
bead_files,
main_channel,
default_transform,
)
}
/// return a new view with transformations applied (channel alignment, drift correction)
#[cfg(feature = "transforms")]
#[expect(unused_variables)]
#[pyo3(signature = (channels = true, drift = false, file = None, bead_files = None, main_channel = None, default_transform = None))]
fn with_transforms<'py>(
&self,
channels: bool,
drift: bool,
file: Option<Bound<'py, PyAny>>,
bead_files: Option<Bound<'py, PyAny>>,
main_channel: Option<usize>,
default_transform: Option<Vec<f64>>,
) -> PyResult<PyView> {
todo!()
}
/// get the transformation
#[cfg(feature = "transforms")]
#[getter]
fn get_transforms(&self) -> PyResult<PyTransforms> {
Ok(PyTransforms {
inner: self.view.get_transforms().clone(),
})
}
/// set the transformation
#[cfg(feature = "transforms")]
#[setter]
fn set_transforms(&mut self, transform: PyTransforms) {
self.view.set_transforms(transform.inner);
}
#[cfg(feature = "transforms")]
fn load_transforms_from_yaml(&mut self, path: PathBuf) -> PyResult<()> {
self.view.load_transforms_from_yaml(path.as_path())?;
Ok(())
}
#[cfg(feature = "transforms")]
fn calculate_channel_transforms_2d(
&self,
py: Python,
main_channel: usize,
) -> PyResult<HashMap<String, PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_channel_transforms_2d(main_channel)?
.into_iter()
.map(|(n, t)| {
(
n,
PyTransform {
inner: t.into_dyn(),
},
)
})
.collect())
})
}
#[cfg(feature = "transforms")]
fn calculate_channel_transforms_3d(
&self,
py: Python,
main_channel: usize,
) -> PyResult<HashMap<String, PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_channel_transforms_3d(main_channel)?
.into_iter()
.map(|(n, t)| {
(
n,
PyTransform {
inner: t.into_dyn(),
},
)
})
.collect())
})
}
#[cfg(feature = "transforms")]
fn calculate_drift_transforms_2d(&self, py: Python) -> PyResult<Vec<PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_drift_transforms_2d()?
.into_iter()
.map(|t| PyTransform {
inner: t.into_dyn(),
})
.collect())
})
}
#[cfg(feature = "transforms")]
fn calculate_drift_transforms_3d(&self, py: Python) -> PyResult<Vec<PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_drift_transforms_3d()?
.into_iter()
.map(|t| PyTransform {
inner: t.into_dyn(),
})
.collect())
})
}
#[gen_stub(override_return_type(type_repr="numpy.ndarray | int | float", imports=("numpy")))]
fn squeeze<'py>(&self, py: Python<'py>) -> PyResult<Bound<'py, PyAny>> {
let view = self.view.squeeze()?;
@@ -1481,6 +1344,7 @@ impl PyView {
))]
axes: Option<Vec<Bound<PyAny>>>,
) -> PyResult<Self> {
// TODO: transpose("xyz"), transpose(*axes)
let view = if let Some(axes) = axes {
let ax = axes
.into_iter()
ndbioimage/src/py/transforms.rs
+133
View File
@@ -1,4 +1,5 @@
use crate::error::Error;
use crate::py::imread::PyView;
use crate::transforms::Transforms;
use image_registration::metric::{SamplingArg, Sigma};
use image_registration::register::{Optimizer, RegistrationStep};
@@ -709,6 +710,138 @@ pub(crate) struct PyTransforms {
pub(crate) inner: Transforms,
}
#[gen_stub_pymethods]
#[pymethods]
impl PyView {
/// return a new view with transformations applied (channel alignment, drift correction)
/// backward compatibility
#[cfg(feature = "transforms")]
#[pyo3(signature = (channels = true, drift = false, file = None, bead_files = None, main_channel = None, default_transform = None))]
fn with_transform<'py>(
&self,
channels: bool,
drift: bool,
file: Option<Bound<'py, PyAny>>,
bead_files: Option<Bound<'py, PyAny>>,
main_channel: Option<usize>,
default_transform: Option<Vec<f64>>,
) -> PyResult<PyView> {
self.with_transforms(
channels,
drift,
file,
bead_files,
main_channel,
default_transform,
)
}
/// return a new view with transformations applied (channel alignment, drift correction)
#[cfg(feature = "transforms")]
#[expect(unused_variables)]
#[pyo3(signature = (channels = true, drift = false, file = None, bead_files = None, main_channel = None, default_transform = None))]
fn with_transforms<'py>(
&self,
channels: bool,
drift: bool,
file: Option<Bound<'py, PyAny>>,
bead_files: Option<Bound<'py, PyAny>>,
main_channel: Option<usize>,
default_transform: Option<Vec<f64>>,
) -> PyResult<PyView> {
todo!()
}
/// get the transformation
#[getter]
fn get_transforms(&self) -> PyResult<PyTransforms> {
Ok(PyTransforms {
inner: self.view.get_transforms().clone(),
})
}
/// set the transformation
#[setter]
fn set_transforms(&mut self, transform: PyTransforms) {
self.view.set_transforms(transform.inner);
}
fn load_transforms_from_yaml(&mut self, path: PathBuf) -> PyResult<()> {
self.view.load_transforms_from_yaml(path.as_path())?;
Ok(())
}
fn calculate_channel_transforms_2d(
&self,
py: Python,
main_channel: usize,
) -> PyResult<HashMap<String, PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_channel_transforms_2d(main_channel)?
.into_iter()
.map(|(n, t)| {
(
n,
PyTransform {
inner: t.into_dyn(),
},
)
})
.collect())
})
}
fn calculate_channel_transforms_3d(
&self,
py: Python,
main_channel: usize,
) -> PyResult<HashMap<String, PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_channel_transforms_3d(main_channel)?
.into_iter()
.map(|(n, t)| {
(
n,
PyTransform {
inner: t.into_dyn(),
},
)
})
.collect())
})
}
fn calculate_drift_transforms_2d(&self, py: Python) -> PyResult<Vec<PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_drift_transforms_2d()?
.into_iter()
.map(|t| PyTransform {
inner: t.into_dyn(),
})
.collect())
})
}
fn calculate_drift_transforms_3d(&self, py: Python) -> PyResult<Vec<PyTransform>> {
py.detach(|| {
Ok(self
.view
.calculate_drift_transforms_3d()?
.into_iter()
.map(|t| PyTransform {
inner: t.into_dyn(),
})
.collect())
})
}
}
#[gen_stub_pymethods]
#[pymethods]
impl PyTransforms {
ndbioimage/src/readers.rs
+14 -11
View File
@@ -2,7 +2,7 @@ use crate::axes::{Axis, Shape};
use crate::error::Error;
use crate::view::View;
use ndarray::{Array, Dimension, Ix2, Ix5, s};
use num::{FromPrimitive, Zero};
use num::{FromPrimitive, One, Zero};
use ome_metadata::Ome;
use regex::Regex;
use serde::{Deserialize, Serialize};
@@ -148,6 +148,8 @@ impl PixelType {
/// array data with a specific pixel type
#[derive(Clone, Debug)]
pub enum ArrayT<D: Dimension> {
/// boolean array
Bool(Array<bool, D>),
/// signed 8-bit integer array
I8(Array<i8, D>),
/// unsigned 8-bit integer array
@@ -197,7 +199,7 @@ pub trait Reader: Clone + Sized + Debug + Send + Hash + Into<DynReader> {
/// get a sliceable view on the image file
fn view(&self) -> View<Ix5, Self> {
let shape = self.shape();
let shape = self.raw_shape();
let slice = s![0..shape.c, 0..shape.z, 0..shape.t, 0..shape.y, 0..shape.x,];
View::new(
self.clone(),
@@ -216,7 +218,7 @@ pub trait Reader: Clone + Sized + Debug + Send + Hash + Into<DynReader> {
/// the position index
fn position(&self) -> usize;
/// the shape of the image data
fn shape(&self) -> &Shape;
fn raw_shape(&self) -> &Shape;
/// the pixel type of the image data
fn pixel_type(&self) -> &PixelType;
/// get all available positions for a given series
@@ -314,6 +316,7 @@ macro_rules! try_into_array_body {
($self:expr, $map:ident) => {{
let mut err = Ok(());
let arr = match $self {
ArrayT::Bool(v) => v.$map(|x| if x { T::one() } else { T::zero() }),
ArrayT::I8(v) => v.$map(|x| {
T::from_i8(x).unwrap_or_else(|| {
err = Err(Error::Cast(x.to_string(), type_name::<T>().to_string()));
@@ -397,7 +400,7 @@ macro_rules! try_into_array_body {
impl<D, T> TryInto<Array<T, D>> for ArrayT<D>
where
D: Dimension,
T: FromPrimitive + Zero + 'static,
T: FromPrimitive + Zero + One + 'static,
{
type Error = Error;
@@ -409,7 +412,7 @@ where
impl<D, T> TryInto<Array<T, D>> for &ArrayT<D>
where
D: Dimension,
T: Clone + FromPrimitive + Zero + 'static,
T: Clone + FromPrimitive + Zero + One + 'static,
{
type Error = Error;
@@ -577,18 +580,18 @@ impl Reader for DynReader {
}
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
match self {
#[cfg(feature = "tiff")]
DynReader::Tiff(r) => r.shape(),
DynReader::Tiff(r) => r.raw_shape(),
#[cfg(feature = "tiffseq")]
DynReader::TiffSeq(r) => r.shape(),
DynReader::TiffSeq(r) => r.raw_shape(),
#[cfg(feature = "czi")]
DynReader::Czi(r) => r.shape(),
DynReader::Czi(r) => r.raw_shape(),
#[cfg(feature = "bioformats_rust")]
DynReader::BioFormatsRust(r) => r.shape(),
DynReader::BioFormatsRust(r) => r.raw_shape(),
#[cfg(feature = "bioformats_java")]
DynReader::BioFormatsJava(r) => r.shape(),
DynReader::BioFormatsJava(r) => r.raw_shape(),
#[expect(unreachable_patterns)]
_ => unreachable!(),
}
ndbioimage/src/readers/bioformats_java.rs
+2 -2
View File
@@ -38,7 +38,7 @@ static JVM_PID: AtomicU32 = AtomicU32::new(0);
/// else -- verified by `JNI_GetCreatedJavaVMs`, which still reports 1 in the child, so a second
/// `JNI_CreateJavaVM` (which is what j4rs' `Jvm::attach_thread` calls) takes the already-exists
/// path. An `exec` replaces the address space, and that does reset the count to 0.
fn check_not_forked() -> Result<(), Error> {
pub fn check_not_forked() -> Result<(), Error> {
let built_in = JVM_PID.load(Ordering::Relaxed);
if built_in != 0 && built_in != std::process::id() {
return Err(Error::ForkedJvm(built_in, std::process::id()));
@@ -698,7 +698,7 @@ impl Reader for BioFormatsJavaReader {
0
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
&self.shape
}
ndbioimage/src/readers/bioformats_rust.rs
+1 -1
View File
@@ -309,7 +309,7 @@ impl Reader for BioFormatsRustReader {
0
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
&self.shape
}
ndbioimage/src/readers/czi.rs
+1 -1
View File
@@ -1149,7 +1149,7 @@ impl Reader for CziReader {
self.position
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
&self.shape
}
ndbioimage/src/readers/tiff.rs
+1 -1
View File
@@ -459,7 +459,7 @@ impl Reader for TiffReader {
fn position(&self) -> usize {
self.position
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
&self.shape
}
fn pixel_type(&self) -> &PixelType {
ndbioimage/src/readers/tiffseq.rs
+1 -1
View File
@@ -572,7 +572,7 @@ impl Reader for TiffSeqReader {
self.position
}
fn shape(&self) -> &Shape {
fn raw_shape(&self) -> &Shape {
&self.shape
}
ndbioimage/src/transforms.rs
+234 -6
View File
@@ -1,13 +1,16 @@
use crate::axes::Axis;
use crate::error::Error;
use crate::readers::{DynReader, Reader};
use crate::readers::{DynReader, Frame, PixelType, Reader};
use crate::view::View;
pub use image_registration::transform::Transform;
use ndarray::{Array2, Dimension, Ix2, Ix3, IxDyn, s};
use ndarray::{Array2, Dimension, Ix2, Ix3, IxDyn, array, s};
use serde::{Deserialize, Serialize};
use std::borrow::Cow;
use std::collections::{HashMap, HashSet};
use std::hash::{Hash, Hasher};
use std::io::Seek;
use std::path::Path;
use std::sync::Arc;
#[derive(Clone, Debug, Default, PartialEq, Serialize, Deserialize)]
struct TransformOld {
@@ -54,6 +57,24 @@ pub struct Transforms {
pub drift: Vec<Transform<IxDyn>>,
}
impl Hash for Transforms {
fn hash<H: Hasher>(&self, state: &mut H) {
let mut keys: Vec<_> = self.channel.keys().collect();
keys.sort();
for key in keys {
key.hash(state);
self.channel[key].hash(state);
}
self.drift.hash(state);
}
}
impl PartialEq<Transforms> for &Transforms {
fn eq(&self, other: &Transforms) -> bool {
**self == *other
}
}
impl Serialize for Transforms {
fn serialize<S>(&self, serializer: S) -> Result<S::Ok, S::Error>
where
@@ -89,7 +110,23 @@ impl Serialize for Transforms {
}
}
pub enum TransformKind {
None,
Channel,
Drift,
Both,
}
impl Transforms {
pub fn kind(&self) -> TransformKind {
match (self.channel.is_empty(), self.drift.is_empty()) {
(true, true) => TransformKind::None,
(false, true) => TransformKind::Channel,
(true, false) => TransformKind::Drift,
(false, false) => TransformKind::Both,
}
}
pub fn load<P: AsRef<Path>>(path: P) -> Result<Self, Error> {
let mut file = std::fs::File::open(&path)?;
Ok(
@@ -210,16 +247,32 @@ impl Transforms {
}
Ok(average_transforms)
}
pub fn get_transform(
&self,
channel: &str,
t: usize,
) -> Result<Option<Cow<'_, Transform<IxDyn>>>, Error> {
match (self.channel.get(channel), self.drift.get(t)) {
(Some(channel), Some(drift)) => match (channel.ndim(), drift.ndim()) {
(2, 2) | (3, 3) => Ok(Some(Cow::Owned(drift * channel))),
(2, 3) => Ok(Some(Cow::Owned(drift.add_dimension(0, 1) * channel))),
(3, 2) => Ok(Some(Cow::Owned(drift * channel.add_dimension(0, 1)))),
(2 | 3, n) | (n, _) => Err(Error::TransformationNotImplemented(n)),
},
(Some(channel), None) => Ok(Some(Cow::Borrowed(channel))),
(None, Some(drift)) => Ok(Some(Cow::Borrowed(drift))),
(None, None) => Ok(None),
}
}
}
impl<D: Dimension, R: Reader> View<D, R> {
#[cfg(feature = "transforms")]
pub(crate) fn get_transforms(&self) -> &Transforms {
pub fn get_transforms(&self) -> &Transforms {
&self.transforms
}
#[cfg(feature = "transforms")]
pub(crate) fn set_transforms(&mut self, transforms: Transforms) {
pub fn set_transforms(&mut self, transforms: Transforms) {
self.transforms = transforms;
}
@@ -318,6 +371,181 @@ impl<D: Dimension, R: Reader> View<D, R> {
pub fn calculate_drift_transforms_3d(&self) -> Result<Vec<Transform<Ix3>>, Error> {
todo!()
}
pub(crate) fn get_transformed_frame(
&self,
c: usize,
z: usize,
t: usize,
transform: &Transform<IxDyn>,
order: usize,
) -> Result<Arc<Frame>, Error> {
match transform.ndim() {
2 => {
let tr = transform.clone().into_dimensionality::<Ix2>()?;
macro_rules! interpolate {
($array:expr ,$var:ident, bool) => {
Ok(Arc::from(Frame::$var(
tr.interpolate_par_order(
&$array.mapv(|i| if i { 1.0 } else { 0.0 }),
order,
)?
.mapv(|i| i > 0.5),
)))
};
($array:expr ,$var:ident, f64) => {
Ok(Arc::from(Frame::$var(
tr.interpolate_par_order($array, order)?,
)))
};
($array:expr ,$var:ident, $t:ty) => {
Ok(Arc::from(Frame::$var(
tr.interpolate_par_order($array, order)?.mapv(|i| i as $t),
)))
};
}
match self.get_cached_frame(c, z, t, None)?.as_ref() {
Frame::Bool(frame) => interpolate!(frame, Bool, bool),
Frame::I8(frame) => interpolate!(frame, I8, i8),
Frame::U8(frame) => interpolate!(frame, U8, u8),
Frame::I16(frame) => interpolate!(frame, I16, i16),
Frame::U16(frame) => interpolate!(frame, U16, u16),
Frame::I32(frame) => interpolate!(frame, I32, i32),
Frame::U32(frame) => interpolate!(frame, U32, u32),
Frame::U64(frame) => interpolate!(frame, U64, u64),
Frame::I64(frame) => interpolate!(frame, I64, i64),
Frame::I128(frame) => interpolate!(frame, I128, i128),
Frame::U128(frame) => interpolate!(frame, U128, u128),
Frame::F32(frame) => interpolate!(frame, F32, f32),
Frame::F64(frame) => interpolate!(frame, F64, f64),
Frame::F128(frame) => interpolate!(frame, F128, f64),
}
}
3 => {
let tr = transform.clone().into_dimensionality::<Ix3>()?;
let shape = self.shape();
let zf = z as f64;
let yf = shape.y as f64 - 1.0;
let xf = shape.x as f64 - 1.0;
let corners = array![[zf, 0.0, 0.0], [zf, 0.0, xf], [zf, yf, 0.0], [zf, yf, xf]];
let extent = tr.transform_points(&corners)?;
let z_extent = extent.slice(s![.., 2]);
let z_min = (z_extent
.iter()
.min_by(|a, b| a.total_cmp(b))
.unwrap()
.floor() as isize
- 2)
.max(0) as usize;
let z_max = (z_extent
.iter()
.max_by(|a, b| a.total_cmp(b))
.unwrap()
.ceil() as isize
+ 2)
.min(shape.z as isize - 1)
.max(0) as usize;
let frames = (z_min..=z_max)
.map(|zi| self.get_cached_frame(c, zi, t, None))
.collect::<Result<Vec<_>, Error>>()?;
// TODO; check these translations
let u0 = tr.as_translation(&[(z_min as f64), 0.0, 0.0]);
let u1 = tr.as_translation(&[-(z_min as f64), 0.0, 0.0]);
let tr = u0 * tr * u1;
match frames.first().map(|f| f.as_ref()) {
Some(first) => {
macro_rules! interpolate {
($var:ident, bool) => {{
let views = frames
.iter()
.map(|f| match f.as_ref() {
Frame::$var(a) => Ok(a.view()),
_ => unreachable!(),
})
.collect::<Result<Vec<_>, Error>>()?;
let stack = ndarray::stack(ndarray::Axis(0), &views)?;
let transformed = tr.interpolate_par_order(
&stack.mapv(|i| if i { 1.0 } else { 0.0 }),
order,
)?;
Ok(Arc::from(Frame::$var(
transformed.slice(s![0, .., ..]).mapv(|i| i > 0.5),
)))
}};
($var:ident, f64) => {{
let views = frames
.iter()
.map(|f| match f.as_ref() {
Frame::$var(a) => Ok(a.view()),
_ => unreachable!(),
})
.collect::<Result<Vec<_>, Error>>()?;
let stack = ndarray::stack(ndarray::Axis(0), &views)?;
let transformed = tr.interpolate_par_order(&stack, order)?;
Ok(Arc::from(Frame::$var(
transformed.slice(s![0, .., ..]).to_owned(),
)))
}};
($var:ident, $t:ty) => {{
let views = frames
.iter()
.map(|f| match f.as_ref() {
Frame::$var(a) => Ok(a.view()),
_ => unreachable!(),
})
.collect::<Result<Vec<_>, Error>>()?;
let stack = ndarray::stack(ndarray::Axis(0), &views)?;
let transformed = tr.interpolate_par_order(&stack, order)?;
Ok(Arc::from(Frame::$var(
transformed.slice(s![0, .., ..]).mapv(|i| i as $t),
)))
}};
}
match first {
Frame::Bool(_) => interpolate!(Bool, bool),
Frame::I8(_) => interpolate!(I8, i8),
Frame::U8(_) => interpolate!(U8, u8),
Frame::I16(_) => interpolate!(I16, i16),
Frame::U16(_) => interpolate!(U16, u16),
Frame::I32(_) => interpolate!(I32, i32),
Frame::U32(_) => interpolate!(U32, u32),
Frame::F32(_) => interpolate!(F32, f32),
Frame::F64(_) => interpolate!(F64, f64),
Frame::I64(_) => interpolate!(I64, i64),
Frame::U64(_) => interpolate!(U64, u64),
Frame::I128(_) => interpolate!(I128, i128),
Frame::U128(_) => interpolate!(U128, u128),
Frame::F128(_) => interpolate!(F128, f64),
}
}
None => {
let shape = self.shape();
Ok(Arc::new(match self.pixel_type() {
PixelType::Bool => Frame::Bool(Array2::default((shape.y, shape.x))),
PixelType::I8 => Frame::I8(Array2::zeros((shape.y, shape.x))),
PixelType::U8 => Frame::U8(Array2::zeros((shape.y, shape.x))),
PixelType::I16 => Frame::I16(Array2::zeros((shape.y, shape.x))),
PixelType::U16 => Frame::U16(Array2::zeros((shape.y, shape.x))),
PixelType::I32 => Frame::I32(Array2::zeros((shape.y, shape.x))),
PixelType::U32 => Frame::U32(Array2::zeros((shape.y, shape.x))),
PixelType::F32 => Frame::F32(Array2::zeros((shape.y, shape.x))),
PixelType::F64 => Frame::F64(Array2::zeros((shape.y, shape.x))),
PixelType::I64 => Frame::I64(Array2::zeros((shape.y, shape.x))),
PixelType::U64 => Frame::U64(Array2::zeros((shape.y, shape.x))),
PixelType::I128 => Frame::I128(Array2::zeros((shape.y, shape.x))),
PixelType::U128 => Frame::U128(Array2::zeros((shape.y, shape.x))),
PixelType::F128 => Frame::F128(Array2::zeros((shape.y, shape.x))),
}))
}
}
}
d => Err(Error::TransformationNotImplemented(d)),
}
}
}
#[cfg(test)]
ndbioimage/src/view.rs
+304 -117
View File
@@ -5,18 +5,23 @@ use crate::metadata::Metadata;
use crate::readers::{Dimensions, DynReader, Frame, Reader};
use crate::stats::MinMax;
#[cfg(feature = "transforms")]
use crate::transforms::Transforms;
use crate::transforms::{TransformKind, Transforms};
#[cfg(feature = "transforms")]
use image_registration::transform::Transform;
use indexmap::IndexMap;
use itertools::{Itertools, iproduct};
use ndarray::{
Array, Array0, Array1, Array2, ArrayD, Dimension, IntoDimension, Ix0, Ix1, Ix2, Ix5, IxDyn,
SliceArg, SliceInfoElem, s,
};
use num::traits::ToBytes;
use num::{Bounded, FromPrimitive, ToPrimitive, Zero};
use num::traits::{AsPrimitive, ToBytes};
use num::{Bounded, FromPrimitive, One, ToPrimitive, Zero};
use serde::{Deserialize, Serialize};
use serde_with::serde_as;
use std::any::type_name;
#[cfg(feature = "transforms")]
use std::borrow::Cow;
use std::cell::OnceCell;
use std::collections::{HashMap, HashSet};
use std::fmt::{Debug, Display, Formatter};
use std::hash::{Hash, Hasher};
@@ -57,27 +62,43 @@ pub trait Number:
+ Sync
+ AddAssign
+ Bounded
+ Clone
+ Copy
+ Div<Self, Output = Self>
+ FromPrimitive
+ PartialOrd
+ Zero
+ One
+ AsPrimitive<f64>
{
const ISFLOAT: bool;
}
impl<T> Number for T where
T: 'static
+ Send
+ Sync
+ AddAssign
+ Bounded
+ Clone
+ Div<Self, Output = Self>
+ FromPrimitive
+ PartialOrd
+ Zero
{
macro_rules! impl_int {
($($t:tt),+ $(,)?) => {
$(
impl Number for $t {
const ISFLOAT: bool = false;
}
)*
};
}
impl_int!(
i8, i16, i32, i64, i128, isize, u8, u16, u32, u64, u128, usize
);
macro_rules! impl_float {
($($t:tt),+ $(,)?) => {
$(
impl Number for $t {
const ISFLOAT: bool = true;
}
)*
};
}
impl_float!(f32, f64);
/// sliceable view on an image file
#[serde_as]
#[derive(Clone, Debug, PartialEq, Eq, Serialize, Deserialize)]
@@ -89,6 +110,8 @@ pub struct View<D: Dimension, R: Reader = DynReader> {
/// always has all of cztyx with possibly some new axes added
axes: Vec<Axis>,
operations: IndexMap<Axis, Operation>,
#[serde(skip)]
channel_names: OnceCell<HashMap<usize, String>>,
dimensionality: PhantomData<D>,
#[cfg(feature = "transforms")]
pub(crate) transforms: Transforms,
@@ -117,6 +140,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
slice,
axes,
operations: IndexMap::new(),
channel_names: OnceCell::new(),
dimensionality: PhantomData,
#[cfg(feature = "transforms")]
transforms: Transforms::default(),
@@ -126,7 +150,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
#[expect(dead_code)]
pub(crate) fn new_with_axes(reader: R, axes: Vec<Axis>) -> Result<Self, Error> {
let mut slice = Vec::new();
let shape = reader.shape();
let shape = reader.raw_shape();
for axis in axes.iter() {
match axis {
Axis::C => slice.push(SliceInfoElem::Slice {
@@ -183,6 +207,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
axes,
operations: IndexMap::new(),
dimensionality: PhantomData,
channel_names: OnceCell::new(),
#[cfg(feature = "transforms")]
transforms: Transforms::default(),
})
@@ -198,6 +223,26 @@ impl<D: Dimension, R: Reader> View<D, R> {
self.reader.series()
}
/// the name of a particular channel
pub fn channel_name(&self, channel: usize) -> Result<&str, Error> {
if self.channel_names.get().is_none() {
let names = self.metadata()?.image[0]
.pixels
.channel
.iter()
.enumerate()
.map(|(i, c)| (i, c.name.clone().unwrap_or_else(|| format!("{}", i))))
.collect::<HashMap<_, _>>();
let _ = self.channel_names.set(names);
}
self.channel_names
.get()
.and_then(|names| names.get(&channel))
.map(String::as_str)
.ok_or(Error::NoChannel(channel))
}
/// set the maximum number of frames the process-wide cache may hold,
/// evicting the least recently used frames
pub fn with_cache_capacity(self, capacity: usize) -> Self {
@@ -234,6 +279,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
slice: self.slice,
axes: self.axes,
operations: self.operations,
channel_names: self.channel_names,
dimensionality: PhantomData,
#[cfg(feature = "transforms")]
transforms: Transforms::default(),
@@ -249,6 +295,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
slice: self.slice,
axes: self.axes,
operations: self.operations,
channel_names: self.channel_names,
dimensionality: PhantomData,
#[cfg(feature = "transforms")]
transforms: self.transforms,
@@ -262,6 +309,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
slice: self.slice,
axes: self.axes,
operations: self.operations,
channel_names: self.channel_names,
dimensionality: PhantomData,
#[cfg(feature = "transforms")]
transforms: self.transforms,
@@ -671,7 +719,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
let slice: Vec<_> = index.iter().map(|s| SliceInfoElem::Index(*s)).collect();
let view = self.clone().into_dyn().slice(slice.as_slice())?;
let arr = view.as_array()?;
Ok(arr.first().unwrap().clone())
Ok(*arr.first().unwrap())
}
/// collect the view into an ndarray
@@ -714,7 +762,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
IndexMap::new()
};
let mut shape_out = Vec::new();
let mut slice = Vec::new();
let mut slice_out = Vec::new();
let mut ax_out = Vec::new();
for (s, a) in self.slice.iter().zip(&self.axes) {
match s {
@@ -722,7 +770,7 @@ impl<D: Dimension, R: Reader> View<D, R> {
let end = end.expect("slice has no end");
if !op_xy.contains_key(a) && !op_czt.contains_key(a) {
shape_out.push(((end - start).max(0) / step) as usize);
slice.push(SliceInfoElem::Slice {
slice_out.push(SliceInfoElem::Slice {
start: 0,
end: None,
step: 1,
@@ -733,59 +781,54 @@ impl<D: Dimension, R: Reader> View<D, R> {
SliceInfoElem::Index(_) => {}
SliceInfoElem::NewAxis => {
shape_out.push(1);
slice.push(SliceInfoElem::Index(0));
slice_out.push(SliceInfoElem::Index(0));
ax_out.push(*a);
}
}
}
let mut slice_reader = vec![Slice::empty(); 5];
let mut xy_dim = 0usize;
let shape_reader = self.reader.shape();
let shape_reader = self.reader.raw_shape();
for (s, &axis) in self.slice.iter().zip(&self.axes) {
match axis {
Axis::New => {}
_ => match s {
SliceInfoElem::Slice { start, end, step } => {
if let Axis::X | Axis::Y = axis
&& !op_xy.contains_key(&axis)
{
xy_dim += 1;
}
slice_reader[axis as usize] = Slice::new(
idx_bnd(*start, shape_reader[axis] as isize)?,
slc_bnd(end.unwrap(), shape_reader[axis] as isize)?,
*step,
);
match (axis, s) {
(Axis::New, _) => {}
(_, SliceInfoElem::Slice { start, end, step }) => {
if let Axis::X | Axis::Y = axis
&& !op_xy.contains_key(&axis)
{
xy_dim += 1;
}
SliceInfoElem::Index(j) => {
slice_reader[axis as usize] = Slice::new(
idx_bnd(*j, shape_reader[axis] as isize)?,
slc_bnd(*j + 1, shape_reader[axis] as isize)?,
1,
);
}
SliceInfoElem::NewAxis => panic!("axis cannot be a new axis"),
},
slice_reader[axis as usize] = Slice::new(
idx_bnd(*start, shape_reader[axis] as isize)?,
slc_bnd(end.unwrap(), shape_reader[axis] as isize)?,
*step,
);
}
(_, SliceInfoElem::Index(j)) => {
slice_reader[axis as usize] = Slice::new(
idx_bnd(*j, shape_reader[axis] as isize)?,
slc_bnd(*j + 1, shape_reader[axis] as isize)?,
1,
);
}
(_, SliceInfoElem::NewAxis) => panic!("axis cannot be a new axis"),
}
}
let xy = [
self.slice[Axis::Y.pos(&self.axes, &self.slice)?],
self.slice[Axis::X.pos(&self.axes, &self.slice)?],
];
let mut array = if let Some((_, op)) = op_czt.first() {
match op {
Operation::Max => {
ArrayD::<T>::from_elem(shape_out.into_dimension(), T::min_value())
}
Operation::Min => {
ArrayD::<T>::from_elem(shape_out.into_dimension(), T::max_value())
}
_ => ArrayD::<T>::zeros(shape_out.into_dimension()),
let mut array = match op_czt.first() {
Some((_, Operation::Max)) => {
ArrayD::<T>::from_elem(shape_out.into_dimension(), T::min_value())
}
} else {
ArrayD::<T>::zeros(shape_out.into_dimension())
Some((_, Operation::Min)) => {
ArrayD::<T>::from_elem(shape_out.into_dimension(), T::max_value())
}
_ => ArrayD::<T>::zeros(shape_out.into_dimension()),
};
let shape = self.reader.shape();
let shape = self.reader.raw_shape();
let mut axes_out_idx = [None; 5];
for (i, ax) in ax_out.iter().enumerate() {
if *ax < Axis::New {
@@ -793,22 +836,56 @@ impl<D: Dimension, R: Reader> View<D, R> {
}
}
#[cfg(feature = "transforms")]
let kind = self.transforms.kind();
#[cfg(feature = "transforms")]
let order = if T::ISFLOAT { 3 } else { 0 };
#[cfg(feature = "transforms")]
let mut transform_done = HashSet::new();
for (c, z, t) in iproduct!(&slice_reader[0], &slice_reader[1], &slice_reader[2]) {
if let Some(i) = axes_out_idx[0] {
slice[i] = SliceInfoElem::Index(c)
slice_out[i] = SliceInfoElem::Index(c)
};
if let Some(i) = axes_out_idx[1] {
slice[i] = SliceInfoElem::Index(z)
slice_out[i] = SliceInfoElem::Index(z)
};
if let Some(i) = axes_out_idx[2] {
slice[i] = SliceInfoElem::Index(t)
slice_out[i] = SliceInfoElem::Index(t)
};
#[cfg(feature = "transforms")]
let transform = self
.transforms
.get_transform(self.channel_name(c as usize)?, t as usize)?;
#[cfg(feature = "transforms")]
let transform_now = if let (TransformKind::Channel, Some((&Axis::C, _)), _)
| (TransformKind::Drift, Some((&Axis::T, _)), _)
| (TransformKind::Both, Some((&Axis::C | &Axis::T, _)), _)
| (
TransformKind::Channel | TransformKind::Drift | TransformKind::Both,
_,
false,
) = (&kind, op_czt.first(), xy_dim == 2)
{
transform_done.insert((c, z, t));
transform.map(|transform| (transform, order))
} else {
None
};
#[cfg(feature = "transforms")]
let frame = self.get_cached_frame(
(c % shape.c as isize) as usize,
(z % shape.z as isize) as usize,
(t % shape.t as isize) as usize,
transform_now,
)?;
#[cfg(not(feature = "transforms"))]
let frame = self.get_cached_frame(
(c % shape.c as isize) as usize,
(z % shape.z as isize) as usize,
(t % shape.t as isize) as usize,
)?;
// 0, 1 or 2d frame depending on x or y operations
let arr_frame: Array2<T> = frame.as_ref().try_into()?;
let arr_frame = match xy_dim {
0 => {
@@ -855,35 +932,28 @@ impl<D: Dimension, R: Reader> View<D, R> {
unreachable!("xy cannot be 3d or more");
}
};
if let Some((_, op)) = op_czt.first() {
match op {
Operation::Max => {
array
.slice_mut(slice.as_slice())
.zip_mut_with(&arr_frame, |x, y| {
*x = if *x >= *y { x.clone() } else { y.clone() }
});
}
Operation::Min => {
array
.slice_mut(slice.as_slice())
.zip_mut_with(&arr_frame, |x, y| {
*x = if *x < *y { x.clone() } else { y.clone() }
});
}
Operation::Sum => {
array
.slice_mut(slice.as_slice())
.zip_mut_with(&arr_frame, |x, y| *x += y.clone());
}
Operation::Mean => {
array
.slice_mut(slice.as_slice())
.zip_mut_with(&arr_frame, |x, y| *x += y.clone());
}
// assign arr_frame to output array, apply first of c, z or t operation
match op_czt.first() {
Some((_, Operation::Max)) => {
array
.slice_mut(slice_out.as_slice())
.zip_mut_with(&arr_frame, |x, y| {
*x = if *x >= *y { *x } else { *y }
});
}
} else {
array.slice_mut(slice.as_slice()).assign(&arr_frame)
Some((_, Operation::Min)) => {
array
.slice_mut(slice_out.as_slice())
.zip_mut_with(&arr_frame, |x, y| {
*x = if *x < *y { *x } else { *y }
});
}
Some((_, Operation::Mean | Operation::Sum)) => {
array
.slice_mut(slice_out.as_slice())
.zip_mut_with(&arr_frame, |x, y| *x += *y);
}
_ => array.slice_mut(slice_out.as_slice()).assign(&arr_frame),
}
}
let mut out = Some(array);
@@ -894,23 +964,61 @@ impl<D: Dimension, R: Reader> View<D, R> {
.collect();
for (ax, op) in self.operations.iter().skip(op_xy.len() + op_czt.len()) {
if let Some(idx) = ax_out.remove(ax) {
#[cfg(feature = "transforms")]
let mut arr = out.take().unwrap();
#[cfg(not(feature = "transforms"))]
let arr = out.take().unwrap();
#[cfg(feature = "transforms")]
if let (Axis::C, TransformKind::Channel | TransformKind::Both)
| (Axis::T, TransformKind::Drift | TransformKind::Both) = (ax, &kind)
{
let mut slice_transform = vec![
SliceInfoElem::Slice {
start: 0,
end: None,
step: 1
};
ax_out.len()
];
for (c, z, t) in iproduct!(&slice_reader[0], &slice_reader[1], &slice_reader[2])
{
if let Some(transform) = self
.transforms
.get_transform(self.channel_name(c as usize)?, t as usize)?
&& !transform_done.contains(&(c, z, t))
{
for (axis, i) in [&Axis::C, &Axis::Z, &Axis::T].iter().zip([c, z, t]) {
if let Some(axis_pos) = ax_out.get(axis) {
slice_transform[*axis_pos] = SliceInfoElem::Index(i);
}
}
let untransformed = arr.slice(slice_transform.as_slice());
let transformed = transform
.interpolate_par_order(untransformed, order)?
.mapv(|i| T::from_f64(i).unwrap());
arr.slice_mut(slice_transform.as_slice())
.assign(&transformed);
}
transform_done.insert((c, z, t));
}
}
let a = op.operate(arr, idx)?;
let _ = out.insert(a);
for i in ax_out.values_mut() {
if *i > idx {
*i -= 1;
}
}
let arr = out.take().unwrap();
let a = op.operate(arr, idx)?;
let _ = out.insert(a);
}
}
let n = if let Some((&ax, op)) = op_czt.first()
&& *op == Operation::Mean
{
let n = if let Some((ax, Operation::Mean)) = op_czt.first() {
self.axes
.iter()
.zip(self.slice.iter())
.find(|(a, _)| **a == ax)
.find(|(a, _)| *a == ax)
.and_then(|(_, s)| match s {
SliceInfoElem::Slice { start, end, step } => {
end.map(|e| (((e - start).max(0) / step) as usize).max(1))
@@ -921,12 +1029,29 @@ impl<D: Dimension, R: Reader> View<D, R> {
} else {
1
};
#[cfg(not(feature = "transforms"))]
let array = if n == 1 {
out.take().unwrap()
} else {
let m = T::from_usize(n).unwrap_or_else(|| T::zero());
out.take().unwrap().mapv(|x| x / m.clone())
out.take().unwrap().mapv(|x| x / m)
};
#[cfg(feature = "transforms")]
let mut array = if n == 1 {
out.take().unwrap()
} else {
let m = T::from_usize(n).unwrap_or_else(|| T::zero());
out.take().unwrap().mapv(|x| x / m)
};
#[cfg(feature = "transforms")]
if transform_done.is_empty() {
let c = slice_reader[0].start;
let t = slice_reader[2].start;
if let Some(transform) = self.transforms.get_transform(self.channel_name(c as usize)?, t as usize)? {
array = transform.interpolate_par_order(&array, order)?.mapv(|i| T::from_f64(i).unwrap());
}
}
ArrayCache::global().insert(key.to_owned(), array.clone());
Ok(array)
}
@@ -957,7 +1082,13 @@ impl<D: Dimension, R: Reader> View<D, R> {
.collect())
}
fn frame_key(&self, c: usize, z: usize, t: usize) -> FrameKeyRef<'_> {
fn frame_key<'a>(
&'a self,
c: usize,
z: usize,
t: usize,
#[cfg(feature = "transforms")] transform: Option<(Cow<'a, Transform<IxDyn>>, usize)>,
) -> FrameKeyRef<'a> {
FrameKeyRef {
reader: ReaderKeyRef {
name: self.reader.reader_name(),
@@ -969,6 +1100,8 @@ impl<D: Dimension, R: Reader> View<D, R> {
c,
z,
t,
#[cfg(feature = "transforms")]
transform,
}
}
@@ -985,14 +1118,37 @@ impl<D: Dimension, R: Reader> View<D, R> {
slice: &self.slice,
axes: &self.axes,
operations: &self.operations,
#[cfg(feature = "transforms")]
transforms: &self.transforms,
}
}
fn get_cached_frame(&self, c: usize, z: usize, t: usize) -> Result<Arc<Frame>, Error> {
let key = self.frame_key(c, z, t);
pub(crate) fn get_cached_frame(
&self,
c: usize,
z: usize,
t: usize,
#[cfg(feature = "transforms")] transform: Option<(Cow<Transform<IxDyn>>, usize)>,
) -> Result<Arc<Frame>, Error> {
let key = self.frame_key(
c,
z,
t,
#[cfg(feature = "transforms")]
transform
.as_ref()
.map(|(transform, order)| (Cow::Borrowed(&**transform), *order)),
);
if let Some(frame) = FrameCache::global().get(&key) {
return Ok(frame);
}
#[cfg(feature = "transforms")]
let frame = if let Some((transform, order)) = transform.as_ref() {
self.get_transformed_frame(c, z, t, transform, *order)?
} else {
Arc::new(self.reader.get_frame(c, z, t)?)
};
#[cfg(not(feature = "transforms"))]
let frame = Arc::new(self.reader.get_frame(c, z, t)?);
FrameCache::global().insert(key.to_owned(), frame.clone());
Ok(frame)
@@ -1187,7 +1343,7 @@ impl<R: Reader> Item for View<Ix0, R> {
Array1<T>: MinMax<Output = Array0<T>>,
Array2<T>: MinMax<Output = Array1<T>>,
{
Ok(self.as_array()?.first().ok_or(Error::EmptyView)?.clone())
Ok(*self.as_array()?.first().ok_or(Error::EmptyView)?)
}
}
@@ -1229,7 +1385,7 @@ to_bytes_vec_impl!(
mod tests {
use crate::axes::{Axis, Operation};
use crate::cache::{
ArrayCache, ArrayKey, ArrayKeyRef, FileStamp, FrameCache, ReaderKey, ReaderKeyRef,
ArrayCache, ArrayKey, ArrayKeyRef, FileStamp, FrameCache, FrameKey, ReaderKey, ReaderKeyRef,
};
use crate::error::Error;
use crate::readers::{DynReader, Frame, Reader};
@@ -1241,6 +1397,8 @@ mod tests {
use postcard::{from_bytes, to_stdvec};
use std::path::{Path, PathBuf};
use std::sync::Arc;
#[cfg(feature = "transforms")]
use crate::transforms::Transforms;
fn open(file: &str) -> Result<DynReader, Error> {
let path = std::env::current_dir()?
@@ -1512,19 +1670,19 @@ mod tests {
#[test]
fn frame_cache() -> Result<(), Error> {
let cache = FrameCache::new(2);
let key = |c: usize, z: usize, t: usize| {
(
ReaderKey {
name: "test".to_string(),
path: PathBuf::from("test.tif"),
series: 0,
position: 0,
stamp: FileStamp::default(),
},
c,
z,
t,
)
let key = |c: usize, z: usize, t: usize| FrameKey {
reader_key: ReaderKey {
name: "test".to_string(),
path: PathBuf::from("test.tif"),
series: 0,
position: 0,
stamp: FileStamp::default(),
},
c,
z,
t,
#[cfg(feature = "transforms")]
transform: None,
};
cache.insert(
key(0, 0, 0),
@@ -1561,6 +1719,8 @@ mod tests {
slice: vec![SliceInfoElem::Index(index)],
axes: axes.to_vec(),
operations: vec![],
#[cfg(feature = "transforms")]
transforms: Transforms::default(),
};
let k0 = key("u16", 0, &[Axis::T]);
let k1 = key("u16", 1, &[Axis::T]);
@@ -1587,6 +1747,8 @@ mod tests {
slice: &[SliceInfoElem::Index(0)],
axes: &[Axis::T],
operations: &ops,
#[cfg(feature = "transforms")]
transforms: &Transforms::default(),
};
assert!(cache.get::<f64, _>(&borrowed).is_none());
assert!(cache.get::<u16, _>(&borrowed).is_some());
@@ -1602,12 +1764,37 @@ mod tests {
let (c, z, t) = (shape[0] - 1, shape[1] - 1, shape[2] - 1);
let cache = FrameCache::global();
let a: Array5<usize> = view.clone().try_into()?;
assert!(cache.get(&view.frame_key(c, z, t)).is_some());
assert!(
cache
.get(&view.frame_key(
c,
z,
t,
#[cfg(feature = "transforms")]
None
))
.is_some()
);
let b: Array5<usize> = view.clone().try_into()?;
assert_eq!(a, b);
// a fresh reader on the same file shares the same cache key
let view2 = open(file)?.view();
assert_eq!(view.frame_key(c, z, t), view2.frame_key(c, z, t));
assert_eq!(
view.frame_key(
c,
z,
t,
#[cfg(feature = "transforms")]
None
),
view2.frame_key(
c,
z,
t,
#[cfg(feature = "transforms")]
None
)
);
let c2: Array5<usize> = view2.clone().try_into()?;
assert_eq!(a, c2);
Ok(())