JWS-112 Bumping version of Mafft to version 7.310.

[jabaws.git] / binaries / src / mafft / core / mafft-distance.c

diff --git a/binaries/src/mafft/core/mafft-distance.c b/binaries/src/mafft/core/mafft-distance.c
index 3ece8a8..cb4e7d5 100644 (file)
--- a/binaries/src/mafft/core/mafft-distance.c
+++ b/binaries/src/mafft/core/mafft-distance.c
@@ -6,9 +6,8 @@
 
 #define END_OF_VEC -1
 
-static int maxl;
-static int tsize;
-static float lenfaca, lenfacb, lenfacc, lenfacd;
+static char outputformat;
+static int nadd;
 #define PLENFACA 0.01
 #define PLENFACB 10000
 #define PLENFACC 10000
@@ -23,9 +22,11 @@ void arguments( int argc, char *argv[] )
        int c;
 
        inputfile = NULL;
+       outputformat = 's';
        scoremtx = 1;
        nblosum = 62;
        dorp = NOTSPECIFIED;
+       nadd = 0;
        alg = 'X';
 
     while( --argc > 0 && (*++argv)[0] == '-' )
@@ -39,6 +40,18 @@ void arguments( int argc, char *argv[] )
                                        fprintf( stderr, "inputfile = %s\n", inputfile );
                                        --argc;
                                        goto nextoption;
+                               case 'I':
+                                       nadd = myatoi(*++argv);
+                                       if( nadd == 0 )
+                                       {
+                                               fprintf( stderr, "nadd = %d?\n", nadd );
+                                               exit( 1 );
+                                       }
+                                       --argc;
+                                       goto nextoption;
+                               case 'p':
+                                       outputformat = 'p';
+                                       break;
                                case 'D':
                                        dorp = 'd';
                                        break;
@@ -117,7 +130,7 @@ void makecompositiontable_p( short *table, int *pointt )
                table[point]++;
 }
 
-int commonsextet_p( short *table, int *pointt )
+static int localcommonsextet_p( short *table, int *pointt )
 {
        int value = 0;
        short tmp;
@@ -217,7 +230,7 @@ void makepointtable( int *pointt, int *n )
 
 int main( int argc, char **argv )
 {
-       int i, j;
+       int i, j, initj;
        FILE *infp;
        char **seq;
        int *grpseq;
@@ -226,14 +239,16 @@ int main( int argc, char **argv )
        static char **name;
        static int *nlen;
        double *mtxself;
-       float score;
+       double score;
        static short *table1;
-       float longer, shorter;
-       float lenfac;
-       float bunbo;
+       double longer, shorter;
+       double lenfac;
+       double bunbo;
+       int norg;
 
        arguments( argc, argv );
 
+
        if( inputfile )
        {
                infp = fopen( inputfile, "r" );
@@ -278,6 +293,10 @@ int main( int argc, char **argv )
 
        constants( njob, seq );
 
+
+       if( nadd ) outputformat = 's';
+       norg = njob - nadd;
+
        if( dorp == 'd' ) tsize = (int)pow( 4, 6 );
        else              tsize = (int)pow( 6, 6 );
 
@@ -325,15 +344,17 @@ int main( int argc, char **argv )
                if( !table1 ) ErrorExit( "Cannot allocate table1\n" );
                makecompositiontable_p( table1, pointt[i] );
 
-               score = commonsextet_p( table1, pointt[i] );
+               score = localcommonsextet_p( table1, pointt[i] );
                mtxself[i] = score;
                free( table1 );
        }
 
        fprintf( stderr, "done.\n" );
        fprintf( stderr, "\nCalculating i-j scores ... \n" );
-       for( i=0; i<njob; i++ )
+       if( outputformat == 'p' ) fprintf( stdout, "%-5d", njob );
+       for( i=0; i<norg; i++ )
        {
+               if( outputformat == 'p' ) fprintf( stdout, "\n%-9d ", i+1 );
                table1 = (short *)calloc( tsize, sizeof( short ) );
                if( !table1 ) ErrorExit( "Cannot allocate table1\n" );
                if( i % 10 == 0 )
@@ -343,30 +364,50 @@ int main( int argc, char **argv )
                makecompositiontable_p( table1, pointt[i] );
 
 
-        for( j=i+1; j<njob; j++ ) 
-        {
+               if( nadd == 0 )
+               {
+                       if( outputformat == 'p' ) initj = 0;
+                       else initj = i+1;
+               }
+               else 
+               {
+                       initj = norg;
+               }
+               for( j=initj; j<njob; j++ ) 
+               {
                        if( nlen[i] > nlen[j] )
                        {
-                               longer=(float)nlen[i];
-                               shorter=(float)nlen[j];
+                               longer=(double)nlen[i];
+                               shorter=(double)nlen[j];
                        }
                        else
                        {
-                               longer=(float)nlen[j];
-                               shorter=(float)nlen[i];
+                               longer=(double)nlen[j];
+                               shorter=(double)nlen[i];
                        }
 //                     lenfac = 3.0 / ( LENFACA + LENFACB / ( longer + LENFACC ) + shorter / longer * LENFACD );
                        lenfac = 1.0 / ( shorter / longer * lenfacd + lenfacb / ( longer + lenfacc ) + lenfaca );
 //                     lenfac = 1.0;
 //                     fprintf( stderr, "lenfac = %f (%.0f,%.0f)\n", lenfac, longer, shorter );
-                       score = commonsextet_p( table1, pointt[j] );
+                       score = localcommonsextet_p( table1, pointt[j] );
                        bunbo = MIN( mtxself[i], mtxself[j] );
-                       if( bunbo == 0.0 )
-                               fprintf( stdout, "%d-%d d=%4.2f l=%d,%d\n", i+1, j+1, 1.0, nlen[i], nlen[j] );
+                       if( outputformat == 'p' )
+                       {
+                               if( bunbo == 0.0 )
+                                       fprintf( stdout, " %8.6f", 1.0 );
+                               else
+                                       fprintf( stdout, " %8.6f", ( 1.0 - score / bunbo ) * lenfac );
+                               if( j % 7 == 6 ) fprintf( stdout, "\n" );
+                       }
                        else
-                               fprintf( stdout, "%d-%d d=%4.2f l=%d,%d\n", i+1, j+1, ( 1.0 - score / bunbo ) * lenfac, nlen[i], nlen[j] );
+                       {
+                               if( bunbo == 0.0 )
+                                       fprintf( stdout, "%d-%d d=%4.2f l=%d,%d\n", i+1, j+1, 1.0, nlen[i], nlen[j] );
+                               else
+                                       fprintf( stdout, "%d-%d d=%4.2f l=%d,%d\n", i+1, j+1, ( 1.0 - score / bunbo ) * lenfac, nlen[i], nlen[j] );
+                       }
 //                     fprintf( stderr, "##### mtx = %f, mtx[i][0]=%f, mtx[j][0]=%f, bunbo=%f\n", mtx[i][j-i], mtx[i][0], mtx[j][0], bunbo );
-//          score = (double)commonsextet_p( table1, pointt[j] );
+//          score = (double)localcommonsextet_p( table1, pointt[j] );
 //                     fprintf( stdout, "%d-%d d=%4.2f l=%d,%d\n", i+1, j+1, ( 1.0 - score / MIN( mtxself[i], mtxself[j] ) ) * 3, nlen[i], nlen[j] );
 
 
@@ -375,6 +416,7 @@ int main( int argc, char **argv )
        }
                
        fprintf( stderr, "\n" );
+       if( outputformat == 'p' ) fprintf( stdout, "\n" );
        SHOWVERSION;
        exit( 0 );
 }