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}