legume-numeric 0.8.11

Numeric and ML foundation for the legume ecosystem (matrix, Leiden, candle, MCMC)
Documentation
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
/// Default `env_logger` filter when verbose mode is on. Promotes everything to
/// info. (The k-NN backend, instant-distance, emits no `log` output, so no
/// per-crate pin is needed here.) Override via `RUST_LOG=...`.
pub const VERBOSE_LOG_FILTER: &str = "info";

/// Default `env_logger` filter when verbose mode is off.
pub const QUIET_LOG_FILTER: &str = "warn";

use flate2::read::MultiGzDecoder;
use rayon::prelude::*;
use std::ffi::OsStr;
use std::fs::File;
use std::io::{BufRead, BufReader, BufWriter, Write};
use std::path::Path;
use tempfile::tempdir;

#[cfg(test)]
mod tests;

/// Define a Delimiter enum to handle both &str and `Vec<char>`
pub enum Delimiter {
    Str(String),
    Chars(Vec<char>),
}

impl From<&str> for Delimiter {
    fn from(s: &str) -> Self {
        Delimiter::Str(s.to_string())
    }
}

impl From<Vec<char>> for Delimiter {
    fn from(chars: Vec<char>) -> Self {
        Delimiter::Chars(chars)
    }
}

impl From<&[char]> for Delimiter {
    fn from(chars: &[char]) -> Self {
        Delimiter::Chars(chars.to_vec())
    }
}

impl<const N: usize> From<&[char; N]> for Delimiter {
    fn from(chars: &[char; N]) -> Self {
        Delimiter::Chars(chars.to_vec())
    }
}

///
/// Read every line of the input_file into memory
///
/// * `input_file` - file name--either gzipped or not
///
pub fn read_lines(input_file_path: &str) -> anyhow::Result<Vec<Box<str>>> {
    let buf: Box<dyn BufRead> = open_buf_reader(input_file_path)?;
    let mut lines = vec![];
    for x in buf.lines() {
        lines.push(x?.into_boxed_str());
    }
    Ok(lines)
}

///
/// Write every line into the output_file
///
/// * `lines` - vector of lines
/// * `output_file` - file name--either gzipped or not
///
pub fn write_lines(lines: &Vec<Box<str>>, output_file_path: &str) -> anyhow::Result<()> {
    write_types(lines, output_file_path)
}

///
/// Write every line into the output_file
///
/// * `lines` - vector of lines
/// * `output_file` - file name--either gzipped or not
///
pub fn write_types<T>(lines: &Vec<T>, output_file_path: &str) -> anyhow::Result<()>
where
    T: std::fmt::Display,
{
    let mut buf = open_buf_writer(output_file_path)?;
    for line in lines {
        if let Err(e) = writeln!(buf, "{}", line) {
            if e.kind() == std::io::ErrorKind::BrokenPipe {
                return Ok(());
            } else {
                return Err(anyhow::anyhow!("unexpected error: {}", e));
            }
        }
    }
    buf.flush()?;
    Ok(())
}

pub struct ReadLinesOut<T: Send> {
    pub lines: Vec<Vec<T>>,
    pub header: Vec<Box<str>>,
}

///
/// Generic function to read lines and parse them into a vector of words or types.
///
/// * `input_file` - file name--either gzipped or not
/// * `hdr_line` - location of a header line (-1 = no header line)
/// * `parse_fn` - function to parse each line into the desired type
///
pub fn read_lines_of_words_generic<T>(
    input_file: &str,
    hdr_line: i64,
    parse_header_fn: impl Fn(&str) -> Vec<Box<str>> + Sync,
    parse_fn: impl Fn(&str) -> Vec<T> + Sync,
) -> anyhow::Result<ReadLinesOut<T>>
where
    T: Send,
{
    let buf_reader: Box<dyn BufRead> = open_buf_reader(input_file)?;

    fn is_not_comment_line(line: &str) -> bool {
        if line.starts_with('#') || line.starts_with('%') {
            return false;
        }
        true
    }

    let lines_raw: Vec<Box<str>> = buf_reader
        .lines()
        .map_while(Result::ok)
        .map(|x| x.into_boxed_str())
        .filter(|x| is_not_comment_line(x.as_ref()))
        .collect();

    let mut header = vec![];

    // Parsing takes more time, so split them into parallel jobs
    let mut lines: Vec<(usize, Vec<T>)> = if hdr_line < 0 {
        lines_raw
            .iter()
            .enumerate()
            .par_bridge()
            .map(|(i, s)| (i, parse_fn(s)))
            .collect()
    } else {
        let n_skip = hdr_line as usize;
        if lines_raw.len() < (n_skip + 1) {
            return Err(anyhow::anyhow!("not enough data"));
        }

        header.extend(parse_header_fn(&lines_raw[n_skip]));

        lines_raw[(n_skip + 1)..]
            .iter()
            .enumerate()
            .par_bridge()
            .map(|(i, s)| (i, parse_fn(s)))
            .collect()
    };

    if lines.len() > 100_000 {
        lines.par_sort_by_key(|&(i, _)| i);
    } else {
        lines.sort_by_key(|&(i, _)| i);
    }

    let lines = lines.into_iter().map(|(_, x)| x).collect();
    Ok(ReadLinesOut { lines, header })
}

///
/// Specialized function to read lines and parse them into a vector of types.
///
/// * `input_file` - file name--either gzipped or not
/// * `delim` - delimiter
/// * `hdr_line` - location of a header line (-1 = no header line)
///
pub fn read_lines_of_types<T>(
    input_file: &str,
    delim: impl Into<Delimiter>,
    hdr_line: i64,
) -> anyhow::Result<ReadLinesOut<T>>
where
    T: Send + std::str::FromStr + std::fmt::Display,
    <T as std::str::FromStr>::Err: std::fmt::Debug,
{
    let delim = delim.into(); // Convert the input delimiter into the Delimiter enum

    let parse_fn = move |line: &str| -> Vec<T> {
        match &delim {
            Delimiter::Str(s) => line
                .split(s.as_str())
                .map(|x| x.parse::<T>().expect("failed to parse"))
                .collect(),
            Delimiter::Chars(chars) => line
                .split(chars.as_slice())
                .map(|x| x.parse::<T>().expect("failed to parse"))
                .collect(),
        }
    };

    let parse_header_fn = |line: &str| -> Vec<Box<str>> {
        line.split_whitespace()
            .map(|x| x.to_owned().into_boxed_str())
            .collect()
    };

    read_lines_of_words_generic(input_file, hdr_line, parse_header_fn, parse_fn)
}

///
/// Specialized function to read lines and parse them into a vector of words.
///
/// * `input_file` - file name--either gzipped or not
/// * `hdr_line` - location of a header line (-1 = no header line)
///
pub fn read_lines_of_words(
    input_file: &str,
    hdr_line: i64,
) -> anyhow::Result<ReadLinesOut<Box<str>>> {
    let parse_fn = |line: &str| -> Vec<Box<str>> {
        line.split_whitespace()
            .map(|x| x.to_owned().into_boxed_str())
            .collect()
    };

    read_lines_of_words_generic(input_file, hdr_line, parse_fn, parse_fn)
}

///
/// Specialized function to read lines and parse them into a vector of words.
///
/// * `input_file` - file name--either gzipped or not
/// * `delim` - delimiter
/// * `hdr_line` - location of a header line (-1 = no header line)
///
/// Trim a field and strip one pair of symmetric quotes.
///
/// Only a field quoted at BOTH ends is a quoted field. Trimming either end
/// on its own corrupts content that merely happens to start or finish with a
/// quote: a GTF attribute column reads `gene_id "X"; gene_name "Y"`, which
/// ends in a quote it needs, and losing it leaves the attribute unparseable.
/// This is the same rule the delimited tokenizer applies to every field, so
/// any caller inspecting raw lines (header sniffing, previews) can match it.
pub fn unquote_field(x: &str) -> &str {
    let t = x.trim();
    for q in ['"', '\''] {
        if t.len() >= 2 && t.starts_with(q) && t.ends_with(q) {
            return &t[q.len_utf8()..t.len() - q.len_utf8()];
        }
    }
    t
}

pub fn read_lines_of_words_delim(
    input_file: &str,
    delim: impl Into<Delimiter>,
    hdr_line: i64,
) -> anyhow::Result<ReadLinesOut<Box<str>>> {
    let delim = delim.into(); // Convert the input delimiter into the Delimiter enum

    // Outer quotes come off here, once, rather than at each consumer.
    //
    // A csv writer that quotes every field yields `"x"` where the caller asked
    // for `x`, so a name match fails and a reader that falls back to reading by
    // POSITION then silently takes whichever columns happen to sit there. The
    // name-list reader in this same file already unquoted; the general reader
    // did not, which is the inconsistency this removes.
    //
    // This splitter cannot honour a quoted field containing the delimiter
    // anyway, so trimming the outer quotes loses nothing it had.
    // Every field is also trimmed. That is wider than unquoting and worth
    // stating: it is what removes the trailing \r on CRLF input, and it applies
    // to every delimited file the workspace reads through here.
    //
    // Only a field quoted at BOTH ends is a quoted field. Trimming either end
    // on its own corrupts content that merely happens to start or finish with a
    // quote: a GTF attribute column reads `gene_id "X"; gene_name "Y"`, which
    // ends in a quote it needs, and losing it leaves the attribute unparseable.
    let parse_fn = |line: &str| -> Vec<Box<str>> {
        match &delim {
            Delimiter::Str(s) => line
                .split(s.as_str())
                .map(|x| unquote_field(x).to_owned().into_boxed_str())
                .collect(),
            Delimiter::Chars(chars) => line
                .split(chars.as_slice())
                .map(|x| unquote_field(x).to_owned().into_boxed_str())
                .collect(),
        }
    };

    read_lines_of_words_generic(input_file, hdr_line, parse_fn, parse_fn)
}

////////////////////////////////////////////////////////////
// Name lists (a column of feature / gene names)           //
////////////////////////////////////////////////////////////

/// Header labels recognized as the name-bearing column of a name list, matched
/// case-insensitively.
const NAME_LIST_HEADERS: [&str; 12] = [
    "gene",
    "genes",
    "gene_name",
    "gene_names",
    "gene_id",
    "gene_symbol",
    "feature",
    "features",
    "feature_name",
    "symbol",
    "name",
    "id",
];

/// Position of the name-bearing column in `header`, if one is labelled.
fn name_list_column(header: &[Box<str>]) -> Option<usize> {
    header.iter().position(|h| {
        let h = h.trim().trim_matches('"').to_ascii_lowercase();
        NAME_LIST_HEADERS.contains(&h.as_str())
    })
}

/// Read a flat list of names (genes / features) from a file, keeping one column
/// and dropping every other.
///
/// The format is inferred from the extension: `.parquet`, else delimited text —
/// tab, comma, or whitespace, optionally gzipped (`.txt`, `.tsv`, `.csv`,
/// `.tsv.gz`, …). A header whose label is gene-like (`gene`, `feature`,
/// `symbol`, …) selects the column to read; without one the first column is
/// used. Every other column is ignored, so a two-column `gene<TAB>celltype`
/// marker table doubles as a plain gene list. Names are de-duplicated, keeping
/// first-seen order.
///
/// Names are returned verbatim — matching them against a data vocabulary
/// (symbol / Ensembl / case) is the caller's job.
pub fn read_name_list(file_path: &str) -> anyhow::Result<Vec<Box<str>>> {
    let is_parquet = Path::new(file_path)
        .extension()
        .and_then(OsStr::to_str)
        .is_some_and(|e| e.eq_ignore_ascii_case("parquet"));

    let names: Vec<Box<str>> = if is_parquet {
        let header = crate::matrix::parquet::peek_parquet_field_names(file_path)?;
        let col = name_list_column(&header).unwrap_or(0);
        crate::matrix::parquet::read_parquet_string_column(file_path, col)?
    } else {
        let raw = read_lines(file_path)?;
        let data_lines: Vec<&str> = raw
            .iter()
            .map(|line| line.trim())
            .filter(|line| !line.is_empty() && !line.starts_with('#') && !line.starts_with('%'))
            .collect();

        // Sniff ONE delimiter from the first row — tab, else comma, else whitespace.
        // Splitting on all three at once would tear a value like `T cell` into two
        // fields and shift every column after it.
        let delim: &[char] = match data_lines.first() {
            Some(first) if first.contains('\t') => &['\t'],
            Some(first) if first.contains(',') => &[','],
            _ => &[' ', '\t'],
        };

        // Blank fields (repeated or leading delimiters) are dropped before the column
        // is indexed, so ragged indentation does not shift it either.
        let rows: Vec<Vec<&str>> = data_lines
            .iter()
            .map(|line| {
                line.split(delim)
                    .map(|w| w.trim().trim_matches('"'))
                    .filter(|w| !w.is_empty())
                    .collect::<Vec<&str>>()
            })
            .filter(|row| !row.is_empty())
            .collect();

        // A gene-like label in the first row means it is a header: read the column it
        // names and skip it. Otherwise the file is headerless and data starts at row 0
        // — treating that first row as a header would silently drop a gene.
        let header: Vec<Box<str>> = rows
            .first()
            .map(|row| row.iter().map(|w| (*w).into()).collect())
            .unwrap_or_default();
        let (col, skip) = match name_list_column(&header) {
            Some(col) => (col, 1),
            None => (0, 0),
        };

        rows.iter()
            .skip(skip)
            .filter_map(|row| row.get(col).map(|w| (*w).into()))
            .collect()
    };

    let mut seen: std::collections::HashSet<Box<str>> = std::collections::HashSet::new();
    let names: Vec<Box<str>> = names
        .into_iter()
        .filter(|n| !n.is_empty())
        .filter(|n| seen.insert(n.clone()))
        .collect();

    if names.is_empty() {
        return Err(anyhow::anyhow!("no names found in {file_path}"));
    }
    Ok(names)
}

///
/// Open a file for reading, and return a buffered reader
/// * `input_file` - file name--either gzipped or not
pub fn open_buf_reader(input_file: &str) -> anyhow::Result<Box<dyn BufRead>> {
    // take a look at the extension
    // return buffered reader accordingly
    let ext = Path::new(input_file).extension().and_then(|x| x.to_str());
    match ext {
        Some("gz") | Some("bgz") | Some("bgzf") => {
            let input_file = File::open(input_file)?;
            // `MultiGzDecoder`, not `GzDecoder`: BGZF and anything else
            // written by `bgzip` is a concatenation of gzip members, and a
            // decoder that stops after the first one drops the rest of the
            // file silently, with no error and no short read.
            let decoder = MultiGzDecoder::new(input_file);
            // A wide buffer over the inflater: a body read one line at a time
            // would otherwise refill through it every few kilobytes.
            Ok(Box::new(BufReader::with_capacity(1 << 20, decoder)))
        }
        _ => {
            // dbg!(input_file);
            let input_file = File::open(input_file)?;
            Ok(Box::new(BufReader::new(input_file)))
        }
    }
}

////////////////////////////////////
// First-line inspection of a table //
////////////////////////////////////

/// The first line of a delimited file, split on `delimiters`, unquoted and
/// trimmed, empty fields dropped. Reads through gzip.
///
/// For error messages and for header detection. It makes no claim about
/// whether that line is a header; it only shows what the file actually holds.
pub fn first_line_fields(path: &str, delimiters: &[char]) -> anyhow::Result<Vec<Box<str>>> {
    let mut first = String::new();
    open_buf_reader(path)?.read_line(&mut first)?;
    Ok(first
        .trim_end_matches(['\n', '\r'])
        .split(delimiters)
        .map(|f| unquote_field(f).to_string().into_boxed_str())
        .filter(|f| !f.is_empty())
        .collect())
}

/// Is the first line a header? Decided by type, not by guessing intent: a line
/// whose fields after column 0 fail to parse as numbers cannot be a data row.
///
/// Column 0 is skipped because a row-name column is non-numeric either way.
/// Fields are unquoted with the tokenizer's own rule first, so a fully quoted
/// numeric field (`"100.5"`) reads as numeric and a quoted headerless file does
/// not lose its first data row to a phantom header. Reads through gzip.
///
/// `Some(0)` when the first line is a header; `None` when every field after
/// column 0 is numeric, or the file cannot be read.
pub fn detect_header_row_numeric(file_path: &str, delimiters: &[char]) -> Option<usize> {
    let mut first = String::new();
    open_buf_reader(file_path)
        .ok()?
        .read_line(&mut first)
        .ok()?;
    let fields: Vec<&str> = first
        .trim_end_matches(['\n', '\r'])
        .split(delimiters)
        .map(unquote_field)
        .collect();
    let any_non_numeric_after_col0 = fields
        .iter()
        .skip(1)
        .any(|t| !t.is_empty() && !is_numeric_or_missing(t));
    // Both outcomes are logged WITH the fields, so a header swallowed as data
    // (all-numeric sample IDs) or a data row taken as a header shows up in the
    // log next to the evidence, not only as a wrong row count later.
    let preview = fields
        .iter()
        .take(6)
        .copied()
        .collect::<Vec<_>>()
        .join(", ");
    if any_non_numeric_after_col0 {
        log::info!("{file_path}: first line treated as a header (non-numeric): [{preview}]");
        Some(0)
    } else {
        log::info!(
            "{file_path}: first line treated as data (all numeric after column 0): [{preview}]"
        );
        None
    }
}

/// A field a count row may legitimately hold: a number, or one of the missing
/// value spellings R and friends write (`NA`, `N/A`; `NaN` already parses).
fn is_numeric_or_missing(t: &str) -> bool {
    t.parse::<f64>().is_ok() || matches!(t, "NA" | "N/A" | "na" | "n/a")
}

///
/// Open a file for writing, and return a buffered writer
/// * `output_file` - file name--either gzipped or not
pub fn open_buf_writer(output_file: &str) -> anyhow::Result<Box<dyn std::io::Write>> {
    // we can simply override with stdout
    if output_file.eq_ignore_ascii_case("stdout") {
        return Ok(Box::new(std::io::BufWriter::new(std::io::stdout())));
    }

    if output_file.eq_ignore_ascii_case("stderr") {
        return Ok(Box::new(std::io::BufWriter::new(std::io::stderr())));
    }

    // take a look at the extension
    let output_file = Path::new(output_file);
    let ext = output_file.extension().and_then(|x| x.to_str());

    match ext {
        Some("gz") => {
            let output_file = File::create(output_file)?;
            let encoder =
                flate2::write::GzEncoder::new(output_file, flate2::Compression::default());
            Ok(Box::new(BufWriter::new(encoder)))
        }
        _ => {
            let output_file = File::create(output_file)?;
            Ok(Box::new(BufWriter::new(output_file)))
        }
    }
}

///
/// Create a directory if needed
/// * `file` - file name
///
pub fn mkdir(file: &str) -> anyhow::Result<()> {
    let path = Path::new(file);
    std::fs::create_dir_all(path)?;
    Ok(())
}

/// Ensure the parent directory of an output path/prefix exists.
///
/// CLI subcommands typically take an `--out` value used as a *file
/// prefix* (e.g. `results/run1` → writes `results/run1.foo.parquet`).
/// Creates the parent directory tree (`results/`) if missing, but
/// never creates a directory named after the prefix itself
/// (`results/run1/`).
///
/// `Path::parent()` returns `Some("")` for a bare filename like
/// `"run1"`; `create_dir_all` would error on that, so skip it.
pub fn mkdir_parent(path: &str) -> anyhow::Result<()> {
    if let Some(parent) = Path::new(path).parent() {
        if !parent.as_os_str().is_empty() {
            std::fs::create_dir_all(parent)?;
        }
    }
    Ok(())
}

pub trait PathOsToStr {
    #[allow(clippy::wrong_self_convention)]
    fn into_boxed_str(&self) -> Box<str>;
}

impl PathOsToStr for Path {
    fn into_boxed_str(&self) -> Box<str> {
        self.to_str()
            .expect("failed to convert to string")
            .to_string()
            .into_boxed_str()
    }
}

impl PathOsToStr for OsStr {
    fn into_boxed_str(&self) -> Box<str> {
        self.to_str()
            .expect("failed to convert to string")
            .to_string()
            .into_boxed_str()
    }
}

/// more general file/directory copy function
pub fn recursive_copy(src_path: &str, dst_path: &str) -> anyhow::Result<()> {
    let src = Path::new(src_path);
    let dst = Path::new(dst_path);

    if src.is_dir() {
        mkdir(dst_path)?;
        for entry in std::fs::read_dir(src)? {
            let entry = entry?;
            if let (Some(src_path), Some(dst_path)) =
                (entry.path().to_str(), dst.join(entry.file_name()).to_str())
            {
                let file_type = entry.file_type()?;
                if file_type.is_dir() {
                    recursive_copy(src_path, dst_path)?;
                } else if file_type.is_file() {
                    std::fs::copy(src_path, dst_path)?;
                }
            }
        }
    } else if src.is_file() {
        if let Some(dir) = dirname(dst_path).as_deref() {
            mkdir(dir)?;
        }
        std::fs::copy(src, dst)?;
    } else if src.is_symlink() {
        if let Ok(abs_src) = std::fs::read_link(src) {
            if let Some(abs_src_path) = abs_src.to_str() {
                recursive_copy(abs_src_path, dst_path)?;
            }
        }
    }

    Ok(())
}

/// Unzip `zip_path` into the `extract_path` If `extract_path` is
/// `None`, just use a current directory.
/// * Returns `extract_path`
pub fn unzip_dir(zip_path: &str, extract_path: Option<&str>) -> anyhow::Result<Box<str>> {
    let zip_file = std::fs::File::open(zip_path)?;
    let mut archive = zip::ZipArchive::new(zip_file)?;

    let extract_path = extract_path
        .map(std::path::PathBuf::from)
        .unwrap_or(std::env::current_dir()?);

    for i in 0..archive.len() {
        let mut file = archive.by_index(i)?;
        let out_path = extract_path.join(file.name());
        // println!("{}", out_path.to_str().unwrap());
        if file.is_dir() {
            std::fs::create_dir_all(&out_path)?;
        } else {
            if let Some(parent) = out_path.parent() {
                std::fs::create_dir_all(parent)?;
            }
            let mut outfile = std::fs::File::create(&out_path)?;
            std::io::copy(&mut file, &mut outfile)?;
        }
    }

    Ok(extract_path.into_boxed_str())
}

/// Zip a directory into a zip archive using `Stored` compression
/// (zarr chunks are already zstd-compressed internally). In-zip entries are
/// prefixed with the source directory's basename.
pub fn zip_dir(source_dir: &str, zip_path: &str) -> anyhow::Result<()> {
    zip_dir_as(source_dir, zip_path, None)
}

/// Like [`zip_dir`] but lets the caller override the in-zip root name, so the
/// archive can use a different prefix than the source directory's basename
/// without renaming or moving the source on disk.
pub fn zip_dir_as(
    source_dir: &str,
    zip_path: &str,
    entry_root: Option<&str>,
) -> anyhow::Result<()> {
    use std::io::Write;
    use zip::write::SimpleFileOptions;
    use zip::ZipWriter;

    let file = std::fs::File::create(zip_path)?;
    let mut zip = ZipWriter::new(file);
    let options = SimpleFileOptions::default().compression_method(zip::CompressionMethod::Stored);
    let source = Path::new(source_dir);

    fn collect_entries(dir: &Path, out: &mut Vec<std::path::PathBuf>) -> std::io::Result<()> {
        for entry in std::fs::read_dir(dir)? {
            let entry = entry?;
            let path = entry.path();
            out.push(path.clone());
            // Use entry.file_type() (does not follow symlinks) to avoid
            // infinite recursion on circular symlinks.
            if entry.file_type()?.is_dir() {
                collect_entries(&path, out)?;
            }
        }
        Ok(())
    }

    let mut entries = vec![];
    collect_entries(source, &mut entries)?;
    entries.sort();

    let root_name =
        entry_root.unwrap_or_else(|| source.file_name().and_then(|s| s.to_str()).unwrap_or(""));

    for path in &entries {
        let rel_under_source = path.strip_prefix(source).unwrap_or(path);
        let rel = if root_name.is_empty() {
            rel_under_source.to_path_buf()
        } else {
            Path::new(root_name).join(rel_under_source)
        };
        if path.symlink_metadata()?.is_dir() {
            zip.add_directory(format!("{}/", rel.display()), options)?;
        } else {
            zip.start_file(rel.display().to_string(), options)?;
            let data = std::fs::read(path)?;
            zip.write_all(&data)?;
        }
    }

    zip.finish()?;
    Ok(())
}

/// just get the directory (parent) name
pub fn dirname(file_path: &str) -> Option<Box<str>> {
    Path::new(file_path).parent().map(|x| x.into_boxed_str())
}

///
/// Take the parent directory, basename, and extension of a file
/// * `file` - file name
///
pub fn dir_base_ext(file_path: &str) -> anyhow::Result<(Box<str>, Box<str>, Box<str>)> {
    let path = Path::new(file_path);

    let dir = path
        .parent()
        .map_or(".".to_string().into_boxed_str(), |x| x.into_boxed_str());

    let ext = path
        .extension()
        .map_or("".to_string().into_boxed_str(), |x| x.into_boxed_str());

    let base = path
        .file_stem()
        .and_then(|x| x.to_str())
        .map(|x| strip_data_ext(x).to_string().into_boxed_str())
        .ok_or(anyhow::anyhow!("failed to find base here: {}", file_path))?;

    Ok((dir, base, ext))
}

///
/// Take the basename of a file
/// * `file` - file name
///
pub fn basename(file: &str) -> anyhow::Result<Box<str>> {
    let path = Path::new(file);
    if let Some(base) = path.file_stem().and_then(|s| s.to_str()) {
        Ok(strip_data_ext(base).to_string().into_boxed_str())
    } else {
        Err(anyhow::anyhow!("no file stem"))
    }
}

/// Strip a leftover sparse-data extension (`.zarr`, `.h5`, `.h5ad`) that
/// `Path::file_stem` leaves behind on doubly-extended paths like
/// `sample.zarr.zip` (file_stem → `sample.zarr`).
fn strip_data_ext(stem: &str) -> &str {
    for sfx in [".zarr", ".h5ad", ".h5"] {
        if let Some(s) = stem.strip_suffix(sfx) {
            return s;
        }
    }
    stem
}

///
/// Take the extension of a file
/// * `file` - file name
///
pub fn file_ext(file: &str) -> anyhow::Result<Box<str>> {
    let path = Path::new(file);
    if let Some(ext) = path.extension() {
        Ok(ext.into_boxed_str())
    } else {
        Err(anyhow::anyhow!("failed to extract extension"))
    }
}

///
/// Create a temporary directory and suggest a file name
/// * `suffix` - suffix of the file name
///
pub fn create_temp_dir_file(suffix: &str) -> anyhow::Result<std::path::PathBuf> {
    let temp_dir = tempdir()?.path().to_path_buf();
    std::fs::create_dir_all(&temp_dir)?;
    let temp_file = tempfile::Builder::new()
        .suffix(suffix)
        .tempfile_in(temp_dir)?
        .path()
        .to_owned();

    Ok(temp_file)
}

///
/// Remove a file if it exists
/// * `file` - file name
///
pub fn remove_file(file: &str) -> anyhow::Result<()> {
    let path = Path::new(file);
    if path.exists() {
        if path.is_file() {
            std::fs::remove_file(path)?;
        } else {
            std::fs::remove_dir_all(path)?;
        }
    }
    Ok(())
}

///
/// Remove a file if it exists
/// * `files` - file name
///
pub fn remove_all_files(files: &Vec<Box<str>>) -> anyhow::Result<()> {
    for file in files {
        remove_file(file)?;
    }
    Ok(())
}

/// Extensions a data file's name is read through, stripped from the end of
/// the name repeatedly by [`file_stem`]; dots inside the name stay.
pub const DATA_FILE_EXTENSIONS: &[&str] = &[
    "gz", "bgz", "bz2", "zst", "tsv", "csv", "txt", "tab", "gaf", "gmt", "obo", "bed", "vcf",
    "parquet", "pq",
];

/// The file name minus its known extensions
/// (`a/b/goa_human.gaf.gz` → `goa_human`, `c2.cp.v1.symbols.gmt` → `c2.cp.v1.symbols`):
/// the name tools derive a relation or a source label from, so they agree.
pub fn file_stem(path: &str) -> String {
    let mut stem = std::path::Path::new(path)
        .file_name()
        .map(|s| s.to_string_lossy().into_owned())
        .unwrap_or_else(|| path.to_string());
    loop {
        let Some((base, ext)) = stem.rsplit_once('.') else {
            break;
        };
        if base.is_empty() || !DATA_FILE_EXTENSIONS.contains(&ext.to_lowercase().as_str()) {
            break;
        }
        stem.truncate(base.len());
    }
    stem
}