#include <stdio.h>
#include <stdlib.h>
#include <string.h>
typedef unsigned int uint;
typedef unsigned char byte;
uint flen( FILE* f ) {
fseek( f, 0, SEEK_END );
uint len = ftell(f);
fseek( f, 0, SEEK_SET );
return len;
}
int obwt_bwts( byte* T, byte* U, int* A, int n );
int obwt_unbwts( byte* T, byte* U, int* LF, int n );
// Usage: BWTS c/d <inpfile> <outfile>
int main( int argc, char** argv ) {
uint mode = (argc<2) ? 0 : (argv[1][0]=='d');
FILE* f = argc<3 ? 0 : fopen( argv[2], "rb" ); if( f==0 ) f=stdin;
FILE* g = argc<4 ? 0 : fopen( argv[3], "wb" ); if( g==0 ) g=stdout;
uint inplen = flen(f);
byte* inpbuf = new byte[inplen]; if( inpbuf==0 ) return 1;
byte* outbuf = new byte[inplen]; if( outbuf==0 ) return 1;
int * ptrtbl = new int [inplen]; if( ptrtbl==0 ) return 1;
fread( inpbuf, 1,inplen, f );
mode ?
obwt_unbwts( inpbuf, outbuf, ptrtbl, inplen ) :
obwt_bwts( inpbuf, outbuf, ptrtbl, inplen );
fwrite( outbuf, 1,inplen, g );
}
#define chr(i) (cs==sizeof(int) ? ((int*)T)[i]:((byte*)T)[i])
#define isLMS(_a) ((bvec_get(F1,(_a))!=0) || ((bvec_get(F2,(_a))!=0) && (bvec_get(F2,(_a)-1)==0)))
int bvec_get( byte* B, int i ) {
return ( B[i >> 3] >> ( i& 7 ) )& 1;
}
void bvec_set( byte* B, int i, int b ) {
if( b ) {
B[i>>3] |= byte( 1U<<(i&7) );
} else {
B[i>>3] &= byte(~(1U<<(i&7)));
}
}
uint bvec_next( byte* B, uint i ) {
uint j;
uint c;
i += 1;
j = i >> 3;
c = B[j] >> ( i& 7 );
if( c==0 ) {
i += 8-(i&7);
j += 1;
for(; (c=B[j])==0; j++, i+=8 );
}
for(; (c&1)==0; i++, c>>=1 );
return i;
}
int bvec_prev( byte* B, int i ) {
int j;
uint c;
i -= 1;
if ( 0 <= i ) {
j = i >> 3;
c = ( B[j] << ( 7-( i& 7 ) ) )& 0xff;
if ( c == 0 ) {
i -= ( i& 7 ) + 1;
j -= 1;
for ( ; ( 0 <= j ) && ( ( c = B[j] ) == 0 ); --j, i -= 8 ){}
if ( c == 0 ) {
c = 128;
}
}
for ( ; ( c& 128 ) == 0; --i, c <<= 1 ){}
}
return i;
}
void getCounts( byte* T, int* C, int n, int k, int cs ) {
int i;
for( i=0; i<k; ++i ) C[i]=0;
for( i=0; i<n; ++i ) ++C[chr(i)];
}
void getBuckets( int* C, int* B, int k, int end ) {
int i, sum = 0;
if ( end ) {
for( i=0; i<k; ++i ) {
sum += C[i];
B[i] = sum;
}
} else {
for( i=0; i<k; ++i ) {
sum += C[i];
B[i] = sum - C[i];
}
}
}
void induceSA( byte* T, int* SA, int* C, int* B, byte* F1, byte* F2, byte* F3, int n, int k, int cs ) {
int i, j, p;
if( C==B ) getCounts( T, C, n, k, cs );
getBuckets( C, B, k, 0 );
p = bvec_prev( F3, n );
for( i=0; i<n; ++i ) {
while( (p>=0) && (B[chr(p)]==i) ) {
j = bvec_next( F1, p ) - 1;
SA[B[chr(j)]++] = j;
p = bvec_prev( F3, p );
}
j = SA[i];
if( j<0 ) continue;
if( bvec_get(F1,j)==0 ) j-=1;
else j=bvec_next(F1,j)-1;
if( bvec_get(F2,j)==0 ) SA[B[chr(j)]++]=j;
}
if( C==B ) getCounts( T, C, n, k, cs );
getBuckets( C, B, k, 1 );
for( i=n-1; i>=0; i-- ) {
j = SA[i];
if( j<0 ) continue;
if( bvec_get(F1,j)==0 ) if( bvec_get(F2,--j)!=0 ) SA[--B[chr(j)]]=j;
}
}
int lw_suffixsort( byte* T, int* SA, byte* F1, int fs, int n, int k, int cs ) {
byte* F2, * F3, * F4;
int* C, * B, * RA;
int i, j, c, m, p, q, plen, qlen, name;
int c0, c1;
int diff;
F2 = ( byte* )calloc( ( n + 7 ) / 8, sizeof(byte));
F3 = ( byte* )calloc( ( n / 8+1 ), sizeof(byte));
if ( ( F2 == NULL ) || ( F3 == NULL ) ) {
free( F2 );
free( F3 );
return - 2;
}
for ( p = 0; p < n; p = q ) {
q = bvec_next( F1, p );
for ( i = q - 2, c = 0, c1 = chr( q - 1 ), m = 0; p <= i; --i, c1 = c0 ) {
bvec_set( F2, i + 1, c );
if ( ( c0 = chr( i ) ) < ( c1 + c ) ) {
c = 1;
} else if ( c != 0 ) {
m += 1;
c = 0;
}
}
bvec_set( F2, p, 1 );
if ( ( m == 0 ) && ( c == 0 ) ) {
bvec_set( F3, p, 1 );
}
}
bvec_set( F3, n, 1 );
if ( k <= fs ) {
C = SA + n;
B = ( k <= ( fs - k ) ) ? C + k: C;
} else if ( ( C = B = ( int* )malloc( k* sizeof( int ) ) ) == NULL ) {
free( F2 );
free( F3 );
return - 2;
}
getCounts( T, C, n, k, cs );
getBuckets( C, B, k, 1 );
for( i=0; i<n; ++i ) {
SA[i] = - 1;
}
for( i=0; i<n; ++i ) {
if ( isLMS( i ) && ( bvec_get( F3, i ) == 0 ) ) {
SA[--B[chr( i )]] = i;
}
}
induceSA( T, SA, C, B, F1, F2, F3, n, k, cs );
if ( fs < k ) {
free( C );
}
for ( i = 0, m = 0; i<n; ++i ) {
p = SA[i];
if ( isLMS( p ) && ( bvec_get( F3, p ) == 0 ) ) {
SA[m++] = p;
}
}
for ( i = m; i<n; ++i ) {
SA[i] = 0;
}
for ( p = 0; p < n; p = q ) {
q = bvec_next( F1, p );
for ( i = q - 2, c = 0, c1 = chr( q - 1 ), j = q; p <= i; --i, c1 = c0 ) {
if ( ( c0 = chr( i ) ) < ( c1 + c ) ) {
c = 1;
} else if ( c != 0 ) {
SA[m + ( ( i + 1 ) >> 1 )] = j - i - 1;
j = i + 1;
c = 0;
}
}
if ( ( j < q ) || ( c != 0 ) ) {
SA[m + ( p >> 1 )] = j - p;
}
}
for ( i = 0, name = 0, q = n, qlen = 0; i<m; ++i ) {
p = SA[i], plen = SA[m + ( p >> 1 )], diff = 1;
if ( plen == qlen ) {
for ( j = 0; ( j < plen ) && ( chr( p + j ) == chr( q + j ) ); j++ ){}
if ( j == plen ) {
diff = 0;
}
}
if ( diff != 0 ) {
++name, q = p, qlen = plen;
}
SA[m + ( p >> 1 )] = name;
}
if ( name < m ) {
F4 = ( byte* )calloc( ( m / 8+1 ), sizeof(byte));
if ( F4 == NULL ) {
free( F2 );
free( F3 );
return - 1;
}
for ( i = 0, j = 0; i<n; ++i ) {
if ( isLMS( i ) && ( bvec_get( F3, i ) == 0 ) ) {
bvec_set( F4, j++, bvec_get( F1, i ) );
}
}
bvec_set( F4, m, 1 );
RA = SA + n + fs - m;
for ( i = n - 1, j = m - 1; m <= i; --i ) {
if ( SA[i] != 0 ) {
RA[j--] = SA[i] - 1;
}
}
if ( lw_suffixsort( ( byte* )RA, SA, F4, fs + n - m * 2, m, name, sizeof( int ) ) != 0 ) {
free( F2 );
free( F3 );
free( F4 );
return - 2;
}
free( F4 );
for ( i = 0, j = 0; i<n; ++i ) {
if ( isLMS( i ) && ( bvec_get( F3, i ) == 0 ) ) {
RA[j++] = i;
}
}
for( i=0; i<m; ++i ) {
SA[i] = RA[SA[i]];
}
}
if( k<=fs ) {
C = SA + n;
B = ( k <= ( fs - k ) ) ? C + k: C;
} else if ( ( C = B = ( int* )malloc( k* sizeof( int ) ) ) == NULL ) {
free( F2 );
free( F3 );
return -2;
}
getCounts( T, C, n, k, cs );
getBuckets( C, B, k, 1 );
for( i=m; i<n; i++ ) SA[i]=-1;
for ( i = m - 1; 0 <= i; --i ) {
j = SA[i], SA[i] = - 1;
SA[--B[chr( j )]] = j;
}
induceSA( T, SA, C, B, F1, F2, F3, n, k, cs );
if( fs<k ) free( C );
free( F2 );
free( F3 );
return 0;
}
int obwt_bwts( byte* T, byte* U, int* A, int n ) {
byte* F;
int i, j, k, p, q;
int c1, c2;
if( (T==0) || (U==0) || (A==0) || (n<0) ) return -1;
if( n<=1 ) { if( n==1 ) U[0] = T[0]; return 0; }
F = (byte*)calloc( (n/8+1), sizeof(byte) ); if( F==NULL ) return -2;
for( i=0, j=k=1; j<=n; ) {
for ( p = i, q = j;; ) {
c1 = ( p < n ) ? T[p]: - 1;
c2 = ( q < n ) ? T[q]: - 1;
++p, ++q;
if ( c1 != c2 ) {
break;
}
if ( p == j ) {
k = q;
p = i;
}
}
if ( c1 < c2 ) {
j = k = q;
} else {
bvec_set( F, i, 1 );
i = k;
j = ( k += 1 );
}
}
bvec_set( F, n, 1 );
if ( lw_suffixsort( T, A, F, 0, n, 256, sizeof(byte)) != 0 ) {
free( A );
free( F );
return - 2;
}
if ( T == U ) {
for( i=0; i<n; ++i ) {
p = A[i];
if ( bvec_get( F, p ) == 0 ) {
p -= 1;
} else {
p = bvec_next( F, p ) - 1;
}
c1 = T[i];
U[i] = byte( ( i <= p ) ? T[p]: A[p] );
A[i] = c1;
}
} else if ( ( ( T < U ) && ( U < ( T + n ) ) ) || ( ( U < T ) && ( T < ( U + n ) ) ) ) {
for( i=0; i<n; ++i ) {
p = A[i];
if ( bvec_get( F, p ) == 0 ) {
p -= 1;
} else {
p = bvec_next( F, p ) - 1;
}
A[i] = T[p];
}
for( i=0; i<n; ++i ) {
U[i] = ( byte )A[i];
}
} else {
for( i=0; i<n; ++i ) {
p = A[i];
if ( bvec_get( F, p ) == 0 ) {
p -= 1;
} else {
p = bvec_next( F, p ) - 1;
}
U[i] = T[p];
}
}
free( F );
return 0;
}
int obwt_unbwts( byte* T, byte* U, int* LF, int n ) {
int C[256];
byte* V;
int i, j, p, t;
if( (T==0) || (U==0) || (LF==0) || (n<0) ) return -1;
if( n<=1 ) { if( n==1 ) U[0]=T[0]; return 0; }
V = new byte[n];
if( LF && V ) {
for( i=0; i<256; ++i ) C[i] = 0;
for( i=0; i<n; ++i ) C[V[i] = T[i]] += 1;
for( i=0, j=0; i<256; ++i ) { t=C[i]; C[i]=j; j+=t; }
for( i=0; i<n; ++i ) LF[i] = C[V[i]]++;
for( i=0, j=n-1; j>=0; ++i ) {
if( LF[i]>=0 ) {
p = i;
do {
U[j--] = V[p];
t = LF[p];
LF[p]=-1;
p = t;
} while( LF[p]>=0 );
}
}
}
if( V ) delete[] V;
return 0;
}
#undef chr
#undef isLMS
I2luY2x1ZGUgPHN0ZGlvLmg+CiNpbmNsdWRlIDxzdGRsaWIuaD4KI2luY2x1ZGUgPHN0cmluZy5oPgoKdHlwZWRlZiB1bnNpZ25lZCBpbnQgIHVpbnQ7CnR5cGVkZWYgdW5zaWduZWQgY2hhciAgYnl0ZTsKCnVpbnQgZmxlbiggRklMRSogZiApIHsKICBmc2VlayggZiwgMCwgU0VFS19FTkQgKTsKICB1aW50IGxlbiA9IGZ0ZWxsKGYpOwogIGZzZWVrKCBmLCAwLCBTRUVLX1NFVCApOwogIHJldHVybiBsZW47Cn0KCmludCBvYnd0X2J3dHMoIGJ5dGUqIFQsIGJ5dGUqIFUsIGludCogQSwgaW50IG4gKTsKaW50IG9id3RfdW5id3RzKCBieXRlKiBULCBieXRlKiBVLCBpbnQqIExGLCBpbnQgbiApOwoKLy8gVXNhZ2U6IEJXVFMgYy9kIDxpbnBmaWxlPiA8b3V0ZmlsZT4gCmludCBtYWluKCBpbnQgYXJnYywgY2hhcioqIGFyZ3YgKSB7CgogIHVpbnQgbW9kZSA9IChhcmdjPDIpID8gMCA6IChhcmd2WzFdWzBdPT0nZCcpOwogIEZJTEUqIGYgPSBhcmdjPDMgPyAwIDogZm9wZW4oIGFyZ3ZbMl0sICJyYiIgKTsgaWYoIGY9PTAgKSBmPXN0ZGluOwogIEZJTEUqIGcgPSBhcmdjPDQgPyAwIDogZm9wZW4oIGFyZ3ZbM10sICJ3YiIgKTsgaWYoIGc9PTAgKSBnPXN0ZG91dDsKCiAgdWludCBpbnBsZW4gPSBmbGVuKGYpOwogIGJ5dGUqIGlucGJ1ZiA9IG5ldyBieXRlW2lucGxlbl07IGlmKCBpbnBidWY9PTAgKSByZXR1cm4gMTsKICBieXRlKiBvdXRidWYgPSBuZXcgYnl0ZVtpbnBsZW5dOyBpZiggb3V0YnVmPT0wICkgcmV0dXJuIDE7CiAgaW50ICogcHRydGJsID0gbmV3IGludCBbaW5wbGVuXTsgaWYoIHB0cnRibD09MCApIHJldHVybiAxOwoKICBmcmVhZCggaW5wYnVmLCAxLGlucGxlbiwgZiApOwoKICBtb2RlID8KICBvYnd0X3VuYnd0cyggaW5wYnVmLCBvdXRidWYsIHB0cnRibCwgaW5wbGVuICkgOgogIG9id3RfYnd0cyggaW5wYnVmLCBvdXRidWYsIHB0cnRibCwgaW5wbGVuICk7CgogIGZ3cml0ZSggb3V0YnVmLCAxLGlucGxlbiwgZyApOwoKfQoKI2RlZmluZSBjaHIoaSkgKGNzPT1zaXplb2YoaW50KSA/ICgoaW50KilUKVtpXTooKGJ5dGUqKVQpW2ldKQojZGVmaW5lIGlzTE1TKF9hKSAoKGJ2ZWNfZ2V0KEYxLChfYSkpIT0wKSB8fCAoKGJ2ZWNfZ2V0KEYyLChfYSkpIT0wKSAmJiAoYnZlY19nZXQoRjIsKF9hKS0xKT09MCkpKQoKaW50IGJ2ZWNfZ2V0KCBieXRlKiBCLCBpbnQgaSApIHsKICByZXR1cm4gKCBCW2kgPj4gM10gPj4gKCBpJiA3ICkgKSYgMTsKfQoKdm9pZCBidmVjX3NldCggYnl0ZSogQiwgaW50IGksIGludCBiICkgewogIGlmKCBiICkgewogICAgQltpPj4zXSB8PSBieXRlKCAgMVU8PChpJjcpICk7CiAgfSBlbHNlIHsKICAgIEJbaT4+M10gJj0gYnl0ZSh+KDFVPDwoaSY3KSkpOwogIH0KfQoKdWludCBidmVjX25leHQoIGJ5dGUqIEIsIHVpbnQgaSApIHsKICB1aW50IGo7CiAgdWludCBjOwogIGkgKz0gMTsKICBqID0gaSA+PiAzOwogIGMgPSBCW2pdID4+ICggaSYgNyApOwogIGlmKCBjPT0wICkgewogICAgaSArPSA4LShpJjcpOwogICAgaiArPSAxOwogICAgZm9yKDsgKGM9QltqXSk9PTA7IGorKywgaSs9OCApOwogIH0KICBmb3IoOyAoYyYxKT09MDsgaSsrLCBjPj49MSApOwogIHJldHVybiBpOwp9CgppbnQgYnZlY19wcmV2KCBieXRlKiBCLCBpbnQgaSApIHsKICBpbnQgajsKICB1aW50IGM7CiAgaSAtPSAxOwogIGlmICggMCA8PSBpICkgewogICAgaiA9IGkgPj4gMzsKICAgIGMgPSAoIEJbal0gPDwgKCA3LSggaSYgNyApICkgKSYgMHhmZjsKICAgIGlmICggYyA9PSAwICkgewogICAgICBpIC09ICggaSYgNyApICsgMTsKICAgICAgaiAtPSAxOwogICAgICBmb3IgKCA7ICggMCA8PSBqICkgJiYgKCAoIGMgPSBCW2pdICkgPT0gMCApOyAtLWosIGkgLT0gOCApe30KICAgICAgaWYgKCBjID09IDAgKSB7CiAgICAgICAgYyA9IDEyODsKICAgICAgfQogICAgfQogICAgZm9yICggOyAoIGMmIDEyOCApID09IDA7IC0taSwgYyA8PD0gMSApe30KICB9CiAgcmV0dXJuIGk7Cn0KCgp2b2lkIGdldENvdW50cyggYnl0ZSogVCwgaW50KiBDLCBpbnQgbiwgaW50IGssIGludCBjcyApIHsKICBpbnQgaTsKICBmb3IoIGk9MDsgaTxrOyArK2kgKSBDW2ldPTA7CiAgZm9yKCBpPTA7IGk8bjsgKytpICkgKytDW2NocihpKV07Cn0KCnZvaWQgZ2V0QnVja2V0cyggaW50KiBDLCBpbnQqIEIsIGludCBrLCBpbnQgZW5kICkgewogIGludCBpLCBzdW0gPSAwOwogIGlmICggZW5kICkgewogICAgZm9yKCBpPTA7IGk8azsgKytpICkgewogICAgICBzdW0gKz0gQ1tpXTsKICAgICAgQltpXSA9IHN1bTsKICAgIH0KICB9IGVsc2UgewogICAgZm9yKCBpPTA7IGk8azsgKytpICkgewogICAgICBzdW0gKz0gQ1tpXTsKICAgICAgQltpXSA9IHN1bSAtIENbaV07CiAgICB9CiAgfQp9CgoKdm9pZCBpbmR1Y2VTQSggYnl0ZSogVCwgaW50KiBTQSwgaW50KiBDLCBpbnQqIEIsIGJ5dGUqIEYxLCBieXRlKiBGMiwgYnl0ZSogRjMsIGludCBuLCBpbnQgaywgaW50IGNzICkgewogIGludCBpLCBqLCBwOwoKICBpZiggQz09QiApIGdldENvdW50cyggVCwgQywgbiwgaywgY3MgKTsKICBnZXRCdWNrZXRzKCBDLCBCLCBrLCAwICk7CgogIHAgPSBidmVjX3ByZXYoIEYzLCBuICk7CgogIGZvciggaT0wOyBpPG47ICsraSApIHsKIAogICAgd2hpbGUoIChwPj0wKSAmJiAoQltjaHIocCldPT1pKSApIHsKICAgICAgaiA9IGJ2ZWNfbmV4dCggRjEsIHAgKSAtIDE7CiAgICAgIFNBW0JbY2hyKGopXSsrXSA9IGo7CiAgICAgIHAgPSBidmVjX3ByZXYoIEYzLCBwICk7CiAgICB9CiAgICBqID0gU0FbaV07CiAgICBpZiggajwwICkgY29udGludWU7CgogICAgaWYoIGJ2ZWNfZ2V0KEYxLGopPT0wICkgai09MTsKICAgIGVsc2Ugaj1idmVjX25leHQoRjEsaiktMTsKCiAgICBpZiggYnZlY19nZXQoRjIsaik9PTAgKSBTQVtCW2NocihqKV0rK109ajsKICB9CgogIGlmKCBDPT1CICkgZ2V0Q291bnRzKCBULCBDLCBuLCBrLCBjcyApOwoKICBnZXRCdWNrZXRzKCBDLCBCLCBrLCAxICk7CgogIGZvciggaT1uLTE7IGk+PTA7IGktLSApIHsKICAgIGogPSBTQVtpXTsKICAgIGlmKCBqPDAgKSBjb250aW51ZTsKICAgIGlmKCBidmVjX2dldChGMSxqKT09MCApIGlmKCBidmVjX2dldChGMiwtLWopIT0wICkgU0FbLS1CW2NocihqKV1dPWo7CiAgfQp9CgppbnQgbHdfc3VmZml4c29ydCggYnl0ZSogVCwgaW50KiBTQSwgYnl0ZSogRjEsIGludCBmcywgaW50IG4sIGludCBrLCBpbnQgY3MgKSB7CgogIGJ5dGUqIEYyLCAqIEYzLCAqIEY0OwogIGludCogQywgKiBCLCAqIFJBOwogIGludCBpLCBqLCBjLCBtLCBwLCBxLCBwbGVuLCBxbGVuLCBuYW1lOwogIGludCBjMCwgYzE7CiAgaW50IGRpZmY7CgogIEYyID0gKCBieXRlKiApY2FsbG9jKCAoIG4gKyA3ICkgLyA4LCBzaXplb2YoYnl0ZSkpOwogIEYzID0gKCBieXRlKiApY2FsbG9jKCAoIG4gLyA4KzEgKSwgc2l6ZW9mKGJ5dGUpKTsKICBpZiAoICggRjIgPT0gTlVMTCApIHx8ICggRjMgPT0gTlVMTCApICkgewogICAgZnJlZSggRjIgKTsKICAgIGZyZWUoIEYzICk7CiAgICByZXR1cm4gIC0gMjsKICB9CiAgZm9yICggcCA9IDA7IHAgPCBuOyBwID0gcSApIHsKICAgIHEgPSBidmVjX25leHQoIEYxLCBwICk7CiAgICBmb3IgKCBpID0gcSAtIDIsIGMgPSAwLCBjMSA9IGNociggcSAtIDEgKSwgbSA9IDA7IHAgPD0gaTsgLS1pLCBjMSA9IGMwICkgewogICAgICBidmVjX3NldCggRjIsIGkgKyAxLCBjICk7CiAgICAgIGlmICggKCBjMCA9IGNociggaSApICkgPCAoIGMxICsgYyApICkgewogICAgICAgIGMgPSAxOwogICAgICB9IGVsc2UgaWYgKCBjICE9IDAgKSB7CiAgICAgICAgbSArPSAxOwogICAgICAgIGMgPSAwOwogICAgICB9CiAgICB9CiAgICBidmVjX3NldCggRjIsIHAsIDEgKTsKICAgIGlmICggKCBtID09IDAgKSAmJiAoIGMgPT0gMCApICkgewogICAgICBidmVjX3NldCggRjMsIHAsIDEgKTsKICAgIH0KICB9CiAgYnZlY19zZXQoIEYzLCBuLCAxICk7CgogIGlmICggayA8PSBmcyApIHsKICAgIEMgPSBTQSArIG47CiAgICBCID0gKCBrIDw9ICggZnMgLSBrICkgKSA/IEMgKyBrOiBDOwogIH0gZWxzZSBpZiAoICggQyA9IEIgPSAoIGludCogKW1hbGxvYyggayogc2l6ZW9mKCBpbnQgKSApICkgPT0gTlVMTCApIHsKICAgIGZyZWUoIEYyICk7CiAgICBmcmVlKCBGMyApOwogICAgcmV0dXJuICAtIDI7CiAgfQogIGdldENvdW50cyggVCwgQywgbiwgaywgY3MgKTsKICBnZXRCdWNrZXRzKCBDLCBCLCBrLCAxICk7CiAgZm9yKCBpPTA7IGk8bjsgKytpICkgewogICAgU0FbaV0gPSAgLSAxOwogIH0KICBmb3IoIGk9MDsgaTxuOyArK2kgKSB7CiAgICBpZiAoIGlzTE1TKCBpICkgJiYgKCBidmVjX2dldCggRjMsIGkgKSA9PSAwICkgKSB7CiAgICAgIFNBWy0tQltjaHIoIGkgKV1dID0gaTsKICAgIH0KICB9CiAgaW5kdWNlU0EoIFQsIFNBLCBDLCBCLCBGMSwgRjIsIEYzLCBuLCBrLCBjcyApOwogIGlmICggZnMgPCBrICkgewogICAgZnJlZSggQyApOwogIH0KCiAgZm9yICggaSA9IDAsIG0gPSAwOyBpPG47ICsraSApIHsKICAgIHAgPSBTQVtpXTsKICAgIGlmICggaXNMTVMoIHAgKSAmJiAoIGJ2ZWNfZ2V0KCBGMywgcCApID09IDAgKSApIHsKICAgICAgU0FbbSsrXSA9IHA7CiAgICB9CiAgfQogIGZvciAoIGkgPSBtOyBpPG47ICsraSApIHsKICAgIFNBW2ldID0gMDsKICB9CgogIGZvciAoIHAgPSAwOyBwIDwgbjsgcCA9IHEgKSB7CiAgICBxID0gYnZlY19uZXh0KCBGMSwgcCApOwogICAgZm9yICggaSA9IHEgLSAyLCBjID0gMCwgYzEgPSBjaHIoIHEgLSAxICksIGogPSBxOyBwIDw9IGk7IC0taSwgYzEgPSBjMCApIHsKICAgICAgaWYgKCAoIGMwID0gY2hyKCBpICkgKSA8ICggYzEgKyBjICkgKSB7CiAgICAgICAgYyA9IDE7CiAgICAgIH0gZWxzZSBpZiAoIGMgIT0gMCApIHsKICAgICAgICBTQVttICsgKCAoIGkgKyAxICkgPj4gMSApXSA9IGogLSBpIC0gMTsKICAgICAgICBqID0gaSArIDE7CiAgICAgICAgYyA9IDA7CiAgICAgIH0KICAgIH0KICAgIGlmICggKCBqIDwgcSApIHx8ICggYyAhPSAwICkgKSB7CiAgICAgIFNBW20gKyAoIHAgPj4gMSApXSA9IGogLSBwOwogICAgfQogIH0KCiAgZm9yICggaSA9IDAsIG5hbWUgPSAwLCBxID0gbiwgcWxlbiA9IDA7IGk8bTsgKytpICkgewogICAgcCA9IFNBW2ldLCBwbGVuID0gU0FbbSArICggcCA+PiAxICldLCBkaWZmID0gMTsKICAgIGlmICggcGxlbiA9PSBxbGVuICkgewogICAgICBmb3IgKCBqID0gMDsgKCBqIDwgcGxlbiApICYmICggY2hyKCBwICsgaiApID09IGNociggcSArIGogKSApOyBqKysgKXt9CiAgICAgIGlmICggaiA9PSBwbGVuICkgewogICAgICAgIGRpZmYgPSAwOwogICAgICB9CiAgICB9CiAgICBpZiAoIGRpZmYgIT0gMCApIHsKICAgICAgKytuYW1lLCBxID0gcCwgcWxlbiA9IHBsZW47CiAgICB9CiAgICBTQVttICsgKCBwID4+IDEgKV0gPSBuYW1lOwogIH0KCiAgaWYgKCBuYW1lIDwgbSApIHsKICAgIEY0ID0gKCBieXRlKiApY2FsbG9jKCAoIG0gLyA4KzEgKSwgc2l6ZW9mKGJ5dGUpKTsKICAgIGlmICggRjQgPT0gTlVMTCApIHsKICAgICAgZnJlZSggRjIgKTsKICAgICAgZnJlZSggRjMgKTsKICAgICAgcmV0dXJuICAtIDE7CiAgICB9CiAgICBmb3IgKCBpID0gMCwgaiA9IDA7IGk8bjsgKytpICkgewogICAgICBpZiAoIGlzTE1TKCBpICkgJiYgKCBidmVjX2dldCggRjMsIGkgKSA9PSAwICkgKSB7CiAgICAgICAgYnZlY19zZXQoIEY0LCBqKyssIGJ2ZWNfZ2V0KCBGMSwgaSApICk7CiAgICAgIH0KICAgIH0KICAgIGJ2ZWNfc2V0KCBGNCwgbSwgMSApOwogICAgUkEgPSBTQSArIG4gKyBmcyAtIG07CiAgICBmb3IgKCBpID0gbiAtIDEsIGogPSBtIC0gMTsgbSA8PSBpOyAtLWkgKSB7CiAgICAgIGlmICggU0FbaV0gIT0gMCApIHsKICAgICAgICBSQVtqLS1dID0gU0FbaV0gLSAxOwogICAgICB9CiAgICB9CiAgICBpZiAoIGx3X3N1ZmZpeHNvcnQoICggYnl0ZSogKVJBLCBTQSwgRjQsIGZzICsgbiAtIG0gKiAyLCBtLCBuYW1lLCBzaXplb2YoIGludCApICkgIT0gMCApIHsKICAgICAgZnJlZSggRjIgKTsKICAgICAgZnJlZSggRjMgKTsKICAgICAgZnJlZSggRjQgKTsKICAgICAgcmV0dXJuICAtIDI7CiAgICB9CiAgICBmcmVlKCBGNCApOwogICAgZm9yICggaSA9IDAsIGogPSAwOyBpPG47ICsraSApIHsKICAgICAgaWYgKCBpc0xNUyggaSApICYmICggYnZlY19nZXQoIEYzLCBpICkgPT0gMCApICkgewogICAgICAgIFJBW2orK10gPSBpOwogICAgICB9CiAgICB9CiAgICBmb3IoIGk9MDsgaTxtOyArK2kgKSB7CiAgICAgIFNBW2ldID0gUkFbU0FbaV1dOwogICAgfQogIH0KCiAgaWYoIGs8PWZzICkgewogICAgQyA9IFNBICsgbjsKICAgIEIgPSAoIGsgPD0gKCBmcyAtIGsgKSApID8gQyArIGs6IEM7CiAgfSBlbHNlIGlmICggKCBDID0gQiA9ICggaW50KiApbWFsbG9jKCBrKiBzaXplb2YoIGludCApICkgKSA9PSBOVUxMICkgewogICAgZnJlZSggRjIgKTsKICAgIGZyZWUoIEYzICk7CiAgICByZXR1cm4gLTI7CiAgfQoKICBnZXRDb3VudHMoIFQsIEMsIG4sIGssIGNzICk7CiAgZ2V0QnVja2V0cyggQywgQiwgaywgMSApOwoKICBmb3IoIGk9bTsgaTxuOyBpKysgKSBTQVtpXT0tMTsKCiAgZm9yICggaSA9IG0gLSAxOyAwIDw9IGk7IC0taSApIHsKICAgIGogPSBTQVtpXSwgU0FbaV0gPSAgLSAxOwogICAgU0FbLS1CW2NociggaiApXV0gPSBqOwogIH0KICBpbmR1Y2VTQSggVCwgU0EsIEMsIEIsIEYxLCBGMiwgRjMsIG4sIGssIGNzICk7CgogIGlmKCBmczxrICkgZnJlZSggQyApOwoKICBmcmVlKCBGMiApOwogIGZyZWUoIEYzICk7CiAgcmV0dXJuIDA7Cn0KCmludCBvYnd0X2J3dHMoIGJ5dGUqIFQsIGJ5dGUqIFUsIGludCogQSwgaW50IG4gKSB7CiAgYnl0ZSogRjsKICBpbnQgaSwgaiwgaywgcCwgcTsKICBpbnQgYzEsIGMyOwoKICBpZiggKFQ9PTApIHx8IChVPT0wKSB8fCAoQT09MCkgfHwgKG48MCkgKSByZXR1cm4gLTE7CgogIGlmKCBuPD0xICkgeyBpZiggbj09MSApIFVbMF0gPSBUWzBdOyByZXR1cm4gMDsgfQoKICBGID0gKGJ5dGUqKWNhbGxvYyggKG4vOCsxKSwgc2l6ZW9mKGJ5dGUpICk7IGlmKCBGPT1OVUxMICkgcmV0dXJuIC0yOwoKICBmb3IoIGk9MCwgaj1rPTE7IGo8PW47ICkgewogICAgZm9yICggcCA9IGksIHEgPSBqOzsgKSB7CiAgICAgIGMxID0gKCBwIDwgbiApID8gVFtwXTogIC0gMTsKICAgICAgYzIgPSAoIHEgPCBuICkgPyBUW3FdOiAgLSAxOwogICAgICArK3AsICsrcTsKICAgICAgaWYgKCBjMSAhPSBjMiApIHsKICAgICAgICBicmVhazsKICAgICAgfQogICAgICBpZiAoIHAgPT0gaiApIHsKICAgICAgICBrID0gcTsKICAgICAgICBwID0gaTsKICAgICAgfQogICAgfQogICAgaWYgKCBjMSA8IGMyICkgewogICAgICBqID0gayA9IHE7CiAgICB9IGVsc2UgewogICAgICBidmVjX3NldCggRiwgaSwgMSApOwogICAgICBpID0gazsKICAgICAgaiA9ICggayArPSAxICk7CiAgICB9CiAgfQogIGJ2ZWNfc2V0KCBGLCBuLCAxICk7CgogIGlmICggbHdfc3VmZml4c29ydCggVCwgQSwgRiwgMCwgbiwgMjU2LCBzaXplb2YoYnl0ZSkpICE9IDAgKSB7CiAgICBmcmVlKCBBICk7CiAgICBmcmVlKCBGICk7CiAgICByZXR1cm4gIC0gMjsKICB9CgogIGlmICggVCA9PSBVICkgewogICAgZm9yKCBpPTA7IGk8bjsgKytpICkgewogICAgICBwID0gQVtpXTsKICAgICAgaWYgKCBidmVjX2dldCggRiwgcCApID09IDAgKSB7CiAgICAgICAgcCAtPSAxOwogICAgICB9IGVsc2UgewogICAgICAgIHAgPSBidmVjX25leHQoIEYsIHAgKSAtIDE7CiAgICAgIH0KICAgICAgYzEgPSBUW2ldOwogICAgICBVW2ldID0gYnl0ZSggKCBpIDw9IHAgKSA/IFRbcF06IEFbcF0gKTsKICAgICAgQVtpXSA9IGMxOwogICAgfQogIH0gZWxzZSBpZiAoICggKCBUIDwgVSApICYmICggVSA8ICggVCArIG4gKSApICkgfHwgKCAoIFUgPCBUICkgJiYgKCBUIDwgKCBVICsgbiApICkgKSApIHsKICAgIGZvciggaT0wOyBpPG47ICsraSApIHsKICAgICAgcCA9IEFbaV07CiAgICAgIGlmICggYnZlY19nZXQoIEYsIHAgKSA9PSAwICkgewogICAgICAgIHAgLT0gMTsKICAgICAgfSBlbHNlIHsKICAgICAgICBwID0gYnZlY19uZXh0KCBGLCBwICkgLSAxOwogICAgICB9CiAgICAgIEFbaV0gPSBUW3BdOwogICAgfQogICAgZm9yKCBpPTA7IGk8bjsgKytpICkgewogICAgICBVW2ldID0gKCBieXRlIClBW2ldOwogICAgfQogIH0gZWxzZSB7CiAgICBmb3IoIGk9MDsgaTxuOyArK2kgKSB7CiAgICAgIHAgPSBBW2ldOwogICAgICBpZiAoIGJ2ZWNfZ2V0KCBGLCBwICkgPT0gMCApIHsKICAgICAgICBwIC09IDE7CiAgICAgIH0gZWxzZSB7CiAgICAgICAgcCA9IGJ2ZWNfbmV4dCggRiwgcCApIC0gMTsKICAgICAgfQogICAgICBVW2ldID0gVFtwXTsKICAgIH0KICB9CgogIGZyZWUoIEYgKTsKCiAgcmV0dXJuIDA7Cn0KCgppbnQgb2J3dF91bmJ3dHMoIGJ5dGUqIFQsIGJ5dGUqIFUsIGludCogTEYsIGludCBuICkgewogIGludCBDWzI1Nl07CiAgYnl0ZSogVjsKICBpbnQgaSwgaiwgcCwgdDsKCiAgaWYoIChUPT0wKSB8fCAoVT09MCkgfHwgKExGPT0wKSB8fCAobjwwKSApIHJldHVybiAtMTsKCiAgaWYoIG48PTEgKSB7IGlmKCBuPT0xICkgVVswXT1UWzBdOyByZXR1cm4gMDsgfQoKICBWICA9IG5ldyBieXRlW25dOwoKICBpZiggTEYgJiYgViApIHsKCiAgICBmb3IoIGk9MDsgaTwyNTY7ICsraSApIENbaV0gPSAwOwoKICAgIGZvciggaT0wOyBpPG47ICsraSApIENbVltpXSA9IFRbaV1dICs9IDE7CgogICAgZm9yKCBpPTAsIGo9MDsgaTwyNTY7ICsraSApIHsgdD1DW2ldOyBDW2ldPWo7IGorPXQ7IH0KCiAgICBmb3IoIGk9MDsgaTxuOyArK2kgKSBMRltpXSA9IENbVltpXV0rKzsKCiAgICBmb3IoIGk9MCwgaj1uLTE7IGo+PTA7ICsraSApIHsKICAgICAgaWYoIExGW2ldPj0wICkgewogICAgICAgIHAgPSBpOwogICAgICAgIGRvIHsKICAgICAgICAgIFVbai0tXSA9IFZbcF07CiAgICAgICAgICB0ID0gTEZbcF07CiAgICAgICAgICBMRltwXT0tMTsKICAgICAgICAgIHAgPSB0OwogICAgICAgIH0gd2hpbGUoIExGW3BdPj0wICk7CiAgICAgIH0KICAgIH0KICB9CgogIGlmKCBWICkgZGVsZXRlW10gVjsKCiAgcmV0dXJuIDA7Cn0KCiN1bmRlZiBjaHIKI3VuZGVmIGlzTE1TCgo=