382 lines
9.6 KiB
C
382 lines
9.6 KiB
C
/*
|
|
* utils.c for libdivsufsort
|
|
* Copyright (c) 2003-2008 Yuta Mori All Rights Reserved.
|
|
*
|
|
* Permission is hereby granted, free of charge, to any person
|
|
* obtaining a copy of this software and associated documentation
|
|
* files (the "Software"), to deal in the Software without
|
|
* restriction, including without limitation the rights to use,
|
|
* copy, modify, merge, publish, distribute, sublicense, and/or sell
|
|
* copies of the Software, and to permit persons to whom the
|
|
* Software is furnished to do so, subject to the following
|
|
* conditions:
|
|
*
|
|
* The above copyright notice and this permission notice shall be
|
|
* included in all copies or substantial portions of the Software.
|
|
*
|
|
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
|
|
* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
|
|
* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
|
|
* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
|
|
* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
|
|
* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
|
|
* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
|
|
* OTHER DEALINGS IN THE SOFTWARE.
|
|
*/
|
|
|
|
#include "divsufsort_private.h"
|
|
|
|
|
|
/*- Private Function -*/
|
|
|
|
/* Binary search for inverse bwt. */
|
|
static
|
|
saidx_t
|
|
binarysearch_lower(const saidx_t *A, saidx_t size, saidx_t value) {
|
|
saidx_t half, i;
|
|
for(i = 0, half = size >> 1;
|
|
0 < size;
|
|
size = half, half >>= 1) {
|
|
if(A[i + half] < value) {
|
|
i += half + 1;
|
|
half -= (size & 1) ^ 1;
|
|
}
|
|
}
|
|
return i;
|
|
}
|
|
|
|
|
|
/*- Functions -*/
|
|
|
|
/* Burrows-Wheeler transform. */
|
|
saint_t
|
|
bw_transform(const sauchar_t *T, sauchar_t *U, saidx_t *SA,
|
|
saidx_t n, saidx_t *idx) {
|
|
saidx_t *A, i, j, p, t;
|
|
saint_t c;
|
|
|
|
/* Check arguments. */
|
|
if((T == NULL) || (U == NULL) || (n < 0) || (idx == NULL)) { return -1; }
|
|
if(n <= 1) {
|
|
if(n == 1) { U[0] = T[0]; }
|
|
*idx = n;
|
|
return 0;
|
|
}
|
|
|
|
if((A = SA) == NULL) {
|
|
i = divbwt(T, U, NULL, n);
|
|
if(0 <= i) { *idx = i; i = 0; }
|
|
return (saint_t)i;
|
|
}
|
|
|
|
/* BW transform. */
|
|
if(T == U) {
|
|
t = n;
|
|
for(i = 0, j = 0; i < n; ++i) {
|
|
p = t - 1;
|
|
t = A[i];
|
|
if(0 <= p) {
|
|
c = T[j];
|
|
U[j] = (j <= p) ? T[p] : (sauchar_t)A[p];
|
|
A[j] = c;
|
|
j++;
|
|
} else {
|
|
*idx = i;
|
|
}
|
|
}
|
|
p = t - 1;
|
|
if(0 <= p) {
|
|
c = T[j];
|
|
U[j] = (j <= p) ? T[p] : (sauchar_t)A[p];
|
|
A[j] = c;
|
|
} else {
|
|
*idx = i;
|
|
}
|
|
} else {
|
|
U[0] = T[n - 1];
|
|
for(i = 0; A[i] != 0; ++i) { U[i + 1] = T[A[i] - 1]; }
|
|
*idx = i + 1;
|
|
for(++i; i < n; ++i) { U[i] = T[A[i] - 1]; }
|
|
}
|
|
|
|
if(SA == NULL) {
|
|
/* Deallocate memory. */
|
|
free(A);
|
|
}
|
|
|
|
return 0;
|
|
}
|
|
|
|
/* Inverse Burrows-Wheeler transform. */
|
|
saint_t
|
|
inverse_bw_transform(const sauchar_t *T, sauchar_t *U, saidx_t *A,
|
|
saidx_t n, saidx_t idx) {
|
|
saidx_t C[ALPHABET_SIZE];
|
|
sauchar_t D[ALPHABET_SIZE];
|
|
saidx_t *B;
|
|
saidx_t i, p;
|
|
saint_t c, d;
|
|
|
|
/* Check arguments. */
|
|
if((T == NULL) || (U == NULL) || (n < 0) || (idx < 0) ||
|
|
(n < idx) || ((0 < n) && (idx == 0))) {
|
|
return -1;
|
|
}
|
|
if(n <= 1) { return 0; }
|
|
|
|
if((B = A) == NULL) {
|
|
/* Allocate n*sizeof(saidx_t) bytes of memory. */
|
|
if((B = (saidx_t *)malloc((size_t)n * sizeof(saidx_t))) == NULL) { return -2; }
|
|
}
|
|
|
|
/* Inverse BW transform. */
|
|
for(c = 0; c < ALPHABET_SIZE; ++c) { C[c] = 0; }
|
|
for(i = 0; i < n; ++i) { ++C[T[i]]; }
|
|
for(c = 0, d = 0, i = 0; c < ALPHABET_SIZE; ++c) {
|
|
p = C[c];
|
|
if(0 < p) {
|
|
C[c] = i;
|
|
D[d++] = (sauchar_t)c;
|
|
i += p;
|
|
}
|
|
}
|
|
for(i = 0; i < idx; ++i) { B[C[T[i]]++] = i; }
|
|
for( ; i < n; ++i) { B[C[T[i]]++] = i + 1; }
|
|
for(c = 0; c < d; ++c) { C[c] = C[D[c]]; }
|
|
for(i = 0, p = idx; i < n; ++i) {
|
|
U[i] = D[binarysearch_lower(C, d, p)];
|
|
p = B[p - 1];
|
|
}
|
|
|
|
if(A == NULL) {
|
|
/* Deallocate memory. */
|
|
free(B);
|
|
}
|
|
|
|
return 0;
|
|
}
|
|
|
|
/* Checks the suffix array SA of the string T. */
|
|
saint_t
|
|
sufcheck(const sauchar_t *T, const saidx_t *SA,
|
|
saidx_t n, saint_t verbose) {
|
|
saidx_t C[ALPHABET_SIZE];
|
|
saidx_t i, p, q, t;
|
|
saint_t c;
|
|
|
|
if(verbose) { fprintf(stderr, "sufcheck: "); }
|
|
|
|
/* Check arguments. */
|
|
if((T == NULL) || (SA == NULL) || (n < 0)) {
|
|
if(verbose) { fprintf(stderr, "Invalid arguments.\n"); }
|
|
return -1;
|
|
}
|
|
if(n == 0) {
|
|
if(verbose) { fprintf(stderr, "Done.\n"); }
|
|
return 0;
|
|
}
|
|
|
|
/* check range: [0..n-1] */
|
|
for(i = 0; i < n; ++i) {
|
|
if((SA[i] < 0) || (n <= SA[i])) {
|
|
if(verbose) {
|
|
fprintf(stderr, "Out of the range [0,%" PRIdSAIDX_T "].\n"
|
|
" SA[%" PRIdSAIDX_T "]=%" PRIdSAIDX_T "\n",
|
|
n - 1, i, SA[i]);
|
|
}
|
|
return -2;
|
|
}
|
|
}
|
|
|
|
/* check first characters. */
|
|
for(i = 1; i < n; ++i) {
|
|
if(T[SA[i - 1]] > T[SA[i]]) {
|
|
if(verbose) {
|
|
fprintf(stderr, "Suffixes in wrong order.\n"
|
|
" T[SA[%" PRIdSAIDX_T "]=%" PRIdSAIDX_T "]=%d"
|
|
" > T[SA[%" PRIdSAIDX_T "]=%" PRIdSAIDX_T "]=%d\n",
|
|
i - 1, SA[i - 1], T[SA[i - 1]], i, SA[i], T[SA[i]]);
|
|
}
|
|
return -3;
|
|
}
|
|
}
|
|
|
|
/* check suffixes. */
|
|
for(i = 0; i < ALPHABET_SIZE; ++i) { C[i] = 0; }
|
|
for(i = 0; i < n; ++i) { ++C[T[i]]; }
|
|
for(i = 0, p = 0; i < ALPHABET_SIZE; ++i) {
|
|
t = C[i];
|
|
C[i] = p;
|
|
p += t;
|
|
}
|
|
|
|
q = C[T[n - 1]];
|
|
C[T[n - 1]] += 1;
|
|
for(i = 0; i < n; ++i) {
|
|
p = SA[i];
|
|
if(0 < p) {
|
|
c = T[--p];
|
|
t = C[c];
|
|
} else {
|
|
c = T[p = n - 1];
|
|
t = q;
|
|
}
|
|
if((t < 0) || (p != SA[t])) {
|
|
if(verbose) {
|
|
fprintf(stderr, "Suffix in wrong position.\n"
|
|
" SA[%" PRIdSAIDX_T "]=%" PRIdSAIDX_T " or\n"
|
|
" SA[%" PRIdSAIDX_T "]=%" PRIdSAIDX_T "\n",
|
|
t, (0 <= t) ? SA[t] : -1, i, SA[i]);
|
|
}
|
|
return -4;
|
|
}
|
|
if(t != q) {
|
|
++C[c];
|
|
if((n <= C[c]) || (T[SA[C[c]]] != c)) { C[c] = -1; }
|
|
}
|
|
}
|
|
|
|
if(1 <= verbose) { fprintf(stderr, "Done.\n"); }
|
|
return 0;
|
|
}
|
|
|
|
|
|
static
|
|
int
|
|
_compare(const sauchar_t *T, saidx_t Tsize,
|
|
const sauchar_t *P, saidx_t Psize,
|
|
saidx_t suf, saidx_t *match) {
|
|
saidx_t i, j;
|
|
saint_t r;
|
|
for(i = suf + *match, j = *match, r = 0;
|
|
(i < Tsize) && (j < Psize) && ((r = T[i] - P[j]) == 0); ++i, ++j) { }
|
|
*match = j;
|
|
return (r == 0) ? -(j != Psize) : r;
|
|
}
|
|
|
|
/* Search for the pattern P in the string T. */
|
|
saidx_t
|
|
sa_search(const sauchar_t *T, saidx_t Tsize,
|
|
const sauchar_t *P, saidx_t Psize,
|
|
const saidx_t *SA, saidx_t SAsize,
|
|
saidx_t *idx) {
|
|
saidx_t size, lsize, rsize, half;
|
|
saidx_t match, lmatch, rmatch;
|
|
saidx_t llmatch, lrmatch, rlmatch, rrmatch;
|
|
saidx_t i, j, k;
|
|
saint_t r;
|
|
|
|
if(idx != NULL) { *idx = -1; }
|
|
if((T == NULL) || (P == NULL) || (SA == NULL) ||
|
|
(Tsize < 0) || (Psize < 0) || (SAsize < 0)) { return -1; }
|
|
if((Tsize == 0) || (SAsize == 0)) { return 0; }
|
|
if(Psize == 0) { if(idx != NULL) { *idx = 0; } return SAsize; }
|
|
|
|
for(i = j = k = 0, lmatch = rmatch = 0, size = SAsize, half = size >> 1;
|
|
0 < size;
|
|
size = half, half >>= 1) {
|
|
match = MIN(lmatch, rmatch);
|
|
r = _compare(T, Tsize, P, Psize, SA[i + half], &match);
|
|
if(r < 0) {
|
|
i += half + 1;
|
|
half -= (size & 1) ^ 1;
|
|
lmatch = match;
|
|
} else if(r > 0) {
|
|
rmatch = match;
|
|
} else {
|
|
lsize = half, j = i, rsize = size - half - 1, k = i + half + 1;
|
|
|
|
/* left part */
|
|
for(llmatch = lmatch, lrmatch = match, half = lsize >> 1;
|
|
0 < lsize;
|
|
lsize = half, half >>= 1) {
|
|
lmatch = MIN(llmatch, lrmatch);
|
|
r = _compare(T, Tsize, P, Psize, SA[j + half], &lmatch);
|
|
if(r < 0) {
|
|
j += half + 1;
|
|
half -= (lsize & 1) ^ 1;
|
|
llmatch = lmatch;
|
|
} else {
|
|
lrmatch = lmatch;
|
|
}
|
|
}
|
|
|
|
/* right part */
|
|
for(rlmatch = match, rrmatch = rmatch, half = rsize >> 1;
|
|
0 < rsize;
|
|
rsize = half, half >>= 1) {
|
|
rmatch = MIN(rlmatch, rrmatch);
|
|
r = _compare(T, Tsize, P, Psize, SA[k + half], &rmatch);
|
|
if(r <= 0) {
|
|
k += half + 1;
|
|
half -= (rsize & 1) ^ 1;
|
|
rlmatch = rmatch;
|
|
} else {
|
|
rrmatch = rmatch;
|
|
}
|
|
}
|
|
|
|
break;
|
|
}
|
|
}
|
|
|
|
if(idx != NULL) { *idx = (0 < (k - j)) ? j : i; }
|
|
return k - j;
|
|
}
|
|
|
|
/* Search for the character c in the string T. */
|
|
saidx_t
|
|
sa_simplesearch(const sauchar_t *T, saidx_t Tsize,
|
|
const saidx_t *SA, saidx_t SAsize,
|
|
saint_t c, saidx_t *idx) {
|
|
saidx_t size, lsize, rsize, half;
|
|
saidx_t i, j, k, p;
|
|
saint_t r;
|
|
|
|
if(idx != NULL) { *idx = -1; }
|
|
if((T == NULL) || (SA == NULL) || (Tsize < 0) || (SAsize < 0)) { return -1; }
|
|
if((Tsize == 0) || (SAsize == 0)) { return 0; }
|
|
|
|
for(i = j = k = 0, size = SAsize, half = size >> 1;
|
|
0 < size;
|
|
size = half, half >>= 1) {
|
|
p = SA[i + half];
|
|
r = (p < Tsize) ? T[p] - c : -1;
|
|
if(r < 0) {
|
|
i += half + 1;
|
|
half -= (size & 1) ^ 1;
|
|
} else if(r == 0) {
|
|
lsize = half, j = i, rsize = size - half - 1, k = i + half + 1;
|
|
|
|
/* left part */
|
|
for(half = lsize >> 1;
|
|
0 < lsize;
|
|
lsize = half, half >>= 1) {
|
|
p = SA[j + half];
|
|
r = (p < Tsize) ? T[p] - c : -1;
|
|
if(r < 0) {
|
|
j += half + 1;
|
|
half -= (lsize & 1) ^ 1;
|
|
}
|
|
}
|
|
|
|
/* right part */
|
|
for(half = rsize >> 1;
|
|
0 < rsize;
|
|
rsize = half, half >>= 1) {
|
|
p = SA[k + half];
|
|
r = (p < Tsize) ? T[p] - c : -1;
|
|
if(r <= 0) {
|
|
k += half + 1;
|
|
half -= (rsize & 1) ^ 1;
|
|
}
|
|
}
|
|
|
|
break;
|
|
}
|
|
}
|
|
|
|
if(idx != NULL) { *idx = (0 < (k - j)) ? j : i; }
|
|
return k - j;
|
|
}
|