Skip to main content

data_beans/
zarr_io.rs

1//! General-purpose zarr I/O for reading numeric arrays and attributes from
2//! zarr directories or `.zarr.zip` archives.
3//!
4//! Also provides a convenience function for reading coordinate-style data
5//! (row names + numeric matrix) that returns `MatWithNames<DMatrix<f32>>`,
6//! matching the parquet reader interface in legume_numeric::matrix.
7
8use crate::sparse_io::SparseIoBackend;
9use legume_numeric::matrix::traits::MatWithNames;
10use log::info;
11use nalgebra::DMatrix;
12use rand_distr::num_traits::FromPrimitive;
13use std::sync::Arc;
14use zarrs::array::{data_type, Array as ZArray};
15use zarrs::config::MetadataRetrieveVersion;
16use zarrs::filesystem::FilesystemStore;
17use zarrs::storage::ReadableListableStorageTraits as ZReadStorageTraits;
18use zarrs_zip::ZipStorageAdapter;
19
20// ── Store management ────────────────────────────────────────────────────
21
22/// A read store paired with an optional write store (present for directories, `None` for zips).
23pub type ZarrStoreRw = (Arc<dyn ZReadStorageTraits>, Option<Arc<FilesystemStore>>);
24
25/// The directory a zarr zip holds its store under, read from the stored
26/// names, which [`ZipStorageAdapter`] matches too: empty when a metadata file
27/// (`zarr.json`, or v2's `.zgroup` / `.zarray`) sits at the archive's root;
28/// else the single top-level directory (as
29/// [`legume_numeric::matrix::common_io::zip_dir`] writes, whatever the
30/// archive is named), ignoring loose files and `__MACOSX`; else, among
31/// several, the one named like the archive (`foo` for `foo.zarr.zip`, or the
32/// older `foo.zarr`); else empty.
33fn zip_zarr_root(zip_path: &std::path::Path) -> anyhow::Result<String> {
34    const METADATA: [&str; 3] = ["zarr.json", ".zgroup", ".zarray"];
35    let file = std::fs::File::open(zip_path)?;
36    let archive = zip::ZipArchive::new(std::io::BufReader::new(file))?;
37    let mut dirs = std::collections::BTreeSet::new();
38    for name in archive.file_names().flatten() {
39        match name.split_once('/') {
40            None if METADATA.contains(&&*name) => return Ok(String::new()),
41            None | Some(("__MACOSX", _)) => {}
42            Some((top, _)) => {
43                dirs.insert(top.to_owned());
44            }
45        }
46    }
47    if dirs.len() == 1 {
48        return Ok(dirs.pop_first().unwrap_or_default());
49    }
50    let file_name = zip_path.file_name().and_then(|f| f.to_str()).unwrap_or("");
51    let stem = file_name.strip_suffix(".zip").unwrap_or(file_name);
52    let named = [stem.strip_suffix(".zarr").unwrap_or(stem), stem];
53    Ok(named
54        .into_iter()
55        .find(|d| dirs.contains(*d))
56        .unwrap_or_default()
57        .to_owned())
58}
59
60/// The path prefix inside a zarr zip: its [`zip_zarr_root`] with a trailing
61/// `/`, or empty.
62fn detect_zip_zarr_prefix(zip_path: &std::path::Path) -> anyhow::Result<Box<str>> {
63    let root = zip_zarr_root(zip_path)?;
64    Ok(if root.is_empty() {
65        root.into()
66    } else {
67        format!("{root}/").into()
68    })
69}
70
71/// Open a zarr store, returning both a read store and an optional write store.
72///
73/// Supports both `.zarr` directories and `.zarr.zip` archives.
74/// For zip archives, uses [`ZipStorageAdapter`] for direct random-access
75/// reads without extracting to a temp directory; the write store is `None`.
76pub fn open_zarr_store_rw(path: &str) -> anyhow::Result<ZarrStoreRw> {
77    let p = std::path::Path::new(path);
78    let is_zip = p.extension().is_some_and(|e| e == "zip");
79
80    if is_zip {
81        let parent = p
82            .parent()
83            .ok_or_else(|| anyhow::anyhow!("no parent directory for {}", path))?;
84        let filename = p
85            .file_name()
86            .ok_or_else(|| anyhow::anyhow!("no filename for {}", path))?
87            .to_str()
88            .ok_or_else(|| anyhow::anyhow!("non-UTF8 filename in {}", path))?;
89        let zarr_prefix = detect_zip_zarr_prefix(p)?;
90        let fs = Arc::new(FilesystemStore::new(parent)?);
91        let key = zarrs::storage::StoreKey::new(filename)?;
92        info!(
93            "Opening zarr zip store: {} (prefix: {:?})",
94            path, zarr_prefix
95        );
96        Ok((
97            Arc::new(ZipStorageAdapter::new_with_path(fs, key, &*zarr_prefix)?),
98            None,
99        ))
100    } else {
101        let fs = Arc::new(FilesystemStore::new(path)?);
102        Ok((fs.clone(), Some(fs)))
103    }
104}
105
106/// Open a zarr store for reading (convenience wrapper around [`open_zarr_store_rw`]).
107pub fn open_zarr_store(path: &str) -> anyhow::Result<Arc<dyn ZReadStorageTraits>> {
108    open_zarr_store_rw(path).map(|(r, _)| r)
109}
110
111/// The file `-o output` ends as, under `--no-zip` (`zip` false) and
112/// `--backend`: `x.zarr.zip`, `x.zarr` or `x.h5`, for a plan to name.
113pub fn output_path(output: &str, zip: bool, backend: &SparseIoBackend) -> String {
114    let effective = apply_zip_flag(output, zip, backend);
115    if effective.ends_with(".zarr.zip") {
116        return effective.to_string();
117    }
118    crate::hdf5_io::resolve_backend_file(&effective, Some(backend.clone()))
119        .map_or_else(|_| effective.to_string(), |(_, f)| f.to_string())
120}
121
122/// If `zip` is true and the target is zarr, ensure the output path ends with
123/// `.zarr.zip`; otherwise pass through unchanged.
124///
125/// The target is treated as HDF5 (and left alone) when `backend` is
126/// [`SparseIoBackend::HDF5`] or the name already carries a `.h5` / `.hdf5`
127/// suffix. Without the `backend` signal, a bare `-o foo` under `--backend hdf5`
128/// would be zarr-ified here and then flipped back to Zarr by
129/// `resolve_backend_file` — silently ignoring the requested backend.
130pub fn apply_zip_flag(output: &str, zip: bool, backend: &SparseIoBackend) -> Box<str> {
131    // The `.zarr.zip` suffix only makes sense for a zarr target. An HDF5 target
132    // — selected by `--backend hdf5` or an explicit `.h5`/`.hdf5` name — is left
133    // untouched, so a bare `-o foo` under `--backend hdf5` resolves to `foo.h5`
134    // rather than being silently zarr-ified by `resolve_backend_file`.
135    let target_is_hdf5 = matches!(backend, SparseIoBackend::HDF5)
136        || output.ends_with(".h5")
137        || output.ends_with(".hdf5");
138    if zip && !target_is_hdf5 && !output.ends_with(".zarr.zip") {
139        let base = crate::hdf5_io::strip_backend_suffix(output);
140        format!("{}.zarr.zip", base).into()
141    } else {
142        output.into()
143    }
144}
145
146/// Where a handler given `-o output` writes, before it writes anything.
147pub struct OutputTarget {
148    /// The output as the handler takes it: `.zarr.zip` added under zip.
149    pub effective: Box<str>,
150    pub backend: SparseIoBackend,
151    /// The file the backend is written to (for a `.zarr.zip` target, the
152    /// un-zipped `.zarr` directory [`staging_zarr_dir`] names).
153    pub working: Box<str>,
154    /// What the user finds once [`finalize_output`] has run: the archive for
155    /// a `.zarr.zip` target, else the working file.
156    pub found: Box<str>,
157}
158
159/// [`OutputTarget`] for `output`, touching nothing.
160pub fn output_target(
161    output: &str,
162    backend: SparseIoBackend,
163    zip: bool,
164) -> anyhow::Result<OutputTarget> {
165    let effective = apply_zip_flag(output, zip, &backend);
166    let (backend, working) = crate::hdf5_io::resolve_backend_file(&effective, Some(backend))?;
167    let found = if effective.ends_with(".zarr.zip") {
168        effective.clone()
169    } else {
170        working.clone()
171    };
172    Ok(OutputTarget {
173        effective,
174        backend,
175        working,
176        found,
177    })
178}
179
180/// Resolve the output target for a handler that writes a fresh backend, and
181/// clear any pre-existing file at that path. Returns
182/// `(effective_output, backend, working_file)`, where `working_file` is the
183/// path to write (the un-zipped `.zarr` directory for a `.zarr.zip` target) —
184/// pair it with [`finalize_output`] once the backend is populated.
185pub fn prepare_output(
186    output: &str,
187    backend: SparseIoBackend,
188    zip: bool,
189) -> anyhow::Result<(Box<str>, SparseIoBackend, Box<str>)> {
190    use legume_numeric::matrix::common_io::remove_file;
191
192    let t = output_target(output, backend, zip)?;
193    if std::path::Path::new(t.working.as_ref()).exists() {
194        remove_file(&t.working)?;
195    }
196    Ok((t.effective, t.backend, t.working))
197}
198
199/// Finalize a handler's output: re-zip the working `.zarr` directory when the
200/// target is `.zarr.zip`, and return the path the user will actually find (the
201/// archive for a `.zarr.zip` target, otherwise the working file).
202pub fn finalize_output<'a>(
203    working_file: &'a str,
204    effective_output: &'a str,
205) -> anyhow::Result<&'a str> {
206    finalize_zarr_output(working_file, effective_output)?;
207    Ok(if effective_output.ends_with(".zarr.zip") {
208        effective_output
209    } else {
210        working_file
211    })
212}
213
214/// Extract a `.zarr.zip` archive into `target_dir` as a flat `.zarr`
215/// directory: the store under its root directory, with legume-numeric's
216/// checks against entries that would leave the target.
217pub fn extract_zarr_zip(zip_path: &str, target_dir: &str) -> anyhow::Result<()> {
218    let root = zip_zarr_root(std::path::Path::new(zip_path))?;
219    legume_numeric::matrix::common_io::unzip_subdir(zip_path, &root, target_dir)
220}
221
222/// Copy input backend to a writable output location. Transparently extracts
223/// `.zarr.zip` archives into a `.zarr` directory so the output can be opened
224/// read/write.
225pub fn materialize_writable_backend(src: &str, dst: &str) -> anyhow::Result<()> {
226    if src.ends_with(".zarr.zip") {
227        info!("extracting {} → {}", src, dst);
228        extract_zarr_zip(src, dst)
229    } else {
230        legume_numeric::matrix::common_io::recursive_copy(src, dst)
231    }
232}
233
234/// The hidden folder, in the output's own folder, that `.zarr.zip` targets
235/// are built in.
236const STAGE: &str = ".data-beans-stage";
237
238/// The working `.zarr` directory for a `.zarr.zip` target `x.zarr.zip`:
239/// `.data-beans-stage/x.zarr` in the same folder, on the same disk. Built as a
240/// plain `x.zarr` beside the archive, every write made a directory there (and
241/// subset a `.subset_tmp` one) and removed it seconds later, and a sync client
242/// watching the folder, such as Dropbox, put the emptied directories back.
243/// [`make_stage`] marks the stage for Dropbox to leave alone, and
244/// [`finalize_zarr_output`] removes it once the archive is in place.
245pub fn staging_zarr_dir(target: &str) -> String {
246    let path = std::path::Path::new(target);
247    let name = path
248        .file_name()
249        .and_then(|f| f.to_str())
250        .and_then(|f| f.strip_suffix(".zip"))
251        .unwrap_or("data.zarr");
252    path.with_file_name(STAGE)
253        .join(name)
254        .to_string_lossy()
255        .into_owned()
256}
257
258/// `name` in the stage of the folder `beside` is in (or of the stage `beside`
259/// is itself in): where scratch kept with an output goes, such as a subset's
260/// `.subset_tmp` or `qc`'s pool. Pair with [`make_stage`] and [`leave_stage`].
261pub fn in_stage(beside: &str, name: &str) -> String {
262    let stage = match stage_of(beside) {
263        Some(stage) => stage.to_path_buf(),
264        None => std::path::Path::new(beside).with_file_name(STAGE),
265    };
266    stage.join(name).to_string_lossy().into_owned()
267}
268
269/// The stage `working` is built in, when it is one ([`staging_zarr_dir`]).
270fn stage_of(working: &str) -> Option<&std::path::Path> {
271    std::path::Path::new(working)
272        .parent()
273        .filter(|p| p.file_name().is_some_and(|n| n == STAGE))
274}
275
276/// Make the stage `working` is built in, if it is in one, and mark it for
277/// Dropbox to ignore before anything is written in it.
278pub fn make_stage(working: &str) -> anyhow::Result<()> {
279    if let Some(stage) = stage_of(working) {
280        std::fs::create_dir_all(stage)?;
281        mark_sync_ignored(stage);
282    }
283    Ok(())
284}
285
286/// Remove the stage `working` was built in once nothing else is being built
287/// there (another run's build keeps it).
288pub fn leave_stage(working: &str) {
289    if let Some(stage) = stage_of(working) {
290        let _ = std::fs::remove_dir(stage);
291    }
292}
293
294/// Dropbox's "ignored" mark: it neither uploads the directory nor puts it
295/// back. Best effort; elsewhere the mark is a harmless attribute.
296fn mark_sync_ignored(dir: &std::path::Path) {
297    #[cfg(target_os = "macos")]
298    let _ = xattr::set(dir, "com.dropbox.ignored", b"1");
299    #[cfg(target_os = "linux")]
300    let _ = xattr::set(dir, "user.com.dropbox.ignored", b"1");
301    #[cfg(not(any(target_os = "macos", target_os = "linux")))]
302    let _ = dir;
303}
304
305/// If `target_path` ends with `.zarr.zip`, zip the `zarr_dir` into it and
306/// remove the directory and its [`staging_zarr_dir`] stage. Otherwise this is
307/// a no-op (the directory IS the target).
308///
309/// To keep in-zip entry names compact, entries are prefixed with the short
310/// stem (`foo/...` rather than `foo.zarr/...`).
311pub fn finalize_zarr_output(zarr_dir: &str, target_path: &str) -> anyhow::Result<()> {
312    if !target_path.ends_with(".zarr.zip") {
313        return Ok(());
314    }
315
316    let zarr_path = std::path::Path::new(zarr_dir);
317    let stem = zarr_path
318        .file_name()
319        .and_then(|s| s.to_str())
320        .map(|s| s.strip_suffix(".zarr").unwrap_or(s))
321        .ok_or_else(|| anyhow::anyhow!("invalid zarr_dir: {}", zarr_dir))?
322        .to_string();
323
324    // Zip beside the source first (in the stage, so nothing comes and goes in
325    // the output folder) so the source .zarr survives any failure (disk full,
326    // etc.). Only after the zip succeeds do we atomically rename into place
327    // and remove the source — so a crash leaves either the source directory
328    // or the final archive, never neither.
329    let tmp_target = format!("{zarr_dir}.zip.tmp");
330    info!("Zipping zarr output: {} → {}", zarr_dir, target_path);
331    if let Err(e) =
332        legume_numeric::matrix::common_io::zip_dir_as(zarr_dir, &tmp_target, Some(&stem))
333    {
334        let _ = std::fs::remove_file(&tmp_target);
335        return Err(e);
336    }
337    if let Err(e) = std::fs::rename(&tmp_target, target_path) {
338        let _ = std::fs::remove_file(&tmp_target);
339        return Err(e.into());
340    }
341    std::fs::remove_dir_all(zarr_dir)?;
342    leave_stage(zarr_dir);
343    Ok(())
344}
345
346// ── Attribute reading ───────────────────────────────────────────────────
347
348/// Read an attribute from a zarr array node.
349///
350/// ```text
351/// read_zarr_array_attr::<Vec<String>>(store, "/cell_summary", "column_names")
352/// ```
353pub fn read_zarr_array_attr<V: serde::de::DeserializeOwned>(
354    store: Arc<dyn ZReadStorageTraits>,
355    array_path: &str,
356    attr_name: &str,
357) -> anyhow::Result<V> {
358    let arr = ZArray::open_opt(store, array_path, &MetadataRetrieveVersion::Default)?;
359    let attr = arr.attributes().get(attr_name).ok_or_else(|| {
360        anyhow::anyhow!(
361            "attribute '{}' not found on array '{}'",
362            attr_name,
363            array_path
364        )
365    })?;
366    Ok(serde_json::from_value(attr.clone())?)
367}
368
369// ── Flat array readers ──────────────────────────────────────────────────
370
371/// Retrieve a zarr array as a flat `Vec<f32>` (row-major) and its shape.
372///
373/// Handles f32, f64, u32, u64 source types with automatic conversion.
374pub fn read_zarr_flat_f32(
375    store: Arc<dyn ZReadStorageTraits>,
376    key: &str,
377) -> anyhow::Result<(Vec<f32>, Vec<u64>)> {
378    let arr = ZArray::open_opt(store, key, &MetadataRetrieveVersion::Default)?;
379    let shape = arr.shape().to_vec();
380    let subset = arr.subset_all();
381
382    let dt = arr.data_type();
383    let data = if *dt == data_type::float32() {
384        arr.retrieve_array_subset::<Vec<f32>>(&subset)?
385    } else if *dt == data_type::float64() {
386        arr.retrieve_array_subset::<Vec<f64>>(&subset)?
387            .into_iter()
388            .map(|x| x as f32)
389            .collect()
390    } else if *dt == data_type::uint32() {
391        arr.retrieve_array_subset::<Vec<u32>>(&subset)?
392            .into_iter()
393            .map(|x| x as f32)
394            .collect()
395    } else if *dt == data_type::uint64() {
396        arr.retrieve_array_subset::<Vec<u64>>(&subset)?
397            .into_iter()
398            .map(|x| x as f32)
399            .collect()
400    } else {
401        anyhow::bail!("unsupported zarr data type: {:?}", dt)
402    };
403
404    Ok((data, shape))
405}
406
407/// Retrieve a zarr array as a flat `Vec<u32>` and its shape.
408pub fn read_zarr_flat_u32(
409    store: Arc<dyn ZReadStorageTraits>,
410    key: &str,
411) -> anyhow::Result<(Vec<u32>, Vec<u64>)> {
412    let arr = ZArray::open_opt(store, key, &MetadataRetrieveVersion::Default)?;
413    let shape = arr.shape().to_vec();
414    let subset = arr.subset_all();
415
416    let dt = arr.data_type();
417    let data = if *dt == data_type::uint32() {
418        arr.retrieve_array_subset::<Vec<u32>>(&subset)?
419    } else if *dt == data_type::uint64() {
420        arr.retrieve_array_subset::<Vec<u64>>(&subset)?
421            .into_iter()
422            .map(|x| x as u32)
423            .collect()
424    } else {
425        anyhow::bail!("unsupported zarr data type for u32: {:?}", dt)
426    };
427
428    Ok((data, shape))
429}
430
431// ── Generic array/attribute readers (moved from misc.rs) ────────────────
432
433#[cfg(feature = "ndarray")]
434/// Read a full ndarray from zarr storage.
435pub fn read_zarr_ndarray<T>(
436    store: Arc<dyn ZReadStorageTraits>,
437    key_name: &str,
438) -> anyhow::Result<ndarray::ArrayD<T>>
439where
440    T: zarrs::array::ElementOwned + FromPrimitive,
441{
442    let arr = ZArray::open_opt(store, key_name, &MetadataRetrieveVersion::Default)?;
443
444    let dt = arr.data_type();
445    if *dt == data_type::float32() {
446        let array: ndarray::ArrayD<f32> =
447            arr.retrieve_array_subset::<ndarray::ArrayD<f32>>(&arr.subset_all())?;
448        Ok(array.mapv(|x| T::from_f32(x).unwrap()))
449    } else if *dt == data_type::float64() {
450        let array: ndarray::ArrayD<f64> =
451            arr.retrieve_array_subset::<ndarray::ArrayD<f64>>(&arr.subset_all())?;
452        Ok(array.mapv(|x| T::from_f64(x).unwrap()))
453    } else if *dt == data_type::uint32() {
454        let array: ndarray::ArrayD<u32> =
455            arr.retrieve_array_subset::<ndarray::ArrayD<u32>>(&arr.subset_all())?;
456        Ok(array.mapv(|x| T::from_u32(x).unwrap()))
457    } else if *dt == data_type::uint64() {
458        let array: ndarray::ArrayD<u64> =
459            arr.retrieve_array_subset::<ndarray::ArrayD<u64>>(&arr.subset_all())?;
460        Ok(array.mapv(|x| T::from_u64(x).unwrap()))
461    } else {
462        anyhow::bail!("unsupported zarr data type: {:?}", dt)
463    }
464}
465
466/// Read a numeric vector from zarr storage.
467pub fn read_zarr_numerics<T>(
468    store: Arc<dyn ZReadStorageTraits>,
469    key_name: &str,
470) -> anyhow::Result<Vec<T>>
471where
472    T: zarrs::array::ElementOwned + FromPrimitive,
473{
474    let arr = ZArray::open_opt(store, key_name, &MetadataRetrieveVersion::Default)?;
475
476    let dt = arr.data_type();
477    let ret = if *dt == data_type::float32() {
478        arr.retrieve_array_subset::<Vec<f32>>(&arr.subset_all())?
479            .into_iter()
480            .map(|x| T::from_f32(x).unwrap())
481            .collect()
482    } else if *dt == data_type::float64() {
483        arr.retrieve_array_subset::<Vec<f64>>(&arr.subset_all())?
484            .into_iter()
485            .map(|x| T::from_f64(x).unwrap())
486            .collect()
487    } else if *dt == data_type::uint32() {
488        arr.retrieve_array_subset::<Vec<u32>>(&arr.subset_all())?
489            .into_iter()
490            .map(|x| T::from_u32(x).unwrap())
491            .collect()
492    } else if *dt == data_type::uint64() {
493        arr.retrieve_array_subset::<Vec<u64>>(&arr.subset_all())?
494            .into_iter()
495            .map(|x| T::from_u64(x).unwrap())
496            .collect()
497    } else {
498        anyhow::bail!("unsupported zarr data type: {:?}", dt);
499    };
500
501    Ok(ret)
502}
503
504/// Extract an attribute from a zarr group node.
505///
506/// `key_name` is parsed as `"group_path/attr_name"`, e.g.
507/// `"/cell_features/features/id"` → group `"/cell_features/features"`, attr `"id"`.
508pub fn read_zarr_group_attr<V>(
509    store: Arc<dyn ZReadStorageTraits>,
510    key_name: &str,
511) -> anyhow::Result<V>
512where
513    V: serde::de::DeserializeOwned,
514{
515    use anyhow::Context;
516
517    fn parse_key_name(key_name: &str) -> (Box<str>, Box<str>) {
518        let trimmed = key_name.strip_prefix('/').unwrap_or(key_name);
519        match trimmed.rsplit_once('/') {
520            Some((left, right)) => (
521                format!("/{}", left).into_boxed_str(),
522                right.to_string().into_boxed_str(),
523            ),
524            None => (
525                "/".to_string().into_boxed_str(),
526                trimmed.to_string().into_boxed_str(),
527            ),
528        }
529    }
530
531    let (group_name, attr_name) = parse_key_name(key_name);
532
533    let group = zarrs::group::Group::open_opt(
534        store,
535        group_name.as_ref(),
536        &MetadataRetrieveVersion::Default,
537    )
538    .with_context(|| format!("Failed to open group '{}'", group_name))?;
539
540    let attr_value = group
541        .attributes()
542        .get(attr_name.as_ref())
543        .with_context(|| {
544            format!(
545                "Attribute '{}' not found in group '{}'",
546                attr_name, group_name
547            )
548        })?;
549
550    Ok(serde_json::from_value(attr_value.clone())?)
551}
552
553/// Read a string array from zarr storage.
554pub fn read_zarr_strings(
555    store: Arc<dyn ZReadStorageTraits>,
556    key_name: &str,
557) -> anyhow::Result<Vec<Box<str>>> {
558    let arr = ZArray::open_opt(store, key_name, &MetadataRetrieveVersion::Default)?;
559
560    Ok(arr
561        .retrieve_array_subset::<Vec<String>>(&arr.subset_all())?
562        .into_iter()
563        .map(|x| x.into_boxed_str())
564        .collect())
565}
566
567// ── 10x cell ID encoding ────────────────────────────────────────────────
568
569/// Hex-digit → shifted-alpha lookup table for 10x cell ID encoding.
570fn hex_to_shifted_lookup() -> [Option<char>; 256] {
571    let mut lookup = [None; 256];
572    for (i, ch) in "0123456789abcdef".chars().enumerate() {
573        lookup[ch as usize] = Some(
574            [
575                'a', 'b', 'c', 'd', 'e', 'f', 'g', 'h', 'i', 'j', 'k', 'l', 'm', 'n', 'o', 'p',
576            ][i],
577        );
578    }
579    lookup
580}
581
582/// Encode a single `(barcode_u32, suffix_u32)` pair into a 10x cell-ID string.
583fn encode_10x_cell_id(
584    lookup: &[Option<char>; 256],
585    barcode: u32,
586    suffix: u32,
587) -> anyhow::Result<Box<str>> {
588    let barcode: String = format!("{:08x}", barcode)
589        .chars()
590        .map(|ch| lookup[ch as usize].ok_or_else(|| anyhow::anyhow!("invalid hex char: {}", ch)))
591        .collect::<anyhow::Result<String>>()?;
592    Ok(format!("{}-{}", barcode, suffix).into_boxed_str())
593}
594
595#[cfg(feature = "ndarray")]
596/// Parse a 10x Xenium `[N, 2]` u32 cell_id ndarray into string barcodes.
597///
598/// See [10x docs](https://www.10xgenomics.com/support/software/xenium-onboard-analysis/3.4/advanced/xoa-output-zarr#cellID).
599pub fn parse_10x_cell_id(
600    input: ndarray::ArrayView<u32, ndarray::IxDyn>,
601) -> anyhow::Result<Vec<Box<str>>> {
602    anyhow::ensure!(
603        input.ndim() == 2 && input.shape()[1] == 2,
604        "Must be 2D with shape [N, 2]"
605    );
606    let lookup = hex_to_shifted_lookup();
607    input
608        .outer_iter()
609        .map(|row| encode_10x_cell_id(&lookup, row[0], row[1]))
610        .collect()
611}
612
613/// Parse a flat row-major `[N, 2]` u32 buffer into 10x Xenium cell-ID strings.
614///
615/// Same encoding as `parse_10x_cell_id` but works on a flat slice.
616pub fn parse_10x_cell_id_flat(data: &[u32], nrows: usize) -> anyhow::Result<Vec<Box<str>>> {
617    anyhow::ensure!(data.len() == nrows * 2, "cell_id buffer size mismatch");
618    let lookup = hex_to_shifted_lookup();
619    (0..nrows)
620        .map(|i| encode_10x_cell_id(&lookup, data[i * 2], data[i * 2 + 1]))
621        .collect()
622}
623
624// ── Matrix reading (generic) ────────────────────────────────────────────
625
626/// Build a column-major `DMatrix<f32>` by selecting columns from a flat
627/// row-major `[nrows, ncols_total]` buffer.
628fn select_columns_to_dmatrix(
629    data: &[f32],
630    nrows: usize,
631    ncols_total: usize,
632    selected: &[usize],
633) -> anyhow::Result<DMatrix<f32>> {
634    let k = selected.len();
635    let mut col_major = vec![0.0f32; nrows * k];
636    for (col_out, &col_in) in selected.iter().enumerate() {
637        anyhow::ensure!(col_in < ncols_total, "column index {} out of range", col_in);
638        for row in 0..nrows {
639            col_major[row + col_out * nrows] = data[row * ncols_total + col_in];
640        }
641    }
642    Ok(DMatrix::from_vec(nrows, k, col_major))
643}
644
645/// Resolve column selection: explicit indices, name lookup, or default.
646fn resolve_columns(
647    column_indices: &[usize],
648    column_names: &[Box<str>],
649    all_col_names: &[Box<str>],
650    default_cols: &[usize],
651) -> anyhow::Result<Vec<usize>> {
652    if !column_indices.is_empty() {
653        Ok(column_indices.to_vec())
654    } else if !column_names.is_empty() {
655        // Keep only names that exist — allows generous defaults
656        // covering multiple platforms (e.g. Visium + Xenium).
657        let matched: Vec<usize> = column_names
658            .iter()
659            .filter_map(|name| all_col_names.iter().position(|c| c == name))
660            .collect();
661        if matched.is_empty() {
662            anyhow::bail!(
663                "none of the requested columns {:?} found (available: {:?})",
664                column_names,
665                all_col_names
666            );
667        }
668        Ok(matched)
669    } else {
670        Ok(default_cols.to_vec())
671    }
672}
673
674/// Read a numeric matrix from a zarr array, selecting columns by index or
675/// name (looked up via an attribute on the array).
676///
677/// * `file_path` — zarr directory or `.zarr.zip`
678/// * `data_array` — path to the numeric `[N, C]` array (e.g. `/cell_summary`)
679/// * `col_names_attr` — attribute name on `data_array` that holds column
680///   names (e.g. `"column_names"`); pass `None` to skip name lookup
681/// * `row_names_array` — path to the row-name array (e.g. `/cell_id`);
682///   `None` generates numeric row names `"0", "1", …`
683/// * `row_names_10x` — if `true`, parse `row_names_array` as 10x u32 `[N,2]`
684///   cell IDs; if `false`, read as string array
685/// * `column_indices` / `column_names` — which columns to select
686pub fn read_zarr_matrix(
687    file_path: &str,
688    data_array: &str,
689    col_names_attr: Option<&str>,
690    row_names_array: Option<&str>,
691    row_names_10x: bool,
692    column_indices: &[usize],
693    column_names: &[Box<str>],
694) -> anyhow::Result<MatWithNames<DMatrix<f32>>> {
695    let store = open_zarr_store(file_path)?;
696
697    // Read the data array
698    let (flat, shape) = read_zarr_flat_f32(store.clone(), data_array)?;
699    let nrows = shape[0] as usize;
700    let ncols_total = if shape.len() > 1 {
701        shape[1] as usize
702    } else {
703        1
704    };
705
706    // Read column names from attribute (if available)
707    let all_col_names: Vec<Box<str>> = if let Some(attr) = col_names_attr {
708        read_zarr_array_attr(store.clone(), data_array, attr).unwrap_or_else(|_| {
709            (0..ncols_total)
710                .map(|i| i.to_string().into_boxed_str())
711                .collect()
712        })
713    } else {
714        (0..ncols_total)
715            .map(|i| i.to_string().into_boxed_str())
716            .collect()
717    };
718
719    // Resolve column selection
720    let selected = resolve_columns(column_indices, column_names, &all_col_names, &[0, 1])?;
721    let sel_names: Vec<Box<str>> = selected
722        .iter()
723        .map(|&i| {
724            all_col_names
725                .get(i)
726                .cloned()
727                .unwrap_or_else(|| i.to_string().into_boxed_str())
728        })
729        .collect();
730
731    // Build the matrix
732    let mat = select_columns_to_dmatrix(&flat, nrows, ncols_total, &selected)?;
733
734    // Read row names
735    let row_names = if let Some(rn_path) = row_names_array {
736        if row_names_10x {
737            let (id_data, id_shape) = read_zarr_flat_u32(store.clone(), rn_path)?;
738            anyhow::ensure!(
739                id_shape[0] as usize == nrows,
740                "row names array rows ({}) != data rows ({})",
741                id_shape[0],
742                nrows
743            );
744            parse_10x_cell_id_flat(&id_data, nrows)?
745        } else {
746            read_zarr_strings(store.clone(), rn_path)?
747        }
748    } else {
749        (0..nrows).map(|i| i.to_string().into_boxed_str()).collect()
750    };
751
752    info!(
753        "Read {} x {} from zarr {}{}: {:?}",
754        nrows,
755        selected.len(),
756        file_path,
757        data_array,
758        sel_names
759    );
760
761    Ok(MatWithNames {
762        rows: row_names,
763        cols: sel_names,
764        mat,
765    })
766}
767
768// ── Xenium convenience ──────────────────────────────────────────────────
769
770/// Read cell coordinates from a Xenium-style zarr file.
771///
772/// Shorthand for [`read_zarr_matrix`] with Xenium defaults:
773///   - data array: `/cell_summary`
774///   - column names attribute: `column_names`
775///   - row names: `/cell_id` (10x encoded)
776pub fn read_zarr_coordinates(
777    file_path: &str,
778    column_indices: &[usize],
779    column_names: &[Box<str>],
780) -> anyhow::Result<MatWithNames<DMatrix<f32>>> {
781    read_zarr_matrix(
782        file_path,
783        "/cell_summary",
784        Some("column_names"),
785        Some("/cell_id"),
786        true,
787        column_indices,
788        column_names,
789    )
790}
791
792#[cfg(test)]
793mod tests {
794    use super::*;
795
796    fn xenium_path() -> Option<std::path::PathBuf> {
797        let p = std::path::Path::new(env!("CARGO_MANIFEST_DIR"))
798            .parent()
799            .unwrap()
800            .join("docs/temp/cells.zarr.zip");
801        p.exists().then_some(p)
802    }
803
804    #[test]
805    fn test_parse_10x_cell_id_flat() {
806        let data = vec![16844u32, 1, 22527, 1];
807        let ids = parse_10x_cell_id_flat(&data, 2).unwrap();
808        assert_eq!(ids.len(), 2);
809        assert!(ids[0].ends_with("-1"), "got: {}", ids[0]);
810        assert!(ids[1].ends_with("-1"), "got: {}", ids[1]);
811        assert_eq!(ids[0].len(), 10);
812    }
813
814    #[test]
815    fn test_parse_10x_cell_id_flat_mismatch() {
816        assert!(parse_10x_cell_id_flat(&[1, 2, 3], 2).is_err());
817    }
818
819    #[test]
820    fn test_read_zarr_coordinates_by_name() {
821        let Some(p) = xenium_path() else { return };
822        let result = read_zarr_coordinates(
823            p.to_str().unwrap(),
824            &[],
825            &["cell_centroid_x".into(), "cell_centroid_y".into()],
826        )
827        .unwrap();
828
829        assert_eq!(result.mat.ncols(), 2);
830        assert!(result.mat.nrows() > 0);
831        assert_eq!(result.rows.len(), result.mat.nrows());
832        assert_eq!(result.cols[0].as_ref(), "cell_centroid_x");
833        assert_eq!(result.cols[1].as_ref(), "cell_centroid_y");
834        assert!(result.mat.min() >= 0.0);
835    }
836
837    #[test]
838    fn test_read_zarr_coordinates_by_index() {
839        let Some(p) = xenium_path() else { return };
840        let result = read_zarr_coordinates(p.to_str().unwrap(), &[0, 1], &[]).unwrap();
841        assert_eq!(result.mat.ncols(), 2);
842        assert_eq!(result.cols[0].as_ref(), "cell_centroid_x");
843    }
844
845    #[test]
846    fn test_read_zarr_coordinates_default() {
847        let Some(p) = xenium_path() else { return };
848        let result = read_zarr_coordinates(p.to_str().unwrap(), &[], &[]).unwrap();
849        assert_eq!(result.mat.ncols(), 2);
850    }
851
852    #[test]
853    fn test_read_zarr_coordinates_bad_column() {
854        let Some(p) = xenium_path() else { return };
855        let result = read_zarr_coordinates(
856            p.to_str().unwrap(),
857            &[],
858            &["nonexistent".to_string().into_boxed_str()],
859        );
860        assert!(result.is_err());
861    }
862
863    #[test]
864    fn test_read_zarr_coordinates_mixed_defaults() {
865        // Simulates pinto's generous defaults: Visium names that won't match
866        // plus Xenium names that will — only matching names should be kept.
867        let Some(p) = xenium_path() else { return };
868        let result = read_zarr_coordinates(
869            p.to_str().unwrap(),
870            &[],
871            &[
872                "pxl_row_in_fullres".into(),
873                "pxl_col_in_fullres".into(),
874                "cell_centroid_x".into(),
875                "cell_centroid_y".into(),
876            ],
877        )
878        .unwrap();
879        assert_eq!(result.mat.ncols(), 2);
880        assert_eq!(result.cols[0].as_ref(), "cell_centroid_x");
881        assert_eq!(result.cols[1].as_ref(), "cell_centroid_y");
882        assert!(result.mat.nrows() > 0);
883    }
884
885    #[test]
886    fn test_read_zarr_array_attr() {
887        let Some(p) = xenium_path() else { return };
888        let store = open_zarr_store(p.to_str().unwrap()).unwrap();
889        let names: Vec<Box<str>> =
890            read_zarr_array_attr(store, "/cell_summary", "column_names").unwrap();
891        assert!(!names.is_empty());
892        assert!(names.contains(&"cell_centroid_x".to_string().into_boxed_str()));
893    }
894
895    #[test]
896    fn test_read_zarr_matrix_generic() {
897        let Some(p) = xenium_path() else { return };
898        // Read cell_summary selecting columns 2,3 (cell_area, nucleus_centroid_x)
899        let result = read_zarr_matrix(
900            p.to_str().unwrap(),
901            "/cell_summary",
902            Some("column_names"),
903            Some("/cell_id"),
904            true,
905            &[2, 3],
906            &[],
907        )
908        .unwrap();
909        assert_eq!(result.mat.ncols(), 2);
910        assert_eq!(result.cols[0].as_ref(), "cell_area");
911        assert_eq!(result.cols[1].as_ref(), "nucleus_centroid_x");
912    }
913
914    /// (entry name, contents) pairs for a test zip.
915    type Entries<'a> = &'a [(&'a str, &'a [u8])];
916
917    /// A zip `name` in `dir` of the given entries, written as is.
918    fn zip_of(dir: &tempfile::TempDir, name: &str, entries: Entries) -> std::path::PathBuf {
919        use std::io::Write;
920        let path = dir.path().join(name);
921        let mut w = zip::ZipWriter::new(std::fs::File::create(&path).unwrap());
922        let opts = zip::write::SimpleFileOptions::default();
923        for (entry, bytes) in entries {
924            w.start_file(*entry, opts).unwrap();
925            w.write_all(bytes).unwrap();
926        }
927        w.finish().unwrap();
928        path
929    }
930
931    fn extract(zip: &std::path::Path, out: &std::path::Path) -> anyhow::Result<()> {
932        extract_zarr_zip(zip.to_str().unwrap(), out.to_str().unwrap())
933    }
934
935    /// An entry that climbs out of the target (`run/../../x`, behind the
936    /// store's own prefix) refuses the whole archive before anything is
937    /// written.
938    #[test]
939    fn extracting_refuses_an_entry_that_leaves_the_target() {
940        let dir = tempfile::tempdir().unwrap();
941        let zip = zip_of(
942            &dir,
943            "run.zarr.zip",
944            &[
945                ("run/zarr.json", b"{}"),
946                ("run/../../escape.txt", b"escaped"),
947            ],
948        );
949        let target = dir.path().join("out").join("run.zarr");
950        let err = extract(&zip, &target).unwrap_err().to_string();
951        assert!(err.contains("no safe path"), "{err}");
952        assert!(
953            !dir.path().join("out").join("escape.txt").exists()
954                && !dir.path().join("escape.txt").exists(),
955            "an entry was written outside the target"
956        );
957        assert!(
958            !target.join("zarr.json").exists(),
959            "a refused archive left files behind"
960        );
961    }
962
963    /// A store zipped by `zip_dir_as` behind its `run/` root, or by `zip_dir`
964    /// behind the older `run.zarr/` one, extracts back to the same files.
965    #[test]
966    fn a_zipped_store_extracts_back_to_its_files() {
967        for root in [Some("run"), None] {
968            let dir = tempfile::tempdir().unwrap();
969            let src = dir.path().join("run.zarr");
970            std::fs::create_dir_all(src.join("X/c/0")).unwrap();
971            std::fs::write(src.join("zarr.json"), b"{}").unwrap();
972            std::fs::write(src.join("X/c/0/1"), b"chunk").unwrap();
973            let zip = dir.path().join("run.zarr.zip");
974            legume_numeric::matrix::common_io::zip_dir_as(
975                src.to_str().unwrap(),
976                zip.to_str().unwrap(),
977                root,
978            )
979            .unwrap();
980            let out = dir.path().join("back.zarr");
981            extract(&zip, &out).unwrap();
982            assert_eq!(
983                std::fs::read(out.join("zarr.json")).unwrap(),
984                b"{}",
985                "{root:?}"
986            );
987            assert_eq!(
988                std::fs::read(out.join("X/c/0/1")).unwrap(),
989                b"chunk",
990                "{root:?}"
991            );
992        }
993    }
994
995    /// The store's prefix is read from the entries: the archive's single
996    /// top-level directory, so a renamed archive keeps it and a store without
997    /// a root group (a nested array) still gets its root, with loose files
998    /// and an OS's metadata beside it ignored and left out of the extracted
999    /// store; beside other directories (a backup), the one named like the
1000    /// archive; and none for a store at the archive's root.
1001    #[test]
1002    fn the_store_prefix_comes_from_the_entries() {
1003        let dir = tempfile::tempdir().unwrap();
1004        let cases: [(&str, Entries, &str); 7] = [
1005            (
1006                "nested.zarr.zip",
1007                &[("run/m/X/zarr.json", b"{}"), ("run/m/X/c/0", b"chunk")],
1008                "run/",
1009            ),
1010            (
1011                "run.zarr.zip",
1012                &[("run/zarr.json", b"{}"), ("run.zarr.bak/zarr.json", b"{}")],
1013                "run/",
1014            ),
1015            (
1016                "renamed.zarr.zip",
1017                &[("run/zarr.json", b"{}"), ("run/X/zarr.json", b"{}")],
1018                "run/",
1019            ),
1020            (
1021                "no_group.zarr.zip",
1022                &[("run/X/zarr.json", b"{}"), ("run/rows/zarr.json", b"{}")],
1023                "run/",
1024            ),
1025            (
1026                "stray.zarr.zip",
1027                &[
1028                    ("README.txt", b""),
1029                    ("__MACOSX/run/._zarr.json", b""),
1030                    ("run/zarr.json", b"{}"),
1031                ],
1032                "run/",
1033            ),
1034            (
1035                "v2.zarr.zip",
1036                &[("old.zarr/.zgroup", b"{}"), ("old.zarr/X/.zarray", b"{}")],
1037                "old.zarr/",
1038            ),
1039            (
1040                "at_root.zarr.zip",
1041                &[("zarr.json", b"{}"), ("X/zarr.json", b"{}")],
1042                "",
1043            ),
1044        ];
1045        for (name, entries, prefix) in cases {
1046            let zip = zip_of(&dir, name, entries);
1047            assert_eq!(&*detect_zip_zarr_prefix(&zip).unwrap(), prefix, "{name}");
1048        }
1049        for (name, kept, dropped) in [
1050            ("stray.zarr.zip", "zarr.json", "README.txt"),
1051            ("no_group.zarr.zip", "rows/zarr.json", "run"),
1052        ] {
1053            let out = dir.path().join(name).with_extension("out");
1054            extract(&dir.path().join(name), &out).unwrap();
1055            assert!(
1056                out.join(kept).exists() && !out.join(dropped).exists(),
1057                "{name}"
1058            );
1059        }
1060    }
1061}