nx_alignment_test.nx source

1// nx_alignment_test.nx -- KAT for AlignmentRecord + mapq + revcomp remap. 2// 3// expect_exit: 0 4// 5// license_tier: ORIGINAL 6 7import "nx_syscalls.nx" 8import "nx_const.nx" 9import "nx_alignment.nx" 10 11func main() -> i64 { 12 13 // ============================================================ 14 // Section A -- AlignmentRecord constructor + field round-trip. 15 // ============================================================ 16 17 let fields: *i64 = sys_mmap(128) as *i64 18 fields[0] = 7 // query_id 19 fields[1] = 3 // ref_id 20 fields[2] = 100 // q_start 21 fields[3] = 150 // q_end 22 fields[4] = 50000 // r_start 23 fields[5] = 50050 // r_end 24 fields[6] = NX_STRAND_FWD 25 fields[7] = 96 // score 26 fields[8] = 45 // mapq 27 fields[9] = 12 // n_seeds 28 29 let rec: *AlignmentRecord = nx_alignment_new(fields) 30 if rec.query_id != 7 { return 1 } 31 if rec.ref_id != 3 { return 2 } 32 if rec.q_start != 100 { return 3 } 33 if rec.q_end != 150 { return 4 } 34 if rec.r_start != 50000 { return 5 } 35 if rec.r_end != 50050 { return 6 } 36 if rec.strand != NX_STRAND_FWD { return 7 } 37 if rec.score != 96 { return 8 } 38 if rec.mapq != 45 { return 9 } 39 if rec.n_seeds != 12 { return 10 } 40 41 // Reverse-strand variant. 42 fields[6] = NX_STRAND_REV 43 let rec_rev: *AlignmentRecord = nx_alignment_new(fields) 44 if rec_rev.strand != NX_STRAND_REV { return 11 } 45 46 // ============================================================ 47 // Section B -- mapq formula edge cases. 48 // ============================================================ 49 50 // Uniquely placed (no secondary): mapq = 60. 51 if nx_alignment_compute_mapq(100, 0) != 60 { return 20 } 52 53 // secondary == primary -> cannot distinguish, mapq = 0. 54 if nx_alignment_compute_mapq(100, 100) != 0 { return 21 } 55 56 // secondary == primary/2 -> 60 * (100-50)/100 = 30. 57 if nx_alignment_compute_mapq(100, 50) != 30 { return 22 } 58 59 // secondary < primary by tiny amount (high confidence boundary). 60 // 60 * (100-99)/100 = 0 (integer truncation). 61 if nx_alignment_compute_mapq(100, 99) != 0 { return 23 } 62 63 // secondary > primary -> impossible-but-defensive, clamps to 0. 64 if nx_alignment_compute_mapq(50, 100) != 0 { return 24 } 65 66 // primary == 0 -> no alignment -> 0. 67 if nx_alignment_compute_mapq(0, 0) != 0 { return 25 } 68 69 // primary < 0 -> defensive 0. 70 if nx_alignment_compute_mapq(-5, 0) != 0 { return 26 } 71 72 // High-quality typical: primary=200 secondary=50 -> 60*150/200 = 45. 73 if nx_alignment_compute_mapq(200, 50) != 45 { return 27 } 74 75 // Mid-quality: primary=60 secondary=20 -> 60*40/60 = 40. 76 if nx_alignment_compute_mapq(60, 20) != 40 { return 28 } 77 78 // Verify cap: large primary, small secondary still clamps at 60. 79 if nx_alignment_compute_mapq(10000, 0) != 60 { return 29 } 80 81 // ============================================================ 82 // Section C -- revcomp coordinate remap. 83 // ============================================================ 84 85 let out_s: *i64 = sys_mmap(16) as *i64 86 let out_e: *i64 = sys_mmap(16) as *i64 87 88 // q_len=10, rc-coords [2..6) -> orig-coords [10-6, 10-2) = [4..8) 89 if nx_remap_revcomp_coords(10, 2, 6, out_s, out_e) != 0 { return 40 } 90 if out_s[0] != 4 { return 41 } 91 if out_e[0] != 8 { return 42 } 92 93 // Involution: remap(remap(x)) == x. 94 // Apply remap again to [4, 8) with q_len=10 -> [10-8, 10-4) = [2, 6). Same as start. 95 let inv_s: *i64 = sys_mmap(16) as *i64 96 let inv_e: *i64 = sys_mmap(16) as *i64 97 if nx_remap_revcomp_coords(10, 4, 8, inv_s, inv_e) != 0 { return 43 } 98 if inv_s[0] != 2 { return 44 } 99 if inv_e[0] != 6 { return 45 } 100 101 // Edge: full range [0..q_len) maps to [0..q_len). 102 if nx_remap_revcomp_coords(10, 0, 10, out_s, out_e) != 0 { return 50 } 103 if out_s[0] != 0 { return 51 } 104 if out_e[0] != 10 { return 52 } 105 106 // Edge: single-base [3..4) with q_len=10 -> [6..7). 107 if nx_remap_revcomp_coords(10, 3, 4, out_s, out_e) != 0 { return 53 } 108 if out_s[0] != 6 { return 54 } 109 if out_e[0] != 7 { return 55 } 110 111 // Bad input rejection. 112 if nx_remap_revcomp_coords(0, 0, 5, out_s, out_e) != -1 { return 60 } // q_len=0 113 if nx_remap_revcomp_coords(10, -1, 5, out_s, out_e) != -1 { return 61 } // negative start 114 if nx_remap_revcomp_coords(10, 5, 3, out_s, out_e) != -1 { return 62 } // end <= start 115 if nx_remap_revcomp_coords(10, 5, 11, out_s, out_e) != -1 { return 63 } // end > q_len 116 117 return 0 118}