src/static_bitsequence/sdarray.cpp

   1 #include <sdarray.h>
   2 #include <stdint.h>
   3 using std::min;
   4 using std::max;
   5 #if 0
   6 typedef unsigned int qword;
   7 #define logD 4
   8 #else
   9 typedef unsigned long long qword;
  10 #define logD 5
  11 #endif
  12 #define PBS (sizeof(uint)*8)
  13 #define D (1<<logD)
  14 #define logM 5
  15 #define M (1<<logM)
  16 #define logP 8
  17 #define P (1<<logP)
  18 #define logLL 16                 // size of word
  19 #define LL (1<<logLL)
  20 //#define logLLL 7
  21 #define logLLL 5
  22 //#define LLL 128
  23 //#define LLL 32
  24 #define LLL (1<<logLLL)
  25 //#define logL 10
  26 //#define logL (logLL-3)
  27 #define logL (logLL-1-5)
  28 #define L (1<<logL)
  29
  30 int __blog(int x) {
  31   int l;
  32   l = 0;
  33   while (x>0) {
  34     x>>=1;
  35     l++;
  36   }
  37   return l;
  38 }
  39
  40
  41 int __setbit(uint *B, int i,int x) {
  42   int j,l;
  43   //printf("%u\n",D);
  44   j = i / D;
  45   l = i % D;
  46   if (x==0) B[j] &= (~(1<<(D-1-l)));
  47   else if (x==1) B[j] |= (1<<(D-1-l));
  48   else {
  49     printf("error __setbit x=%d\n",x);
  50     exit(1);
  51   }
  52   return x;
  53 }
  54
  55
  56 int __setbit2(uchar *B, int i,int x) {
  57   int j,l;
  58
  59   j = i / 8;
  60   l = i % 8;
  61   if (x==0) B[j] &= (~(1<<(8-1-l)));
  62   else if (x==1) B[j] |= (1<<(8-1-l));
  63   else {
  64     printf("error __setbit2 x=%d\n",x);
  65     exit(1);
  66   }
  67   return x;
  68 }
  69
  70
  71 int __setbits(uint *B, int i, int d, int x) {
  72   int j;
  73
  74   for (j=0; j<d; j++) {
  75     __setbit(B,i+j,(x>>(d-j-1))&1);
  76   }
  77   return x;
  78 }
  79
  80
  81 int __getbit(uint *B, int i) {
  82   int j,l;
  83
  84   //j = i / D;
  85   //l = i % D;
  86   j = i >> logD;
  87   l = i & (D-1);
  88   return (B[j] >> (D-1-l)) & 1;
  89 }
  90
  91
  92 static int __getbit2(uchar *B, int i) {
  93   int j,l;
  94
  95   //j = i / D;
  96   //l = i % D;
  97   j = i >> 3;
  98   l = i & (8-1);
  99   return (B[j] >> (8-1-l)) & 1;
 100 }
 101
 102
 103 #if 1
 104 uint __getbits_aux(uint *B, int i, int d) {
 105   qword x,z;
 106   int j = i;
 107   B += (i >> logD);
 108   i &= (D-1);
 109   if (d==8 && (j & 7) == 0) {
 110     i = (24 - i) >> 3;
 111     x = (uint) (((unsigned char*) B)[i]);
 112   } else if (i+d <= 2*D) {
 113     x = (((qword)B[0]) << D) + B[1];
 114     x <<= i;
 115     x >>= (D*2-1-d);
 116     x >>= 1;
 117   }
 118   else {
 119     fprintf(stderr, "Warning: %d, %d\n", D, d);
 120     x = (((qword)B[0])<<D)+B[1];
 121     z = (x<<D)+B[2];
 122     x <<= i;
 123     x &= (((qword)1L<<D)-1)<<D;
 124     z <<= i;
 125     z >>= D;
 126     x += z;
 127     x >>= (2*D-d);
 128   }
 129
 130   return x;
 131 }
 132
 133 static uint __getbits(uint *B, int i, int d)
 134 {
 135   uint64_t x;
 136   B += (i >> logD);
 137   i &= (D-1);
 138   x = ((uint64_t *) B)[0];
 139   x = (x << 32)|(x >> 32);
 140   x = (x << i) >> (2*D - d);
 141   return x;
 142 }
 143
 144 #endif
 145
 146 #if HAS_NATIVE_POPCOUNT
 147 static inline unsigned int popcount(unsigned int n){
 148   asm ("popcnt %1, %0" : "=r" (n) : "0" (n));
 149   return n;
 150 }
 151
 152 // static inline unsigned int popcount8(unsigned int n) {
 153 //   return __builtin_popcount(n & 0xff);
 154 // }
 155
 156 static inline unsigned int _fast_popcount(uchar x)
 157 {
 158   return __builtin_popcount(x);
 159 }
 160
 161
 162 #else
 163 static const unsigned int _popCount[] = {
 164   0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,
 165   1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,
 166   1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,
 167   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 168   1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,
 169   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 170   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 171   3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,
 172   1,2,2,3,2,3,3,4,2,3,3,4,3,4,4,5,
 173   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 174   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 175   3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,
 176   2,3,3,4,3,4,4,5,3,4,4,5,4,5,5,6,
 177   3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,
 178   3,4,4,5,4,5,5,6,4,5,5,6,5,6,6,7,
 179   4,5,5,6,5,6,6,7,5,6,6,7,6,7,7,8
 180 };
 181
 182 static inline unsigned int
 183 _fast_popcount(int x) {
 184   return _popCount[x & 0xff];
 185 }
 186 #endif
 187
 188
 189
 190 //static unsigned int __selecttbl[8*256];
 191 static int built = 0;
 192
 193 void make___selecttbl(void) {
 194   if(built) return;
 195   built = 1;
 196   int i,x,r;
 197   uint buf[1] = { 0 };
 198
 199   for (x = 0; x < 256; x++) {
 200     __setbits(buf,0,8,x);
 201 //    for (r=0; r<8; r++) __selecttbl[(r<<8)+x] = -1;
 202     r = 0;
 203     for (i=0; i<8; i++) {
 204       if (__getbit(buf,i)) {
 205         //      __selecttbl[(r<<8)+x] = i;
 206         r++;
 207       }
 208     }
 209   }
 210 }
 211
 212 int selectd2_save(selectd2 * s, FILE * fp) {
 213   uint wr = 0;
 214   wr += fwrite(&s->n,sizeof(uint),1,fp);
 215   wr += fwrite(&s->m,sizeof(uint),1,fp);
 216   wr += fwrite(&s->size,sizeof(uint),1,fp);
 217   wr += fwrite(&s->ss_len,sizeof(uint),1,fp);
 218   wr += fwrite(&s->sl_len,sizeof(uint),1,fp);
 219   wr += fwrite(s->buf,sizeof(uchar),(s->n+7)/8+1,fp);
 220   uint nl = (s->m-1) / L + 1;
 221   wr += fwrite(s->lp,sizeof(uint),nl+1,fp);
 222   wr += fwrite(s->p,sizeof(uint),nl+1,fp);
 223   wr += fwrite(s->ss,sizeof(ushort),s->ss_len,fp);
 224   wr += fwrite(s->sl,sizeof(uint),s->sl_len,fp);
 225   if(wr!=s->sl_len+s->ss_len+2*(nl+1)+(s->n+7)/8+1+5)
 226     return 1;
 227   return 0;
 228 }
 229
 230
 231 int selectd2_load(selectd2 * s, FILE * fp) {
 232   uint rd = 0;
 233   rd += fread(&s->n,sizeof(uint),1,fp);
 234   rd += fread(&s->m,sizeof(uint),1,fp);
 235   rd += fread(&s->size,sizeof(uint),1,fp);
 236   rd += fread(&s->ss_len,sizeof(uint),1,fp);
 237   rd += fread(&s->sl_len,sizeof(uint),1,fp);
 238   s->buf = new uchar[(s->n+7)/8+1];
 239   rd += fread(s->buf,sizeof(uchar),(s->n+7)/8+1,fp);
 240   uint nl = (s->m-1) / L + 1;
 241   s->lp = new uint[nl+1];
 242   rd += fread(s->lp,sizeof(uint),nl+1,fp);
 243   s->p = new uint[nl+1];
 244   rd += fread(s->p,sizeof(uint),nl+1,fp);
 245   s->ss = new ushort[s->ss_len];
 246   rd += fread(s->ss,sizeof(ushort),s->ss_len,fp);
 247   s->sl = new uint[s->sl_len];
 248   rd += fread(s->sl,sizeof(uint),s->sl_len,fp);
 249   if(rd!=s->sl_len+s->ss_len+2*(nl+1)+(s->n+7)/8+1+5)
 250     return 1;
 251   return 0;
 252 }
 253
 254
 255 void selectd2_free(selectd2 * s) {
 256   //delete [] s->buf;
 257   delete [] s->lp;
 258   delete [] s->p;
 259   delete [] s->ss;
 260   delete [] s->sl;
 261 }
 262
 263
 264 int selectd2_construct(selectd2 *select, int n, uchar *buf) {
 265   int i,m;
 266   int nl;
 267   int p,pp;
 268   int il,is,ms,ml;
 269   int r;
 270   uint *s;
 271
 272   //make___selecttbl();
 273
 274   if (L/LLL == 0) {
 275     printf("ERROR: L=%d LLL=%d\n",L,LLL);
 276     exit(1);
 277   }
 278
 279   m = 0;
 280   for (i=0; i<n; i++) m += __getbit2(buf,i);
 281   select->n = n;
 282   select->m = m;
 283   //printf("n=%d m=%d\n",n,m);
 284
 285   select->buf = buf;
 286
 287   s = new uint[m];
 288   m = 0;
 289   for (i=0; i<n; i++) {
 290     if (__getbit2(buf,i)) {
 291       m++;
 292       s[m-1] = i;
 293     }
 294   }
 295
 296   nl = (m-1) / L + 1;
 297   select->size = 0;              //ignoring buf, shared with selects3
 298   select->lp = new uint[nl+1];
 299   for(int k=0;k<nl+1;k++) select->lp[k]=0;
 300   select->size += (nl+1)*sizeof(uint);
 301   select->p = new uint[nl+1];
 302   for(int k=0;k<nl+1;k++) select->p[k]=0;
 303   select->size += (nl+1)*sizeof(uint);
 304
 305
 306   for (r = 0; r < 2; r++) {
 307     ml = ms = 0;
 308     for (il = 0; il < nl; il++) {
 309       pp = s[il*L];
 310       select->lp[il] = pp;
 311       i = min((il+1)*L-1,m-1);
 312       p = s[i];
 313       // printf("p-pp=%d, LL=%d\n",p-pp, LL);
 314
 315       if (p - pp >= LL) {
 316         if (r == 1) {
 317           for (is = 0; is < L; is++) {
 318             if (il*L+is >= m) break;
 319             select->sl[ml*L+is] = s[il*L+is];
 320           }
 321         }
 322
 323         select->p[il] = -((ml<<logL)+1);
 324         ml++;
 325       }
 326       else {
 327         if (r == 1) {
 328           for (is = 0; is < L/LLL; is++) {
 329             if (il*L+is*LLL >= m) break;
 330             select->ss[ms*(L/LLL)+is] = s[il*L+is*LLL] - pp;
 331           }
 332         }
 333         select->p[il] = ms << (logL-logLLL);
 334         ms++;
 335       }
 336     }
 337
 338     if (r == 0) {
 339       select->sl = new uint[ml*L+1];
 340       for(int k=0;k<ml*L+1;k++) select->sl[k]=0;
 341       select->size += sizeof(uint)*(ml*L+1);
 342       select->sl_len = ml*L+1;
 343       select->ss = new ushort[ms*(L/LLL)+1];
 344       for(int k=0;k<ms*(L/LLL)+1;k++) select->ss[k]=0;
 345       select->ss_len = ms*(L/LLL)+1;
 346       select->size += sizeof(ushort)*(ms*(L/LLL)+1);
 347     }
 348   }
 349   delete [] s;
 350   return 0;
 351 }
 352
 353
 354 int selectd2_select1(selectd2 *select, int i) {
 355   int p,r;
 356   int il;
 357   int rr;
 358   uchar *q;
 359   if (i <= 0) return -1;
 360   i--;
 361
 362   il = select->p[i>>logL];
 363   if (il < 0) {
 364     il = -il-1;
 365     //p = select->sl[(il<<logL)+(i & (L-1))];
 366     p = select->sl[il+(i & (L-1))];
 367   }
 368   else {
 369     p = select->lp[i>>logL];
 370     p += select->ss[il+((i & (L-1))>>logLLL)];
 371     r = i - (i & (LLL-1));
 372
 373     q = &(select->buf[p>>3]);
 374
 375     rr = p & (8-1);
 376     r -= _fast_popcount(*q >> (8-1-rr));
 377
 378     while (1) {
 379         //rr = _popCount[*q];
 380       rr = _fast_popcount(*q);
 381       if (r + rr >= i) break;
 382       r += rr;
 383       //p += 8;
 384       q++;
 385     }
 386       p = (q - select->buf) << 3;
 387       p += selecttbl[((i-r-1)<<8)+(*q)];
 388   }
 389   return p;
 390 }
 391
 392 int selectd2_select0(selectd2 *select, int i) {
 393   int p,r;
 394   int il;
 395   int rr;
 396   uchar *q;
 397
 398   if (i <= 0) return -1;
 399   i--;
 400
 401   il = select->p[i>>logL];
 402   if (il < 0) {
 403     il = -il-1;
 404     return select->sl[il+(i & (L-1))];
 405
 406   } else {
 407     p = select->lp[i>>logL];
 408
 409     p += select->ss[il+((i & (L-1))>>logLLL)];
 410     r = i - (i & (LLL-1));
 411
 412     q = &(select->buf[p>>3]);
 413
 414     rr = p & 7;
 415
 416     uint masked_q = *q ^ 0xff;
 417
 418     r -= _fast_popcount(masked_q >> (7-rr));
 419
 420     for(;;) {
 421       rr = _fast_popcount(masked_q);
 422       if (rr >= i-r) {
 423         p = (q - select->buf) << 3;
 424         p += selecttbl[((i-r-1)<<8)+masked_q];
 425         return p;
 426       };
 427       r += rr;
 428       q++;
 429       masked_q = *q ^ 0xff;
 430     };
 431   }
 432 }
 433
 434 int selectd2_select(selectd2 *select, int i,int f) {
 435   return f ? selectd2_select1(select, i) :
 436     selectd2_select0(select, i);
 437 }
 438
 439
 440 int selectd2_select2(selectd2 *select, int i,int f, int *st, int *en) {
 441   int p,r,p2;
 442   int il;
 443   int rr;
 444   uchar *q;
 445
 446   if (i == 0) {
 447     *st = -1;
 448     return -1;
 449   }
 450
 451   #if 0
 452   if (i > select->m) {
 453     printf("ERROR: m=%d i=%d\n",select->m,i);
 454     exit(1);
 455   }
 456   #endif
 457
 458   i--;
 459
 460   il = select->p[i>>logL];
 461   if (il < 0) {
 462     il = -il-1;
 463     //p = select->sl[(il<<logL)+(i & (L-1))];
 464     p = select->sl[il+(i & (L-1))];
 465
 466     if ((i>>logL) == ((i+1)>>logL)) {
 467       p2 = select->sl[il+((i+1) & (L-1))];
 468     }
 469     else {
 470       p2 = selectd2_select(select,i+2,f);
 471     }
 472   }
 473   else {
 474     p = select->lp[i>>logL];
 475     //p += select->ss[(il<<(logL-logLLL))+(i & (L-1))/LLL];
 476     p += select->ss[il+((i & (L-1))>>logLLL)];
 477     r = i - (i & (LLL-1));
 478
 479     q = &(select->buf[p>>3]);
 480
 481     if (f == 1) {
 482       rr = p & (8-1);
 483       //r -= _popCount[*q >> (8-1-rr)];
 484       r -= _fast_popcount(*q >> (8-1-rr));
 485       //p = p - rr;
 486
 487       while (1) {
 488         //rr = _popCount[*q];
 489         rr = _fast_popcount(*q);
 490         if (r + rr >= i) break;
 491         r += rr;
 492         //p += 8;
 493         q++;
 494       }
 495       p = (q - select->buf) << 3;
 496       p += selecttbl[((i-r-1)<<8)+(*q)];
 497
 498       if ((i>>logL) == ((i+1)>>logL)) {
 499         i++;
 500         while (1) {
 501           //rr = _popCount[*q];
 502           r = _fast_popcount(*q);
 503           if (r + rr >= i) break;
 504           r += rr;
 505           q++;
 506         }
 507         p2 = (q - select->buf) << 3;
 508         p2 += selecttbl[((i-r-1)<<8)+(*q)];
 509       }
 510       else {
 511         p2 = selectd2_select(select,i+2,f);
 512       }
 513
 514     }
 515     else {
 516       rr = p & (8-1);
 517       //r -= _popCount[(*q ^ 0xff) >> (8-1-rr)];
 518       r -= _fast_popcount((*q ^ 0xff) >> (8-1-rr));
 519       //p = p - rr;
 520
 521       while (1) {
 522         //rr = _popCount[*q ^ 0xff];
 523         rr = _fast_popcount(*q ^ 0xff);
 524         if (r + rr >= i) break;
 525         r += rr;
 526         //p += 8;
 527         q++;
 528       }
 529       p = (q - select->buf) << 3;
 530       p += selecttbl[((i-r-1)<<8)+(*q ^ 0xff)];
 531
 532       if ((i>>logL) == ((i+1)>>logL)) {
 533         i++;
 534         while (1) {
 535           //rr = _popCount[*q ^ 0xff];
 536           rr = _fast_popcount(*q ^ 0xff);
 537           if (r + rr >= i) break;
 538           r += rr;
 539           q++;
 540         }
 541         p2 = (q - select->buf) << 3;
 542         p2 += selecttbl[((i-r-1)<<8)+(*q ^ 0xff)];
 543       }
 544       else {
 545         p2 = selectd2_select(select,i+2,f);
 546       }
 547     }
 548   }
 549   *st = p;
 550   *en = p2;
 551   return p;
 552 }
 553
 554
 555 int selects3_save(selects3 * s, FILE * fp) {
 556   uint wr = 0;
 557   wr += fwrite(&s->n,sizeof(uint),1,fp);
 558   wr += fwrite(&s->m,sizeof(uint),1,fp);
 559   wr += fwrite(&s->size,sizeof(uint),1,fp);
 560   wr += fwrite(&s->d,sizeof(uint),1,fp);
 561   wr += fwrite(&s->hi_len,sizeof(uint),1,fp);
 562   wr += fwrite(&s->low_len,sizeof(uint),1,fp);
 563   wr += fwrite(s->hi,sizeof(uchar),s->hi_len,fp);
 564   wr += fwrite(s->low,sizeof(uint),s->low_len,fp);
 565   if(wr!=(6+s->hi_len+s->low_len))
 566     return 1;
 567   if(selectd2_save(s->sd0,fp)) return 2;
 568   if(selectd2_save(s->sd1,fp)) return 3;
 569   return 0;
 570 }
 571
 572
 573 int selects3_load(selects3 * s, FILE * fp) {
 574   uint rd = 0;
 575   rd += fread(&s->n,sizeof(uint),1,fp);
 576   rd += fread(&s->m,sizeof(uint),1,fp);
 577   rd += fread(&s->size,sizeof(uint),1,fp);
 578   rd += fread(&s->d,sizeof(uint),1,fp);
 579   rd += fread(&s->hi_len,sizeof(uint),1,fp);
 580   rd += fread(&s->low_len,sizeof(uint),1,fp);
 581   s->hi = new uchar[s->hi_len];
 582   rd += fread(s->hi,sizeof(uchar),s->hi_len,fp);
 583   s->low = new uint[s->low_len];
 584   rd += fread(s->low,sizeof(uint),s->low_len,fp);
 585   if(rd!=(6+s->hi_len+s->low_len))
 586     return 1;
 587   s->sd0 = new selectd2;
 588   if(selectd2_load(s->sd0,fp)) return 2;
 589   s->sd1 = new selectd2;
 590   if(selectd2_load(s->sd1,fp)) return 3;
 591   delete [] s->sd0->buf;
 592   delete [] s->sd1->buf;
 593   s->sd0->buf = s->hi;
 594   s->sd1->buf = s->hi;
 595   return 0;
 596 }
 597
 598
 599 void selects3_free(selects3 * s) {
 600   delete [] s->hi;
 601   delete [] s->low;
 602   //delete [] s->sd0->buf;
 603   selectd2_free(s->sd0);
 604   delete s->sd0;
 605   selectd2_free(s->sd1);
 606   delete s->sd1;
 607 }
 608
 609
 610 int selects3_construct(selects3 *select, int n, uint *buf) {
 611   int i,m;
 612   int d,mm;
 613   uint *low;
 614   uchar *buf2;
 615   selectd2 *sd0,*sd1;
 616
 617   m = 0;
 618   for (i=0; i<n; i++) m += __getbit(buf,i);
 619   select->n = n;
 620   select->m = m;
 621
 622   if (m == 0) return 0;
 623
 624   mm = m;
 625   d = 0;
 626   while (mm < n) {
 627     mm <<= 1;
 628     d++;
 629   }
 630
 631   select->d = d;
 632
 633   buf2 = new uchar[(2*m+8-1)/8+1];
 634   for(int k=0;k<(2*m+8-1)/8+1;k++) buf2[k]=0;
 635   select->hi_len = (2*m+8-1)/8+1;
 636   low = new uint[(d*m+PBS-1)/PBS+1];
 637   for(uint k=0;k<(d*m+PBS-1)/PBS+1;k++) low[k]=0;
 638   select->low_len = (d*m+PBS-1)/PBS+1;
 639
 640   select->hi = buf2;
 641   select->low = low;
 642   select->size = sizeof(uchar)*((2*m+8-1)/8+1) + sizeof(uint)*((d*m+PBS-1)/PBS+1);
 643
 644   for (i=0; i<m*2; i++) __setbit2(buf2,i,0);
 645
 646   m = 0;
 647   for (i=0; i<n; i++) {
 648     if (__getbit(buf,i)) {
 649       __setbit2(buf2,(i>>d)+m,1);
 650       __setbits(low,m*d,d,i & ((1<<d)-1));
 651       m++;
 652     }
 653   }
 654
 655   sd1 = new selectd2;
 656   sd0 = new selectd2;
 657   select->size += 2*sizeof(selectd2);
 658
 659   selectd2_construct(sd1,m*2,buf2);
 660   select->sd1 = sd1;
 661
 662   for (i=0; i<m*2; i++) __setbit2(buf2,i,1-__getbit2(buf2,i));
 663   selectd2_construct(sd0,m*2,buf2);
 664   select->sd0 = sd0;
 665
 666   for (i=0; i<m*2; i++) __setbit2(buf2,i,1-__getbit2(buf2,i));
 667   return 0;
 668 }
 669
 670 //selects3 * lasts3=NULL;
 671 //int lasti=0;
 672 //int lasts=0;
 673
 674 int selects3_select(selects3 *select, int i) {
 675   int d,x;
 676
 677   #if 0
 678   if (i > select->m) {
 679     printf("ERROR: m=%d i=%d\n",select->m,i);
 680     exit(1);
 681   }
 682   #endif
 683
 684   if (i == 0) return -1;
 685
 686   d = select->d;
 687         /*if(select->lasti==(uint)i-1) {
 688                 while(!__getbit2(select->sd1->buf,++select->lasts));
 689         }
 690         else {
 691           select->lasts = selectd2_select(select->sd1,i,1);
 692         }
 693         select->lasti = i;
 694         //lasts3 = select; */
 695         x = selectd2_select1(select->sd1,i) - (i-1);
 696   //x = (select->lasts-(i-1)) << d;
 697   x <<= d;
 698   x += __getbits(select->low,(i-1)*d,d);
 699   return x;
 700 }
 701
 702 void pr_byte(FILE* fp, uchar b)
 703 {
 704   uchar * buff = &b;
 705   for(int i = 0; i < 8; i++){
 706     fprintf(stderr, "%i", __getbit2(buff, i));
 707   };
 708 }
 709
 710 int selects3_selectnext(selects3 *select, int i) {
 711   //return selects3_select(select,selects3_rank(select,i)+1);
 712   int d,x,w,y;
 713   int r,j;
 714   int z,ii;
 715   int xoffset;
 716   uint *q;
 717   d = select->d;
 718   q = select->low;
 719   ii = i>>d;
 720   y = selectd2_select0(select->sd0,ii)+1;
 721   int k2=y-ii;
 722   x = y - ii;
 723   int x_orig = x;
 724   j = i - (ii<<d);
 725   r = y & 7;
 726   y >>= 3;
 727   z = select->hi[y];
 728   xoffset = x * d;
 729   while (1) {
 730     if (((z << r) & 0x80) == 0) {
 731       k2 += (x!=x_orig);
 732       break;
 733     };
 734
 735     w = __getbits(q,xoffset,d);
 736     if (w >= j) {
 737       if (w == j) {
 738         if (__getbit2(select->hi,((y << 3)+r)))  k2 ++;
 739         x++;
 740         r++;
 741         xoffset += d;
 742       };
 743       break;
 744     };
 745
 746     x++;
 747     r++;
 748     xoffset += d;
 749     if(__getbit2(select->hi,( (y << 3)+r))) k2++;
 750     if (r == 8) {
 751       r = 0;
 752       y++;
 753       z = select->hi[y];
 754     }
 755   };
 756
 757   if(x==select->m) return (uint)-1;
 758
 759
 760   int c = (y << 3)+r;
 761
 762   uint mask = ~(0xffu << (8*sizeof(uint) - 8));
 763   uint tmp = (select->hi[y] << ((8*sizeof(uint) - 8)+r)) | mask;
 764   uint c_incr = __builtin_clz(tmp);
 765   c += (c_incr & 7);
 766   if (c_incr == 8) {
 767     c += (8-r);
 768
 769     uchar * pp = &(select->hi[c >> 3]);
 770
 771     while(*pp == 0) {
 772       pp++;
 773       c += 8;
 774     };
 775
 776     c += __builtin_clz(*pp) - 24;
 777
 778   };
 779   c -= k2;
 780
 781   return __getbits(q,xoffset,d)+((c)<<d);
 782 }
 783
 784 int selects3_rank(selects3 *select, int i) {
 785   int d,x,w,y;
 786   int r,j;
 787   int z,ii;
 788   int xoffset;
 789   uint *q;
 790
 791   d = select->d;
 792   q = select->low;
 793
 794   ii = i >> d;
 795
 796   y = selectd2_select0(select->sd0, ii)+1;
 797
 798   x = y - ii;
 799
 800   j = i - (ii << d);
 801
 802   r = y & 7;
 803   y >>= 3;
 804   z = select->hi[y];
 805   xoffset = x * d;
 806
 807   while (1) {
 808     if (((z << r) & 0x80) == 0) return x;
 809
 810     w = __getbits(q, xoffset, d);
 811
 812     if (w >= j) return x + (w == j);
 813
 814     x++;
 815     r++;
 816     xoffset += d;
 817     if (r == 8) {
 818       r = 0;
 819       y++;
 820       z = select->hi[y];
 821     }
 822   }
 823 }
 824