| 154 | |
| 155 | template <class indexT> |
| 156 | sequence<indexT> suffix_array(sequence<uchar> const &ss) { |
| 157 | if (verbose) sa_timer.start(); |
| 158 | size_t n = ss.size(); |
| 159 | |
| 160 | // renumber characters densely |
| 161 | // start numbering at 1 leaving 0 to indicate end-of-string |
| 162 | size_t pad = 48; |
| 163 | sequence<indexT> flags(256, (indexT) 0); |
| 164 | parallel_for (0, n, [&] (size_t i) { |
| 165 | if (!flags[ss[i]]) flags[ss[i]] = 1;}, 1000); |
| 166 | auto add = [&] (indexT a, indexT b) {return a + b;}; |
| 167 | indexT m; |
| 168 | std::tie(flags, m) = scan(flags, make_monoid(add,(indexT) 1)); |
| 169 | |
| 170 | // pad the end of string with 0s |
| 171 | sequence<uchar> s(n + pad, [&] (size_t i) { |
| 172 | return (i < n) ? flags[ss[i]] : 0;}); |
| 173 | |
| 174 | if (verbose) std::cout << "distinct characters = " << m-1 << std::endl; |
| 175 | |
| 176 | // pack characters into 128-bit word, along with the location i |
| 177 | // 96 bits for characters, and 32 for location |
| 178 | double logm = log2((double) m); |
| 179 | indexT nchars = floor(96.0/logm); |
| 180 | |
| 181 | sequence<uint128> Cl(n, [&] (size_t i) { |
| 182 | uint128 r = s[i]; |
| 183 | for (indexT j=1; j < nchars; j++) r = r*m + s[i+j]; |
| 184 | return (r << 32) + i; |
| 185 | }); |
| 186 | sa_timer.next("copy into 128bit int"); |
| 187 | |
| 188 | // sort based on packed words |
| 189 | sample_sort_inplace(Cl.slice(), std::less<uint128>()); |
| 190 | sa_timer.next("sort"); |
| 191 | |
| 192 | // identify segments of equal values |
| 193 | sequence<indexT> ranks(n); |
| 194 | sequence<seg<indexT>> seg_outs(n); |
| 195 | sequence<ipair<indexT>> C = split_segment_top(seg_outs, ranks, Cl); |
| 196 | Cl.clear(); |
| 197 | sa_timer.next("split"); |
| 198 | |
| 199 | indexT offset = nchars; |
| 200 | uint round =0; |
| 201 | indexT nKeys = n; |
| 202 | |
| 203 | // offset is how many characters for each suffix have already been sorted |
| 204 | // each round doubles offset so there should be at most log n rounds |
| 205 | // The segments keep regions that have not yet been fully sorted |
| 206 | while (1) { |
| 207 | if (round++ > 40) { |
| 208 | cout << "Suffix Array: Too many rounds" << std::endl; |
| 209 | abort(); |
| 210 | } |
| 211 | |
| 212 | auto is_seg = [&] (seg<indexT> s) {return s.length > 1;}; |
| 213 | // only keep segments that are longer than 1 (otherwise already sorted) |
nothing calls this directly
no test coverage detected