wincfg: Add swapfile.sys
[wimlib] / src / divsufsort / divsufsort.c
1 /*
2  * divsufsort.c for libdivsufsort
3  * Copyright (c) 2003-2008 Yuta Mori All Rights Reserved.
4  *
5  * Permission is hereby granted, free of charge, to any person
6  * obtaining a copy of this software and associated documentation
7  * files (the "Software"), to deal in the Software without
8  * restriction, including without limitation the rights to use,
9  * copy, modify, merge, publish, distribute, sublicense, and/or sell
10  * copies of the Software, and to permit persons to whom the
11  * Software is furnished to do so, subject to the following
12  * conditions:
13  *
14  * The above copyright notice and this permission notice shall be
15  * included in all copies or substantial portions of the Software.
16  *
17  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
18  * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
19  * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
20  * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
21  * HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
22  * WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
23  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
24  * OTHER DEALINGS IN THE SOFTWARE.
25  */
26
27 #include "divsufsort_private.h"
28 #ifdef _OPENMP
29 # include <omp.h>
30 #endif
31
32 /*- Private Functions -*/
33
34 /* Sorts suffixes of type B*. */
35 static
36 saidx_t
37 sort_typeBstar(const sauchar_t *T, saidx_t *SA,
38                saidx_t *bucket_A, saidx_t *bucket_B,
39                saidx_t n) {
40   saidx_t *PAb, *ISAb, *buf;
41   saidx_t i, j, k, t, m, bufsize;
42   saint_t c0, c1;
43
44   /* Initialize bucket arrays. */
45   for(i = 0; i < BUCKET_A_SIZE; ++i) { bucket_A[i] = 0; }
46   for(uint32_t i = 0; i < BUCKET_B_SIZE; ++i) { bucket_B[i] = 0; }
47
48   /* Count the number of occurrences of the first one or two characters of each
49      type A, B and B* suffix. Moreover, store the beginning position of all
50      type B* suffixes into the array SA. */
51
52   for(i = n - 1, m = n, c0 = T[n - 1]; 0 <= i;) {
53     /* type A suffix. */
54     do { ++BUCKET_A(c1 = c0); } while((0 <= --i) && ((c0 = T[i]) >= c1));
55     if(0 <= i) {
56       /* type B* suffix. */
57       ++BUCKET_BSTAR(c0, c1);
58       SA[--m] = i;
59       /* type B suffix. */
60       for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) <= c1); --i, c1 = c0) {
61         ++BUCKET_B(c0, c1);
62       }
63     }
64   }
65   m = n - m;
66 /*
67 note:
68   A type B* suffix is lexicographically smaller than a type B suffix that
69   begins with the same first two characters.
70 */
71
72   /* Calculate the index of start/end point of each bucket. */
73   for(c0 = 0, i = 0, j = 0; c0 < ALPHABET_SIZE; ++c0) {
74     t = i + BUCKET_A(c0);
75     BUCKET_A(c0) = i + j; /* start point */
76     i = t + BUCKET_B(c0, c0);
77     for(c1 = c0 + 1; c1 < ALPHABET_SIZE; ++c1) {
78       j += BUCKET_BSTAR(c0, c1);
79       BUCKET_BSTAR(c0, c1) = j; /* end point */
80       i += BUCKET_B(c0, c1);
81     }
82   }
83
84   if(0 < m) {
85     /* Sort the type B* suffixes by their first two characters. */
86     PAb = SA + n - m; ISAb = SA + m;
87     for(i = m - 2; 0 <= i; --i) {
88       t = PAb[i], c0 = T[t], c1 = T[t + 1];
89       SA[--BUCKET_BSTAR(c0, c1)] = i;
90     }
91     t = PAb[m - 1], c0 = T[t], c1 = T[t + 1];
92     SA[--BUCKET_BSTAR(c0, c1)] = m - 1;
93
94     /* Sort the type B* substrings using sssort. */
95     buf = SA + m, bufsize = n - (2 * m);
96     for(c0 = ALPHABET_SIZE - 2, j = m; 0 < j; --c0) {
97       for(c1 = ALPHABET_SIZE - 1; c0 < c1; j = i, --c1) {
98         i = BUCKET_BSTAR(c0, c1);
99         if(1 < (j - i)) {
100           sssort(T, PAb, SA + i, SA + j,
101                  buf, bufsize, 2, n, *(SA + i) == (m - 1));
102         }
103       }
104     }
105
106     /* Compute ranks of type B* substrings. */
107     for(i = m - 1; 0 <= i; --i) {
108       if(0 <= SA[i]) {
109         j = i;
110         do { ISAb[SA[i]] = i; } while((0 <= --i) && (0 <= SA[i]));
111         SA[i + 1] = i - j;
112         if(i <= 0) { break; }
113       }
114       j = i;
115       do { ISAb[SA[i] = ~SA[i]] = j; } while(SA[--i] < 0);
116       ISAb[SA[i]] = j;
117     }
118
119     /* Construct the inverse suffix array of type B* suffixes using trsort. */
120     trsort(ISAb, SA, m, 1);
121
122     /* Set the sorted order of tyoe B* suffixes. */
123     for(i = n - 1, j = m, c0 = T[n - 1]; 0 <= i;) {
124       for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) >= c1); --i, c1 = c0) { }
125       if(0 <= i) {
126         t = i;
127         for(--i, c1 = c0; (0 <= i) && ((c0 = T[i]) <= c1); --i, c1 = c0) { }
128         SA[ISAb[--j]] = ((t == 0) || (1 < (t - i))) ? t : ~t;
129       }
130     }
131
132     /* Calculate the index of start/end point of each bucket. */
133     BUCKET_B(ALPHABET_SIZE - 1, ALPHABET_SIZE - 1) = n; /* end point */
134     for(c0 = ALPHABET_SIZE - 2, k = m - 1; 0 <= c0; --c0) {
135       i = BUCKET_A(c0 + 1) - 1;
136       for(c1 = ALPHABET_SIZE - 1; c0 < c1; --c1) {
137         t = i - BUCKET_B(c0, c1);
138         BUCKET_B(c0, c1) = i; /* end point */
139
140         /* Move all type B* suffixes to the correct position. */
141         for(i = t, j = BUCKET_BSTAR(c0, c1);
142             j <= k;
143             --i, --k) { SA[i] = SA[k]; }
144       }
145       BUCKET_BSTAR(c0, c0 + 1) = i - BUCKET_B(c0, c0) + 1; /* start point */
146       BUCKET_B(c0, c0) = i; /* end point */
147     }
148   }
149
150   return m;
151 }
152
153 /* Constructs the suffix array by using the sorted order of type B* suffixes. */
154 static
155 void
156 construct_SA(const sauchar_t *T, saidx_t *SA,
157              saidx_t *bucket_A, saidx_t *bucket_B,
158              saidx_t n, saidx_t m) {
159   saidx_t *i, *j, *k;
160   saidx_t s;
161   saint_t c0, c1, c2;
162
163   if(0 < m) {
164     /* Construct the sorted order of type B suffixes by using
165        the sorted order of type B* suffixes. */
166     for(c1 = ALPHABET_SIZE - 2; 0 <= c1; --c1) {
167       /* Scan the suffix array from right to left. */
168       for(i = SA + BUCKET_BSTAR(c1, c1 + 1),
169           j = SA + BUCKET_A(c1 + 1) - 1, k = NULL, c2 = -1;
170           i <= j;
171           --j) {
172         if(0 < (s = *j)) {
173           assert(T[s] == c1);
174           assert(((s + 1) < n) && (T[s] <= T[s + 1]));
175           assert(T[s - 1] <= T[s]);
176           *j = ~s;
177           c0 = T[--s];
178           if((0 < s) && (T[s - 1] > c0)) { s = ~s; }
179           if(c0 != c2) {
180             if(0 <= c2) { BUCKET_B(c2, c1) = k - SA; }
181             k = SA + BUCKET_B(c2 = c0, c1);
182           }
183           assert(k < j);
184           *k-- = s;
185         } else {
186           assert(((s == 0) && (T[s] == c1)) || (s < 0));
187           *j = ~s;
188         }
189       }
190     }
191   }
192
193   /* Construct the suffix array by using
194      the sorted order of type B suffixes. */
195   k = SA + BUCKET_A(c2 = T[n - 1]);
196   *k++ = (T[n - 2] < c2) ? ~(n - 1) : (n - 1);
197   /* Scan the suffix array from left to right. */
198   for(i = SA, j = SA + n; i < j; ++i) {
199     if(0 < (s = *i)) {
200       assert(T[s - 1] >= T[s]);
201       c0 = T[--s];
202       if((s == 0) || (T[s - 1] < c0)) { s = ~s; }
203       if(c0 != c2) {
204         BUCKET_A(c2) = k - SA;
205         k = SA + BUCKET_A(c2 = c0);
206       }
207       assert(i < k);
208       *k++ = s;
209     } else {
210       assert(s < 0);
211       *i = ~s;
212     }
213   }
214 }
215
216 #if 0
217 /* Constructs the burrows-wheeler transformed string directly
218    by using the sorted order of type B* suffixes. */
219 static
220 saidx_t
221 construct_BWT(const sauchar_t *T, saidx_t *SA,
222               saidx_t *bucket_A, saidx_t *bucket_B,
223               saidx_t n, saidx_t m) {
224   saidx_t *i, *j, *k, *orig;
225   saidx_t s;
226   saint_t c0, c1, c2;
227
228   if(0 < m) {
229     /* Construct the sorted order of type B suffixes by using
230        the sorted order of type B* suffixes. */
231     for(c1 = ALPHABET_SIZE - 2; 0 <= c1; --c1) {
232       /* Scan the suffix array from right to left. */
233       for(i = SA + BUCKET_BSTAR(c1, c1 + 1),
234           j = SA + BUCKET_A(c1 + 1) - 1, k = NULL, c2 = -1;
235           i <= j;
236           --j) {
237         if(0 < (s = *j)) {
238           assert(T[s] == c1);
239           assert(((s + 1) < n) && (T[s] <= T[s + 1]));
240           assert(T[s - 1] <= T[s]);
241           c0 = T[--s];
242           *j = ~((saidx_t)c0);
243           if((0 < s) && (T[s - 1] > c0)) { s = ~s; }
244           if(c0 != c2) {
245             if(0 <= c2) { BUCKET_B(c2, c1) = k - SA; }
246             k = SA + BUCKET_B(c2 = c0, c1);
247           }
248           assert(k < j);
249           *k-- = s;
250         } else if(s != 0) {
251           *j = ~s;
252 #ifndef NDEBUG
253         } else {
254           assert(T[s] == c1);
255 #endif
256         }
257       }
258     }
259   }
260
261   /* Construct the BWTed string by using
262      the sorted order of type B suffixes. */
263   k = SA + BUCKET_A(c2 = T[n - 1]);
264   *k++ = (T[n - 2] < c2) ? ~((saidx_t)T[n - 2]) : (n - 1);
265   /* Scan the suffix array from left to right. */
266   for(i = SA, j = SA + n, orig = SA; i < j; ++i) {
267     if(0 < (s = *i)) {
268       assert(T[s - 1] >= T[s]);
269       c0 = T[--s];
270       *i = c0;
271       if((0 < s) && (T[s - 1] < c0)) { s = ~((saidx_t)T[s - 1]); }
272       if(c0 != c2) {
273         BUCKET_A(c2) = k - SA;
274         k = SA + BUCKET_A(c2 = c0);
275       }
276       assert(i < k);
277       *k++ = s;
278     } else if(s != 0) {
279       *i = ~s;
280     } else {
281       orig = i;
282     }
283   }
284
285   return orig - SA;
286 }
287 #endif
288
289
290 /*---------------------------------------------------------------------------*/
291
292 /*- Function -*/
293
294 /* XXX Modified from original: use provided temporary space instead of
295  * allocating it.  */
296 saint_t
297 divsufsort(const sauchar_t *T, saidx_t *SA, saidx_t n,
298      saidx_t *bucket_A, saidx_t *bucket_B) {
299   saidx_t m;
300
301   switch (n) {
302     case 0:
303       break;
304
305     case 1:
306       SA[0] = 0;
307       break;
308
309     case 2:
310       m = (T[0] < T[1]);
311       SA[m ^ 1] = 0;
312       SA[m] = 1;
313       break;
314
315     default:
316       m = sort_typeBstar(T, SA, bucket_A, bucket_B, n);
317       construct_SA(T, SA, bucket_A, bucket_B, n, m);
318       break;
319   }
320   return 0;
321 }
322
323 #if 0
324 saidx_t
325 divbwt(const sauchar_t *T, sauchar_t *U, saidx_t *A, saidx_t n) {
326   saidx_t *B;
327   saidx_t *bucket_A, *bucket_B;
328   saidx_t m, pidx, i;
329
330   /* Check arguments. */
331   if((T == NULL) || (U == NULL) || (n < 0)) { return -1; }
332   else if(n <= 1) { if(n == 1) { U[0] = T[0]; } return n; }
333
334   if((B = A) == NULL) { B = (saidx_t *)malloc((size_t)(n + 1) * sizeof(saidx_t)); }
335   bucket_A = (saidx_t *)malloc(BUCKET_A_SIZE * sizeof(saidx_t));
336   bucket_B = (saidx_t *)malloc(BUCKET_B_SIZE * sizeof(saidx_t));
337
338   /* Burrows-Wheeler Transform. */
339   if((B != NULL) && (bucket_A != NULL) && (bucket_B != NULL)) {
340     m = sort_typeBstar(T, B, bucket_A, bucket_B, n);
341     pidx = construct_BWT(T, B, bucket_A, bucket_B, n, m);
342
343     /* Copy to output string. */
344     U[0] = T[n - 1];
345     for(i = 0; i < pidx; ++i) { U[i + 1] = (sauchar_t)B[i]; }
346     for(i += 1; i < n; ++i) { U[i] = (sauchar_t)B[i]; }
347     pidx += 1;
348   } else {
349     pidx = -2;
350   }
351
352   free(bucket_B);
353   free(bucket_A);
354   if(A == NULL) { free(B); }
355
356   return pidx;
357 }
358 #endif
359
360 #if 0
361 const char *
362 divsufsort_version(void) {
363   return PROJECT_VERSION_FULL;
364 }
365 #endif