[med-svn] [Git][med-team/mafft][master] 4 commits: New upstream version 7.464
Steffen Möller
gitlab at salsa.debian.org
Sun May 3 18:23:28 BST 2020
Steffen Möller pushed to branch master at Debian Med / mafft
Commits:
a9b2d97a by Steffen Moeller at 2020-05-03T19:21:44+02:00
New upstream version 7.464
- - - - -
cb61b047 by Steffen Moeller at 2020-05-03T19:21:44+02:00
routine-update: New upstream version
- - - - -
3d0ae048 by Steffen Moeller at 2020-05-03T19:21:46+02:00
Update upstream source from tag 'upstream/7.464'
Update to upstream version '7.464'
with Debian dir 52aa8526d15ce73fa73c65bafa092e3cd0525b14
- - - - -
44912116 by Steffen Moeller at 2020-05-03T19:22:24+02:00
routine-update: Ready to upload to unstable
- - - - -
17 changed files:
- MPI/mltaln9_mpi.c
- core/Falign.c
- core/Galign11.c
- core/Makefile
- − core/README
- core/Salignmm.c
- core/addsingle.c
- core/blosum.c
- core/disttbfast.c
- core/functions.h
- core/io.c
- core/mafft.tmpl
- core/mltaln.h
- core/mltaln9.c
- core/tbfast.c
- debian/changelog
- readme
Changes:
=====================================
MPI/mltaln9_mpi.c
=====================================
@@ -1450,6 +1450,56 @@ static void msaresetnearest( int nseq, Bchain *acpt, double **distfrompt, double
}
#endif
+static int getdensest( int *m, double *d )
+{
+ int i;
+ double dmax = -100.0;
+ int pmax = -1;
+ for( i=0; m[i]>-1; i++ )
+ {
+ if( d[m[i]] > dmax )
+ {
+ dmax = d[m[i]];
+ pmax = m[i];
+ }
+ }
+ return( pmax );
+}
+
+
+static void setdensity( int nseq, Bchain *acpt, double **eff, double *density, int pos )
+{
+ int j;
+ double tmpdouble;
+// double **effptpt;
+ Bchain *acptj;
+
+// printf( "[%d], %f, dist=%d ->", pos, *mindisfrompt, *nearestpt );
+
+// if( (acpt+pos)->next ) effpt = eff[pos]+(acpt+pos)->next->pos-pos;
+
+ tmpdouble = 0.0;
+// for( j=pos+1; j<nseq; j++ )
+ for( acptj=(acpt+pos)->next; acptj!=NULL; acptj=acptj->next )
+ {
+ j = acptj->pos;
+ if( eff[pos][j-pos] < 1.0 )
+ tmpdouble += (2.0-eff[pos][j-pos]);
+ }
+// effptpt = eff;
+// for( j=0; j<pos; j++ )
+ for( acptj=acpt; (acptj&&acptj->pos!=pos); acptj=acptj->next )
+ {
+ j = acptj->pos;
+ if( eff[j][pos-j] < 1.0 )
+ tmpdouble += (2.0-eff[j][pos-j]);
+ }
+
+ *density = tmpdouble;
+// printf( "p=%d, d=%f \n", pos, *density );
+}
+
+
static void setnearest( int nseq, Bchain *acpt, double **eff, double *mindisfrompt, int *nearestpt, int pos )
{
int j;
@@ -7540,7 +7590,7 @@ void compacttree_memsaveselectable( int nseq, double **partmtx, int *nearest, do
if( result ) free( result );
}
-void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree )
+void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree, int treeout )
{
int i, j, k, miniim, maxiim, minijm, maxijm;
@@ -7566,6 +7616,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
FILE *fp;
double (*clusterfuncpt[1])(double,double);
char namec;
+ double *density;
sueff1 = 1 - (double)sueff_global;
@@ -7595,6 +7646,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
nametmp = AllocateCharVec( 1000 ); // nagasugi
// tree = AllocateCharMtx( njob, njob*600 );
tree = AllocateCharMtx( njob, 0 );
+ if( treeout == 2 ) density = AllocateDoubleVec( njob );
}
@@ -7637,6 +7689,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
ac[nseq-1].next = NULL;
for( i=0; i<nseq; i++ ) setnearest( nseq, ac, eff, mindisfrom+i, nearest+i, i ); // muscle
+ if( treeout == 2 ) for( i=0; i<nseq; i++ ) setdensity( nseq, ac, eff, density+i, i );
for( i=0; i<nseq; i++ ) tmptmplen[i] = 0.0;
@@ -7896,18 +7949,36 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
#endif
-#if 0
- printf( "\nooSTEP-%03d:\n", k+1 );
- printf( "len0 = %f\n", len[k][0] );
- for( i=0; topol[k][0][i]>-1; i++ ) printf( " %03d", topol[k][0][i]+1 );
- printf( "\n" );
- printf( "len1 = %f\n", len[k][1] );
- for( i=0; topol[k][1][i]>-1; i++ ) printf( " %03d", topol[k][1][i]+1 );
- printf( "\n" );
-#endif
}
fp = fopen( "infile.tree", "w" );
fprintf( fp, "%s;\n", treetmp );
+ if( treeout == 2 )
+ {
+ int *mem = calloc( sizeof( int ), nseq );
+ fprintf( fp, "\nDensity:" );
+ for( k=0; k<nseq; k++ ) fprintf( fp, "\nSequence %d, %7.4f", k+1, density[k] );
+
+ fprintf( fp, "\n\nNode info:" );
+ for( k=0; k<nseq-1; k++ )
+ {
+ if( dep )
+ fprintf( fp, "\nNode %d, Height=%f\n", k+1, dep[k].distfromtip );
+// fprintf( fp, "len0 = %f\n", len[k][0] );
+ topolorderz( mem, topol, dep, k, 0 );
+// for( i=0; topol[k][0][i]>-1; i++ ) fprintf( fp, " %03d", topol[k][0][i]+1 );
+ fprintf( fp, "%d:", getdensest( mem, density )+1 );
+ for( i=0; mem[i]>-1; i++ ) fprintf( fp, " %d", mem[i]+1 );
+ fprintf( fp, "\n" );
+
+ topolorderz( mem, topol, dep, k, 1 );
+// fprintf( fp, "len1 = %f\n", len[k][1] );
+// for( i=0; topol[k][1][i]>-1; i++ ) fprintf( fp, " %03d", topol[k][1][i]+1 );
+ fprintf( fp, "%d:", getdensest( mem, density )+1 );
+ for( i=0; mem[i]>-1; i++ ) fprintf( fp, " %d", mem[i]+1 );
+ fprintf( fp, "\n" );
+ }
+ free( mem );
+ }
fclose( fp );
free( tree[0] );
@@ -7920,8 +7991,10 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
free( (void *)nmemar ); nmemar = NULL;
free( mindisfrom );
free( nearest );
+ if( treeout == 2 ) free( density );
}
+
void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree )
{
int i, j, k, miniim, maxiim, minijm, maxijm;
=====================================
core/Falign.c
=====================================
@@ -774,19 +774,18 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
{
int i, j, k, l, m, maxk;
int nlen, nlen2, nlen4;
- static TLS int prevalloclen = 0;
static TLS int crossscoresize = 0;
- static TLS char **tmpseq1 = NULL;
- static TLS char **tmpseq2 = NULL;
- static TLS char **tmpptr1 = NULL;
- static TLS char **tmpptr2 = NULL;
- static TLS char **tmpres1 = NULL;
- static TLS char **tmpres2 = NULL;
- static TLS char **result1 = NULL;
- static TLS char **result2 = NULL;
+ char **tmpseq1 = NULL;
+ char **tmpseq2 = NULL;
+ char **tmpptr1 = NULL;
+ char **tmpptr2 = NULL;
+ char **tmpres1 = NULL;
+ char **tmpres2 = NULL;
+ char **result1 = NULL;
+ char **result2 = NULL;
#if RND
- static TLS char **rndseq1 = NULL;
- static TLS char **rndseq2 = NULL;
+ char **rndseq1 = NULL;
+ char **rndseq2 = NULL;
#endif
static TLS Fukusosuu **seqVector1 = NULL;
static TLS Fukusosuu **seqVector2 = NULL;
@@ -803,7 +802,7 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
static TLS Segment **sortedseg2 = NULL;
static TLS int *cut1 = NULL;
static TLS int *cut2 = NULL;
- static TLS char *sgap1, *egap1, *sgap2, *egap2;
+ char *sgap1, *egap1, *sgap2, *egap2;
static TLS int localalloclen = 0;
int lag;
int tmpint;
@@ -817,11 +816,10 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
if( seq1 == NULL )
{
- if( result1 )
+ if( kouho )
{
// fprintf( stderr, "Freeing localarrays in Falign\n" );
localalloclen = 0;
- prevalloclen = 0;
crossscoresize = 0;
mymergesort( 0, 0, NULL );
alignableReagion( 0, 0, NULL, NULL, NULL, NULL, NULL );
@@ -834,21 +832,10 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
blockAlign2( NULL, NULL, NULL, NULL, NULL, NULL );
if( crossscore ) FreeDoubleMtx( crossscore );
crossscore = NULL;
- FreeCharMtx( result1 ); result1 = NULL;
- FreeCharMtx( result2 );
- FreeCharMtx( tmpres1 );
- FreeCharMtx( tmpres2 );
- FreeCharMtx( tmpseq1 );
- FreeCharMtx( tmpseq2 );
- free( sgap1 );
- free( egap1 );
- free( sgap2 );
- free( egap2 );
free( kouho );
+ kouho = NULL;
free( cut1 );
free( cut2 );
- free( tmpptr1 );
- free( tmpptr2 );
free( segment );
free( segment1 );
free( segment2 );
@@ -888,33 +875,33 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
fprintf( stderr, "nlentmp = %d, nlen = %d\n", nlentmp, nlen );
#endif
- if( prevalloclen != alloclen ) // Falign_noudp mo kaeru
- {
- if( prevalloclen )
- {
- FreeCharMtx( result1 );
- FreeCharMtx( result2 );
- FreeCharMtx( tmpres1 );
- FreeCharMtx( tmpres2 );
- }
-// fprintf( stderr, "\n\n\nreallocating ...\n" );
- result1 = AllocateCharMtx( njob, alloclen );
- result2 = AllocateCharMtx( njob, alloclen );
- tmpres1 = AllocateCharMtx( njob, alloclen );
- tmpres2 = AllocateCharMtx( njob, alloclen );
- prevalloclen = alloclen;
- }
+ result1 = AllocateCharMtx( clus1, alloclen );
+ result2 = AllocateCharMtx( clus2, alloclen );
+ tmpres1 = AllocateCharMtx( clus1, alloclen );
+ tmpres2 = AllocateCharMtx( clus2, alloclen );
+ sgap1 = AllocateCharVec( clus1 );
+ egap1 = AllocateCharVec( clus1 );
+ sgap2 = AllocateCharVec( clus2 );
+ egap2 = AllocateCharVec( clus2 );
+ tmpptr1 = calloc( clus1, sizeof( char * ) );
+ tmpptr2 = calloc( clus2, sizeof( char * ) );
+ tmpseq1 = AllocateCharMtx( clus1, nlen );
+ tmpseq2 = AllocateCharMtx( clus2, nlen );
+#if RND
+ rndseq1 = AllocateCharMtx( clus1, nlen );
+ rndseq2 = AllocateCharMtx( clus2, nlen );
+ for( i=0; i<clus1; i++ )
+ generateRndSeq( rndseq1[i], nlen );
+ for( i=0; i<clus2; i++ )
+ generateRndSeq( rndseq2[i], nlen );
+#endif
+
+
if( !localalloclen )
{
- sgap1 = AllocateCharVec( njob );
- egap1 = AllocateCharVec( njob );
- sgap2 = AllocateCharVec( njob );
- egap2 = AllocateCharVec( njob );
kouho = AllocateIntVec( NKOUHO );
cut1 = AllocateIntVec( MAXSEG );
cut2 = AllocateIntVec( MAXSEG );
- tmpptr1 = AllocateCharMtx( njob, 0 );
- tmpptr2 = AllocateCharMtx( njob, 0 );
// crossscore = AllocateDoubleMtx( MAXSEG, MAXSEG );
segment = (Segment *)calloc( MAXSEG, sizeof( Segment ) );
segment1 = (Segment *)calloc( MAXSEG, sizeof( Segment ) );
@@ -941,17 +928,9 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
FreeFukusosuuMtx( naiseki );
FreeDoubleVec( soukan );
}
- FreeCharMtx( tmpseq1 );
- FreeCharMtx( tmpseq2 );
-#endif
-#if RND
- FreeCharMtx( rndseq1 );
- FreeCharMtx( rndseq2 );
#endif
}
- tmpseq1 = AllocateCharMtx( njob, nlen );
- tmpseq2 = AllocateCharMtx( njob, nlen );
if( !kobetsubunkatsu )
{
naisekiNoWa = AllocateFukusosuuVec( nlen );
@@ -960,18 +939,9 @@ double Falign( int **whichmtx, double ***scoringmatrices, double **n_dynamicmtx,
seqVector2 = AllocateFukusosuuMtx( n20or4or2+1, nlen+1 );
soukan = AllocateDoubleVec( nlen+1 );
}
-#if RND
- rndseq1 = AllocateCharMtx( njob, nlen );
- rndseq2 = AllocateCharMtx( njob, nlen );
- for( i=0; i<njob; i++ )
- {
- generateRndSeq( rndseq1[i], nlen );
- generateRndSeq( rndseq2[i], nlen );
- }
-#endif
localalloclen = nlen;
}
-
+
for( j=0; j<clus1; j++ ) strcpy( tmpseq1[j], seq1[j] );
for( j=0; j<clus2; j++ ) strcpy( tmpseq2[j], seq2[j] );
@@ -1610,6 +1580,7 @@ system( "less seqVec2 < /dev/tty > /dev/tty" );
}
#endif
}
+
#if KEIKA
fprintf( stderr, "DP ... done \n" );
#endif
@@ -1627,6 +1598,24 @@ system( "less seqVec2 < /dev/tty > /dev/tty" );
fprintf( stderr, "in Falign, %s\n", result2[j] );
}
#endif
+
+ FreeCharMtx( result1 );
+ FreeCharMtx( result2 );
+ FreeCharMtx( tmpres1 );
+ FreeCharMtx( tmpres2 );
+ FreeCharMtx( tmpseq1 );
+ FreeCharMtx( tmpseq2 );
+ free( sgap1 );
+ free( egap1 );
+ free( sgap2 );
+ free( egap2 );
+ free( tmpptr1 );
+ free( tmpptr2 );
+#if RND
+ FreeCharMtx( rndseq1 );
+ FreeCharMtx( rndseq2 );
+#endif
+
return( totalscore );
}
@@ -1671,6 +1660,7 @@ static int terminalmargin( int lshorter, double groupsizefac )
// return ( lshorter * 1.1 + 10 ) * groupsizefac;
}
+#if 0
static int estimatenogaplen( int n, char **s, int start, int end )
{
int i, j, l, d, o;
@@ -1732,6 +1722,7 @@ static int nogapmargin( int n, char **s, int start, int end, int m )
// reporterr( "minl=%d, so returning %d\n", minl, start+minl*d );
return( start + minl*d );
}
+#endif
double Falign_givenanchors( ExtAnch *pairanch,
int **whichmtx, double ***scoringmatrices,
@@ -1876,7 +1867,7 @@ double Falign_givenanchors( ExtAnch *pairanch,
FreeCharMtx( tmpres2 );
}
// fprintf( stderr, "\n\n\nreallocating ...\n" );
- result1 = AllocateCharMtx( njob, alloclen );
+ result1 = AllocateCharMtx( njob, alloclen ); // ato de loca nseq ni kakihaosu
result2 = AllocateCharMtx( njob, alloclen );
tmpres1 = AllocateCharMtx( njob, alloclen );
tmpres2 = AllocateCharMtx( njob, alloclen );
@@ -2388,8 +2379,9 @@ double Falign_givenanchors( ExtAnch *pairanch,
cut1[count+1] = pairanch[count0].endi+1;
cut2[count+1] = pairanch[count0].endj+1;
alignorcopy[count+1] = 'a';
-// reporterr( "\n###cut at %d and %d / %d\n", cut1[count+1], cut1[count+2], len1 );
-// reporterr( "\n###cut at %d and %d / %d\n", cut2[count+1], cut2[count+2], len2 );
+// reporterr( "\n###cut1 at %d / %d\n", cut1[count+1], len1 );
+// reporterr( "###cut2 at %d / %d\n", cut2[count+1], len2 );
+// reporterr( "sa1=%d, sa2=%d\n", cut1[count+1]-cut1[count], cut2[count+1]-cut2[count] );
count += 1;
count0++;
}
@@ -2455,6 +2447,7 @@ double Falign_givenanchors( ExtAnch *pairanch,
count += 1;
+
#if 0
for( i=0; i<count; i++ )
{
@@ -2477,6 +2470,7 @@ double Falign_givenanchors( ExtAnch *pairanch,
count += 2;
#endif
+
totallen = 0;
for( j=0; j<clus1; j++ ) result1[j][0] = 0;
for( j=0; j<clus2; j++ ) result2[j][0] = 0;
@@ -2492,8 +2486,8 @@ double Falign_givenanchors( ExtAnch *pairanch,
//reporterr( "i=%d, headgp=%d\n", i, headgp );
//reporterr( "i=%d, tailgp=%d\n", i, tailgp );
- //for( j=0; j<clus1; j++ ) reporterr( "Cut seq 0-%d at %d\n", j, cut1[i] );
- //for( j=0; j<clus2; j++ ) reporterr( "Cut seq 1-%d at %d\n", j, cut2[i] );
+// for( j=0; j<clus1; j++ ) reporterr( "Cut seq 0-%d at %d\n", j, cut1[i] );
+// for( j=0; j<clus2; j++ ) reporterr( "Cut seq 1-%d at %d\n", j, cut2[i] );
#if 0 // hoka mo atode henkou suru
if( cut1[i] )
@@ -2684,13 +2678,13 @@ double Falign_givenanchors( ExtAnch *pairanch,
int tmplen;
if( endtermcut1 )
{
- tmplen = strlen( tmpres2[j] );
+ tmplen = strlen( tmpres2[0] );
for( j=0; j<clus2; j++ ) if( tmpres2[j][tmplen-1] != '-' ) break;
if( j<clus2 ) reporterr( "There may be a problem in a hard-coded parameter (3). Please contact katoh at ifrec.osaka-u.ac.jp\n" );
}
else if( endtermcut2 )
{
- tmplen = strlen( tmpres1[j] );
+ tmplen = strlen( tmpres1[0] );
for( j=0; j<clus1; j++ ) if( tmpres1[j][tmplen-1] != '-' ) break;
if( j<clus1 ) reporterr( "There may be a problem in a hard-coded parameter (4). Please contact katoh at ifrec.osaka-u.ac.jp\n" );
}
@@ -2733,14 +2727,14 @@ double Falign_givenanchors( ExtAnch *pairanch,
fprintf( stderr, "keika \n\n" );
for( j=0; j<clus1; j++ )
{
- fprintf( stderr, ">group1-%d\n%s\n", j, result1[j] );
+ fprintf( stderr, ">group1-%d\n%100.100s\n", j, result1[j] );
}
fprintf( stderr, "- - - - - - - - - - -\n" );
for( j=0; j<clus2; j++ )
{
- fprintf( stderr, ">group2-%d\n%s\n", j, result2[j] );
+ fprintf( stderr, ">group2-%d\n%100.100s\n", j, result2[j] );
}
- if( clus1 == 2 ) exit( 1 );
+// if( clus1 == 1 && clus2 == 5 ) exit( 1 );
#endif
return( totalscore );
}
@@ -2752,7 +2746,7 @@ double Falign_givenanchors( ExtAnch *pairanch,
/*
sakujo wo kentou (2010/10/05)
*/
-double Falign_udpari_long(
+double Falign_udpari_long(
int **whichmtx, double ***scoringmatrices,
double **n_dynamicmtx,
char **seq1, char **seq2,
@@ -2763,19 +2757,18 @@ double Falign_udpari_long(
{
int i, j, k, l, m, maxk;
int nlen, nlen2, nlen4;
- static TLS int prevalloclen = 0;
static TLS int crossscoresize = 0;
- static TLS char **tmpseq1 = NULL;
- static TLS char **tmpseq2 = NULL;
- static TLS char **tmpptr1 = NULL;
- static TLS char **tmpptr2 = NULL;
- static TLS char **tmpres1 = NULL;
- static TLS char **tmpres2 = NULL;
- static TLS char **result1 = NULL;
- static TLS char **result2 = NULL;
+ char **tmpseq1 = NULL;
+ char **tmpseq2 = NULL;
+ char **tmpptr1 = NULL;
+ char **tmpptr2 = NULL;
+ char **tmpres1 = NULL;
+ char **tmpres2 = NULL;
+ char **result1 = NULL;
+ char **result2 = NULL;
#if RND
- static TLS char **rndseq1 = NULL;
- static TLS char **rndseq2 = NULL;
+ char **rndseq1 = NULL;
+ char **rndseq2 = NULL;
#endif
static TLS Fukusosuu **seqVector1 = NULL;
static TLS Fukusosuu **seqVector2 = NULL;
@@ -2792,7 +2785,7 @@ double Falign_udpari_long(
static TLS Segment **sortedseg2 = NULL;
static TLS int *cut1 = NULL;
static TLS int *cut2 = NULL;
- static TLS char *sgap1, *egap1, *sgap2, *egap2;
+ char *sgap1, *egap1, *sgap2, *egap2;
static TLS int localalloclen = 0;
int lag;
int tmpint;
@@ -2806,11 +2799,10 @@ double Falign_udpari_long(
if( seq1 == NULL )
{
- if( result1 )
+ if( kouho )
{
// fprintf( stderr, "### Freeing localarrays in Falign\n" );
localalloclen = 0;
- prevalloclen = 0;
crossscoresize = 0;
mymergesort( 0, 0, NULL );
alignableReagion( 0, 0, NULL, NULL, NULL, NULL, NULL );
@@ -2822,21 +2814,10 @@ double Falign_udpari_long(
blockAlign2( NULL, NULL, NULL, NULL, NULL, NULL );
if( crossscore ) FreeDoubleMtx( crossscore );
crossscore = NULL; // reallocate sareru kanousei ga arunode.
- FreeCharMtx( result1 ); result1 = NULL;
- FreeCharMtx( result2 );
- FreeCharMtx( tmpres1 );
- FreeCharMtx( tmpres2 );
- FreeCharMtx( tmpseq1 );
- FreeCharMtx( tmpseq2 );
- free( sgap1 );
- free( egap1 );
- free( sgap2 );
- free( egap2 );
free( kouho );
+ kouho = NULL;
free( cut1 );
free( cut2 );
- free( tmpptr1 );
- free( tmpptr2 );
free( segment );
free( segment1 );
free( segment2 );
@@ -2859,8 +2840,6 @@ double Falign_udpari_long(
return( 0.0 );
}
-
-
len1 = strlen( seq1[0] );
len2 = strlen( seq2[0] );
nlentmp = MAX( len1, len2 );
@@ -2878,34 +2857,34 @@ double Falign_udpari_long(
fprintf( stderr, "nlentmp = %d, nlen = %d\n", nlentmp, nlen );
#endif
- if( prevalloclen != alloclen ) // Falign_noudp mo kaeru
- {
- if( prevalloclen )
- {
- FreeCharMtx( result1 );
- FreeCharMtx( result2 );
- FreeCharMtx( tmpres1 );
- FreeCharMtx( tmpres2 );
- }
-// fprintf( stderr, "\n\n\nreallocating ...\n" );
- result1 = AllocateCharMtx( njob, alloclen );
- result2 = AllocateCharMtx( njob, alloclen );
- tmpres1 = AllocateCharMtx( njob, alloclen );
- tmpres2 = AllocateCharMtx( njob, alloclen );
- prevalloclen = alloclen;
- }
+ result1 = AllocateCharMtx( clus1, alloclen );
+ result2 = AllocateCharMtx( clus2, alloclen );
+ tmpres1 = AllocateCharMtx( clus1, alloclen );
+ tmpres2 = AllocateCharMtx( clus2, alloclen );
+ sgap1 = AllocateCharVec( clus1 );
+ egap1 = AllocateCharVec( clus1 );
+ sgap2 = AllocateCharVec( clus2 );
+ egap2 = AllocateCharVec( clus2 );
+
+ tmpseq1 = AllocateCharMtx( clus1, nlen );
+ tmpseq2 = AllocateCharMtx( clus2, nlen );
+ tmpptr1 = calloc( clus1, sizeof(char*) );
+ tmpptr2 = calloc( clus2, sizeof(char*) );
+
+#if RND
+ rndseq1 = AllocateCharMtx( clus1, nlen );
+ rndseq2 = AllocateCharMtx( clus2, nlen );
+ for( i=0; i<clus1; i++ )
+ generateRndSeq( rndseq1[i], nlen );
+ for( i=0; i<clus2; i++ )
+ generateRndSeq( rndseq2[i], nlen );
+#endif
if( !localalloclen )
{
- sgap1 = AllocateCharVec( njob );
- egap1 = AllocateCharVec( njob );
- sgap2 = AllocateCharVec( njob );
- egap2 = AllocateCharVec( njob );
kouho = AllocateIntVec( NKOUHO_LONG );
cut1 = AllocateIntVec( MAXSEG );
cut2 = AllocateIntVec( MAXSEG );
- tmpptr1 = AllocateCharMtx( njob, 0 );
- tmpptr2 = AllocateCharMtx( njob, 0 );
segment = (Segment *)calloc( MAXSEG, sizeof( Segment ) );
segment1 = (Segment *)calloc( MAXSEG, sizeof( Segment ) );
segment2 = (Segment *)calloc( MAXSEG, sizeof( Segment ) );
@@ -2918,6 +2897,7 @@ double Falign_udpari_long(
else if( fftscore ) n20or4or2 = 1;
else n20or4or2 = 20;
}
+
if( localalloclen < nlen )
{
if( localalloclen )
@@ -2931,18 +2911,10 @@ double Falign_udpari_long(
FreeFukusosuuMtx( naiseki );
FreeDoubleVec( soukan );
}
- FreeCharMtx( tmpseq1 );
- FreeCharMtx( tmpseq2 );
-#endif
-#if RND
- FreeCharMtx( rndseq1 );
- FreeCharMtx( rndseq2 );
#endif
}
- tmpseq1 = AllocateCharMtx( njob, nlen );
- tmpseq2 = AllocateCharMtx( njob, nlen );
if( !kobetsubunkatsu )
{
naisekiNoWa = AllocateFukusosuuVec( nlen );
@@ -2951,15 +2923,6 @@ double Falign_udpari_long(
seqVector2 = AllocateFukusosuuMtx( n20or4or2, nlen+1 );
soukan = AllocateDoubleVec( nlen+1 );
}
-#if RND
- rndseq1 = AllocateCharMtx( njob, nlen );
- rndseq2 = AllocateCharMtx( njob, nlen );
- for( i=0; i<njob; i++ )
- {
- generateRndSeq( rndseq1[i], nlen );
- generateRndSeq( rndseq2[i], nlen );
- }
-#endif
localalloclen = nlen;
}
@@ -3554,6 +3517,7 @@ system( "less seqVec2 < /dev/tty > /dev/tty" );
totalscore += MSalignmm_variousdist( NULL, scoringmatrices, NULL, tmpres1, tmpres2, eff1, eff2, eff1s, eff2s, clus1, clus2, alloclen, sgap1, sgap2, egap1, egap2, NULL, 0, NULL, headgp, tailgp );
else
totalscore += MSalignmm( n_dynamicmtx, tmpres1, tmpres2, eff1, eff2, clus1, clus2, alloclen, sgap1, sgap2, egap1, egap2, NULL, 0, NULL, headgp, tailgp, NULL, NULL, NULL, 0.0, 0.0 );
+// totalscore += G__align11( n_dynamicmtx, tmpres1, tmpres2, alloclen, headgp, tailgp ); // CHUUI!!!
break;
default:
fprintf( stderr, "alg = %c\n", alg );
@@ -3602,5 +3566,22 @@ system( "less seqVec2 < /dev/tty > /dev/tty" );
fprintf( stderr, "%s\n", result2[j] );
}
#endif
+
+ FreeCharMtx( result1 );
+ FreeCharMtx( result2 );
+ FreeCharMtx( tmpres1 );
+ FreeCharMtx( tmpres2 );
+ free( sgap1 );
+ free( egap1 );
+ free( sgap2 );
+ free( egap2 );
+ FreeCharMtx( tmpseq1 );
+ FreeCharMtx( tmpseq2 );
+ free( tmpptr1 );
+ free( tmpptr2 );
+#if RND
+ FreeCharMtx( rndseq1 );
+ FreeCharMtx( rndseq2 );
+#endif
return( totalscore );
}
=====================================
core/Galign11.c
=====================================
@@ -352,8 +352,8 @@ double G__align11( double **n_dynamicmtx, char **seq1, char **seq2, int alloclen
if( orlgth1 == 0 )
{
- mseq1 = AllocateCharMtx( njob, 0 );
- mseq2 = AllocateCharMtx( njob, 0 );
+ mseq1 = AllocateCharMtx( 2, 0 ); // 2020/Apr
+ mseq2 = AllocateCharMtx( 2, 0 ); // 2020/Apr
}
@@ -399,7 +399,7 @@ double G__align11( double **n_dynamicmtx, char **seq1, char **seq2, int alloclen
m = AllocateFloatVec( ll2+2 );
mp = AllocateIntVec( ll2+2 );
- mseq = AllocateCharMtx( njob, ll1+ll2 );
+ mseq = AllocateCharMtx( 2, ll1+ll2 ); // 2020/Apr
doublework = AllocateFloatMtx( nalphabets, MAX( ll1, ll2 )+2 );
=====================================
core/Makefile
=====================================
@@ -51,6 +51,9 @@ MYCFLAGS = $(MNO_CYGWIN) $(ENABLE_MULTITHREAD) $(ENABLE_ATOMIC) $(STDF) $(CFLAGS
INSTALL = install
+STRIP = strip
+#STRIP = true # to disable strip
+
PROGS = dvtditr dndfast7 dndblast sextet5 mafft-distance pairlocalalign \
multi2hat3s pairash addsingle maffttext2hex hex2maffttext \
splittbfast disttbfast tbfast nodepair mafft-profile f2cl mccaskillwrap contrafoldwrap countlen \
@@ -159,6 +162,7 @@ dlls : $(DLLS)
$(DASH_CLIENT): dash_client.go
go build dash_client.go
+# go build --ldflags '-extldflags "-static"' dash_client.go
univscript: univscript.tmpl Makefile
sed "s:_PROGS:$(PROGS):" univscript.tmpl > univscript
@@ -556,7 +560,9 @@ install : all
chmod 755 $(SCRIPTS)
$(INSTALL) $(SCRIPTS) $(DESTDIR)$(BINDIR)
chmod 755 $(PROGS) ||: # in MinGW, it's ok if this fails
- $(INSTALL) -s $(PROGS) $(DESTDIR)$(LIBDIR)
+# $(INSTALL) -s $(PROGS) $(DESTDIR)$(LIBDIR)
+ $(STRIP) $(PROGS) ||: # may fail for dash_client on mac.
+ $(INSTALL) $(PROGS) $(DESTDIR)$(LIBDIR)
$(INSTALL) $(PERLPROGS) $(DESTDIR)$(LIBDIR)
$(INSTALL) -m 644 $(MANPAGES) $(DESTDIR)$(LIBDIR)
=====================================
core/README deleted
=====================================
@@ -1 +0,0 @@
-SUPPORT environments without Pthreads
=====================================
core/Salignmm.c
=====================================
@@ -1103,9 +1103,9 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
static TLS double *match;
static TLS double *initverticalw; /* kufuu sureba iranai */
static TLS double *lastverticalw; /* kufuu sureba iranai */
- static TLS char **mseq1;
- static TLS char **mseq2;
- static TLS char **mseq;
+ char **mseq1;
+ char **mseq2;
+ char **mseq;
static TLS double *ogcp1, *ogcp1o;
static TLS double *ogcp2, *ogcp2o;
static TLS double *fgcp1, *fgcp1o;
@@ -1164,8 +1164,6 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
imp_match_init_strict( NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, NULL, NULL, NULL, NULL, 0, 0 );
- free( mseq1 );
- free( mseq2 );
FreeFloatVec( w1 );
FreeFloatVec( w2 );
FreeFloatVec( match );
@@ -1175,7 +1173,6 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
FreeFloatVec( m );
FreeIntVec( mp );
- FreeCharMtx( mseq );
FreeFloatVec( ogcp1 );
FreeFloatVec( ogcp1o );
@@ -1279,12 +1276,10 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
fprintf( stderr, "#### strlen( seq2[0] ) = %d\n", strlen( seq2[0] ) );
for( i=0; i<jcyc; i++ ) fprintf( stderr, "eff2[%d] = %f\n", i, eff2[i] );
#endif
- if( orlgth1 == 0 )
- {
- mseq1 = AllocateCharMtx( njob, 0 );
- mseq2 = AllocateCharMtx( njob, 0 );
- }
+ mseq1 = AllocateCharMtx( icyc, 0 );
+ mseq2 = AllocateCharMtx( jcyc, 0 );
+ mseq = AllocateCharMtx( icyc+jcyc, lgth1+lgth2+100 );
if( lgth1 > orlgth1 || lgth2 > orlgth2 )
@@ -1303,8 +1298,6 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
FreeFloatVec( m );
FreeIntVec( mp );
- FreeCharMtx( mseq );
-
FreeFloatVec( ogcp1 );
FreeFloatVec( ogcp1o );
FreeFloatVec( ogcp2 );
@@ -1342,7 +1335,6 @@ double A__align( double **n_dynamicmtx, int penalty_l, int penalty_ex_l, char **
m = AllocateFloatVec( ll2+2 );
mp = AllocateIntVec( ll2+2 );
- mseq = AllocateCharMtx( njob, ll1+ll2 );
ogcp1 = AllocateFloatVec( ll1+2 );
ogcp1o = AllocateFloatVec( ll1+2 );
@@ -2184,6 +2176,10 @@ fprintf( stderr, "\n" );
previousicyc = icyc;
previouscall = calledbyfulltreebase;
+ free( mseq1 );
+ free( mseq2 );
+ FreeCharMtx( mseq );
+
return( wm );
}
=====================================
core/addsingle.c
=====================================
@@ -1368,7 +1368,7 @@ static double treebase( int nseq, int *nlen, char **aseq, int nadd, char *mergeo
*/
- if( !nevermemsave && ( constraint != 2 && alg != 'M' && ( len1 > 30000 || len2 > 30000 ) ) )
+ if( !nevermemsave && ( constraint != 2 && alg != 'M' && ( len1 > 50000 || len2 > 50000 ) ) )
{
fprintf( stderr, "\nlen1=%d, len2=%d, Switching to the memsave mode.\n", len1, len2 );
alg = 'M';
@@ -1379,12 +1379,22 @@ static double treebase( int nseq, int *nlen, char **aseq, int nadd, char *mergeo
}
+ if( alg == 'M' ) // hoka no thread ga M ni shitakamo shirenainode
+ {
+ if( commonIP ) FreeIntMtx( commonIP );
+ commonIP = NULL;
+ commonAlloc1 = 0;
+ commonAlloc2 = 0;
+ }
+
+
// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000 );
else ffttry = 0;
// ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 5000 && clus2 < 5000 ); // v6.708
// fprintf( stderr, "f=%d, len1/fftlog[m1]=%f, clus1=%d, len2/fftlog[m2]=%f, clus2=%d\n", ffttry, (double)len1/fftlog[m1], clus1, (double)len2/fftlog[m2], clus2 );
// fprintf( stderr, "f=%d, clus1=%d, fftlog[m1]=%d, clus2=%d, fftlog[m2]=%d\n", ffttry, clus1, fftlog[m1], clus2, fftlog[m2] );
+
if( constraint == 2 )
{
if( alg == 'M' )
@@ -1598,6 +1608,7 @@ static void *addsinglethread( void *arg )
norg = njob - nadd;
njobc = norg+1;
+#if 0
alnleninnode = AllocateIntVec( norg );
addmem = AllocateIntVec( nadd+1 );
depc = (Treedep *)calloc( njobc, sizeof( Treedep ) );
@@ -1609,6 +1620,19 @@ static void *addsinglethread( void *arg )
mergeoralign = AllocateCharVec( njob );
nogaplenjusttodecideaddhereornot = AllocateIntVec( njobc );
tmpseq = calloc( alloclen, sizeof( char ) );
+#else
+ alnleninnode = AllocateIntVec( norg );
+ addmem = AllocateIntVec( nadd+1 );
+ depc = (Treedep *)calloc( njobc, sizeof( Treedep ) );
+ mseq1 = AllocateCharMtx( njobc, 0 );
+ mseq2 = AllocateCharMtx( njobc, 0 );
+ bseq = AllocateCharMtx( njobc, alloclen );
+ namec = AllocateCharMtx( njobc, 0 );
+ nlenc = AllocateIntVec( njobc );
+ mergeoralign = AllocateCharVec( njobc );
+ nogaplenjusttodecideaddhereornot = AllocateIntVec( njobc );
+ tmpseq = calloc( alloclen, sizeof( char ) );
+#endif
if( allowlongadds ) // hontou ha iranai.
{
@@ -1701,7 +1725,10 @@ static void *addsinglethread( void *arg )
pthread_mutex_unlock( targ->mutex_counter );
break;
}
- fprintf( stderr, "\r%d / %d (thread %d) \r", iadd, nadd, thread_no );
+ if( iadd < 500 )
+ fprintf( stderr, "\rSTEP %d / %d (thread %d) \r", iadd, nadd, thread_no );
+ else if( iadd % 100 == 0 )
+ fprintf( stderr, "\nSTEP %d / %d (thread %d) \n", iadd, nadd, thread_no );
++(*iaddshare);
targetseq = seq[norg+iadd];
pthread_mutex_unlock( targ->mutex_counter );
@@ -1712,7 +1739,10 @@ static void *addsinglethread( void *arg )
iadd++;
if( iadd == nadd ) break;
targetseq = seq[norg+iadd];
- fprintf( stderr, "\r%d / %d \r", iadd, nadd );
+ if( iadd < 500 )
+ fprintf( stderr, "\rSTEP %d / %d \r", iadd, nadd );
+ else if( iadd % 100 == 0 )
+ fprintf( stderr, "\nSTEP %d / %d \n", iadd, nadd );
}
for( i=0; i<norg; i++ ) strcpy( bseq[i], seq[i] );
@@ -1974,6 +2004,7 @@ static void *addsinglethread( void *arg )
if( mseq1 ) free( mseq1 ); mseq1 = NULL;
if( mseq2 ) free( mseq2 ); mseq2 = NULL;
Falign( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL, NULL, 0, NULL );
+ Falign_udpari_long( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
A__align( NULL, 0, 0, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0, -1, -1, NULL, NULL, NULL, 0.0, 0.0 );
if( commonIP ) FreeIntMtx( commonIP );
commonIP = NULL;
@@ -2929,11 +2960,11 @@ int main( int argc, char *argv[] )
njobc = norg+1;
fprintf( stderr, "norg = %d\n", norg );
fprintf( stderr, "njobc = %d\n", njobc );
- if( norg > 1000 || nadd > 1000 ) use_fft = 0;
+// if( norg > 1000 || nadd > 1000 ) use_fft = 0;
+ if( norg > 1000 ) use_fft = 0;
fullseqlen = alloclen = nlenmax*4+1; //chuui!
-// fullseqlen = alloclen = nlenmax*2+1; //chuui!
seq = AllocateCharMtx( njob, alloclen );
name = AllocateCharMtx( njob, B+1 );
@@ -3392,7 +3423,15 @@ int main( int argc, char *argv[] )
}
free( dep );
FreeFloatMtx( len );
- if( multidist || tuplesize > 0 ) FreeFloatMtx( nscore );
+ if( multidist || tuplesize > 0 )
+ {
+ FreeFloatHalfMtx( iscore, norg );
+ FreeFloatMtx( nscore );
+ }
+ else
+ {
+ FreeFloatHalfMtx( iscore, njob );
+ }
// for( i=0; i<nadd; i++ ) fprintf( stdout, ">%s (%d) \n%s\n", name[norg+i], norg+i, seq[norg+i] );
@@ -3699,6 +3738,11 @@ int main( int argc, char *argv[] )
if( constraint ) FreeLocalHomTable( localhomtable, njob );
}
+ FreeCharMtx( seq );
+ FreeCharMtx( name );
+ closeFiles();
+ commonsextet_p( NULL, NULL );
+
SHOWVERSION;
if( ndeleted > 0 )
{
=====================================
core/blosum.c
=====================================
@@ -290,6 +290,7 @@ void BLOSUMmtx( int n, double **matrix, double *freq, unsigned char *amino, char
else
for( i=0; i<20; i++ ) freq[i] = freqd[i];
+ if( n == -1 ) free( tmpmtx );
#if 0
av = 0.0;
for( i=0; i<20; i++ )
=====================================
core/disttbfast.c
=====================================
@@ -8,6 +8,8 @@
#define SCOREOUT 0
#define SKIP 1
+#define ITERATIVECYCLE 1
+
#define END_OF_VEC -1
static int treein;
@@ -27,6 +29,8 @@ static int smoothing;
static double maxdistmtxsize;
static int nthreadtb;
static int useexternalanchors;
+static int oneiteration;
+static double maxanchorseparation;
#if 0
#define PLENFACA 0.0123
@@ -168,6 +172,7 @@ void arguments( int argc, char *argv[] )
devide = 0;
use_fft = 0;
useexternalanchors = 0;
+ oneiteration = 0;
force_fft = 0;
fftscore = 1;
fftRepeatStop = 0;
@@ -217,6 +222,7 @@ void arguments( int argc, char *argv[] )
mapout = 0;
smoothing = 0;
nwildcard = 0;
+ maxanchorseparation = 1000.0;
while( --argc > 0 && (*++argv)[0] == '-' )
{
@@ -313,14 +319,15 @@ void arguments( int argc, char *argv[] )
case 't':
treeout = 1;
break;
+ case '^':
+ treeout = 2;
+ break;
case 'T':
noalign = 1;
break;
-#if 0
case 'r':
- fmodel = -1;
+ oneiteration = 1;
break;
-#endif
case 'D':
dorp = 'd';
break;
@@ -333,6 +340,10 @@ void arguments( int argc, char *argv[] )
case 'e':
fftscore = 0;
break;
+ case 'x':
+ maxanchorseparation = myatof( *++argv );
+ --argc;
+ goto nextoption;
case 'H':
subalignment = 1;
subalignmentoffset = myatoi( *++argv );
@@ -1346,7 +1357,7 @@ static void recountpositions( ExtAnch *pairanch, int n1, int n2, char **seq1, ch
len = strlen( seq1[0] )+1;
map = calloc( sizeof( int ), len );
-
+
for( k=0; k<n1; k++ )
{
pos = 0;
@@ -1441,6 +1452,8 @@ static void indexanchors( ExtAnch *a, int **idx )
#endif
}
+
+#if 0
static void checkanchors_internal( ExtAnch *a )
{
int p, q, r, s;
@@ -1467,11 +1480,12 @@ static void checkanchors_internal( ExtAnch *a )
qsort( a+s, r-s, sizeof( ExtAnch ), anchscorecomp );
#if 0
- reporterr( "after sortscore\n" );
+ reporterr( "after sortscore, r=%d\n", r );
for( p=s; p<r; p++ )
{
reporterr( "a[%d].starti,j=%d,%d, score=%d\n", p, a[p].starti, a[p].startj, a[p].score );
}
+ exit( 1 );
#endif
for( p=s; p<r; p++ )
@@ -1551,10 +1565,14 @@ static void checkanchors_internal( ExtAnch *a )
exit( 1 );
#endif
}
+#endif
-static void checkanchors_strongestfirst( ExtAnch *a, int s )
+static void checkanchors_strongestfirst( ExtAnch *a, int s, double gapratio1, double gapratio2 )
{
int p, q;
+ double zureij;
+ double nogaplenestimation1;
+ double nogaplenestimation2;
#if 0
reporterr( "before sortscore\n" );
for( p=0; a[p].i>-1; p++ )
@@ -1563,14 +1581,69 @@ static void checkanchors_strongestfirst( ExtAnch *a, int s )
}
#endif
qsort( a, s, sizeof( ExtAnch ), anchscorecomp );
+
+ nogaplenestimation1 = (double)a[0].starti / (1.0+gapratio1);
+ nogaplenestimation2 = (double)a[0].startj / (1.0+gapratio2);
+ zureij = nogaplenestimation1 - nogaplenestimation2;
+
for( p=0; a[p].i>-1; p++ )
{
if( a[p].starti == -1 ) continue;
+
+#if 0
+ nogaplenestimation1 = (double)a[p].starti / (1.0+gapratio1);
+ nogaplenestimation2 = (double)a[p].startj / (1.0+gapratio2);
+ if( fabs( zureij - ( nogaplenestimation1 - nogaplenestimation2 ) ) > maxanchorseparation )
+ {
+// reporterr( "warning: long internal gaps in %d-%d, |%5.2f-%5.2f - %5.2f| = %5.2f > %5.2f\n", a[p].i, a[p].j, nogaplenestimation1, nogaplenestimation2, zureij, fabs( zureij - ( nogaplenestimation1, nogaplenestimation2 ) ), maxanchorseparation );
+ a[p].starti = a[p].startj = a[p].startj = a[p].endj = -1;
+ continue;
+ }
+#else
+ int nearest, mindist;
+ double zurei, zurej;
+ if( p )
+ {
+ mindist = 999999999;
+ for( q=0; q<p; q++ )
+ {
+ if( a[q].starti == -1 ) continue;
+ if( abs( a[p].starti - a[q].starti ) < mindist )
+ {
+ nearest = q;
+ mindist = abs( a[p].starti - a[q].starti );
+ }
+ }
+ //reporterr( "nearest=%d\n", nearest );
+ if( a[nearest].starti < a[p].starti )
+ {
+ zurei = (double)( a[p].starti - a[nearest].endi )/(1.0+gapratio1);
+ zurej = (double)( a[p].startj - a[nearest].endj )/(1.0+gapratio2);
+ }
+ else
+ {
+ zurei = (double)( a[nearest].starti - a[p].endi )/(1.0+gapratio1);
+ zurej = (double)( a[nearest].startj - a[p].endj )/(1.0+gapratio2);
+ }
+ }
+ else
+ zurei = zurej = 0.0;
+ if( fabs( zurei - zurej ) > maxanchorseparation )
+// if( fabs( zurei - zurej ) > maxanchorseparation || zurei > maxanchorseparation || zurej > maxanchorseparation ) // test
+ {
+// reporterr( "warning: long internal gaps in %d-%d, |%5.2f-%5.2f - %5.2f| = %5.2f > %5.2f\n", a[p].i, a[p].j, nogaplenestimation1, nogaplenestimation2, zureij, fabs( zureij - ( nogaplenestimation1, nogaplenestimation2 ) ), maxanchorseparation );
+ a[p].starti = a[p].startj = a[p].startj = a[p].endj = -1;
+ continue;
+ }
+#endif
+
// reporterr( "P score=%d, %d-%d, %d-%d\n", a[p].score, a[p].starti, a[p].endi, a[p].startj, a[p].endj );
for( q=p+1; a[q].i>-1; q++ )
{
if( a[q].starti == -1 ) continue;
// reporterr( "Q score=%d, %d-%d, %d-%d\n", a[q].score, a[q].starti, a[q].endi, a[q].startj, a[q].endj );
+
+
if( a[p].endi < a[q].starti && a[p].endj < a[q].startj )
{
// reporterr( "consistent\n" );
@@ -1604,19 +1677,52 @@ static void checkanchors_strongestfirst( ExtAnch *a, int s )
}
}
}
+
+ qsort( a, s, sizeof( ExtAnch ), anchcomp );
#if 0
- reporterr( "after filtering\n" );
+ reporterr( "after filtering and sorting\n" );
for( p=0; a[p].i>-1; p++ )
{
reporterr( "a[%d].starti,j=%d,%d, score=%d\n", p, a[p].starti, a[p].startj, a[p].score );
}
-// exit( 1 );
#endif
-
- qsort( a, s, sizeof( ExtAnch ), anchcomp );
}
+static double gapnongapratio( int n, char **s )
+{
+ int i, j, len;
+ char *seq, *pt1, *pt2;
+ double fv, ng;
+
+ len = strlen( s[0] );
+ seq = calloc( len+1, sizeof( char ) );
+
+ fv = 0.0;
+ ng = 0.0;
+ for( i=0; i<n; i++ )
+ {
+ pt1 = s[i];
+ while( *pt1 == '-' ) pt1++;
+ pt2 = seq;
+ while( *pt1 != 0 ) *pt2++ = *pt1++;
+ *pt2 = *pt1; // 0
+ pt1 = pt2-1;
+ while( *pt1 == '-' ) pt1--;
+ *(pt1+1) = 0;
+// reporterr( "seq[i]=%s\n", s[i] );
+// reporterr( "seq=%s\n", seq );
+ len = pt1-seq+1;
+ for( j=0; j<len; j++ )
+ if( seq[j] == '-' )
+ fv+=1.0;
+ else
+ ng+=1.0;
+ }
+ free( seq );
+ return( fv/ng );
+}
+
static void pickpairanch( ExtAnch **pairanch, ExtAnch *extanch, int **anchindex, int n1, int n2, int *m1, int *m2, char **seq1, char **seq2 ) // loop no junban wo kaeta hou ga iikamo
{
int i, j, k, s;
@@ -1754,8 +1860,19 @@ static void pickpairanch( ExtAnch **pairanch, ExtAnch *extanch, int **anchindex,
}
#endif
+#if 0
+ reporterr( "\ngroup1=\n" );
+ for( i=0; m1[i]>-1; i++ )
+ reporterr( "%d ", m1[i] );
+ reporterr( "\n" );
+ reporterr( "\ngroup2=\n" );
+ for( i=0; m2[i]>-1; i++ )
+ reporterr( "%d ", m2[i] );
+ reporterr( "\n" );
+#endif
+
+ checkanchors_strongestfirst( *pairanch, s, gapnongapratio( n1, seq1 ), gapnongapratio( n2, seq2 ) );
- checkanchors_strongestfirst( *pairanch, s );
// qsort( *pairanch, s, sizeof( ExtAnch ), anchcomp );
// checkanchors_new( *pairanch );
@@ -2101,10 +2218,13 @@ static void *treebasethread( void *arg )
// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
- if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000);
- else ffttry = 0;
+// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000);
+// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
+// else ffttry = 0;
+ ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
// ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 5000 && clus2 < 5000); // v6.708
// reporterr( "f=%d, len1/fftlog[m1]=%f, clus1=%d, len2/fftlog[m2]=%f, clus2=%d\n", ffttry, (double)len1/fftlog[m1], clus1, (double)len2/fftlog[m2], clus2 );
+// reporterr( "fftlog=%d,%d, ffttry=%d\n", fftlog[m1], fftlog[m2], ffttry );
if( useexternalanchors )
{
@@ -2284,6 +2404,316 @@ static void *treebasethread( void *arg )
}
#endif
+static int dooneiteration( int *nlen, char **aseq, int nadd, char *mergeoralign, char **mseq1, char **mseq2, int ***topol, Treedep *dep, int **memhist, double ***cpmxhist, double *effarr, double **newdistmtx, int *selfscore, ExtAnch *extanch, int **anchindex, int *alloclen, int (*callback)(int, int, char*) )
+{
+ int l, ll, len1, len2, i, j;
+ int clus1, clus2;
+ double pscore, tscore;
+ char *indication1 = NULL, *indication2 = NULL;
+ double *effarr1 = NULL;
+ double *effarr2 = NULL;
+ int *fftlog = NULL; // fixed at 2006/07/26
+// double dumfl = 0.0;
+ double dumdb = 0.0;
+ int ffttry;
+ int m1, m2;
+ int *alreadyaligned = NULL;
+ double **dynamicmtx = NULL;
+ int **localmem = NULL;
+ double **cpmxchild0, **cpmxchild1;
+ double orieff1, orieff2;
+ ExtAnch *pairanch;
+#if SKIP
+ int **skiptable1 = NULL, **skiptable2 = NULL;
+#endif
+#if 0
+ int i, j;
+#endif
+
+
+ if( effarr1 == NULL )
+ {
+ effarr1 = AllocateDoubleVec( njob );
+ effarr2 = AllocateDoubleVec( njob );
+ indication1 = AllocateCharVec( 150 );
+ indication2 = AllocateCharVec( 150 );
+ fftlog = AllocateIntVec( njob );
+ alreadyaligned = AllocateIntVec( njob );
+ if( specificityconsideration )
+ dynamicmtx = AllocateDoubleMtx( nalphabets, nalphabets );
+ localmem = calloc( sizeof( int * ), 2 );
+ }
+ for( i=0; i<njob-nadd; i++ ) alreadyaligned[i] = 1;
+ for( i=njob-nadd; i<njob; i++ ) alreadyaligned[i] = 0;
+
+ if( callback && callback( 0, 50, "Progressive alignment" ) ) goto chudan_tbfast;
+
+ for( l=0; l<njob; l++ ) fftlog[l] = 1;
+
+#if 0 // chain you
+ localmem[0][0] = -1;
+ localmem[1][0] = -1;
+ clus1 = 1;// chain ni hitsuyou
+#endif
+
+#if 0
+ reporterr( "##### fftwinsize = %d, fftthreshold = %d\n", fftWinSize, fftThreshold );
+#endif
+
+#if 0
+ for( i=0; i<njob; i++ )
+ reporterr( "TBFAST effarr[%d] = %f\n", i, effarr[i] );
+#endif
+
+// for( i=0; i<njob; i++ ) strcpy( aseq[i], seq[i] );
+
+
+// writePre( njob, name, nlen, aseq, 0 );
+
+ tscore = 0.0;
+ for( ll=0; ll<njob*ITERATIVECYCLE; ll++ )
+ {
+ l = ll % njob;
+ cpmxchild0 = NULL;
+ cpmxchild1 = NULL;
+
+ localmem[0] = calloc( sizeof( int ), 2 );
+ localmem[0][0] = l;
+ localmem[0][1] = -1;
+ clus1 = 1;
+ m1 = localmem[0][0];
+
+ localmem[1] = calloc( sizeof( int ), njob );
+ for( i=0,j=0; i<njob; i++ )
+ if( i != l ) localmem[1][j++] = i;
+ localmem[1][j] = -1;
+ clus2 = njob-1;
+ m2 = localmem[1][0];
+
+// reporterr( "\ndistfromtip = %f\n", dep[l].distfromtip );
+ if( specificityconsideration )
+ makedynamicmtx( dynamicmtx, n_dis_consweight_multi, dep[l].distfromtip );
+ else
+ dynamicmtx = n_dis_consweight_multi;
+// makedynamicmtx( dynamicmtx, n_dis_consweight_multi, ( dep[l].distfromtip - 0.2 ) * 3 );
+
+ len1 = strlen( aseq[m1] );
+ len2 = strlen( aseq[m2] );
+ if( *alloclen < len1 + len2 )
+ {
+ reporterr( "\nReallocating.." );
+ *alloclen = ( len1 + len2 ) + 1000;
+ ReallocateCharMtx( aseq, njob, *alloclen + 10 );
+ }
+
+#if 1 // CHUUI@@@@
+ clus1 = fastconjuction_noname( localmem[0], aseq, mseq1, effarr1, effarr, indication1, 0.0, &orieff1 );
+ clus2 = fastconjuction_noname( localmem[1], aseq, mseq2, effarr2, effarr, indication2, 0.0, &orieff2 );
+#else
+ clus1 = fastconjuction_noname( topol[l][0], aseq, mseq1, effarr1, effarr, indication1, 0.0 );
+ clus2 = fastconjuction_noname( topol[l][1], aseq, mseq2, effarr2, effarr, indication2, 0.0 );
+// clus1 = fastconjuction_noweight( topol[l][0], aseq, mseq1, effarr1, indication1 );
+// clus2 = fastconjuction_noweight( topol[l][1], aseq, mseq2, effarr2, indication2 );
+#endif
+
+ newgapstr = "-";
+
+ commongappick( clus2, mseq2 );
+ commongappick( clus1, mseq1 );
+
+ if( l < 500 || l % 100 == 0 ) reporterr( "\rIteration % 5d / %d ", ll+1, njob*ITERATIVECYCLE );
+ if( callback && callback( 0, 50+50*l/(njob-1), "Progressive alignment" ) ) goto chudan_tbfast;
+#if 0
+ reporterr( "\nclus1=%d, clus2=%d\n", clus1, clus2 );
+ for( i=0; i<clus1; i++ ) reporterr( "effarr1[%d]=%f\n", i, effarr1[i] );
+ for( i=0; i<clus2; i++ ) reporterr( "effarr2[%d]=%f\n", i, effarr2[i] );
+#endif
+
+#if 0
+ reporterr( "STEP %d /%d\n", l+1, njob-1 );
+ reporterr( "group1 = %.66s", indication1 );
+ if( strlen( indication1 ) > 66 ) reporterr( "..." );
+ reporterr( "\n" );
+ reporterr( "group2 = %.66s", indication2 );
+ if( strlen( indication2 ) > 66 ) reporterr( "..." );
+ reporterr( "\n" );
+#endif
+
+/*
+ reporterr( "before align all\n" );
+ display( aseq, njob );
+ reporterr( "\n" );
+ reporterr( "before align 1 %s \n", indication1 );
+ display( mseq1, clus1 );
+ reporterr( "\n" );
+ reporterr( "before align 2 %s \n", indication2 );
+ display( mseq2, clus2 );
+ reporterr( "\n" );
+*/
+
+ if( !nevermemsave && ( alg != 'M' && ( len1 > 30000 || len2 > 30000 ) ) )
+ {
+ reporterr( "\nlen1=%d, len2=%d, Switching to the memsave mode\n", len1, len2 );
+ alg = 'M';
+ if( commonIP ) FreeIntMtx( commonIP );
+ commonIP = NULL;
+ commonAlloc1 = 0;
+ commonAlloc2 = 0;
+ }
+
+// ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000);
+ ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
+// ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 5000 && clus2 < 5000); // v6.708
+// reporterr( "f=%d, len1/fftlog[m1]=%f, clus1=%d, len2/fftlog[m2]=%f, clus2=%d\n", ffttry, (double)len1/fftlog[m1], clus1, (double)len2/fftlog[m2], clus2 );
+
+ if( useexternalanchors )
+ {
+ pickpairanch( &pairanch, extanch, anchindex, clus1, clus2, localmem[0], localmem[1], mseq1, mseq2 );
+// reporterr( "pairanch: %d:%d\n", pairanch[0].starti, pairanch[0].startj );
+ pscore = Falign_givenanchors( pairanch, NULL, NULL, dynamicmtx, mseq1, mseq2, effarr1, effarr2, NULL, NULL, clus1, clus2, *alloclen, fftlog+m1 );
+ free( pairanch );
+ pairanch = NULL;
+ }
+ else if( force_fft || ( use_fft && ffttry ) )
+ {
+ if( l < 500 || l % 100 == 0 ) reporterr( " f\b\b" );
+ if( alg == 'M' )
+ {
+ if( l < 500 || l % 100 == 0 ) reporterr( "m" );
+ pscore = Falign_udpari_long( NULL, NULL, dynamicmtx, mseq1, mseq2, effarr1, effarr2, NULL, NULL, clus1, clus2, *alloclen, fftlog+m1 );
+ }
+ else
+ {
+ pscore = Falign( NULL, NULL, dynamicmtx, mseq1, mseq2, effarr1, effarr2, NULL, NULL, clus1, clus2, *alloclen, fftlog+m1, NULL, 0, NULL );
+// reporterr( "######### mseq1[0] = %s\n", mseq1[0] );
+ }
+ }
+ else
+ {
+ if( l < 500 || l % 100 == 0 ) reporterr( " d\b\b" );
+ fftlog[m1] = 0;
+ switch( alg )
+ {
+ case( 'a' ):
+ pscore = Aalign( mseq1, mseq2, effarr1, effarr2, clus1, clus2, *alloclen );
+ break;
+ case( 'M' ):
+ if( l < 500 || l % 100 == 0 ) reporterr( "m" );
+ if( l < 500 || l % 100 == 0 ) if( cpmxchild1 || cpmxchild0 ) reporterr( " h" );
+// reporterr( "%d-%d", clus1, clus2 );
+ pscore = MSalignmm( dynamicmtx, mseq1, mseq2, effarr1, effarr2, clus1, clus2, *alloclen, NULL, NULL, NULL, NULL, NULL, 0, NULL, outgap, outgap, cpmxchild0, cpmxchild1, cpmxhist+l, orieff1, orieff2 );
+ break;
+ case( 'd' ):
+ if( 1 && clus1 == 1 && clus2 == 1 )
+ {
+// reporterr( "%d-%d", clus1, clus2 );
+ pscore = G__align11( dynamicmtx, mseq1, mseq2, *alloclen, outgap, outgap );
+ }
+ else
+ {
+// reporterr( "%d-%d", clus1, clus2 );
+ pscore = D__align_ls( dynamicmtx, mseq1, mseq2, effarr1, effarr2, clus1, clus2, *alloclen, 0, &dumdb, NULL, NULL, NULL, NULL, NULL, 0, NULL, outgap, outgap );
+ }
+ break;
+ case( 'A' ):
+ if( clus1 == 1 && clus2 == 1 )
+ {
+// reporterr( "%d-%d", clus1, clus2 );
+ pscore = G__align11( dynamicmtx, mseq1, mseq2, *alloclen, outgap, outgap );
+ }
+ else
+ {
+ if( l < 500 || l % 100 == 0 ) if( cpmxchild1 || cpmxchild0 ) reporterr( " h" );
+// reporterr( "\n\n %d - %d (%d x %d) : \n", topol[l][0][0], topol[l][1][0], clus1, clus2 );
+ pscore = A__align( dynamicmtx, penalty, penalty_ex, mseq1, mseq2, effarr1, effarr2, clus1, clus2, *alloclen, 0, &dumdb, NULL, NULL, NULL, NULL, NULL, 0, NULL, outgap, outgap, localmem[0][0], 1, cpmxchild0, cpmxchild1, cpmxhist+l, orieff1, orieff2 );
+ }
+
+ break;
+ default:
+ ErrorExit( "ERROR IN SOURCE FILE" );
+ }
+ }
+#if SCOREOUT
+ reporterr( "score = %10.2f\n", pscore );
+#endif
+ tscore += pscore;
+ nlen[m1] = 0.5 * ( nlen[m1] + nlen[m2] );
+
+// writePre( njob, name, nlen, aseq, 0 );
+
+ if( disp ) display( aseq, njob );
+// reporterr( "\n" );
+
+
+
+#if 0
+ if( localmem[1][0] == 13 )
+ {
+ reporterr( "OUTPUT!\n" );
+ for( i=0; i<clus1; i++ ) reporterr( ">g1\n%s\n", mseq1[i] );
+ for( i=0; i<clus2; i++ ) reporterr( ">g2\n%s\n", mseq2[i] );
+ exit( 1 );
+ }
+#endif
+
+// free( topol[l][0] ); topol[l][0] = NULL;
+// free( topol[l][1] ); topol[l][1] = NULL;
+// free( topol[l] ); topol[l] = NULL;
+
+
+// reporterr( ">514\n%s\n", aseq[514] );
+ free( localmem[0] );
+ free( localmem[1] );
+
+ }
+
+#if SCOREOUT
+ reporterr( "totalscore = %10.2f\n\n", tscore );
+#endif
+ Falign( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL, NULL, 0, NULL );
+ Falign_udpari_long( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
+ Falign_givenanchors( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
+ A__align( NULL, 0, 0, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0, -1, -1, NULL, NULL, NULL, 0.0, 0.0 );
+ D__align( NULL, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0 );
+ G__align11( NULL, NULL, NULL, 0, 0, 0 ); // iru?
+ free( effarr1 );
+ free( effarr2 );
+ free( indication1 );
+ free( indication2 );
+ free( fftlog );
+ if( specificityconsideration )
+ FreeDoubleMtx( dynamicmtx );
+ free( alreadyaligned );
+ free( localmem );
+ effarr1 = NULL;
+ return( 0 );
+
+ chudan_tbfast:
+
+ Falign( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL, NULL, 0, NULL );
+ Falign_udpari_long( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
+ Falign_givenanchors( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
+ A__align( NULL, 0, 0, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0, -1, -1, NULL, NULL, NULL, 0.0, 0.0 );
+ D__align( NULL, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0 );
+ G__align11( NULL, NULL, NULL, 0, 0, 0 ); // iru?
+ if( effarr1 ) free( effarr1 ); effarr1 = NULL;
+ if( effarr2 ) free( effarr2 ); effarr2 = NULL;
+ if( indication1 ) free( indication1 ); indication1 = NULL;
+ if( indication2 ) free( indication2 ); indication2 = NULL;
+ if( fftlog ) free( fftlog ); fftlog = NULL;
+ if( alreadyaligned ) free( alreadyaligned ); alreadyaligned = NULL;
+ if( specificityconsideration )
+ {
+ if( dynamicmtx ) FreeDoubleMtx( dynamicmtx ); dynamicmtx = NULL;
+ }
+ if( localmem ) free( localmem ); localmem = NULL;
+#if SKIP
+ if( skiptable1 ) FreeIntMtx( skiptable1 ); skiptable1 = NULL;
+ if( skiptable2 ) FreeIntMtx( skiptable2 ); skiptable2 = NULL;
+#endif
+
+ return( 1 );
+}
static int treebase( int *nlen, char **aseq, int nadd, char *mergeoralign, char **mseq1, char **mseq2, int ***topol, Treedep *dep, int **memhist, double ***cpmxhist, double *effarr, double **newdistmtx, int *selfscore, ExtAnch *extanch, int **anchindex, int *alloclen, int (*callback)(int, int, char*) )
{
int l, len1, len2, i, m, immin, immax;
@@ -2593,10 +3023,12 @@ static int treebase( int *nlen, char **aseq, int nadd, char *mergeoralign, char
}
// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
- if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000);
- else ffttry = 0;
-// ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 5000 && clus2 < 5000); // v6.708
+// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 && clus1 < 1000 && clus2 < 1000);
+// if( fftlog[m1] && fftlog[m2] ) ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
+// else ffttry = 0;
+ ffttry = ( nlen[m1] > clus1 && nlen[m2] > clus2 );
// reporterr( "f=%d, len1/fftlog[m1]=%f, clus1=%d, len2/fftlog[m2]=%f, clus2=%d\n", ffttry, (double)len1/fftlog[m1], clus1, (double)len2/fftlog[m2], clus2 );
+// reporterr( "fftlog=%d,%d, ffttry=%d\n", fftlog[m1], fftlog[m2], ffttry );
if( useexternalanchors )
{
@@ -3247,7 +3679,7 @@ int disttbfast( int ngui, int lgui, char **namegui, char **seqgui, int argc, cha
{
compacttree = 1;
treein = 0;
- use_fft = 0; // kankeinai?
+// use_fft = 0; // kankeinai?
// maxdistmtxsize = 5 * 1000 * 1000; // 5GB. ato de kahen ni suru.
// maxdistmtxsize = 1.0 * 1000 * 1000 * 1000; // 5GB. ato de kahen ni suru.
}
@@ -3255,13 +3687,13 @@ int disttbfast( int ngui, int lgui, char **namegui, char **seqgui, int argc, cha
{
compacttree = 4; // youngest linkage, 3 ha tbfast de tsukaunode ichiou sakeru
treein = 0;
- use_fft = 0; // kankeinai?
+// use_fft = 0; // kankeinai?
}
else if( treein == 'S' || treein == 'C' )
{
compacttree = 2; // 3 ha tbfast de tsukaunode ichiou sakeru
treein = 0;
- use_fft = 0; // kankeinai?
+// use_fft = 0; // kankeinai?
}
else if( treein == 'a' )
{
@@ -3297,7 +3729,8 @@ int disttbfast( int ngui, int lgui, char **namegui, char **seqgui, int argc, cha
nthreadtb = 0;
}
#endif
- if( njob > 10000 ) nthreadtb = 0;
+// if( njob > 10000 ) nthreadtb = 0;
+ if( njob > 20000 ) nthreadtb = 0;
// 2018/Jan. Hairetsu ga ooi toki
// 1. topolorder_lessargs no stack ga tarinakunaru
// 2. localcopy no tame kouritsu warui
@@ -3775,7 +4208,7 @@ int disttbfast( int ngui, int lgui, char **namegui, char **seqgui, int argc, cha
else if( treeout )
{
reporterr( "Constructing a UPGMA tree (treeout, efffree=%d) ... ", !calcpairdists );
- fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( njob, mtx, topol, len, name, nogaplen, dep, !calcpairdists );
+ fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( njob, mtx, topol, len, name, nogaplen, dep, !calcpairdists, treeout );
if( !calcpairdists )
{
FreeFloatHalfMtx( mtx, njob ); mtx = NULL;
@@ -4552,6 +4985,14 @@ int disttbfast( int ngui, int lgui, char **namegui, char **seqgui, int argc, cha
// use_getrusage();
}
+
+
+ if( oneiteration )
+ {
+ reporterr( "Iterative refinement (one vs others)\n" );
+ dooneiteration( nlen, bseq, nadd, mergeoralign, mseq1, mseq2, topol, dep, memhist, cpmxhist, eff, NULL, NULL, extanch, anchindex, &alloclen, callback );
+ }
+
#if DEBUG
reporterr( "closing trap_g\n" );
#endif
=====================================
core/functions.h
=====================================
@@ -326,7 +326,7 @@ extern int check_guidetreefile( int *seed, int *npick, double *limitram );
extern void createchain( int nseq, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int treeout, int shuffle, int seed );
//extern void loadtop( int nseq, double **eff, int ***topol, double **len );
extern void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *, int efffree ); // KESU
-extern void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *, int efffree );
+extern void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *, int efffree, int treeout );
extern void fixed_supg_double_realloc_nobk_halfmtx_treeout_constrained( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *, int ncons, int **constraints, int efffree );
extern void fixed_musclesupg_double_treeout( int nseq, double **eff, int ***topol, double **len, char **name );
extern void fixed_supg_double_treeout_constrained( int nseq, double **eff, int ***topol, double **len, char **name, int ncons, int **constraints );
=====================================
core/io.c
=====================================
@@ -70,7 +70,6 @@ static int countc( char *s, char q )
static void ttou( char *s )
{
- int v = 0;
while( *s )
{
if( *s == 't' ) *s = 'u';
=====================================
core/mafft.tmpl
=====================================
@@ -1,7 +1,7 @@
#! /bin/bash
er=0;
myself=`dirname "$0"`/`basename "$0"`; export myself
-version="v7.453 (2019/Nov/8)"; export version
+version="v7.464 (2020/Apr/21)"; export version
LANG=C; export LANG
os=`uname`
progname=`basename "$0"`
@@ -215,6 +215,8 @@ mccaskill=$defaultmccaskill
contrafold=$defaultcontrafold
progressfile="/dev/stderr"
anchorfile="/dev/null"
+anchoropt=""
+maxanchorseparation=1000
debug=0
sw=0
algopt=$defaultalgopt
@@ -237,6 +239,7 @@ partsize=50
partdist="ktuples"
partorderopt=" -x "
treeout=0
+nodeout=0
distout=0
treein=0
topin=0
@@ -379,6 +382,9 @@ if [ $# -gt 0 ]; then
partdist="fasta"
elif [ "$1" = "--treeout" ]; then
treeout=1
+ elif [ "$1" = "--nodeout" ]; then
+ nodeout=1
+ treeout=1
elif [ "$1" = "--distout" ]; then
distout=1
elif [ "$1" = "--fastswpair" ]; then
@@ -521,6 +527,8 @@ if [ $# -gt 0 ]; then
adjustdirection=1
elif [ "$1" = "--adjustdirectionaccurately" ]; then
adjustdirection=2
+ elif [ "$1" = "--oneiteration" ]; then
+ oneiterationopt=" -r "
elif [ "$1" = "--progress" ]; then
shift
progressfile="$1"
@@ -531,6 +539,13 @@ if [ $# -gt 0 ]; then
elif [ "$1" = "--out" ]; then
shift
outputfile="$1"
+ elif [ "$1" = "--skipanchorsremoterthan" ]; then
+ shift
+ if ! expr "$1" : "[0-9]" > /dev/null ; then
+ echo "Specify maximum gap length between anchors." 1>&2
+ exit
+ fi
+ maxanchorseparation=`expr "$1" - 0`
elif [ "$1" = "--anchors" ]; then
shift
anchorfile="$1"
@@ -1664,7 +1679,11 @@ function removetmpfile() { # for MPI
#new
if [ $cycle -eq 0 ]; then
- treeoutopt="-t -T"
+ if [ $nodeout -eq 1 ]; then
+ treeoutopt="-^ -T"
+ else
+ treeoutopt="-t -T"
+ fi
iterate=0
weighti="0.0" # 2016Jul31, tbfast.c kara idou
# if [ $distance = "global" -o $distance = "local" -o $distance = "localgenaf" -o $distance = "globalgenaf" ]; then # 2012/04, localpair --> local alignment distance
@@ -1684,7 +1703,14 @@ function removetmpfile() { # for MPI
fi
fi
else
- if [ $treeout -eq 1 ]; then
+ if [ $nodeout -eq 1 ]; then
+ if [ $iterate -gt 0 ]; then
+ echo "The --nodeout option supports only progressive method (--maxiterate 0) for now." 1>&2
+ exit 1
+ fi
+ parttreeoutopt="-t"
+ treeoutopt="-^"
+ elif [ $treeout -eq 1 ]; then
parttreeoutopt="-t"
treeoutopt="-t"
else
@@ -2547,7 +2573,7 @@ function removetmpfile() { # for MPI
if [ $fragment -ne 0 ]; then
"$prefix/addsingle" -Q 100 $legacygapopt -W $tuplesize -O $outnum $addsinglearg $addarg $add2ndhalfarg -C $numthreads $memopt $weightopt $treeinopt $treeoutopt $distoutopt $seqtype $model -f "-"$gop -h $aof $param_fft $localparam $algopt $treealg $scoreoutarg < infile > /dev/null 2>>"$progressfile" || exit 1
else
- "$prefix/disttbfast" -q $npickup -E $cycledisttbfast -V "-"$gopdist -s $unalignlevel $legacygapopt $mergearg -W $tuplesize $termgapopt $outnum $addarg $add2ndhalfarg -C $numthreads-$numthreadstb $memopt $weightopt $treeinopt $treeoutopt $distoutopt $seqtype $model -g $gexp -f "-"$gop -Q $spfactor -h $aof $param_fft $algopt $treealg $scoreoutarg $anchoropt < infile > pre 2>>"$progressfile" || exit 1
+ "$prefix/disttbfast" -q $npickup -E $cycledisttbfast -V "-"$gopdist -s $unalignlevel $legacygapopt $mergearg -W $tuplesize $termgapopt $outnum $addarg $add2ndhalfarg -C $numthreads-$numthreadstb $memopt $weightopt $treeinopt $treeoutopt $distoutopt $seqtype $model -g $gexp -f "-"$gop -Q $spfactor -h $aof $param_fft $algopt $treealg $scoreoutarg $anchoropt -x $maxanchorseparation $oneiterationopt < infile > pre 2>>"$progressfile" || exit 1
mv hat3.seed hat3
fi
fi
=====================================
core/mltaln.h
=====================================
@@ -37,7 +37,7 @@
-#define VERSION "7.453"
+#define VERSION "7.464"
#define SHOWVERSION reporterr( "%s (%s) Version " VERSION "\nalg=%c, model=%s, amax=%3.1f\n%d thread(s)\n\n", progName( argv[0] ), (dorp=='d')?"nuc":((nblosum==-2)?"text":"aa"), alg, modelname, specificityconsideration, nthread )
#define FFT_THRESHOLD 80
=====================================
core/mltaln9.c
=====================================
@@ -1347,6 +1347,54 @@ static void msaresetnearest( int nseq, Bchain *acpt, double **distfrompt, double
}
#endif
+static int getdensest( int *m, double *d )
+{
+ int i;
+ double dmax = -100.0;
+ int pmax = -1;
+ for( i=0; m[i]>-1; i++ )
+ {
+ if( d[m[i]] > dmax )
+ {
+ dmax = d[m[i]];
+ pmax = m[i];
+ }
+ }
+ return( pmax );
+}
+
+static void setdensity( int nseq, Bchain *acpt, double **eff, double *density, int pos )
+{
+ int j;
+ double tmpdouble;
+// double **effptpt;
+ Bchain *acptj;
+
+// printf( "[%d], %f, dist=%d ->", pos, *mindisfrompt, *nearestpt );
+
+// if( (acpt+pos)->next ) effpt = eff[pos]+(acpt+pos)->next->pos-pos;
+
+ tmpdouble = 0.0;
+// for( j=pos+1; j<nseq; j++ )
+ for( acptj=(acpt+pos)->next; acptj!=NULL; acptj=acptj->next )
+ {
+ j = acptj->pos;
+ if( eff[pos][j-pos] < 1.0 )
+ tmpdouble += (2.0-eff[pos][j-pos]);
+ }
+// effptpt = eff;
+// for( j=0; j<pos; j++ )
+ for( acptj=acpt; (acptj&&acptj->pos!=pos); acptj=acptj->next )
+ {
+ j = acptj->pos;
+ if( eff[j][pos-j] < 1.0 )
+ tmpdouble += (2.0-eff[j][pos-j]);
+ }
+
+ *density = tmpdouble;
+// printf( "p=%d, d=%f \n", pos, *density );
+}
+
static void setnearest( int nseq, Bchain *acpt, double **eff, double *mindisfrompt, int *nearestpt, int pos )
{
int j;
@@ -6079,7 +6127,7 @@ void compacttree_memsaveselectable( int nseq, double **partmtx, int *nearest, do
if( result ) free( result );
}
-void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree )
+void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree, int treeout )
{
int i, j, k, miniim, maxiim, minijm, maxijm;
@@ -6105,6 +6153,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
FILE *fp;
double (*clusterfuncpt[1])(double,double);
char namec;
+ double *density;
sueff1 = 1 - (double)sueff_global;
@@ -6134,6 +6183,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
nametmp = AllocateCharVec( 1000 ); // nagasugi
// tree = AllocateCharMtx( njob, njob*600 );
tree = AllocateCharMtx( njob, 0 );
+ if( treeout == 2 ) density = AllocateDoubleVec( njob );
}
@@ -6176,6 +6226,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
ac[nseq-1].next = NULL;
for( i=0; i<nseq; i++ ) setnearest( nseq, ac, eff, mindisfrom+i, nearest+i, i ); // muscle
+ if( treeout == 2 ) for( i=0; i<nseq; i++ ) setdensity( nseq, ac, eff, density+i, i );
for( i=0; i<nseq; i++ ) tmptmplen[i] = 0.0;
@@ -6435,18 +6486,36 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
#endif
-#if 0
- printf( "\nooSTEP-%03d:\n", k+1 );
- printf( "len0 = %f\n", len[k][0] );
- for( i=0; topol[k][0][i]>-1; i++ ) printf( " %03d", topol[k][0][i]+1 );
- printf( "\n" );
- printf( "len1 = %f\n", len[k][1] );
- for( i=0; topol[k][1][i]>-1; i++ ) printf( " %03d", topol[k][1][i]+1 );
- printf( "\n" );
-#endif
}
fp = fopen( "infile.tree", "w" );
fprintf( fp, "%s;\n", treetmp );
+ if( treeout == 2 )
+ {
+ int *mem = calloc( sizeof( int ), nseq );
+ fprintf( fp, "\nDensity:" );
+ for( k=0; k<nseq; k++ ) fprintf( fp, "\nSequence %d, %7.4f", k+1, density[k] );
+
+ fprintf( fp, "\n\nNode info:" );
+ for( k=0; k<nseq-1; k++ )
+ {
+ if( dep )
+ fprintf( fp, "\nNode %d, Height=%f\n", k+1, dep[k].distfromtip );
+// fprintf( fp, "len0 = %f\n", len[k][0] );
+ topolorderz( mem, topol, dep, k, 0 );
+// for( i=0; topol[k][0][i]>-1; i++ ) fprintf( fp, " %03d", topol[k][0][i]+1 );
+ fprintf( fp, "%d:", getdensest( mem, density )+1 );
+ for( i=0; mem[i]>-1; i++ ) fprintf( fp, " %d", mem[i]+1 );
+ fprintf( fp, "\n" );
+
+ topolorderz( mem, topol, dep, k, 1 );
+// fprintf( fp, "len1 = %f\n", len[k][1] );
+// for( i=0; topol[k][1][i]>-1; i++ ) fprintf( fp, " %03d", topol[k][1][i]+1 );
+ fprintf( fp, "%d:", getdensest( mem, density )+1 );
+ for( i=0; mem[i]>-1; i++ ) fprintf( fp, " %d", mem[i]+1 );
+ fprintf( fp, "\n" );
+ }
+ free( mem );
+ }
fclose( fp );
free( tree[0] );
@@ -6459,6 +6528,7 @@ void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( int nseq, dou
free( (void *)nmemar ); nmemar = NULL;
free( mindisfrom );
free( nearest );
+ if( treeout == 2 ) free( density );
}
void fixed_musclesupg_double_realloc_nobk_halfmtx_treeout( int nseq, double **eff, int ***topol, double **len, char **name, int *nlen, Treedep *dep, int efffree )
=====================================
core/tbfast.c
=====================================
@@ -348,6 +348,9 @@ static void arguments( int argc, char *argv[], int *pac, char **pav, int *tac, c
case 't':
treeout = 1;
break;
+ case '^':
+ treeout = 2;
+ break;
case 'T':
noalign = 1;
break;
@@ -1765,6 +1768,7 @@ void treebase( int *nlen, char **aseq, int nadd, char *mergeoralign, char **mseq
free( effarr1_kozo );
free( effarr2_kozo );
Falign( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL, NULL, 0, NULL );
+ Falign_udpari_long( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, NULL );
D__align( NULL, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0 );
A__align( NULL, 0, 0, NULL, NULL, NULL, NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, 0, 0, -1, -1, NULL, NULL, NULL, 0.0, 0.0 );
imp_match_init_strictD( NULL, 0, 0, 0, 0, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, NULL, NULL, NULL, NULL, NULL, 0, 0 );
@@ -2832,7 +2836,7 @@ int main( int argc, char *argv[] )
else if( treeout )
{
fprintf( stderr, "Constructing a UPGMA tree ... " );
- fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( njob, iscore, topol, len, name, nlen, dep, 1 ); // _memsave demo iihazu
+ fixed_musclesupg_double_realloc_nobk_halfmtx_treeout_memsave( njob, iscore, topol, len, name, nlen, dep, 1, treeout ); // _memsave demo iihazu
}
else
{
=====================================
debian/changelog
=====================================
@@ -1,3 +1,10 @@
+mafft (7.464-1) unstable; urgency=medium
+
+ * Team upload.
+ * New upstream version
+
+ -- Steffen Moeller <moeller at debian.org> Sun, 03 May 2020 19:21:55 +0200
+
mafft (7.453-2) unstable; urgency=medium
* Team upload
=====================================
readme
=====================================
@@ -1,6 +1,6 @@
-----------------------------------------------------------------------
MAFFT: a multiple sequence alignment program
- version 7.453, 2019/Nov/8
+ version 7.464, 2020/Apr/21
http://mafft.cbrc.jp/alignment/software/
katoh at ifrec.osaka-u.ac.jp
View it on GitLab: https://salsa.debian.org/med-team/mafft/-/compare/3412a690ee609bd85a394882a27f315fa8ea0e83...449121169168d0756221667a6cf67fabb0001c35
--
View it on GitLab: https://salsa.debian.org/med-team/mafft/-/compare/3412a690ee609bd85a394882a27f315fa8ea0e83...449121169168d0756221667a6cf67fabb0001c35
You're receiving this email because of your account on salsa.debian.org.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/debian-med-commit/attachments/20200503/6c4d6dc7/attachment-0001.html>
More information about the debian-med-commit
mailing list