code wiki / (root) / nx_bam_test.nx

nx_bam_test.nx source

↩ module page · 306 lines · 12407 B

1// nx_bam_test.nx -- KAT for BAM 4-bit seq + cigar packing. 2// 3// Section A -- nx_bam_iupac_to_nibble: every spec letter + lowercase 4// + unknown -> N fallback. 5// Section B -- nx_bam_pack_seq: even-length, odd-length, mixed-case, 6// capacity overflow. 7// Section C -- nx_bam_pack_cigar: single-op, multi-op, validate 8// little-endian byte order, reject invalid op codes. 9// 10// expect_exit: 0 11// 12// license_tier: ORIGINAL 13 14import "nx_syscalls.nx" 15import "nx_const.nx" 16import "nx_bam.nx" 17 18func main() -> i64 { 19 20 // ============================================================ 21 // Section A -- IUPAC ASCII -> 4-bit code lookup. 22 // Spec values per SAMv1 §4.2 Table 1. 23 // ============================================================ 24 25 if nx_bam_iupac_to_nibble(0x3D) != 0 { return 1 } // = 26 if nx_bam_iupac_to_nibble(0x41) != 1 { return 2 } // A 27 if nx_bam_iupac_to_nibble(0x43) != 2 { return 3 } // C 28 if nx_bam_iupac_to_nibble(0x4D) != 3 { return 4 } // M 29 if nx_bam_iupac_to_nibble(0x47) != 4 { return 5 } // G 30 if nx_bam_iupac_to_nibble(0x52) != 5 { return 6 } // R 31 if nx_bam_iupac_to_nibble(0x53) != 6 { return 7 } // S 32 if nx_bam_iupac_to_nibble(0x56) != 7 { return 8 } // V 33 if nx_bam_iupac_to_nibble(0x54) != 8 { return 9 } // T 34 if nx_bam_iupac_to_nibble(0x57) != 9 { return 10 } // W 35 if nx_bam_iupac_to_nibble(0x59) != 10 { return 11 } // Y 36 if nx_bam_iupac_to_nibble(0x48) != 11 { return 12 } // H 37 if nx_bam_iupac_to_nibble(0x4B) != 12 { return 13 } // K 38 if nx_bam_iupac_to_nibble(0x44) != 13 { return 14 } // D 39 if nx_bam_iupac_to_nibble(0x42) != 14 { return 15 } // B 40 if nx_bam_iupac_to_nibble(0x4E) != 15 { return 16 } // N 41 42 // Lowercase folds. 43 if nx_bam_iupac_to_nibble(0x61) != 1 { return 17 } // a 44 if nx_bam_iupac_to_nibble(0x67) != 4 { return 18 } // g 45 if nx_bam_iupac_to_nibble(0x74) != 8 { return 19 } // t 46 if nx_bam_iupac_to_nibble(0x6E) != 15 { return 20 } // n 47 48 // U/u -> T per RNA-collapse convention. 49 if nx_bam_iupac_to_nibble(0x55) != 8 { return 21 } // U 50 if nx_bam_iupac_to_nibble(0x75) != 8 { return 22 } // u 51 52 // Unknown -> N fallback. 53 if nx_bam_iupac_to_nibble(0x21) != 15 { return 23 } // '!' 54 if nx_bam_iupac_to_nibble(0x5A) != 15 { return 24 } // 'Z' 55 56 // ============================================================ 57 // Section B -- nx_bam_pack_seq packing semantics. 58 // ============================================================ 59 60 // Canonical even-length: "ACGT" = [A=1, C=2, G=4, T=8] 61 // Packed: byte0 = (1<<4)|2 = 0x12, byte1 = (4<<4)|8 = 0x48 62 let seq1: *u8 = sys_mmap(8) 63 seq1[0]=0x41; seq1[1]=0x43; seq1[2]=0x47; seq1[3]=0x54 // ACGT 64 let out1: *u8 = sys_mmap(8) 65 let n1: i64 = nx_bam_pack_seq(seq1, 4, out1, 8) 66 if n1 != 2 { return 30 } 67 if (out1[0] & 0xff) != 0x12 { return 31 } 68 if (out1[1] & 0xff) != 0x48 { return 32 } 69 70 // Odd-length: "ACG" = [A=1, C=2, G=4, _=0] 71 // Packed: byte0 = (1<<4)|2 = 0x12, byte1 = (4<<4)|0 = 0x40 72 seq1[0]=0x41; seq1[1]=0x43; seq1[2]=0x47 // ACG 73 let n2: i64 = nx_bam_pack_seq(seq1, 3, out1, 8) 74 if n2 != 2 { return 40 } 75 if (out1[0] & 0xff) != 0x12 { return 41 } 76 if (out1[1] & 0xff) != 0x40 { return 42 } 77 78 // Mixed case + IUPAC ambiguity: "aNcT" = [a=1, N=15, c=2, T=8] 79 // Packed: byte0 = (1<<4)|15 = 0x1F, byte1 = (2<<4)|8 = 0x28 80 seq1[0]=0x61; seq1[1]=0x4E; seq1[2]=0x63; seq1[3]=0x54 // aNcT 81 let n3: i64 = nx_bam_pack_seq(seq1, 4, out1, 8) 82 if n3 != 2 { return 50 } 83 if (out1[0] & 0xff) != 0x1F { return 51 } 84 if (out1[1] & 0xff) != 0x28 { return 52 } 85 86 // Length-1: "G" = [G=4, _=0] -> 0x40 87 seq1[0]=0x47 88 let n4: i64 = nx_bam_pack_seq(seq1, 1, out1, 8) 89 if n4 != 1 { return 60 } 90 if (out1[0] & 0xff) != 0x40 { return 61 } 91 92 // Length-0 edge: 0 bytes written. 93 let n5: i64 = nx_bam_pack_seq(seq1, 0, out1, 8) 94 if n5 != 0 { return 70 } 95 96 // Capacity overflow: 4 bases -> 2 bytes, max_out=1. 97 seq1[0]=0x41; seq1[1]=0x43; seq1[2]=0x47; seq1[3]=0x54 98 if nx_bam_pack_seq(seq1, 4, out1, 1) != -1 { return 80 } 99 100 // Negative seq_len rejected. 101 if nx_bam_pack_seq(seq1, -1, out1, 8) != -1 { return 81 } 102 103 // ============================================================ 104 // Section C -- nx_bam_pack_cigar packing semantics. 105 // ============================================================ 106 107 // Single op: 100M = (100<<4)|0 = 1600 = 0x640. LE bytes: 40 06 00 00 108 let ops1: *i64 = sys_mmap(8 * 4) as *i64 109 let lens1: *i64 = sys_mmap(8 * 4) as *i64 110 ops1[0] = NX_CIGAR_M 111 lens1[0] = 100 112 let cig_out: *u8 = sys_mmap(64) 113 let cn1: i64 = nx_bam_pack_cigar(ops1, lens1, 1, cig_out, 64) 114 if cn1 != 4 { return 100 } 115 if (cig_out[0] & 0xff) != 0x40 { return 101 } 116 if (cig_out[1] & 0xff) != 0x06 { return 102 } 117 if (cig_out[2] & 0xff) != 0x00 { return 103 } 118 if (cig_out[3] & 0xff) != 0x00 { return 104 } 119 120 // Multi-op: 5M2I3D = [(5<<4)|0=80=0x50, (2<<4)|1=33=0x21, (3<<4)|2=50=0x32] 121 // LE: 50 00 00 00 21 00 00 00 32 00 00 00 122 ops1[0] = NX_CIGAR_M; lens1[0] = 5 123 ops1[1] = NX_CIGAR_I; lens1[1] = 2 124 ops1[2] = NX_CIGAR_D; lens1[2] = 3 125 let cn2: i64 = nx_bam_pack_cigar(ops1, lens1, 3, cig_out, 64) 126 if cn2 != 12 { return 110 } 127 if (cig_out[0] & 0xff) != 0x50 { return 111 } 128 if (cig_out[1] & 0xff) != 0x00 { return 112 } 129 if (cig_out[4] & 0xff) != 0x21 { return 113 } 130 if (cig_out[5] & 0xff) != 0x00 { return 114 } 131 if (cig_out[8] & 0xff) != 0x32 { return 115 } 132 if (cig_out[9] & 0xff) != 0x00 { return 116 } 133 if (cig_out[10] & 0xff) != 0x00 { return 117 } 134 if (cig_out[11] & 0xff) != 0x00 { return 118 } 135 136 // Length-0: 0 bytes written. 137 let cn3: i64 = nx_bam_pack_cigar(ops1, lens1, 0, cig_out, 64) 138 if cn3 != 0 { return 120 } 139 140 // Invalid op code rejected. 141 ops1[0] = 9; lens1[0] = 5 142 if nx_bam_pack_cigar(ops1, lens1, 1, cig_out, 64) != -1 { return 130 } 143 ops1[0] = -1 144 if nx_bam_pack_cigar(ops1, lens1, 1, cig_out, 64) != -1 { return 131 } 145 146 // Negative length rejected. 147 ops1[0] = NX_CIGAR_M; lens1[0] = -5 148 if nx_bam_pack_cigar(ops1, lens1, 1, cig_out, 64) != -1 { return 132 } 149 150 // Capacity overflow. 151 ops1[0] = NX_CIGAR_M; lens1[0] = 100 152 if nx_bam_pack_cigar(ops1, lens1, 1, cig_out, 3) != -1 { return 133 } 153 154 // Negative n_op rejected. 155 if nx_bam_pack_cigar(ops1, lens1, -1, cig_out, 64) != -1 { return 134 } 156 157 // ============================================================ 158 // Section D -- nx_bam_write_record full binary record. 159 // 160 // Canonical mapped read: 161 // read_name = "r1" (l_read_name=3 including NUL) 162 // refID=0, pos=99 (0-indexed; SAM POS 100), 163 // mapq=60, bin=4681 (UCSC bin for pos 99-99+3), 164 // n_cigar_op=1 (3M), flag=0, l_seq=3 (ACG), 165 // next_refID=-1, next_pos=-1, tlen=0, no_qual=0 166 // qual = [37, 37, 37] (raw Phred, not +33) 167 // 168 // Expected layout (LE bytes): 169 // block_size = total - 4 170 // total = 4 + 32 + 3 (name+NUL) + 4 (cigar) + 2 (seq) + 3 (qual) = 48 171 // block_size = 44 = 0x2C 172 // ============================================================ 173 174 let rname: *u8 = sys_mmap(8) 175 rname[0]=0x72; rname[1]=0x31 // "r1" 176 let ops_d: *i64 = sys_mmap(8) as *i64 177 let lens_d: *i64 = sys_mmap(8) as *i64 178 ops_d[0] = NX_CIGAR_M 179 lens_d[0] = 3 180 let seq_d: *u8 = sys_mmap(8) 181 seq_d[0]=0x41; seq_d[1]=0x43; seq_d[2]=0x47 // ACG 182 let qual_d: *u8 = sys_mmap(8) 183 qual_d[0]=37; qual_d[1]=37; qual_d[2]=37 184 185 let fields_d: *i64 = sys_mmap(8 * 12) as *i64 186 fields_d[NX_BAM_REC_REF_ID] = 0 187 fields_d[NX_BAM_REC_POS] = 99 188 fields_d[NX_BAM_REC_MAPQ] = 60 189 fields_d[NX_BAM_REC_BIN] = 4681 190 fields_d[NX_BAM_REC_FLAG] = 0 191 fields_d[NX_BAM_REC_NEXT_REF] = -1 192 fields_d[NX_BAM_REC_NEXT_POS] = -1 193 fields_d[NX_BAM_REC_TLEN] = 0 194 fields_d[NX_BAM_REC_NO_QUAL] = 0 195 fields_d[NX_BAM_REC_READ_NAME_LEN] = 2 196 fields_d[NX_BAM_REC_N_CIGAR_OP] = 1 197 fields_d[NX_BAM_REC_L_SEQ] = 3 198 199 let rec_out: *u8 = sys_mmap(128) 200 let rn: i64 = nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 201 fields_d, rec_out, 128) 202 if rn != 48 { return rn + 200 } 203 204 // block_size = 44 LE = 2C 00 00 00 205 if (rec_out[0] & 0xff) != 0x2C { return 140 } 206 if (rec_out[1] & 0xff) != 0x00 { return 141 } 207 if (rec_out[2] & 0xff) != 0x00 { return 142 } 208 if (rec_out[3] & 0xff) != 0x00 { return 143 } 209 210 // refID = 0 211 if (rec_out[4] & 0xff) != 0x00 { return 144 } 212 if (rec_out[5] & 0xff) != 0x00 { return 145 } 213 if (rec_out[6] & 0xff) != 0x00 { return 146 } 214 if (rec_out[7] & 0xff) != 0x00 { return 147 } 215 216 // pos = 99 = 0x63 217 if (rec_out[8] & 0xff) != 0x63 { return 148 } 218 if (rec_out[9] & 0xff) != 0x00 { return 149 } 219 if (rec_out[10] & 0xff) != 0x00 { return 150 } 220 if (rec_out[11] & 0xff) != 0x00 { return 151 } 221 222 // l_read_name = 3 (incl. NUL) 223 if (rec_out[12] & 0xff) != 0x03 { return 152 } 224 // mapq = 60 = 0x3C 225 if (rec_out[13] & 0xff) != 0x3C { return 153 } 226 227 // bin = 4681 = 0x1249 → 49 12 LE 228 if (rec_out[14] & 0xff) != 0x49 { return 154 } 229 if (rec_out[15] & 0xff) != 0x12 { return 155 } 230 231 // n_cigar_op = 1 → 01 00 232 if (rec_out[16] & 0xff) != 0x01 { return 156 } 233 if (rec_out[17] & 0xff) != 0x00 { return 157 } 234 235 // flag = 0 236 if (rec_out[18] & 0xff) != 0x00 { return 158 } 237 if (rec_out[19] & 0xff) != 0x00 { return 159 } 238 239 // l_seq = 3 → 03 00 00 00 240 if (rec_out[20] & 0xff) != 0x03 { return 160 } 241 if (rec_out[21] & 0xff) != 0x00 { return 161 } 242 243 // next_refID = -1 → FF FF FF FF 244 if (rec_out[24] & 0xff) != 0xFF { return 162 } 245 if (rec_out[27] & 0xff) != 0xFF { return 163 } 246 247 // next_pos = -1 → FF FF FF FF 248 if (rec_out[28] & 0xff) != 0xFF { return 164 } 249 if (rec_out[31] & 0xff) != 0xFF { return 165 } 250 251 // tlen = 0 252 if (rec_out[32] & 0xff) != 0x00 { return 166 } 253 if (rec_out[35] & 0xff) != 0x00 { return 167 } 254 255 // read_name "r1\0" 256 if (rec_out[36] & 0xff) != 0x72 { return 168 } // r 257 if (rec_out[37] & 0xff) != 0x31 { return 169 } // 1 258 if (rec_out[38] & 0xff) != 0x00 { return 170 } // NUL 259 260 // cigar 3M = (3<<4)|0 = 0x30 → LE: 30 00 00 00 261 if (rec_out[39] & 0xff) != 0x30 { return 171 } 262 if (rec_out[40] & 0xff) != 0x00 { return 172 } 263 if (rec_out[41] & 0xff) != 0x00 { return 173 } 264 if (rec_out[42] & 0xff) != 0x00 { return 174 } 265 266 // seq "ACG" packed: (A=1<<4)|C=2 = 0x12, (G=4<<4)|0 = 0x40 267 if (rec_out[43] & 0xff) != 0x12 { return 175 } 268 if (rec_out[44] & 0xff) != 0x40 { return 176 } 269 270 // qual = [37, 37, 37] raw 271 if (rec_out[45] & 0xff) != 37 { return 177 } 272 if (rec_out[46] & 0xff) != 37 { return 178 } 273 if (rec_out[47] & 0xff) != 37 { return 179 } 274 275 // NO_QUAL path: rewrite with NO_QUAL=1; qual section becomes all 0xFF. 276 fields_d[NX_BAM_REC_NO_QUAL] = 1 277 let rn_nq: i64 = nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 278 fields_d, rec_out, 128) 279 if rn_nq != 48 { return rn_nq + 250 } 280 if (rec_out[45] & 0xff) != 0xFF { return 180 } 281 if (rec_out[46] & 0xff) != 0xFF { return 181 } 282 if (rec_out[47] & 0xff) != 0xFF { return 182 } 283 284 // Reject paths. 285 fields_d[NX_BAM_REC_NO_QUAL] = 0 286 fields_d[NX_BAM_REC_READ_NAME_LEN] = 0 287 if nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 288 fields_d, rec_out, 128) != -1 { return 190 } 289 fields_d[NX_BAM_REC_READ_NAME_LEN] = 2 290 291 fields_d[NX_BAM_REC_MAPQ] = 256 292 if nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 293 fields_d, rec_out, 128) != -1 { return 191 } 294 fields_d[NX_BAM_REC_MAPQ] = 60 295 296 fields_d[NX_BAM_REC_BIN] = -1 297 if nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 298 fields_d, rec_out, 128) != -1 { return 192 } 299 fields_d[NX_BAM_REC_BIN] = 4681 300 301 // Capacity overflow. 302 if nx_bam_write_record(rname, ops_d, lens_d, seq_d, qual_d, 303 fields_d, rec_out, 20) != -1 { return 193 } 304 305 return 0 306}