nx_sequence_fm_test.nx source
↩ module page · 168 lines · 5350 B
1// nx_sequence_fm_test.nx -- KAT for FM-index construction + backward search.
2//
3// Canonical "BANANA$" reference (Ferragina-Manzini 2000, textbook
4// across every IR / bioinformatics intro):
5//
6// Suffixes (sorted lexicographically, $=0x24 < A < B < N):
7// 6 $
8// 5 A$
9// 3 ANA$
10// 1 ANANA$
11// 0 BANANA$
12// 4 NA$
13// 2 NANA$
14// SA = [6, 5, 3, 1, 0, 4, 2]
15// BWT = [A, N, N, B, $, A, A] = "ANNB$AA"
16// C[$] = 0 ; C[A] = 1 ; C[B] = 4 ; C[N] = 5 ; C[anything > N] = 7
17// fm_count("ANA") = 2 (positions 1, 3)
18// fm_count("BAN") = 1 (position 0)
19// fm_count("NAB") = 0 (not in text)
20// fm_count("A") = 3
21//
22// Plus DNA "ACGT$" KAT and a not-in-text refusal.
23//
24// expect_exit: 0
25//
26// license_tier: ORIGINAL
27
28import "nx_syscalls.nx"
29import "nx_sequence_fm.nx"
30
31func main() -> i64 {
32
33 // ============================================================
34 // Section A -- "BANANA$" (length 7)
35 // ============================================================
36
37 let t: *u8 = sys_mmap(16)
38 t[0]=0x42; t[1]=0x41; t[2]=0x4E; t[3]=0x41; t[4]=0x4E; t[5]=0x41; t[6]=0x24 // B A N A N A $
39
40 if fm_check_sentinel(t, 7) != 0 { return 1 }
41
42 let sa: *i64 = sys_mmap(128) as *i64
43 if fm_build_sa(t, 7, sa) != 0 { return 2 }
44
45 // Expected SA = [6, 5, 3, 1, 0, 4, 2]
46 if sa[0] != 6 { return 10 }
47 if sa[1] != 5 { return 11 }
48 if sa[2] != 3 { return 12 }
49 if sa[3] != 1 { return 13 }
50 if sa[4] != 0 { return 14 }
51 if sa[5] != 4 { return 15 }
52 if sa[6] != 2 { return 16 }
53
54 let bwt: *u8 = sys_mmap(16)
55 if fm_build_bwt(t, 7, sa, bwt) != 0 { return 3 }
56
57 // Expected BWT = "ANNB$AA" = [A, N, N, B, $, A, A]
58 if (bwt[0] & 0xff) != 0x41 { return 20 }
59 if (bwt[1] & 0xff) != 0x4E { return 21 }
60 if (bwt[2] & 0xff) != 0x4E { return 22 }
61 if (bwt[3] & 0xff) != 0x42 { return 23 }
62 if (bwt[4] & 0xff) != 0x24 { return 24 }
63 if (bwt[5] & 0xff) != 0x41 { return 25 }
64 if (bwt[6] & 0xff) != 0x41 { return 26 }
65
66 let c_arr: *i64 = sys_mmap(4096) as *i64
67 if fm_build_c(bwt, 7, c_arr) != 0 { return 4 }
68
69 // Expected: C[$]=0, C[A]=1, C[B]=4, C[N]=5, C[O]=7
70 if c_arr[0x24] != 0 { return 30 }
71 if c_arr[0x41] != 1 { return 31 }
72 if c_arr[0x42] != 4 { return 32 }
73 if c_arr[0x4E] != 5 { return 33 }
74 if c_arr[0x4F] != 7 { return 34 }
75
76 // Backward-search count tests.
77 let pat: *u8 = sys_mmap(16)
78
79 // "ANA" -- expect 2
80 pat[0]=0x41; pat[1]=0x4E; pat[2]=0x41
81 if fm_count(bwt, 7, c_arr, pat, 3) != 2 { return 40 }
82
83 // "BAN" -- expect 1
84 pat[0]=0x42; pat[1]=0x41; pat[2]=0x4E
85 if fm_count(bwt, 7, c_arr, pat, 3) != 1 { return 41 }
86
87 // "NAB" -- expect 0
88 pat[0]=0x4E; pat[1]=0x41; pat[2]=0x42
89 if fm_count(bwt, 7, c_arr, pat, 3) != 0 { return 42 }
90
91 // "A" -- expect 3
92 pat[0]=0x41
93 if fm_count(bwt, 7, c_arr, pat, 1) != 3 { return 43 }
94
95 // "N" -- expect 2
96 pat[0]=0x4E
97 if fm_count(bwt, 7, c_arr, pat, 1) != 2 { return 44 }
98
99 // "B" -- expect 1
100 pat[0]=0x42
101 if fm_count(bwt, 7, c_arr, pat, 1) != 1 { return 45 }
102
103 // "X" (absent) -- expect 0
104 pat[0]=0x58
105 if fm_count(bwt, 7, c_arr, pat, 1) != 0 { return 46 }
106
107 // "NANA" -- expect 1 (position 2)
108 pat[0]=0x4E; pat[1]=0x41; pat[2]=0x4E; pat[3]=0x41
109 if fm_count(bwt, 7, c_arr, pat, 4) != 1 { return 47 }
110
111 // ============================================================
112 // Section B -- DNA "ACGT$" (length 5)
113 // SA = [4, 0, 1, 2, 3]
114 // BWT = [T, $, A, C, G] = "T$ACG"
115 // C[$]=0, C[A]=1, C[C]=2, C[G]=3, C[T]=4, C[U]=5
116 // ============================================================
117
118 let t2: *u8 = sys_mmap(8)
119 t2[0]=0x41; t2[1]=0x43; t2[2]=0x47; t2[3]=0x54; t2[4]=0x24
120
121 if fm_check_sentinel(t2, 5) != 0 { return 50 }
122
123 let sa2: *i64 = sys_mmap(64) as *i64
124 let bwt2: *u8 = sys_mmap(8)
125 let c2: *i64 = sys_mmap(4096) as *i64
126
127 fm_build_sa(t2, 5, sa2)
128 fm_build_bwt(t2, 5, sa2, bwt2)
129 fm_build_c(bwt2, 5, c2)
130
131 if sa2[0] != 4 { return 51 }
132 if sa2[1] != 0 { return 52 }
133 if sa2[2] != 1 { return 53 }
134 if sa2[3] != 2 { return 54 }
135 if sa2[4] != 3 { return 55 }
136
137 if (bwt2[0] & 0xff) != 0x54 { return 60 }
138 if (bwt2[1] & 0xff) != 0x24 { return 61 }
139 if (bwt2[2] & 0xff) != 0x41 { return 62 }
140 if (bwt2[3] & 0xff) != 0x43 { return 63 }
141 if (bwt2[4] & 0xff) != 0x47 { return 64 }
142
143 // "CG" -- expect 1 (position 1)
144 pat[0]=0x43; pat[1]=0x47
145 if fm_count(bwt2, 5, c2, pat, 2) != 1 { return 70 }
146
147 // "ACGT" -- expect 1 (position 0)
148 pat[0]=0x41; pat[1]=0x43; pat[2]=0x47; pat[3]=0x54
149 if fm_count(bwt2, 5, c2, pat, 4) != 1 { return 71 }
150
151 // "GA" (absent) -- expect 0
152 pat[0]=0x47; pat[1]=0x41
153 if fm_count(bwt2, 5, c2, pat, 2) != 0 { return 72 }
154
155 // ============================================================
156 // Section C -- sentinel-invariant violations are caught.
157 // ============================================================
158
159 let bad: *u8 = sys_mmap(8)
160 bad[0]=0x41; bad[1]=0x43; bad[2]=0x47; bad[3]=0x54 // no $
161 if fm_check_sentinel(bad, 4) != -1 { return 80 }
162
163 let bad2: *u8 = sys_mmap(8)
164 bad2[0]=0x41; bad2[1]=0x24; bad2[2]=0x47; bad2[3]=0x24 // $ in middle
165 if fm_check_sentinel(bad2, 4) != -1 { return 81 }
166
167 return 0
168}