diff --git a/bwt.c b/bwt.c index 7b37fe5..47b06e2 100644 --- a/bwt.c +++ b/bwt.c @@ -161,7 +161,7 @@ void bwt_2occ(const bwt_t *bwt, bwtint_t k, bwtint_t l, ubyte_t c, bwtint_t *ok, void bwt_occ4(const bwt_t *bwt, bwtint_t k, bwtint_t cnt[4]) { bwtint_t l, j, x; - uint32_t *p; + uint32_t *p, tmp; if (k == (bwtint_t)(-1)) { memset(cnt, 0, 4 * sizeof(bwtint_t)); return; @@ -171,9 +171,10 @@ void bwt_occ4(const bwt_t *bwt, bwtint_t k, bwtint_t cnt[4]) memcpy(cnt, p, 4 * sizeof(bwtint_t)); p += sizeof(bwtint_t); j = k >> 4 << 4; - for (l = k / OCC_INTERVAL * OCC_INTERVAL, x = 0; l < j; l += 16, ++p) + for (l = k & ~OCC_INTV_MASK, x = 0; l < j; l += 16, ++p) x += __occ_aux4(bwt, *p); - x += __occ_aux4(bwt, *p & ~((1U<<((~k&15)<<1)) - 1)) - (~k&15); + tmp = *p & ~((1U<<((~k&15)<<1)) - 1); + x += __occ_aux4(bwt, tmp) - (~k&15); cnt[0] += x&0xff; cnt[1] += x>>8&0xff; cnt[2] += x>>16&0xff; cnt[3] += x>>24; } @@ -188,7 +189,7 @@ void bwt_2occ4(const bwt_t *bwt, bwtint_t k, bwtint_t l, bwtint_t cntk[4], bwtin bwt_occ4(bwt, l, cntl); } else { bwtint_t i, j, x, y; - uint32_t *p; + uint32_t *p, tmp; if (k >= bwt->primary) --k; // because $ is not in bwt if (l >= bwt->primary) --l; p = bwt_occ_intv(bwt, k); @@ -196,14 +197,16 @@ void bwt_2occ4(const bwt_t *bwt, bwtint_t k, bwtint_t l, bwtint_t cntk[4], bwtin p += sizeof(bwtint_t); // prepare cntk[] j = k >> 4 << 4; - for (i = k / OCC_INTERVAL * OCC_INTERVAL, x = 0; i < j; i += 16, ++p) + for (i = k & ~OCC_INTV_MASK, x = 0; i < j; i += 16, ++p) x += __occ_aux4(bwt, *p); y = x; - x += __occ_aux4(bwt, *p & ~((1U<<((~k&15)<<1)) - 1)) - (~k&15); + tmp = *p & ~((1U<<((~k&15)<<1)) - 1); + x += __occ_aux4(bwt, tmp) - (~k&15); // calculate cntl[] and finalize cntk[] j = l >> 4 << 4; for (; i < j; i += 16, ++p) y += __occ_aux4(bwt, *p); - y += __occ_aux4(bwt, *p & ~((1U<<((~l&15)<<1)) - 1)) - (~l&15); + tmp = *p & ~((1U<<((~l&15)<<1)) - 1); + y += __occ_aux4(bwt, tmp) - (~l&15); memcpy(cntl, cntk, 4 * sizeof(bwtint_t)); cntk[0] += x&0xff; cntk[1] += x>>8&0xff; cntk[2] += x>>16&0xff; cntk[3] += x>>24; cntl[0] += y&0xff; cntl[1] += y>>8&0xff; cntl[2] += y>>16&0xff; cntl[3] += y>>24; diff --git a/bwt.h b/bwt.h index e06329a..ab5aecd 100644 --- a/bwt.h +++ b/bwt.h @@ -30,8 +30,10 @@ #include -// requirement: (OCC_INTERVAL%16 == 0); please DO NOT change this line -#define OCC_INTERVAL 0x80 +// requirement: (OCC_INTERVAL%16 == 0); please DO NOT change this line because some part of the code assume OCC_INTERVAL=0x80 +#define OCC_INTV_SHIFT 7 +#define OCC_INTERVAL (1LL<