code wiki / (root) / nx_cigar.nx

nx_cigar.nx source

↩ module page · 225 lines · 8571 B

1// nx_cigar.nx -- CIGAR run-length compression + ASCII serialiser. 2// 3// license_tier: INDEPENDENT_REDERIVE 4// genealogy_id: international-research-sources/sam-spec-v1.4.6-section-1.4.6 5// 6// G1.4c of NISHI_GENOMICS_SUBSTRATE_ROADMAP.md. Closes the gap 7// from raw op-stream (nx_align_backtrace output) to SAM/BAM-spec 8// CIGAR string. Three primitives compose: 9// 10// cigar_compress : MMMIIDD -> [(M,3), (I,2), (D,2)] 11// cigar_serialize : [(M,3), (I,2), (D,2)] -> "3M2I2D" ASCII 12// nx_i64_to_ascii_digits : positive int -> ASCII digits (no float) 13// 14// The two-stage compress + serialise split (rather than a one-shot 15// ops -> ASCII string) is deliberate: downstream consumers may want 16// the (op, count) array WITHOUT serialisation (for variant-calling 17// pileup, for binary BAM writer where ops are length-prefixed 18// uint32 not ASCII), and pre-serialised strings are inefficient 19// to walk programmatically. 20// 21// SAM-spec ASCII mapping (from NX_CIGAR_M/I/D codes): 22// NX_CIGAR_M (0) -> 'M' 0x4D 23// NX_CIGAR_I (1) -> 'I' 0x49 24// NX_CIGAR_D (2) -> 'D' 0x44 25// NX_CIGAR_N (3) -> 'N' 0x4E 26// NX_CIGAR_S (4) -> 'S' 0x53 27// NX_CIGAR_H (5) -> 'H' 0x48 28// NX_CIGAR_P (6) -> 'P' 0x50 29// NX_CIGAR_EQ(7) -> '=' 0x3D 30// NX_CIGAR_X (8) -> 'X' 0x58 31// 32// What G1.4c does NOT do (deferred): 33// - BAM binary uint32 packing (op_len << 4 | op_code) -- nx_bam_writer (G6) 34// - Soft-clipping op insertion at alignment ends -- G1.4c.2 35// - = / X explicit match/mismatch distinction -- G1.4c.2 36// - CIGAR parsing (ASCII -> ops) for reading SAM input -- G1.4d 37// 38// API: 39// cigar_compress(ops, n_ops, out_codes, out_lens, max_out) -> i64 40// cigar_serialize(codes, lens, n_groups, out_ascii, max_ascii) -> i64 41// nx_i64_to_ascii_digits(n, out_buf) -> i64 42// nx_cigar_code_to_ascii(code) -> i64 43// 44// nx_safety_envelope: (schema: nishi-library/seeds/safety-critical-standards.toml) 45// intended_use: "Convert SW backtrace op-stream to compressed 46// CIGAR form and SAM-spec ASCII; the bridge 47// between alignment primitives and the SAM/BAM 48// interchange format" 49// sil_target: SIL2 50// asil_target: QM 51// dal_target: DAL C 52// iec_62304_class: B 53// evidence: [no_floating_point, deterministic, 54// bit_equal_reproducible, 55// composes_nx_const_NX_CIGAR_codes, 56// composes_nx_const_NX_ASCII_DIGIT_0, 57// multi_digit_run_length_KAT, 58// end_to_end_backtrace_to_string_KAT, 59// license_tier_INDEPENDENT_REDERIVE] 60// hazard_register: [bug-tape-cigar-compress-runs-cross-op-boundary, 61// bug-tape-cigar-serialize-leading-zero-digit, 62// bug-tape-int-to-ascii-reverses-wrong-direction, 63// bug-tape-cigar-ascii-mapping-uses-wrong-letter, 64// bug-tape-empty-ops-emits-junk-instead-of-zero] 65// residual_risk: "Run-length max is i64; SAM-spec stores uint32 66// per run. BAM writer (G6) must validate runs 67// fit in uint32 before binary encoding. No 68// check here because the runs originate from 69// op_count which fits in i64 by construction." 70// verdict: NOT_YET_EVALUATED 71 72import "nx_syscalls.nx" 73import "nx_const.nx" 74 75// Per SAM-spec consumption flags: 76// M (0): consumes query + ref 77// I (1): consumes query only 78// D (2): consumes ref only 79// N (3): consumes ref only (intron-style skip) 80// S (4): consumes query only (soft-clip, but bases stored) 81// H (5): consumes neither (hard-clip, bases discarded) 82// P (6): consumes neither (padding vs padded ref) 83// = (7): consumes query + ref (explicit match) 84// X (8): consumes query + ref (explicit mismatch) 85 86// Returns 1 if op consumes query bases, 0 otherwise. 87func nx_cigar_consumes_query(op: i64) -> i64 { 88 let c: i64 = op & 0xff 89 if c == NX_CIGAR_M { return 1 } 90 if c == NX_CIGAR_I { return 1 } 91 if c == NX_CIGAR_S { return 1 } 92 if c == NX_CIGAR_EQ { return 1 } 93 if c == NX_CIGAR_X { return 1 } 94 return 0 95} 96 97// Returns 1 if op consumes reference bases, 0 otherwise. 98func nx_cigar_consumes_ref(op: i64) -> i64 { 99 let c: i64 = op & 0xff 100 if c == NX_CIGAR_M { return 1 } 101 if c == NX_CIGAR_D { return 1 } 102 if c == NX_CIGAR_N { return 1 } 103 if c == NX_CIGAR_EQ { return 1 } 104 if c == NX_CIGAR_X { return 1 } 105 return 0 106} 107 108// Map a NX_CIGAR_* code to its SAM-spec ASCII letter. 109// Returns ASCII byte value, or '?' (0x3F) for out-of-range codes. 110func nx_cigar_code_to_ascii(code: i64) -> i64 { 111 let c: i64 = code & 0xff 112 if c == NX_CIGAR_M { return 0x4D } // M 113 if c == NX_CIGAR_I { return 0x49 } // I 114 if c == NX_CIGAR_D { return 0x44 } // D 115 if c == NX_CIGAR_N { return 0x4E } // N 116 if c == NX_CIGAR_S { return 0x53 } // S 117 if c == NX_CIGAR_H { return 0x48 } // H 118 if c == NX_CIGAR_P { return 0x50 } // P 119 if c == NX_CIGAR_EQ { return 0x3D } // = 120 if c == NX_CIGAR_X { return 0x58 } // X 121 return 0x3F // '?' fallback 122} 123 124// Convert a non-negative integer to ASCII decimal digits. 125// out_buf must have capacity for the digits (n digits where n = 1 126// for 0..9, 2 for 10..99, etc; bounded by 19 for i64 max). 127// Returns the digit count written. Negative inputs return 0 (refused). 128func nx_i64_to_ascii_digits(n: i64, out_buf: *u8) -> i64 { 129 if n < 0 { return 0 } 130 if n == 0 { 131 out_buf[0] = NX_ASCII_DIGIT_0 & 0xff 132 return 1 133 } 134 135 let tmp: *u8 = sys_mmap(32) 136 var v: i64 = n 137 var i: i64 = 0 138 while v > 0 { 139 let d: i64 = v % 10 140 tmp[i] = (NX_ASCII_DIGIT_0 + d) & 0xff 141 v = v / 10 142 i = i + 1 143 } 144 145 // tmp now holds digits in reverse order; copy reversed into out_buf. 146 var j: i64 = 0 147 while j < i { 148 out_buf[j] = tmp[i - 1 - j] 149 j = j + 1 150 } 151 return i 152} 153 154// Compress raw op stream into (code, length) run-length groups. 155// Walks ops left to right, starts a new group whenever the op changes. 156// out_codes[k] gets the op-code byte; out_lens[k] gets the count. 157// Returns number of groups (>= 0), or -1 if max_out is too small. 158func cigar_compress(ops: *u8, n_ops: i64, 159 out_codes: *u8, out_lens: *i64, 160 max_out: i64) -> i64 { 161 if n_ops < 0 { return -1 } 162 if n_ops == 0 { return 0 } 163 if max_out <= 0 { return -1 } 164 165 var group_count: i64 = 0 166 var cur_code: i64 = ops[0] & 0xff 167 var cur_len: i64 = 1 168 169 var i: i64 = 1 170 while i < n_ops { 171 let next: i64 = ops[i] & 0xff 172 if next == cur_code { 173 cur_len = cur_len + 1 174 } else { 175 // Flush current group. 176 if group_count >= max_out { return -1 } 177 out_codes[group_count] = cur_code & 0xff 178 out_lens[group_count] = cur_len 179 group_count = group_count + 1 180 cur_code = next 181 cur_len = 1 182 } 183 i = i + 1 184 } 185 186 // Flush trailing group. 187 if group_count >= max_out { return -1 } 188 out_codes[group_count] = cur_code & 0xff 189 out_lens[group_count] = cur_len 190 group_count = group_count + 1 191 return group_count 192} 193 194// Serialise compressed CIGAR groups to SAM-spec ASCII string. 195// Format: <len1><letter1><len2><letter2>... e.g., "3M2I2D". 196// out_ascii must have capacity for the full string (no terminator). 197// Returns the byte length written, or -1 if capacity exceeded. 198func cigar_serialize(codes: *u8, lens: *i64, n_groups: i64, 199 out_ascii: *u8, max_ascii: i64) -> i64 { 200 if n_groups < 0 { return -1 } 201 if n_groups == 0 { return 0 } 202 if max_ascii <= 0 { return -1 } 203 204 let digit_buf: *u8 = sys_mmap(32) 205 206 var written: i64 = 0 207 var g: i64 = 0 208 while g < n_groups { 209 let len: i64 = lens[g] 210 let n_digits: i64 = nx_i64_to_ascii_digits(len, digit_buf) 211 // Need n_digits + 1 (the letter) more bytes. 212 if written + n_digits + 1 > max_ascii { return -1 } 213 var d: i64 = 0 214 while d < n_digits { 215 out_ascii[written] = digit_buf[d] & 0xff 216 written = written + 1 217 d = d + 1 218 } 219 let letter: i64 = nx_cigar_code_to_ascii(codes[g] & 0xff) 220 out_ascii[written] = letter & 0xff 221 written = written + 1 222 g = g + 1 223 } 224 return written 225}