psltools 0.0.1

work with .psl files in Rust
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

// Copyright (c) 2026 Alejandro Gonzales-Irribarren <alejandrxgzi@gmail.com>
// Distributed under the terms of the Apache License, Version 2.0.

//! Structural validation of a PSL record, mirroring UCSC `pslCheck`.
//!
//! Verifies coordinate bounds, list-length agreement, block monotonicity, and
//! that the block spans reconcile with the alignment span and base-insert
//! counts (applying `sizeMul` on the reference side).

use crate::model::psl::PslRecord;
use crate::ops::score::size_mul;

/// The result of [`check`]: a list of human-readable invariant violations.
#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
#[derive(Debug, Clone, Default)]
pub struct CheckReport {
    pub violations: Vec<String>,
}

impl CheckReport {
    /// `true` when no invariants were violated.
    pub fn is_ok(&self) -> bool {
        self.violations.is_empty()
    }

    fn fail(&mut self, msg: impl Into<String>) {
        self.violations.push(msg.into());
    }
}

/// Validates the structural invariants of a record.
// `Coord as u64` is a real widening for the default `u32` build; it is only a
// no-op under `bigcoords`, so the cast is conditionally necessary.
#[allow(clippy::unnecessary_cast)]
pub fn check<P: PslRecord>(p: &P) -> CheckReport {
    let mut report = CheckReport::default();

    let sizes = p.block_sizes();
    let q_starts = p.query_starts();
    let t_starts = p.reference_starts();

    // List lengths must agree with each other and with blockCount.
    let n = sizes.len();
    if q_starts.len() != n || t_starts.len() != n {
        report.fail(format!(
            "block list lengths disagree: blockSizes={}, qStarts={}, tStarts={}",
            n,
            q_starts.len(),
            t_starts.len()
        ));
        // Without equal lengths the remaining per-block checks are unreliable.
        return report;
    }
    if p.block_count() != n {
        report.fail(format!(
            "blockCount {} disagrees with list length {n}",
            p.block_count()
        ));
    }

    // Coordinate bounds.
    if p.query_start() > p.query_end() {
        report.fail(format!(
            "qStart ({}) > qEnd ({})",
            p.query_start(),
            p.query_end()
        ));
    }
    if p.query_end() > p.query_size() {
        report.fail(format!(
            "qEnd ({}) > qSize ({})",
            p.query_end(),
            p.query_size()
        ));
    }
    if p.reference_start() > p.reference_end() {
        report.fail(format!(
            "tStart ({}) > tEnd ({})",
            p.reference_start(),
            p.reference_end()
        ));
    }
    if p.reference_end() > p.reference_size() {
        report.fail(format!(
            "tEnd ({}) > tSize ({})",
            p.reference_end(),
            p.reference_size()
        ));
    }

    let mul = size_mul(p) as u64;

    // Block monotonicity (raw coordinates increase along the alignment).
    for i in 1..n {
        if (q_starts[i - 1] as u64) + (sizes[i - 1] as u64) > q_starts[i] as u64 {
            report.fail(format!("qStarts not monotonic at block {i}"));
        }
        if (t_starts[i - 1] as u64) + (sizes[i - 1] as u64) * mul > t_starts[i] as u64 {
            report.fail(format!("tStarts not monotonic at block {i}"));
        }
    }

    // Span consistency: span == sum(blockSizes [* sizeMul]) + baseInsert.
    let sum_sizes: u64 = sizes.iter().map(|&s| s as u64).sum();
    let query_span = (p.query_end() as u64).saturating_sub(p.query_start() as u64);
    let expected_query = sum_sizes + p.query_base_insert() as u64;
    if query_span != expected_query {
        report.fail(format!(
            "qEnd-qStart ({query_span}) != sum(blockSizes) + qBaseInsert ({expected_query})"
        ));
    }
    let reference_span = (p.reference_end() as u64).saturating_sub(p.reference_start() as u64);
    let expected_reference = sum_sizes * mul + p.reference_base_insert() as u64;
    if reference_span != expected_reference {
        report.fail(format!(
            "tEnd-tStart ({reference_span}) != sum(blockSizes)*sizeMul + tBaseInsert ({expected_reference})"
        ));
    }

    report
}