binaries/src/mafft/core/io.c.fast

   1 #include "mltaln.h"
   2
   3 static int upperCase = 0;
   4
   5 #define DEBUG   0
   6 #define IODEBUG 0
   7
   8
   9 int addlocalhom2( char *al1, char *al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa, int skip )
  10 {
  11         int pos1, pos2, start1, start2, end1, end2;
  12         char *pt1, *pt2;
  13         double score;
  14         double sumscore;
  15         int sumoverlap;
  16         LocalHom *tmppt;
  17         int st;
  18         int nlocalhom = 0;
  19         pt1 = al1; pt2 = al2;
  20         pos1 = off1; pos2 = off2;
  21
  22         sumscore = 0.0;
  23         sumoverlap = 0;
  24
  25 #if 0
  26         fprintf( stderr, "nlocalhom = %d in addlocalhom\n", nlocalhom );
  27         fprintf( stderr, "al1 = %s, al2 = %s\n", al1, al2 );
  28         fprintf( stderr, "off1 = %d, off2 = %d\n", off1, off2 );
  29         fprintf( stderr, "localhopt = %p, skip = %d\n", localhompt, skip );
  30         fprintf( stderr, "pt1 = \n%s\n, pt2 = \n%s\n", pt1, pt2 );
  31 #endif
  32
  33         if( skip )
  34         {
  35                 while( --skip > 0 ) localhompt = localhompt->next;
  36                 localhompt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
  37                 localhompt = localhompt->next;
  38 //              fprintf( stderr, "tmppt = %p, localhompt = %p\n", tmppt, localhompt );
  39         }
  40         tmppt = localhompt;
  41
  42         st = 0;
  43         score = 0.0;
  44         while( *pt1 != 0 )
  45         {
  46 //              fprintf( stderr, "In in while loop\n" );
  47 //              fprintf( stderr, "pt = %c, %c, st=%d\n", *pt1, *pt2, st );
  48                 if( st == 1 && ( *pt1 == '-' || *pt2 == '-' ) )
  49                 {
  50                         end1 = pos1 - 1;
  51                         end2 = pos2 - 1;
  52
  53                         if( nlocalhom++ > 0 )
  54                         {
  55 //                              fprintf( stderr, "reallocating ...\n" );
  56                                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
  57 //                              fprintf( stderr, "done\n" );
  58                                 tmppt = tmppt->next;
  59                                 tmppt->next = NULL;
  60                         }
  61                         tmppt->start1 = start1;
  62                         tmppt->start2 = start2;
  63                         tmppt->end1   = end1  ;
  64                         tmppt->end2   = end2  ;
  65
  66 #if 1
  67                         sumscore += score;
  68                         sumoverlap += end2-start2+1;
  69 #else
  70                         tmppt->overlapaa   = end2-start2+1;
  71                         tmppt->opt = score * 5.8 / 600;
  72                         tmppt->overlapaa   = overlapaa;
  73                         tmppt->opt = (double)opt;
  74 #endif
  75
  76 #if 0
  77                         fprintf( stderr, "score (1)= %f\n", score );
  78                         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
  79                         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
  80 #endif
  81                         score = 0.0;
  82                         st = 0;
  83                 }
  84                 else if( *pt1 != '-' && *pt2 != '-' )
  85                 {
  86                         if( st == 0 )
  87                         {
  88                                 start1 = pos1; start2 = pos2;
  89                                 st = 1;
  90                         }
  91                         score += n_dis[(int)amino_n[(int)*pt1]][(int)amino_n[(int)*pt2]];
  92 //                      fprintf( stderr, "%c-%c, score(0) = %f\n", *pt1, *pt2, score );
  93                 }
  94                 if( *pt1++ != '-' ) pos1++;
  95                 if( *pt2++ != '-' ) pos2++;
  96         }
  97
  98         if( st )
  99         {
 100                 if( nlocalhom++ > 0 )
 101                 {
 102 //                      fprintf( stderr, "reallocating ...\n" );
 103                         tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 104 //                      fprintf( stderr, "done\n" );
 105                         tmppt = tmppt->next;
 106                         tmppt->next = NULL;
 107                 }
 108                 end1 = pos1 - 1;
 109                 end2 = pos2 - 1;
 110                 tmppt->start1 = start1;
 111                 tmppt->start2 = start2;
 112                 tmppt->end1   = end1  ;
 113                 tmppt->end2   = end2  ;
 114
 115 #if 1
 116                 sumscore += score;
 117                 sumoverlap += end2-start2+1;
 118 #else
 119                 tmppt->overlapaa   = end2-start2+1;
 120                 tmppt->opt = score * 5.8 / 600;
 121                 tmppt->overlapaa   = overlapaa;
 122                 tmppt->opt = (double)opt;
 123 #endif
 124 #if 0
 125                 fprintf( stderr, "score (2)= %f\n", score );
 126                 fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 127                 fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 128 #endif
 129         }
 130
 131         for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
 132         {
 133                 tmppt->overlapaa = sumoverlap;
 134                 tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
 135         }
 136         return( nlocalhom );
 137 }
 138
 139
 140
 141 int addlocalhom( char *al1, char *al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa, int skip )
 142 {
 143         int pos1, pos2, start1, start2, end1, end2;
 144         char *pt1, *pt2;
 145         double score;
 146         double sumscore;
 147         int sumoverlap;
 148         LocalHom *tmppt;
 149         int st;
 150         int nlocalhom = 0;
 151         pt1 = al1; pt2 = al2;
 152         pos1 = off1; pos2 = off2;
 153
 154         sumscore = 0.0;
 155         sumoverlap = 0;
 156
 157 #if 1
 158         fprintf( stderr, "nlocalhom = %d in addlocalhom\n", nlocalhom );
 159         fprintf( stderr, "al1 = %s, al2 = %s\n", al1, al2 );
 160         fprintf( stderr, "off1 = %d, off2 = %d\n", off1, off2 );
 161         fprintf( stderr, "localhopt = %p, skip = %d\n", localhompt, skip );
 162 #endif
 163         fprintf( stderr, "pt1 = \n%s\n, pt2 = \n%s\n", pt1, pt2 );
 164
 165         if( skip )
 166         {
 167                 while( --skip > 0 ) localhompt = localhompt->next;
 168                 localhompt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 169                 localhompt = localhompt->next;
 170                 fprintf( stderr, "tmppt = %p, localhompt = %p\n", tmppt, localhompt );
 171         }
 172         tmppt = localhompt;
 173
 174         st = 0;
 175         score = 0.0;
 176         while( *pt1 != 0 )
 177         {
 178                 fprintf( stderr, "In in while loop\n" );
 179                 fprintf( stderr, "pt = %c, %c, st=%d\n", *pt1, *pt2, st );
 180                 if( st == 1 && ( *pt1 == '-' || *pt2 == '-' ) )
 181                 {
 182                         end1 = pos1 - 1;
 183                         end2 = pos2 - 1;
 184
 185                         if( nlocalhom++ > 0 )
 186                         {
 187 //                              fprintf( stderr, "reallocating ...\n" );
 188                                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 189 //                              fprintf( stderr, "done\n" );
 190                                 tmppt = tmppt->next;
 191                                 tmppt->next = NULL;
 192                         }
 193                         tmppt->start1 = start1;
 194                         tmppt->start2 = start2;
 195                         tmppt->end1   = end1  ;
 196                         tmppt->end2   = end2  ;
 197
 198 #if 1
 199                         sumscore += score;
 200                         sumoverlap += end2-start2+1;
 201 #else
 202                         tmppt->overlapaa   = end2-start2+1;
 203                         tmppt->opt = score * 5.8 / 600;
 204                         tmppt->overlapaa   = overlapaa;
 205                         tmppt->opt = (double)opt;
 206 #endif
 207
 208                         fprintf( stderr, "score (1)= %f\n", score );
 209                         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 210                         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 211                         score = 0.0;
 212                         st = 0;
 213                 }
 214                 else if( *pt1 != '-' && *pt2 != '-' )
 215                 {
 216                         if( st == 0 )
 217                         {
 218                                 start1 = pos1; start2 = pos2;
 219                                 st = 1;
 220                         }
 221                         score += n_dis[(int)amino_n[(int)*pt1]][(int)amino_n[(int)*pt2]];
 222 //                      fprintf( stderr, "%c-%c, score(0) = %f\n", *pt1, *pt2, score );
 223                 }
 224                 if( *pt1++ != '-' ) pos1++;
 225                 if( *pt2++ != '-' ) pos2++;
 226         }
 227         if( nlocalhom++ > 0 )
 228         {
 229 //              fprintf( stderr, "reallocating ...\n" );
 230                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 231 //              fprintf( stderr, "done\n" );
 232                 tmppt = tmppt->next;
 233                 tmppt->next = NULL;
 234         }
 235         end1 = pos1 - 1;
 236         end2 = pos2 - 1;
 237         tmppt->start1 = start1;
 238         tmppt->start2 = start2;
 239         tmppt->end1   = end1  ;
 240         tmppt->end2   = end2  ;
 241
 242 #if 1
 243         sumscore += score;
 244         sumoverlap += end2-start2+1;
 245 #else
 246         tmppt->overlapaa   = end2-start2+1;
 247         tmppt->opt = score * 5.8 / 600;
 248         tmppt->overlapaa   = overlapaa;
 249         tmppt->opt = (double)opt;
 250 #endif
 251
 252         fprintf( stderr, "score (2)= %f\n", score );
 253         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 254         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 255
 256         for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
 257         {
 258                 tmppt->overlapaa = sumoverlap;
 259                 tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
 260         }
 261         return( nlocalhom );
 262 }
 263 void putlocalhom2( char *al1, char *al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
 264 {
 265         int pos1, pos2, start1, start2, end1, end2;
 266         char *pt1, *pt2;
 267         double score;
 268         double sumscore;
 269         int sumoverlap;
 270         LocalHom *tmppt = localhompt;
 271         int nlocalhom = 0;
 272         int st;
 273         pt1 = al1; pt2 = al2;
 274         pos1 = off1; pos2 = off2;
 275
 276
 277         sumscore = 0.0;
 278         sumoverlap = 0;
 279
 280         st = 0;
 281         score = 0.0;
 282         while( *pt1 != 0 )
 283         {
 284 //              fprintf( stderr, "pt = %c, %c, st=%d\n", *pt1, *pt2, st );
 285                 if( st == 1 && ( *pt1 == '-' || *pt2 == '-' ) )
 286                 {
 287                         end1 = pos1 - 1;
 288                         end2 = pos2 - 1;
 289
 290                         if( nlocalhom++ > 0 )
 291                         {
 292 //                              fprintf( stderr, "reallocating ...\n" );
 293                                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 294 //                              fprintf( stderr, "done\n" );
 295                                 tmppt = tmppt->next;
 296                                 tmppt->next = NULL;
 297                         }
 298                         tmppt->start1 = start1;
 299                         tmppt->start2 = start2;
 300                         tmppt->end1   = end1  ;
 301                         tmppt->end2   = end2  ;
 302
 303 #if 1
 304                         if( divpairscore )
 305                         {
 306                                 tmppt->overlapaa   = end2-start2+1;
 307                                 tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
 308                         }
 309                         else
 310                         {
 311                                 sumscore += score;
 312                                 sumoverlap += end2-start2+1;
 313                         }
 314 #else
 315                         tmppt->overlapaa   = overlapaa;
 316                         tmppt->opt = (double)opt;
 317 #endif
 318
 319 #if 0
 320                         fprintf( stderr, "score (1)= %f\n", score );
 321                         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 322                         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 323 #endif
 324                         score = 0.0;
 325                         st = 0;
 326                 }
 327                 else if( *pt1 != '-' && *pt2 != '-' )
 328                 {
 329                         if( st == 0 )
 330                         {
 331                                 start1 = pos1; start2 = pos2;
 332                                 st = 1;
 333                         }
 334                         score += n_dis[(int)amino_n[(int)*pt1]][(int)amino_n[(int)*pt2]]; // - offset ¤Ï¤¤¤é¤Ê¤¤¤«¤â
 335 //                      fprintf( stderr, "%c-%c, score(0) = %f\n", *pt1, *pt2, score );
 336                 }
 337                 if( *pt1++ != '-' ) pos1++;
 338                 if( *pt2++ != '-' ) pos2++;
 339         }
 340         if( *(pt1-1) != '-' && *(pt2-1) != '-'  )
 341         {
 342                 if( nlocalhom++ > 0 )
 343                 {
 344 //                      fprintf( stderr, "reallocating ...\n" );
 345                         tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 346 //                      fprintf( stderr, "done\n" );
 347                         tmppt = tmppt->next;
 348                         tmppt->next = NULL;
 349                 }
 350                 end1 = pos1 - 1;
 351                 end2 = pos2 - 1;
 352                 tmppt->start1 = start1;
 353                 tmppt->start2 = start2;
 354                 tmppt->end1   = end1  ;
 355                 tmppt->end2   = end2  ;
 356
 357 #if 1
 358                 if( divpairscore )
 359                 {
 360                         tmppt->overlapaa   = end2-start2+1;
 361                         tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
 362                 }
 363                 else
 364                 {
 365                         sumscore += score;
 366                         sumoverlap += end2-start2+1;
 367                 }
 368 #else
 369                 tmppt->overlapaa   = overlapaa;
 370                 tmppt->opt = (double)opt;
 371 #endif
 372
 373 #if 0
 374                 fprintf( stderr, "score (2)= %f\n", score );
 375                 fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 376                 fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 377 #endif
 378         }
 379
 380         if( !divpairscore )
 381         {
 382                 for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
 383                 {
 384                         tmppt->overlapaa = sumoverlap;
 385                         tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
 386 //                      fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
 387                 }
 388         }
 389 }
 390 void putlocalhom( char *al1, char *al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
 391 {
 392         int pos1, pos2, start1, start2, end1, end2;
 393         char *pt1, *pt2;
 394         double score;
 395         double sumscore;
 396         int sumoverlap;
 397         LocalHom *tmppt = localhompt;
 398         int nlocalhom = 0;
 399         int st;
 400         pt1 = al1; pt2 = al2;
 401         pos1 = off1; pos2 = off2;
 402
 403
 404         sumscore = 0.0;
 405         sumoverlap = 0;
 406
 407         st = 0;
 408         score = 0.0;
 409         while( *pt1 != 0 )
 410         {
 411 //              fprintf( stderr, "pt = %c, %c, st=%d\n", *pt1, *pt2, st );
 412                 if( st == 1 && ( *pt1 == '-' || *pt2 == '-' ) )
 413                 {
 414                         end1 = pos1 - 1;
 415                         end2 = pos2 - 1;
 416
 417                         if( nlocalhom++ > 0 )
 418                         {
 419 //                              fprintf( stderr, "reallocating ...\n" );
 420                                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 421 //                              fprintf( stderr, "done\n" );
 422                                 tmppt = tmppt->next;
 423                                 tmppt->next = NULL;
 424                         }
 425                         tmppt->start1 = start1;
 426                         tmppt->start2 = start2;
 427                         tmppt->end1   = end1  ;
 428                         tmppt->end2   = end2  ;
 429
 430 #if 1
 431                         if( divpairscore )
 432                         {
 433                                 tmppt->overlapaa   = end2-start2+1;
 434                                 tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
 435                         }
 436                         else
 437                         {
 438                                 sumscore += score;
 439                                 sumoverlap += end2-start2+1;
 440                         }
 441 #else
 442                         tmppt->overlapaa   = overlapaa;
 443                         tmppt->opt = (double)opt;
 444 #endif
 445
 446 #if 0
 447                         fprintf( stderr, "score (1)= %f\n", score );
 448                         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 449                         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 450 #endif
 451                         score = 0.0;
 452                         st = 0;
 453                 }
 454                 else if( *pt1 != '-' && *pt2 != '-' )
 455                 {
 456                         if( st == 0 )
 457                         {
 458                                 start1 = pos1; start2 = pos2;
 459                                 st = 1;
 460                         }
 461                         score += n_dis[(int)amino_n[(int)*pt1]][(int)amino_n[(int)*pt2]]; // - offset ¤Ï¤¤¤é¤Ê¤¤¤«¤â
 462 //                      fprintf( stderr, "%c-%c, score(0) = %f\n", *pt1, *pt2, score );
 463                 }
 464                 if( *pt1++ != '-' ) pos1++;
 465                 if( *pt2++ != '-' ) pos2++;
 466         }
 467         if( nlocalhom++ > 0 )
 468         {
 469 //              fprintf( stderr, "reallocating ...\n" );
 470                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
 471 //              fprintf( stderr, "done\n" );
 472                 tmppt = tmppt->next;
 473                 tmppt->next = NULL;
 474         }
 475         end1 = pos1 - 1;
 476         end2 = pos2 - 1;
 477         tmppt->start1 = start1;
 478         tmppt->start2 = start2;
 479         tmppt->end1   = end1  ;
 480         tmppt->end2   = end2  ;
 481
 482 #if 1
 483         if( divpairscore )
 484         {
 485                 tmppt->overlapaa   = end2-start2+1;
 486                 tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
 487         }
 488         else
 489         {
 490                 sumscore += score;
 491                 sumoverlap += end2-start2+1;
 492         }
 493 #else
 494         tmppt->overlapaa   = overlapaa;
 495         tmppt->opt = (double)opt;
 496 #endif
 497
 498 #if 0
 499         fprintf( stderr, "score (2)= %f\n", score );
 500         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
 501         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
 502 #endif
 503
 504         if( !divpairscore )
 505         {
 506                 for( tmppt=localhompt; tmppt; tmppt=tmppt->next )
 507                 {
 508                         tmppt->overlapaa = sumoverlap;
 509                         tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
 510 //                      fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
 511                 }
 512         }
 513 }
 514
 515 char *cutal( char *al, int al_display_start, int start, int end )
 516 {
 517         int pos;
 518         char *pt = al;
 519         char *val;
 520
 521         pos = al_display_start;
 522         do
 523         {
 524                 if( start == pos ) val = pt;
 525                 if( end == pos ) break;
 526 //              fprintf( stderr, "pos=%d, *pt=%c, val=%p\n", pos, *pt, val );
 527                 if( *pt != '-' ) pos++;
 528         } while( *pt++ != 0 );
 529         *(pt+1) = 0;
 530         return( val );
 531 }
 532
 533 void ErrorExit( char *message )
 534 {
 535         fprintf( stderr, "%s\n", message );
 536         exit( 1 );
 537 }
 538
 539 void strncpy_caseC( char *str1, char *str2, int len )
 540 {
 541         if( dorp == 'd' && upperCase > 0 )
 542         {
 543                 while( len-- )
 544                         *str1++ = toupper( *str2++ );
 545         }
 546         else strncpy( str1, str2, len );
 547 }
 548
 549 void seqUpper( int nseq, char **seq ) /* not used */
 550 {
 551         int i, j, len;
 552         for( i=0; i<nseq; i++ )
 553         {
 554                 len = strlen( seq[i] );
 555                 for( j=0; j<len; j++ )
 556                         seq[i][j] = toupper( seq[i][j] );
 557         }
 558 }
 559
 560 void seqLower( int nseq, char **seq )
 561 {
 562         int i, j, len;
 563         for( i=0; i<nseq; i++ )
 564         {
 565                 len = strlen( seq[i] );
 566                 for( j=0; j<len; j++ )
 567                         seq[i][j] = tolower( seq[i][j] );
 568         }
 569 }
 570
 571 int getaline_fp_eof( char *s, int l, FILE *fp )  /* end of file -> return 1 */
 572 {
 573     int c, i = 0 ;
 574     int noteofflag = 0;
 575     for( i=0; i<l && ( noteofflag = ( (c=getc(fp)) != EOF ) ) && c != '\n'; i++ )
 576         *s++ = c;
 577     *s = '\0' ;
 578      return( !noteofflag );
 579 }
 580
 581 int getaline_fp_eof_new(s, l, fp)  /* end of file -> return 1 */
 582 char    s[] ; int l ; FILE *fp ;
 583 {
 584         int c = 0, i = 0 ;
 585                 int noteofflag = 0;
 586
 587                 if( feof( fp ) ) return( 1 );
 588
 589                 for( i=0; i<l && ( noteofflag = ( (c=getc(fp)) != EOF ) ) && c != '\n'; i++ )
 590         { *s++ = c ; }
 591         *s = '\0' ;
 592                 if( c != '\n' && c != EOF ) while( getc(fp) != '\n' )
 593                         ;
 594                 return( !noteofflag );
 595 }
 596
 597 int myfgets(s, l, fp)  /* l°Ê¾å¤Ï¡¢¹ÔËö¤Þ¤ÇÆÉ¤ßÈô¤Ð¤¹ */
 598 char    s[] ; int l ; FILE *fp ;
 599 {
 600         int     c = 0, i = 0 ;
 601
 602                 if( feof( fp ) ) return( 1 );
 603
 604                 for( i=0; i<l && ( c=getc( fp ) ) != '\n'; i++ )
 605                 *s++ = c;
 606         *s = '\0' ;
 607                 if( c != '\n' )
 608                         while( getc(fp) != '\n' )
 609                                 ;
 610                 return( 0 );
 611 }
 612
 613 float input_new( FILE *fp, int d )
 614 {
 615         char mojiretsu[10];
 616         int i, c;
 617
 618         c = getc( fp );
 619         if( c != '\n' )
 620                 ungetc( c, fp );
 621
 622         for( i=0; i<d; i++ )
 623                 mojiretsu[i] = getc( fp );
 624         mojiretsu[i] = 0;
 625
 626         return( atof( mojiretsu ) );
 627 }
 628
 629
 630 void PreRead( FILE *fp, int *locnjob, int *locnlenmax )
 631 {
 632         int i, nleni;
 633         char b[B];
 634
 635         fgets( b, B-1, fp ); *locnjob = atoi( b );
 636         *locnlenmax = 0;
 637         i=0;
 638         while( i<*locnjob )
 639         {
 640                 fgets( b, B-1, fp );
 641                 if( !strncmp( b, "=", 1 ) )
 642                 {
 643                         i++;
 644                         fgets( b, B-1, fp ); nleni = atoi( b );
 645                         if( nleni > *locnlenmax ) *locnlenmax = nleni;
 646                 }
 647         }
 648         if( *locnlenmax > N )
 649         {
 650                 fprintf( stderr, "TOO LONG SEQUENCE!\n" );
 651                 exit( 1 );
 652         }
 653         if( njob > M  )
 654         {
 655                 fprintf( stderr, "TOO MANY SEQUENCE!\n" );
 656                 fprintf( stderr, "%d > %d\n", njob, M );
 657                 exit( 1 );
 658         }
 659 }
 660
 661 int allSpace( char *str )
 662 {
 663         int value = 1;
 664         while( *str ) value *= ( !isdigit( *str++ ) );
 665         return( value );
 666 }
 667
 668 void Read( char name[M][B], int nlen[M], char **seq )
 669 {
 670         extern void FRead( FILE *x, char y[M][B], int z[M], char **w );
 671         FRead( stdin, name, nlen, seq );
 672 }
 673
 674 void FRead( FILE *fp, char name[][B], int nlen[], char **seq )
 675 {
 676         int i, j;
 677         char b[B];
 678
 679         fgets( b, B-1, fp );
 680 #if DEBUG
 681         fprintf( stderr, "b = %s\n", b );
 682 #endif
 683
 684     if( strstr( b, "onnet" ) ) scoremtx = 1;
 685     else if( strstr( b, "DnA" ) )
 686         {
 687                 scoremtx = -1;
 688                 upperCase = -1;
 689         }
 690     else if( strstr( b, "dna" ) )
 691         {
 692                 scoremtx = -1;
 693                 upperCase = 0;
 694         }
 695         else if( strstr( b, "DNA" ) )
 696         {
 697                 scoremtx = -1;
 698                 upperCase = 1;
 699         }
 700     else if( strstr( b, "M-Y" ) || strstr( b, "iyata" ) ) scoremtx = 2;
 701     else scoremtx = 0;
 702 #if DEBUG
 703         fprintf( stderr, " %s->scoremtx = %d\n", b, scoremtx );
 704 #endif
 705
 706         geta2 = GETA2;
 707
 708 #if 0
 709         if( strlen( b ) >=25 )
 710         {
 711                 b[25] = 0;
 712         #if DEBUG
 713                 fprintf( stderr, "kimuraR = %s\n", b+20 );
 714         #endif
 715                 kimuraR = atoi( b+20 );
 716
 717                 if( kimuraR < 0 || 20 < kimuraR ) ErrorExit( "Illeagal kimuraR value.\n" );
 718                 if( allSpace( b+20 ) ) kimuraR = NOTSPECIFIED;
 719         }
 720         else kimuraR = NOTSPECIFIED;
 721         #if DEBUG
 722         fprintf( stderr, "kimuraR = %d\n", kimuraR );
 723         #endif
 724
 725         if( strlen( b ) >=20 )
 726         {
 727                 b[20] = 0;
 728         #if DEBUG
 729                 fprintf( stderr, "pamN = %s\n", b+15 );
 730         #endif
 731                 pamN = atoi( b+15 );
 732                 if( pamN < 0 || 400 < pamN ) ErrorExit( "Illeagal pam value.\n" );
 733                 if( allSpace( b+15 ) ) pamN = NOTSPECIFIED;
 734         }
 735         else pamN = NOTSPECIFIED;
 736
 737         if( strlen( b ) >= 15 )
 738         {
 739                 b[15] = 0;
 740         #if DEBUG
 741                 fprintf( stderr, "poffset = %s\n", b+10 );
 742         #endif
 743                 poffset = atoi( b+10 );
 744                 if( poffset > 500 ) ErrorExit( "Illegal extending gap ppenalty\n" );
 745                 if( allSpace( b+10 ) ) poffset = NOTSPECIFIED;
 746         }
 747         else poffset = NOTSPECIFIED;
 748
 749         if( strlen( b ) >= 10 )
 750         {
 751                 b[10] = 0;
 752         #if DEBUG
 753                 fprintf( stderr, "ppenalty = %s\n", b+5 );
 754         #endif
 755                 ppenalty = atoi( b+5 );
 756                 if( ppenalty > 0 ) ErrorExit( "Illegal opening gap ppenalty\n" );
 757                 if( allSpace( b+5 ) ) ppenalty = NOTSPECIFIED;
 758         }
 759         else ppenalty = NOTSPECIFIED;
 760 #endif
 761
 762         for( i=0; i<njob; i++ )
 763         {
 764                 getaline_fp_eof_new( b, B-1, fp );
 765                 strcpy( name[i], b );
 766 #if DEBUG
 767                 fprintf( stderr, "name[%d] = %s\n", i, name[i] );
 768 #endif
 769                 fgets( b, B-1, fp ); nlen[i] = atoi( b );      /* seq i no nagasa */
 770                 seq[i][0] = 0;
 771                 if( nlen[i] ) for( j=0; j <= (nlen[i]-1)/C; j++ )
 772                 {
 773                         getaline_fp_eof_new( b, B-1, fp );
 774                         /*      b[C] = 0;  */
 775                         strcat( seq[i], b );
 776                 }
 777                 seq[i][nlen[i]] = 0;
 778         }
 779         if( scoremtx == -1 && upperCase != -1 ) seqLower( njob, seq );
 780 }
 781
 782
 783 static int countKUorWA( FILE *fp )
 784 {
 785         int value;
 786         int c, b;
 787
 788         value= 0;
 789         b = '\n';
 790         while( ( c = getc( fp ) ) != EOF )
 791         {
 792                 if( b == '\n' && ( c == '=' || c == '>' ) )
 793                         value++;
 794                 b = c;
 795         }
 796         rewind( fp );
 797         return( value );
 798 }
 799
 800 static void searchKUorWA( FILE *fp )
 801 {
 802         int c, b;
 803         b = '\n';
 804         while( !( ( ( c = getc( fp ) ) == '>' || c == '=' || c == EOF ) && b == '\n' ) )
 805                 b = c;
 806         ungetc( c, fp );
 807 }
 808
 809 static int onlyAlpha_lower( char *str )
 810 {
 811         char tmp;
 812         char *res = str;
 813         char *bk = str;
 814
 815         while( (tmp=*str++) )
 816                 if( isalpha( tmp ) || tmp == '-' || tmp == '*' || tmp == '.' )
 817                         *res++ = tolower( tmp );
 818         *res = 0;
 819         return( res - bk );
 820 }
 821 static int onlyAlpha_upper( char *str )
 822 {
 823         char tmp;
 824         char *res = str;
 825         char *bk = str;
 826
 827         while( (tmp=*str++) )
 828                 if( isalpha( tmp ) || tmp == '-' || tmp == '*' || tmp == '.' )
 829                         *res++ = toupper( tmp );
 830         *res = 0;
 831         return( res - bk );
 832 }
 833
 834 void kake2hiku( char *str )
 835 {
 836         do
 837                 if( *str == '*' ) *str = '-';
 838         while( *str++ );
 839 }
 840
 841 int load1SeqWithoutName_new( FILE *fpp, char *cbuf )
 842 {
 843         int c, b;
 844         char *bk = cbuf;
 845
 846         b = '\n';
 847         while( ( c = getc( fpp ) ) != EOF &&                    /* by T. Nishiyama */
 848           !( ( c == '>' || c == '=' || c == '(' || c == EOF ) && b == '\n' ) )
 849         {
 850                 *cbuf++ = (char)c;  /* Ä¹¤¹¤®¤Æ¤â¤·¤é¤Ê¤¤ */
 851                 b = c;
 852         }
 853         ungetc( c, fpp );
 854         *cbuf = 0;
 855         if( dorp == 'd' )
 856                 onlyAlpha_lower( bk );
 857         else
 858                 onlyAlpha_upper( bk );
 859         kake2hiku( bk );
 860         return( 0 );
 861 }
 862
 863
 864 void readDataforgaln( FILE *fp, char **name, int *nlen, char **seq )
 865 {
 866         int i;
 867         static char *tmpseq = NULL;
 868
 869         if( !tmpseq )
 870         {
 871                 tmpseq = AllocateCharVec( N );
 872         }
 873
 874         rewind( fp );
 875         searchKUorWA( fp );
 876
 877         for( i=0; i<njob; i++ )
 878         {
 879                 name[i][0] = '='; getc( fp );
 880 #if 0
 881                 fgets( name[i]+1, B-2, fp );
 882                 j = strlen( name[i] );
 883                 if( name[i][j-1] != '\n' )
 884                         ErrorExit( "Too long name\n" );
 885                 name[i][j-1] = 0;
 886 #else
 887                 myfgets( name[i]+1, B-2, fp );
 888 #endif
 889 #if 0
 890                 fprintf( stderr, "name[%d] = %s\n", i, name[i] );
 891 #endif
 892                 load1SeqWithoutName_new( fp, tmpseq );
 893                 strcpy( seq[i], tmpseq );
 894                 nlen[i] = strlen( seq[i] );
 895         }
 896         if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
 897 #if 0
 898         free( tmpseq );
 899 #endif
 900 }
 901
 902 void readData( FILE *fp, char name[][B], int nlen[], char **seq )
 903 {
 904         int i;
 905         static char *tmpseq = NULL;
 906
 907         if( !tmpseq )
 908         {
 909                 tmpseq = AllocateCharVec( N );
 910         }
 911
 912         rewind( fp );
 913         searchKUorWA( fp );
 914
 915         for( i=0; i<njob; i++ )
 916         {
 917                 name[i][0] = '='; getc( fp );
 918 #if 0
 919                 fgets( name[i]+1, B-2, fp );
 920                 j = strlen( name[i] );
 921                 if( name[i][j-1] != '\n' )
 922                         ErrorExit( "Too long name\n" );
 923                 name[i][j-1] = 0;
 924 #else
 925                 myfgets( name[i]+1, B-2, fp );
 926 #endif
 927 #if 0
 928                 fprintf( stderr, "name[%d] = %s\n", i, name[i] );
 929 #endif
 930                 load1SeqWithoutName_new( fp, tmpseq );
 931                 strcpy( seq[i], tmpseq );
 932                 nlen[i] = strlen( seq[i] );
 933         }
 934         if( dorp == 'd' && upperCase != -1 ) seqLower( njob, seq );
 935 #if 0
 936         free( tmpseq );
 937 #endif
 938 }
 939
 940
 941 double countATGC( char *s )
 942 {
 943         int nATGC;
 944         int nChar;
 945         char c;
 946         nATGC = nChar = 0;
 947
 948         do
 949         {
 950                 c = tolower( *s );
 951                 if( isalpha( c ) )
 952                 {
 953                         nChar++;
 954                         if( c == 'a' || c == 't' || c == 'g' || c == 'c' || c == 'u' || c == 'n' )
 955                                 nATGC++;
 956                 }
 957         }
 958         while( *++s );
 959         return( (double)nATGC / nChar );
 960 }
 961
 962
 963 void getnumlen( FILE *fp )
 964 {
 965         int i, tmp;
 966         char *tmpseq;
 967         double atgcfreq;
 968         tmpseq = AllocateCharVec( N );
 969         njob = countKUorWA( fp );
 970         searchKUorWA( fp );
 971         nlenmax = 0;
 972         atgcfreq = 0.0;
 973         for( i=0; i<njob; i++ )
 974         {
 975                 fgets( tmpseq, N-1, fp );
 976                 load1SeqWithoutName_new( fp, tmpseq );
 977                 tmp = strlen( tmpseq );
 978                 if( tmp > nlenmax ) nlenmax  = tmp;
 979                 atgcfreq += countATGC( tmpseq );
 980         }
 981         atgcfreq /= (double)njob;
 982         if( dorp == NOTSPECIFIED )
 983         {
 984                 if( atgcfreq > 0.75 )
 985                 {
 986                         dorp = 'd';
 987                         upperCase = -1;
 988                 }
 989                 else
 990                 {
 991                         dorp = 'p';
 992                         upperCase = 0;
 993                 }
 994         }
 995         free( tmpseq );
 996 }
 997
 998
 999
1000 void WriteGapFill( FILE *fp, int locnjob, char name[][B], int nlen[M], char **aseq )
1001 {
1002         static char b[N];
1003         int i, j;
1004         int nalen[M];
1005         static char gap[N];
1006         static char buff[N];
1007
1008 #if IODEBUG
1009         fprintf( stderr, "IMAKARA KAKU\n" );
1010 #endif
1011         nlenmax = 0;
1012         for( i=0; i<locnjob; i++ )
1013         {
1014                 int len = strlen( aseq[i] );
1015                 if( nlenmax < len ) nlenmax = len;
1016         }
1017
1018         for( i=0; i<nlenmax; i++ ) gap[i] = '-';
1019         gap[nlenmax] = 0;
1020
1021         fprintf( fp, "%5d", locnjob );
1022         fprintf( fp, "\n" );
1023
1024         for( i=0; i<locnjob; i++ )
1025         {
1026                 strcpy( buff, aseq[i] );
1027                 strncat( buff, gap, nlenmax-strlen( aseq[i] ) );
1028                 buff[nlenmax] = 0;
1029                 nalen[i] = strlen( buff );
1030                 fprintf( fp, "%s\n", name[i] );
1031                 fprintf( fp, "%5d\n", nalen[i] );
1032                 for( j=0; j<nalen[i]; j=j+C )
1033                 {
1034                         strncpy_caseC( b, buff+j, C ); b[C] = 0;
1035                         fprintf( fp, "%s\n",b );
1036                 }
1037         }
1038 #if DEBUG
1039         fprintf( stderr, "nalen[0] = %d\n", nalen[0] );
1040 #endif
1041 #if IODEBUG
1042         fprintf( stderr, "KAKIOWATTA\n" );
1043 #endif
1044 }
1045
1046 void writeDataforgaln( FILE *fp, int locnjob, char **name, int *nlen, char **aseq )
1047 {
1048         int i, j;
1049         int nalen;
1050
1051         for( i=0; i<locnjob; i++ )
1052         {
1053                 nalen = strlen( aseq[i] );
1054                 fprintf( fp, ">%s\n", name[i]+1 );
1055                 for( j=0; j<nalen; j=j+C )
1056                 {
1057 #if 0
1058                         strncpy( b, aseq[i]+j, C ); b[C] = 0;
1059                         fprintf( fp, "%s\n",b );
1060 #else
1061                         fprintf( fp, "%.*s\n", C, aseq[i]+j );
1062 #endif
1063                 }
1064         }
1065 }
1066
1067 void writeData( FILE *fp, int locnjob, char name[][B], int nlen[], char **aseq )
1068 {
1069         int i, j;
1070         int nalen;
1071
1072         for( i=0; i<locnjob; i++ )
1073         {
1074 #if DEBUG
1075                 fprintf( stderr, "i = %d in writeData\n", i );
1076 #endif
1077                 nalen = strlen( aseq[i] );
1078                 fprintf( fp, ">%s\n", name[i]+1 );
1079                 for( j=0; j<nalen; j=j+C )
1080                 {
1081 #if 0
1082                         strncpy( b, aseq[i]+j, C ); b[C] = 0;
1083                         fprintf( fp, "%s\n",b );
1084 #else
1085                         fprintf( fp, "%.*s\n", C, aseq[i]+j );
1086 #endif
1087                 }
1088         }
1089 }
1090
1091
1092
1093
1094
1095 void readhat2_int( FILE *fp, int nseq, char name[M][B], int **mtx )
1096 {
1097     int i, j, nseq0;
1098     char b[B];
1099
1100     fgets( b, B, fp );
1101     fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
1102     fgets( b, B, fp );
1103     for( i=0; i<nseq; i++ )
1104     {
1105 #if 0
1106         getaline_fp_eof( b, B, fp );
1107 #else
1108                 myfgets( b, B-2, fp );
1109 #endif
1110 #if 0
1111                 j = MIN( strlen( b+6 ), 10 );
1112         if( strncmp( name[i], b+6 , j ) )
1113                 {
1114                         fprintf( stderr, "Error in hat2\n" );
1115                         fprintf( stderr, "%s != %s\n", b, name[i] );
1116                         exit( 1 );
1117                 }
1118 #endif
1119     }
1120     for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
1121     {
1122         mtx[i][j] = (int)( input_new( fp, D ) * INTMTXSCALE );
1123     }
1124 }
1125 void readhat2( FILE *fp, int nseq, char name[M][B], double **mtx )
1126 {
1127     int i, j, nseq0;
1128     char b[B];
1129
1130     fgets( b, B, fp );
1131     fgets( b, B, fp ); b[5] = 0; nseq0 = atoi( b ); if( nseq != nseq0 ) ErrorExit( "hat2 is wrong." );
1132     fgets( b, B, fp );
1133     for( i=0; i<nseq; i++ )
1134     {
1135 #if 0
1136         getaline_fp_eof( b, B, fp );
1137 #else
1138                 myfgets( b, B-2, fp );
1139 #endif
1140 #if 0
1141                 j = MIN( strlen( b+6 ), 10 );
1142         if( strncmp( name[i], b+6 , j ) )
1143                 {
1144                         fprintf( stderr, "Error in hat2\n" );
1145                         fprintf( stderr, "%s != %s\n", b, name[i] );
1146                         exit( 1 );
1147                 }
1148 #endif
1149     }
1150     for( i=0; i<nseq-1; i++ ) for( j=i+1; j<nseq; j++ )
1151     {
1152         mtx[i][j] = (double)input_new( fp, D);
1153     }
1154 }
1155
1156 void WriteFloatHat2( FILE *hat2p, int locnjob, char name[M][B], float **mtx )
1157 {
1158         int i, j;
1159         double max = 0.0;
1160         for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
1161
1162         fprintf( hat2p, "%5d\n", 1 );
1163         fprintf( hat2p, "%5d\n", locnjob );
1164         fprintf( hat2p, " %#6.3f\n", max * 2.5 );
1165
1166         for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
1167         for( i=0; i<locnjob-1; i++ )
1168         {
1169                 for( j=i+1; j<locnjob; j++ )
1170                 {
1171                         fprintf( hat2p, "%#6.3f", mtx[i][j] );
1172                         if( (j-i) % 12 == 0 || j == locnjob-1 ) fprintf( hat2p, "\n" );
1173                 }
1174         }
1175 }
1176
1177 void WriteHat2_int( FILE *hat2p, int locnjob, char name[M][B], int **mtx )
1178 {
1179         int i, j;
1180         double max = 0.0;
1181         for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
1182         max /= INTMTXSCALE;
1183
1184         fprintf( hat2p, "%5d\n", 1 );
1185         fprintf( hat2p, "%5d\n", locnjob );
1186         fprintf( hat2p, " %#6.3f\n", max * 2.5 );
1187
1188         for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
1189         for( i=0; i<locnjob-1; i++ )
1190         {
1191                 for( j=i+1; j<locnjob; j++ )
1192                 {
1193                         fprintf( hat2p, "%#6.3f", (float)mtx[i][j] / INTMTXSCALE );
1194                         if( (j-i) % 12 == 0 || j == locnjob-1 ) fprintf( hat2p, "\n" );
1195                 }
1196         }
1197 }
1198 void WriteHat2( FILE *hat2p, int locnjob, char name[M][B], double **mtx )
1199 {
1200         int i, j;
1201         double max = 0.0;
1202         for( i=0; i<locnjob-1; i++ ) for( j=i+1; j<locnjob; j++ ) if( mtx[i][j] > max ) max = mtx[i][j];
1203
1204         fprintf( hat2p, "%5d\n", 1 );
1205         fprintf( hat2p, "%5d\n", locnjob );
1206         fprintf( hat2p, " %#6.3f\n", max * 2.5 );
1207
1208         for( i=0; i<locnjob; i++ ) fprintf( hat2p, "%4d. %s\n", i+1, name[i] );
1209         for( i=0; i<locnjob-1; i++ )
1210         {
1211                 for( j=i+1; j<locnjob; j++ )
1212                 {
1213                         fprintf( hat2p, "%#6.3f", mtx[i][j] );
1214                         if( (j-i) % 12 == 0 || j == locnjob-1 ) fprintf( hat2p, "\n" );
1215                 }
1216         }
1217 }
1218
1219 int ReadFasta_sub( FILE *fp, double *dis, int nseq, char name[M][B] )
1220 {
1221     int i, count=0;
1222     char b[B];
1223     int junban[M];
1224
1225     count = 0;
1226     for( i=0; i<10000000 && count<nseq; i++ )
1227     {
1228         fgets( b, B-1, fp );
1229         if( !strncmp( "+==========+", b, 12 ) )
1230         {
1231             junban[count] = atoi( b+12 );
1232             count++;
1233         }
1234     }
1235
1236         for( i=0; i<nseq; i++ ) dis[i] = 0.0;
1237     count = 0;
1238     for( i=0; i<100000 && count<nseq; i++ )
1239     {
1240                 if( fgets( b, B-1, fp ) ) break;
1241         if( !strncmp( name[junban[count]], b, 20  ) )
1242         {
1243             fgets( b, B-1, fp );
1244             dis[junban[count]] = atof( b );
1245             count++;
1246         }
1247     }
1248     return 0;
1249 }
1250
1251
1252 int ReadSsearch( FILE *fp, double *dis, int nseq, char name[M][B] )
1253 {
1254     int i, count=0;
1255     char b[B];
1256     int junban[M];
1257         int opt;
1258
1259     count = 0;
1260     for( i=0; i<10000000 && count<nseq; i++ )
1261     {
1262         fgets( b, B-1, fp );
1263         if( !strncmp( "+==========+", b, 12 ) )
1264         {
1265             junban[count] = atoi( b+12 );
1266                         sscanf( b+75, "%d", &opt );
1267             dis[junban[count]] = (double)opt;
1268             count++;
1269         }
1270     }
1271
1272 /*
1273     count = 0;
1274     for( i=0; i<100000 && count<nseq; i++ )
1275     {
1276         fgets( b, B-1, fp );
1277         if( !strncmp( name[junban[count]], b, 20  ) )
1278         {
1279             dis[junban[count]] = atof( b+65 );
1280             count++;
1281         }
1282     }
1283 */
1284     return 0;
1285 }
1286
1287 int ReadBlastm7( FILE *fp, double *dis, int qmem, char name[M][B], LocalHom *localhomlist )
1288 {
1289     int count=0;
1290     char b[B];
1291         char *pt;
1292     static int junban[M];
1293         int overlapaa;
1294         double score, sumscore;
1295         int qstart, qend, tstart, tend;
1296         double z, bits;
1297         int qal_display_start, tal_display_start;
1298         static char qal[N], tal[N], al[N];
1299         char *qal2, *tal2;
1300         int c, i, nlocalhom;
1301
1302
1303
1304         count = 0;
1305         sumscore = 0.0;
1306         score = 0.0;
1307     while( 1 )
1308         {
1309
1310                 if( feof( fp ) ) break;
1311
1312                 while( fgets( b, B-1, fp ) )
1313                 {
1314                         if( !strncmp( "          <Hit_def>", b, 19 ) || !strncmp( "              <Hsp_num>", b, 23 ) ) break;
1315                 }
1316
1317                 if( !strncmp( "          <Hit_def>", b, 19 ) )
1318                 {
1319                         junban[count] = atoi( b+31 );
1320                         nlocalhom = 0;
1321                 }
1322
1323
1324                 while( fgets( b, B-1, fp ) )
1325                         if( !strncmp( "              <Hsp_score>", b, 25 ) ) break;
1326                 pt = b + 25;
1327                 score = atof( pt );
1328                 sumscore += score;
1329
1330
1331                 while( fgets( b, B-1, fp ) )
1332                         if( !strncmp( "              <Hsp_query-from>", b, 30 ) ) break;
1333                 pt = b + 30;
1334                 qstart = atoi( pt ) - 1;
1335
1336
1337                 while( fgets( b, B-1, fp ) )
1338                         if( !strncmp( "              <Hsp_query-to>", b, 28 ) ) break;
1339                 pt = b + 28;
1340                 qend = atoi( pt ) - 1;
1341
1342
1343                 while( fgets( b, B-1, fp ) )
1344                         if( !strncmp( "              <Hsp_hit-from>", b, 28 ) ) break;
1345                 pt = b + 28;
1346                 tstart = atoi( pt ) - 1;
1347
1348
1349                 while( fgets( b, B-1, fp ) )
1350                         if( !strncmp( "              <Hsp_hit-to>", b, 26 ) ) break;
1351                 pt = b + 26;
1352                 tend = atoi( pt ) - 1;
1353
1354
1355                 while( fgets( b, B-1, fp ) )
1356                         if( !strncmp( "              <Hsp_align-len>", b, 29 ) ) break;
1357                 pt = b + 29;
1358                 overlapaa = atoi( pt );
1359
1360
1361                 while( fgets( al, N-100, fp ) )
1362                         if( !strncmp( "              <Hsp_qseq>", al, 24 ) ) break;
1363
1364                 strcpy( qal, al+24 );
1365                 pt = qal;
1366                 while( *++pt != '<' )
1367                         ;
1368                 *pt = 0;
1369
1370
1371                 while( fgets( al, N-100, fp ) )
1372                         if( !strncmp( "              <Hsp_hseq>", al, 24 ) ) break;
1373
1374                 strcpy( tal, al+24 );
1375                 pt = tal;
1376                 while( *++pt != '<' )
1377                         ;
1378                 *pt = 0;
1379
1380
1381 //              fprintf( stderr, "t=%d, score = %f, qstart=%d, qend=%d, tstart=%d, tend=%d, overlapaa=%d\n", junban[count], score, qstart, qend, tstart, tend, overlapaa );
1382
1383                 nlocalhom += addlocalhom( qal, tal, localhomlist+junban[count], qstart, tstart, score, overlapaa, nlocalhom );
1384
1385                 while( fgets( b, B-1, fp ) )
1386                         if( !strncmp( "            </Hsp>:", b, 18 ) ) break;
1387
1388
1389                 fgets( b, B-1, fp );
1390
1391
1392                 if( !strncmp( "          </Hit_hsps>", b, 21 ) )
1393                 {
1394                         dis[junban[count++]] = sumscore;
1395                         sumscore = 0.0;
1396                         fgets( b, B-1, fp );
1397                         fgets( b, B-1, fp );
1398                         if( !strncmp( "      </Iteration_hits>", b, 23 ) ) break;
1399                 }
1400         }
1401     return count;
1402 }
1403
1404 int ReadFasta34noalign( FILE *fp, double *dis, int qmem, char name[M][B], LocalHom *localhomlist )
1405 {
1406     int count=0;
1407     char b[B];
1408         char *pt;
1409     static int junban[M];
1410         int overlapaa;
1411         int opt, qstart, qend, tstart, tend;
1412         double z, bits;
1413         int qal_display_start, tal_display_start;
1414         static char qal[N], tal[N];
1415         char *qal2, *tal2;
1416         int c;
1417
1418
1419     count = 0;
1420 #if 0
1421     for( i=0; i<10000000 && count<nseq; i++ )
1422 #else
1423     while( !feof( fp ) )
1424 #endif
1425     {
1426         fgets( b, B-1, fp );
1427         if( !strncmp( "+==========+", b, 12 ) )
1428         {
1429             junban[count] = atoi( b+12 );
1430
1431                         pt = strchr( b, ')' ) + 1;
1432                         sscanf( pt, "%d %lf %lf",  &opt, &bits, &z );
1433             dis[junban[count]] = (double)opt;
1434             count++;
1435
1436         }
1437     }
1438
1439     return count;
1440 }
1441 int ReadFasta34m10_nuc( FILE *fp, double *dis, int qmem, char name[M][B], LocalHom *localhomlist )
1442 {
1443     int count=0;
1444     char b[B];
1445         char *pt;
1446     static int junban[M];
1447         int overlapaa;
1448         int opt, qstart, qend, tstart, tend;
1449         double z, bits;
1450         int qal_display_start, tal_display_start;
1451         static char qal[N], tal[N];
1452         char *qal2, *tal2;
1453         int c;
1454
1455
1456     count = 0;
1457 #if 0
1458     for( i=0; i<10000000 && count<nseq; i++ )
1459 #else
1460     while( !feof( fp ) )
1461 #endif
1462     {
1463         fgets( b, B-1, fp );
1464         if( !strncmp( "+==========+", b, 12 ) )
1465         {
1466             junban[count] = atoi( b+12 );
1467
1468                         if( strchr( b, 'r' ) ) continue;
1469
1470                         pt = strchr( b, ']' ) + 1;
1471                         sscanf( pt, "%d %lf %lf",  &opt, &bits, &z );
1472             dis[junban[count]] = (double)opt;
1473             count++;
1474
1475         }
1476                 else if( 0 == strncmp( ">>+==========+", b, 14 ) )
1477                 {
1478                         break;
1479                 }
1480
1481     }
1482         if( !count ) return -1;
1483
1484         count = 0;
1485     while( 1 )
1486         {
1487                 if( strncmp( ">>+==========+", b, 14 ) )
1488                 {
1489                         fgets( b, B-1, fp );
1490                         if( feof( fp ) ) break;
1491                         continue;
1492                 }
1493                 junban[count++] = atoi( b+14 );
1494 //              fprintf( stderr, "t = %d\n", atoi( b+14 ) );
1495                 while( fgets( b, B-1, fp ) )
1496                         if( !strncmp( "; fa_opt:", b, 9 ) || !strncmp( "; sw_s-w opt:", b, 13 ) ) break;
1497                 pt = strstr( b, ":" ) +1;
1498                 opt = atoi( pt );
1499
1500
1501                 while( fgets( b, B-1, fp ) )
1502                         if( !strncmp( "_overlap:", b+4, 9 ) ) break;
1503                 pt = strstr( b, ":" ) +1;
1504                 overlapaa = atoi( pt );
1505
1506                 while( fgets( b, B-1, fp ) )
1507                         if( !strncmp( "_start:", b+4, 7 ) ) break;
1508                 pt = strstr( b, ":" ) +1;
1509                 qstart = atoi( pt ) - 1;
1510
1511                 while( fgets( b, B-1, fp ) )
1512                         if( !strncmp( "_stop:", b+4, 6 ) ) break;
1513                 pt = strstr( b, ":" ) +1;
1514                 qend = atoi( pt ) - 1;
1515
1516                 while( fgets( b, B-1, fp ) )
1517                         if( !strncmp( "_display_start:", b+4, 15 ) ) break;
1518                 pt = strstr( b, ":" ) +1;
1519                 qal_display_start = atoi( pt ) - 1;
1520
1521                 pt = qal;
1522                 while( c = fgetc( fp ) )
1523                 {
1524                         if( c == '>' )
1525                         {
1526                                 ungetc( c, fp );
1527                                 break;
1528                         }
1529                         if( isalpha( c ) || c == '-' )
1530                         *pt++ = c;
1531                 }
1532                 *pt = 0;
1533
1534                 while( fgets( b, B-1, fp ) )
1535                         if( !strncmp( "_start:", b+4, 7 ) ) break;
1536                 pt = strstr( b, ":" ) + 1;
1537                 tstart = atoi( pt ) - 1;
1538
1539                 while( fgets( b, B-1, fp ) )
1540                         if( !strncmp( "_stop:", b+4, 6 ) ) break;
1541                 pt = strstr( b, ":" ) + 1;
1542                 tend = atoi( pt ) - 1;
1543
1544                 while( fgets( b, B-1, fp ) )
1545                         if( !strncmp( "_display_start:", b+4, 15 ) ) break;
1546                 pt = strstr( b, ":" ) + 1;
1547                 tal_display_start = atoi( pt ) - 1;
1548
1549                 pt = tal;
1550                 while( c = fgetc( fp ) )
1551                 {
1552                         if( c == '>' )
1553                         {
1554                                 ungetc( c, fp );
1555                                 break;
1556                         }
1557                         if( isalpha( c ) || c == '-' )
1558                         *pt++ = c;
1559                 }
1560                 *pt = 0;
1561
1562 //              fprintf( stderr, "(%d-%d:%d-%d)\n", qstart, qend, tstart, tend );
1563 //              fprintf( stderr, "qal_display_start = %d, tal_display_start = %d\n", qal_display_start, tal_display_start );
1564
1565 //              fprintf( stderr, "qal = %s\n", qal );
1566 //              fprintf( stderr, "tal = %s\n", tal );
1567
1568                 qal2 = cutal( qal, qal_display_start, qstart, qend );
1569                 tal2 = cutal( tal, tal_display_start, tstart, tend );
1570
1571 //              fprintf( stderr, "qal2 = %s\n", qal2 );
1572 //              fprintf( stderr, "tal2 = %s\n", tal2 );
1573
1574 //              fprintf( stderr, "putting   %d - %d, opt = %d\n", qmem, junban[count-1], opt );
1575                 putlocalhom( qal2, tal2, localhomlist+junban[count-1], qstart, tstart, opt, overlapaa );
1576         }
1577 //      fprintf( stderr, "count = %d\n", count );
1578     return count;
1579 }
1580 int ReadFasta34m10( FILE *fp, double *dis, int qmem, char name[M][B], LocalHom *localhomlist )
1581 {
1582     int count=0;
1583     char b[B];
1584         char *pt;
1585     static int junban[M];
1586         int overlapaa;
1587         int opt, qstart, qend, tstart, tend;
1588         double z, bits;
1589         int qal_display_start, tal_display_start;
1590         static char qal[N], tal[N];
1591         char *qal2, *tal2;
1592         int c;
1593
1594
1595     count = 0;
1596 #if 0
1597     for( i=0; i<10000000 && count<nseq; i++ )
1598 #else
1599     while( !feof( fp ) )
1600 #endif
1601     {
1602         fgets( b, B-1, fp );
1603         if( !strncmp( "+==========+", b, 12 ) )
1604         {
1605             junban[count] = atoi( b+12 );
1606
1607                         pt = strchr( b, ')' ) + 1;
1608                         sscanf( pt, "%d %lf %lf",  &opt, &bits, &z );
1609             dis[junban[count]] = (double)opt;
1610             count++;
1611
1612         }
1613                 else if( 0 == strncmp( ">>+==========+", b, 14 ) )
1614                 {
1615                         break;
1616                 }
1617
1618     }
1619         if( !count ) return -1;
1620
1621         count = 0;
1622     while( 1 )
1623         {
1624                 if( strncmp( ">>+==========+", b, 14 ) )
1625                 {
1626                         fgets( b, B-1, fp );
1627                         if( feof( fp ) ) break;
1628                         continue;
1629                 }
1630                 junban[count++] = atoi( b+14 );
1631 //              fprintf( stderr, "t = %d\n", atoi( b+14 ) );
1632                 while( fgets( b, B-1, fp ) )
1633                         if( !strncmp( "; fa_opt:", b, 9 ) || !strncmp( "; sw_s-w opt:", b, 13 ) ) break;
1634                 pt = strstr( b, ":" ) +1;
1635                 opt = atoi( pt );
1636
1637
1638                 while( fgets( b, B-1, fp ) )
1639                         if( !strncmp( "_overlap:", b+4, 9 ) ) break;
1640                 pt = strstr( b, ":" ) +1;
1641                 overlapaa = atoi( pt );
1642
1643                 while( fgets( b, B-1, fp ) )
1644                         if( !strncmp( "_start:", b+4, 7 ) ) break;
1645                 pt = strstr( b, ":" ) +1;
1646                 qstart = atoi( pt ) - 1;
1647
1648                 while( fgets( b, B-1, fp ) )
1649                         if( !strncmp( "_stop:", b+4, 6 ) ) break;
1650                 pt = strstr( b, ":" ) +1;
1651                 qend = atoi( pt ) - 1;
1652
1653                 while( fgets( b, B-1, fp ) )
1654                         if( !strncmp( "_display_start:", b+4, 15 ) ) break;
1655                 pt = strstr( b, ":" ) +1;
1656                 qal_display_start = atoi( pt ) - 1;
1657
1658                 pt = qal;
1659                 while( c = fgetc( fp ) )
1660                 {
1661                         if( c == '>' )
1662                         {
1663                                 ungetc( c, fp );
1664                                 break;
1665                         }
1666                         if( isalpha( c ) || c == '-' )
1667                         *pt++ = c;
1668                 }
1669                 *pt = 0;
1670
1671                 while( fgets( b, B-1, fp ) )
1672                         if( !strncmp( "_start:", b+4, 7 ) ) break;
1673                 pt = strstr( b, ":" ) + 1;
1674                 tstart = atoi( pt ) - 1;
1675
1676                 while( fgets( b, B-1, fp ) )
1677                         if( !strncmp( "_stop:", b+4, 6 ) ) break;
1678                 pt = strstr( b, ":" ) + 1;
1679                 tend = atoi( pt ) - 1;
1680
1681                 while( fgets( b, B-1, fp ) )
1682                         if( !strncmp( "_display_start:", b+4, 15 ) ) break;
1683                 pt = strstr( b, ":" ) + 1;
1684                 tal_display_start = atoi( pt ) - 1;
1685
1686                 pt = tal;
1687                 while( c = fgetc( fp ) )
1688                 {
1689                         if( c == '>' )
1690                         {
1691                                 ungetc( c, fp );
1692                                 break;
1693                         }
1694                         if( isalpha( c ) || c == '-' )
1695                         *pt++ = c;
1696                 }
1697                 *pt = 0;
1698
1699 //              fprintf( stderr, "(%d-%d:%d-%d)\n", qstart, qend, tstart, tend );
1700 //              fprintf( stderr, "qal_display_start = %d, tal_display_start = %d\n", qal_display_start, tal_display_start );
1701
1702 //              fprintf( stderr, "qal = %s\n", qal );
1703 //              fprintf( stderr, "tal = %s\n", tal );
1704
1705                 qal2 = cutal( qal, qal_display_start, qstart, qend );
1706                 tal2 = cutal( tal, tal_display_start, tstart, tend );
1707
1708 //              fprintf( stderr, "qal2 = %s\n", qal2 );
1709 //              fprintf( stderr, "tal2 = %s\n", tal2 );
1710
1711 //              fprintf( stderr, "putting   %d - %d, opt = %d\n", qmem, junban[count-1], opt );
1712                 putlocalhom( qal2, tal2, localhomlist+junban[count-1], qstart, tstart, opt, overlapaa );
1713         }
1714 //      fprintf( stderr, "count = %d\n", count );
1715     return count;
1716 }
1717 int ReadFasta34( FILE *fp, double *dis, int nseq, char name[M][B], LocalHom *localhomlist )
1718 {
1719     int count=0;
1720     char b[B];
1721         char *pt;
1722     static int junban[M];
1723         int overlapaa;
1724         int opt, qstart, qend, tstart, tend;
1725         double z, bits;
1726
1727
1728     count = 0;
1729 #if 0
1730     for( i=0; i<10000000 && count<nseq; i++ )
1731 #else
1732     while( !feof( fp ) )
1733 #endif
1734     {
1735         fgets( b, B-1, fp );
1736         if( !strncmp( "+==========+", b, 12 ) )
1737         {
1738             junban[count] = atoi( b+12 );
1739
1740                         pt = strchr( b, ')' ) + 1;
1741                         sscanf( pt, "%d %lf %lf",  &opt, &bits, &z );
1742             dis[junban[count]] = (double)opt;
1743             count++;
1744
1745         }
1746                 else if( 0 == strncmp( ">>+==========+", b, 14 ) )
1747                 {
1748                         break;
1749                 }
1750
1751     }
1752         if( !count ) return -1;
1753
1754         count = 0;
1755     while( !feof( fp ) )
1756         {
1757                 if( !strncmp(">>+==========+", b, 14 ) )
1758                 {
1759             junban[count] = atoi( b+14 );
1760             count++;
1761                 fgets( b, B-1, fp ); // initn:
1762                         pt = strstr( b, "opt: " ) + 5;
1763                         localhomlist[junban[count-1]].opt = atof( pt );
1764                 fgets( b, B-1, fp ); // Smith-Waterman score
1765                         pt = strstr( b, "ungapped) in " ) + 13;
1766                         sscanf( pt, "%d", &overlapaa );
1767                         fprintf( stderr, "pt = %s, overlapaa = %d\n", pt, overlapaa );
1768                         pt = strstr( b, "overlap (" ) + 8;
1769                         sscanf( pt, "(%d-%d:%d-%d)", &qstart, &qend, &tstart, &tend );
1770                         localhomlist[junban[count-1]].overlapaa = overlapaa;
1771                         localhomlist[junban[count-1]].start1 = qstart-1;
1772                         localhomlist[junban[count-1]].end1   = qend-1;
1773                         localhomlist[junban[count-1]].start2 = tstart-1;
1774                         localhomlist[junban[count-1]].end2   = tend-1;
1775                 }
1776         fgets( b, B-1, fp );
1777         }
1778         fprintf( stderr, "count = %d\n", count );
1779     return count;
1780 }
1781
1782 int ReadFasta3( FILE *fp, double *dis, int nseq, char name[M][B] )
1783 {
1784     int count=0;
1785     char b[B];
1786         char *pt;
1787     int junban[M];
1788         int initn, init1, opt;
1789         double z;
1790
1791     count = 0;
1792 #if 0
1793     for( i=0; i<10000000 && count<nseq; i++ )
1794 #else
1795     while( !feof( fp ) )
1796 #endif
1797     {
1798         fgets( b, B-1, fp );
1799         if( !strncmp( "+==========+", b, 12 ) )
1800         {
1801             junban[count] = atoi( b+12 );
1802
1803                         pt = strchr( b, ')' ) + 1;
1804                         sscanf( pt, "%d %d %d %lf", &initn, &init1, &opt, &z );
1805             dis[junban[count]] = (double)opt;
1806             count++;
1807         }
1808     }
1809     return 0;
1810 }
1811
1812 int ReadFasta( FILE *fp, double *dis, int nseq, char name[M][B] )
1813 {
1814     int i, count=0;
1815     char b[B];
1816     int junban[M];
1817         int initn, init1, opt;
1818
1819     count = 0;
1820         for( i=0; i<nseq; i++ ) dis[i] = 0.0;
1821     for( i=0; !feof( fp ) && count<nseq; i++ )
1822     {
1823         fgets( b, B-1, fp );
1824         if( !strncmp( "+==========+", b, 12 ) )
1825         {
1826             junban[count] = atoi( b+12 );
1827
1828                         sscanf( b+50, "%d %d %d", &initn, &init1, &opt );
1829             dis[junban[count]] = (double)opt;
1830             count++;
1831         }
1832     }
1833
1834 /*
1835     count = 0;
1836     for( i=0; i<100000 && count<nseq; i++ )
1837     {
1838         fgets( b, B-1, fp );
1839         if( !strncmp( name[junban[count]], b, 20  ) )
1840         {
1841             dis[junban[count]] = atof( b+65 );
1842             count++;
1843         }
1844     }
1845 */
1846     return 0;
1847 }
1848
1849
1850 int ReadOpt( FILE *fp, int opt[M], int nseq, char name[M][B] )
1851 {
1852     int i, count=0;
1853     char b[B];
1854     int junban[M];
1855         int optt, initn, init1;
1856
1857     count = 0;
1858     for( i=0; i<10000000 && count<nseq; i++ )
1859     {
1860         fgets( b, B-1, fp );
1861         if( !strncmp( "+==========+", b, 12 ) )
1862         {
1863             junban[count] = atoi( b+12 );
1864                         sscanf( b+50, "%d %d %d", &initn, &init1, &optt );
1865             opt[junban[count]] = (double)optt;
1866             count++;
1867         }
1868     }
1869     return 0;
1870 }
1871
1872 int ReadOpt2( FILE *fp, int opt[M], int nseq, char name[M][B] )
1873 {
1874     int i, count=0;
1875     char b[B];
1876     int junban[M];
1877
1878     count = 0;
1879     for( i=0; i<10000000 && count<nseq; i++ )
1880     {
1881         fgets( b, B-1, fp );
1882         if( !strncmp( "+==========+", b, 12 ) )
1883         {
1884             junban[count] = atoi( b+12 );
1885             opt[junban[count]] = atoi( b+65 );
1886             count++;
1887         }
1888     }
1889     return 0;
1890 }
1891
1892
1893
1894 int writePre( int nseq, char name[][B], int nlen[M], char **aseq, int force )
1895 {
1896 #if USE_XCED
1897         int i, value;
1898         if( !signalSM )
1899         {
1900                 if( force )
1901                 {
1902                         rewind( prep_g );
1903                         if( devide ) dvWrite( prep_g, nseq, name, nlen, aseq );
1904 #if 0
1905                         else    WriteGapFill( prep_g, nseq, name, nlen, aseq );
1906 #else
1907                         else    writeData( prep_g, nseq, name, nlen, aseq );
1908 #endif
1909                 }
1910                 return( 0 );
1911         }
1912         for( i=0; i<10; i++ )
1913         {
1914 #if IODEBUG
1915                 fprintf( stderr, "SEMAPHORE = %d\n", signalSM[SEMAPHORE] );
1916 #endif
1917                 if( signalSM[SEMAPHORE]-- > 0 )
1918                 {
1919 #if 0 /* /tmp/pre ¤Î´Ø·¸¤Ç¤Ï¤º¤·¤¿ */
1920                         if( ferror( prep_g ) ) prep_g = fopen( "pre", "w" );
1921                         if( !prep_g ) ErrorExit( "Cannot re-open pre." );
1922 #endif
1923                         rewind( prep_g );
1924                         signalSM[STATUS] = IMA_KAITERU;
1925 #if IODEBUG
1926                         if( force ) fprintf( stderr, "FINAL " );
1927 #endif
1928                         if( devide ) dvWrite( prep_g, nseq, name, nlen, aseq );
1929                         else    WriteGapFill( prep_g, nseq, name, nlen, aseq );
1930                         /*
1931                         fprintf( prep_g, '\EOF' );
1932                         */
1933                         fflush( prep_g );
1934                         if( force ) signalSM[STATUS] = OSHIMAI;
1935                         else        signalSM[STATUS] = KAKIOWATTA;
1936                         value = 1;
1937                         signalSM[SEMAPHORE]++;
1938 #if IODEBUG
1939                         fprintf( stderr, "signalSM[STATUS] = %c\n", signalSM[STATUS] );
1940 #endif
1941                         break;
1942                 }
1943                 else
1944                 {
1945 #if IODEBUG
1946                         fprintf( stderr, "YONDERUKARA_AKIRAMERU\n" );
1947 #endif
1948                         value = 0;
1949                         signalSM[SEMAPHORE]++;
1950                         if( !force ) break;
1951 #if IODEBUG
1952                         fprintf( stderr, "MATSU\n" );
1953 #endif
1954                         sleep( 1 );
1955                 }
1956         }
1957         if( force && !value ) ErrorExit( "xced ga pre wo hanasanai \n" );
1958         return( value );
1959 #else
1960         if( force )
1961         {
1962                 rewind( prep_g );
1963                         writeData( prep_g, nseq, name, nlen, aseq );
1964         }
1965 #endif
1966         return( 0 );
1967 }
1968
1969
1970 void readOtherOptions( int *ppidptr, int *fftThresholdptr, int *fftWinSizeptr )
1971 {
1972         if( calledByXced )
1973         {
1974                 FILE *fp = fopen( "pre", "r" );
1975                 char b[B];
1976                 if( !fp ) ErrorExit( "Cannot open pre.\n" );
1977                 fgets( b, B-1, fp );
1978                 sscanf( b, "%d %d %d", ppidptr, fftThresholdptr, fftWinSizeptr );
1979                 fclose( fp );
1980 #if IODEBUG
1981         fprintf( stderr, "b = %s\n", b );
1982         fprintf( stderr, "ppid = %d\n", ppid );
1983         fprintf( stderr, "fftThreshold = %d\n", fftThreshold );
1984         fprintf( stderr, "fftWinSize = %d\n", fftWinSize );
1985 #endif
1986         }
1987         else
1988         {
1989                 *ppidptr = 0;
1990                 *fftThresholdptr = FFT_THRESHOLD;
1991                 if( dorp == 'd' )
1992                         *fftWinSizeptr = FFT_WINSIZE_D;
1993                 else
1994                         *fftWinSizeptr = FFT_WINSIZE_P;
1995         }
1996 #if 0
1997         fprintf( stderr, "fftThresholdptr=%d\n", *fftThresholdptr );
1998         fprintf( stderr, "fftWinSizeptr=%d\n", *fftWinSizeptr );
1999 #endif
2000 }
2001
2002 void initSignalSM( void )
2003 {
2004 //      int signalsmid;
2005
2006 #if IODEBUG
2007         if( ppid ) fprintf( stderr, "PID of xced = %d\n", ppid );
2008 #endif
2009         if( !ppid )
2010         {
2011                 signalSM = NULL;
2012                 return;
2013         }
2014
2015 #if 0
2016         signalsmid = shmget( (key_t)ppid, 3, IPC_ALLOC | 0666 );
2017         if( signalsmid == -1 ) ErrorExit( "Cannot get Shared memory for signal.\n" );
2018         signalSM = shmat( signalsmid, 0, 0 );
2019         if( (int)signalSM == -1 ) ErrorExit( "Cannot attatch Shared Memory for signal!\n" );
2020         signalSM[STATUS] = IMA_KAITERU;
2021         signalSM[SEMAPHORE] = 1;
2022 #endif
2023 }
2024
2025 void initFiles( void )
2026 {
2027         char pname[100];
2028         if( ppid )
2029                 sprintf( pname, "/tmp/pre.%d", ppid );
2030         else
2031                 sprintf( pname, "pre" );
2032         prep_g = fopen( pname, "w" );
2033         if( !prep_g ) ErrorExit( "Cannot open pre" );
2034
2035         trap_g = fopen( "trace", "w" );
2036         if( !trap_g ) ErrorExit( "cannot open trace" );
2037         fprintf( trap_g, "PID = %d\n", getpid() );
2038         fflush( trap_g );
2039 }
2040
2041
2042 void WriteForFasta( FILE *fp, int locnjob, char name[][B], int nlen[M], char **aseq )
2043 {
2044     static char b[N];
2045     int i, j;
2046     int nalen[M];
2047
2048     for( i=0; i<locnjob; i++ )
2049     {
2050         nalen[i] = strlen( aseq[i] );
2051         fprintf( fp, ">%s\n", name[i] );
2052         for( j=0; j<nalen[i]; j=j+C )
2053         {
2054             strncpy( b, aseq[i]+j, C ); b[C] = 0;
2055             fprintf( fp, "%s\n",b );
2056         }
2057     }
2058 }
2059
2060 void readlocalhomtable( FILE*fp, int njob, LocalHom **localhomtable )
2061 {
2062         double opt;
2063         static char buff[B];
2064         int i, j, overlapaa, start1, end1, start2, end2;
2065         int **nlocalhom = NULL;
2066         LocalHom *tmpptr1, *tmpptr2;
2067
2068         nlocalhom = AllocateIntMtx( njob, njob );
2069         for( i=0; i<njob; i++ ) for( j=0; j<njob; j++ ) nlocalhom[i][j] = 0;
2070
2071         while ( NULL != fgets( buff, B-1, fp ) )
2072         {
2073 //              fprintf( stderr, "\n" );
2074                 sscanf( buff, "%d %d %d %lf %d %d %d %d",  &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2 );
2075
2076 #if 0
2077                 if( start1 == end1 || start2 == end2 ) continue; //mondai ari
2078 #endif
2079
2080
2081                 if( nlocalhom[i][j]++ > 0 )
2082                 {
2083 //                      fprintf( stderr, "reallocating, nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
2084                         tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2085                         tmpptr1 = tmpptr1->next;
2086                         tmpptr1->next = NULL;
2087                 }
2088                 else
2089                 {
2090                         tmpptr1 = localhomtable[i]+j;
2091 //                      fprintf( stderr, "nlocalhom[%d][%d] = %d\n", i, j, nlocalhom[i][j] );
2092                 }
2093
2094                 tmpptr1->start1 = start1;
2095                 tmpptr1->start2 = start2;
2096                 tmpptr1->end1 = end1;
2097                 tmpptr1->end2 = end2;
2098 //              tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8  * 600;
2099 //              tmpptr1->opt = opt;
2100                 tmpptr1->opt = ( opt + 0.00 ) / 5.8  * 600;
2101                 tmpptr1->overlapaa = overlapaa;
2102
2103 //              fprintf( stderr, "i=%d, j=%d, opt = %f\n", i, j, opt );
2104
2105                 if( nlocalhom[j][i]++ > 0 )
2106                 {
2107                         tmpptr2->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2108                         tmpptr2 = tmpptr2->next;
2109                         tmpptr2->next = NULL;
2110                 }
2111                 else
2112                         tmpptr2 = localhomtable[j]+i;
2113
2114                 tmpptr2->start2 = start1;
2115                 tmpptr2->start1 = start2;
2116                 tmpptr2->end2 = end1;
2117                 tmpptr2->end1 = end2;
2118 //              tmpptr2->opt = ( opt / overlapaa + 0.00 ) / 5.8  * 600;
2119 //              tmpptr2->opt = opt;
2120                 tmpptr2->opt = ( opt + 0.00 ) / 5.8  * 600;
2121                 tmpptr2->overlapaa = overlapaa;
2122
2123         }
2124         FreeIntMtx( nlocalhom );
2125 }
2126
2127 void outlocalhom( LocalHom **localhom, int nseq )
2128 {
2129         int i, j;
2130         LocalHom *tmpptr;
2131         for( i=0; i<nseq; i++ ) for( j=0; j<nseq; j++ )
2132         {
2133                 tmpptr = localhom[i]+j;
2134                 fprintf( stderr, "%d-%d\n", i, j );
2135                 do
2136                 {
2137                         fprintf( stderr, "reg1=%d-%d, reg2=%d-%d, imp=%f, opt=%f\n", tmpptr->start1, tmpptr->end1, tmpptr->start2, tmpptr->end2, tmpptr->importance, tmpptr->opt );
2138                 }
2139                 while( tmpptr=tmpptr->next );
2140         }
2141 }
2142
2143 void outlocalhompt( LocalHom ***localhom, int n1, int n2 )
2144 {
2145         int i, j;
2146         LocalHom *tmpptr;
2147         for( i=0; i<n1; i++ ) for( j=0; j<n2; j++ )
2148         {
2149                 tmpptr = localhom[i][j];
2150                 fprintf( stderr, "%d-%d\n", i, j );
2151                 do
2152                 {
2153                         fprintf( stderr, "reg1=%d-%d, reg2=%d-%d, imp=%f, opt=%f, wimp=%f\n", tmpptr->start1, tmpptr->end1, tmpptr->start2, tmpptr->end2, tmpptr->importance, tmpptr->opt, tmpptr->wimportance );
2154                 }
2155                 while( tmpptr=tmpptr->next );
2156         }
2157 }
2158
2159 void FreeLocalHomTable( LocalHom **localhomtable, int n )
2160 {
2161         int i, j;
2162         LocalHom *ppp, *tmpptr;
2163         fprintf( stderr, "freeing localhom\n" );
2164         for( i=0; i<n; i++ )
2165         {
2166                 for( j=0; j<n; j++ )
2167                 {
2168                         tmpptr=localhomtable[i]+j;
2169                         ppp = tmpptr->next;
2170                         for( ; tmpptr; tmpptr=ppp )
2171                         {
2172 #if DEBUG
2173                                 fprintf( stderr, "i=%d, j=%d\n", i, j );
2174 #endif
2175                                 ppp = tmpptr->next;
2176                                 if( tmpptr!=localhomtable[i]+j )
2177                                 {
2178 #if DEBUG
2179                                         fprintf( stderr, "freeing %p\n", tmpptr );
2180 #endif
2181                                         free( tmpptr );
2182                                 }
2183                         }
2184                 }
2185 #if DEBUG
2186                 fprintf( stderr, "freeing localhomtable[%d]\n", i );
2187 #endif
2188                 free( localhomtable[i] );
2189         }
2190 #if DEBUG
2191         fprintf( stderr, "freeing localhomtable\n" );
2192 #endif
2193         free( localhomtable );
2194 #if DEBUG
2195         fprintf( stderr, "freed\n" );
2196 #endif
2197 }
2198
2199 char *progName( char *str )
2200 {
2201     char *value;
2202     if( ( value = strrchr( str, '/' ) ) != NULL )
2203         return( value+1 );
2204     else
2205         return( str );
2206 }
2207
2208 static char *extractfirstword( char *str )
2209 {
2210         char *val = str;
2211
2212         while( *str )
2213         {
2214                 if( val == str && *str == ' ' )
2215                 {
2216                         val++; str++;
2217                 }
2218                 else if( *str != ' ' )
2219                 {
2220                         str++;
2221                 }
2222                 else if( *str == ' ' )
2223                 {
2224                         *str = 0;
2225                 }
2226         }
2227         return( val );
2228 }
2229
2230
2231 void clustalout( FILE *fp, int nseq, int maxlen, char **seq, char name[][B], char *mark, char *comment, int *order )
2232 {
2233         int pos, j;
2234         pos = 0;
2235         if( comment == NULL )
2236                 fprintf( fp, "CLUSTAL (-like) formatted alignment by MAFFT (v%s)\n\n", VERSION );
2237         else
2238                 fprintf( fp, "CLUSTAL (-like) formatted alignment by MAFFT %s (v%s)\n\n", comment, VERSION );
2239
2240         while( pos < maxlen )
2241         {
2242                 fprintf( fp, "\n" );
2243                 for( j=0; j<nseq; j++ )
2244                 {
2245                         fprintf( fp, "%-15.15s ", extractfirstword( name[order[j]]+1 ) );
2246                         fprintf( fp, "%.60s\n", seq[order[j]]+pos ); // Ä¹¤µ¤¬°ã¤¦¤È¤À¤á¡£
2247                 }
2248                 fprintf( fp, "%-15.15s ", "" );
2249                 fprintf( fp, "%.60s\n", mark + pos ); // Ä¹¤µ¤¬°ã¤¦¤È¤À¤á¡£
2250                 pos += 60;
2251         }
2252 }
2253
2254 void writeData_reorder( FILE *fp, int locnjob, char name[][B], int nlen[], char **aseq, int *order )
2255 {
2256         int i, j, k;
2257         int nalen;
2258
2259         for( i=0; i<locnjob; i++ )
2260         {
2261                 k = order[i];
2262 #if DEBUG
2263                 fprintf( stderr, "i = %d in writeData\n", i );
2264 #endif
2265                 nalen = strlen( aseq[k] );
2266                 fprintf( fp, ">%s\n", name[k]+1 );
2267                 for( j=0; j<nalen; j=j+C )
2268                 {
2269 #if 0
2270                         strncpy( b, aseq[k]+j, C ); b[C] = 0;
2271                         fprintf( fp, "%s\n",b );
2272 #else
2273                         fprintf( fp, "%.*s\n", C, aseq[k]+j );
2274 #endif
2275                 }
2276         }
2277 }
2278
2279 void putlocalhom3( char *al1, char *al2, LocalHom *localhompt, int off1, int off2, int opt, int overlapaa )
2280 {
2281         int pos1, pos2, start1, start2, end1, end2;
2282         char *pt1, *pt2;
2283         double score;
2284         double sumscore;
2285         int sumoverlap;
2286         LocalHom *tmppt;
2287         LocalHom *subnosento;
2288         int st;
2289         int saisho;
2290
2291         pt1 = al1; pt2 = al2;
2292         pos1 = off1; pos2 = off2;
2293
2294         sumscore = 0.0;
2295         sumoverlap = 0;
2296
2297         subnosento = localhompt;
2298         while( subnosento->next ) subnosento = subnosento->next;
2299         tmppt = subnosento;
2300
2301         saisho = ( localhompt->nokori == 0 );
2302
2303         fprintf( stderr, "localhompt = %p\n", localhompt );
2304         fprintf( stderr, "tmppt = %p\n", tmppt );
2305         fprintf( stderr, "subnosento = %p\n", subnosento );
2306
2307         st = 0;
2308         score = 0.0;
2309         while( *pt1 != 0 )
2310         {
2311 //              fprintf( stderr, "pt = %c, %c, st=%d\n", *pt1, *pt2, st );
2312                 if( st == 1 && ( *pt1 == '-' || *pt2 == '-' ) )
2313                 {
2314                         end1 = pos1 - 1;
2315                         end2 = pos2 - 1;
2316
2317                         if( localhompt->nokori++ > 0 )
2318                         {
2319 //                              fprintf( stderr, "reallocating ...\n" );
2320                                 tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2321 //                              fprintf( stderr, "done\n" );
2322                                 tmppt = tmppt->next;
2323                                 tmppt->next = NULL;
2324                         }
2325                         tmppt->start1 = start1;
2326                         tmppt->start2 = start2;
2327                         tmppt->end1   = end1  ;
2328                         tmppt->end2   = end2  ;
2329
2330 #if 1
2331                         if( divpairscore )
2332                         {
2333                                 tmppt->overlapaa   = end2-start2+1;
2334                                 tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
2335                         }
2336                         else
2337                         {
2338                                 sumscore += score;
2339                                 sumoverlap += end2-start2+1;
2340                         }
2341 #else
2342                         tmppt->overlapaa   = overlapaa;
2343                         tmppt->opt = (double)opt;
2344 #endif
2345
2346 #if 0
2347                         fprintf( stderr, "score (1)= %f\n", score );
2348                         fprintf( stderr, "al1: %d - %d\n", start1, end1 );
2349                         fprintf( stderr, "al2: %d - %d\n", start2, end2 );
2350 #endif
2351                         score = 0.0;
2352                         st = 0;
2353                 }
2354                 else if( *pt1 != '-' && *pt2 != '-' )
2355                 {
2356                         if( st == 0 )
2357                         {
2358                                 start1 = pos1; start2 = pos2;
2359                                 st = 1;
2360                         }
2361                         score += n_dis[(int)amino_n[(int)*pt1]][(int)amino_n[(int)*pt2]]; // - offset ¤Ï¤¤¤é¤Ê¤¤¤«¤â
2362 //                      fprintf( stderr, "%c-%c, score(0) = %f\n", *pt1, *pt2, score );
2363                 }
2364                 if( *pt1++ != '-' ) pos1++;
2365                 if( *pt2++ != '-' ) pos2++;
2366         }
2367         if( *(pt1-1) != '-' && *(pt2-1) != '-'  )
2368         {
2369                 if( localhompt->nokori++ > 0 )
2370                 {
2371 //                      fprintf( stderr, "reallocating ...\n" );
2372                         tmppt->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2373 //                      fprintf( stderr, "done\n" );
2374                         tmppt = tmppt->next;
2375                         tmppt->next = NULL;
2376                 }
2377
2378                 end1 = pos1 - 1;
2379                 end2 = pos2 - 1;
2380                 tmppt->start1 = start1;
2381                 tmppt->start2 = start2;
2382                 tmppt->end1   = end1  ;
2383                 tmppt->end2   = end2  ;
2384
2385
2386 #if 1
2387                 if( divpairscore )
2388                 {
2389                         tmppt->overlapaa   = end2-start2+1;
2390                         tmppt->opt = score / tmppt->overlapaa * 5.8 / 600;
2391                 }
2392                 else
2393                 {
2394                         sumscore += score;
2395                         sumoverlap += end2-start2+1;
2396                 }
2397 #else
2398                 tmppt->overlapaa   = overlapaa;
2399                 tmppt->opt = (double)opt;
2400 #endif
2401
2402 #if 0
2403                 fprintf( stderr, "score (2)= %f\n", score );
2404                 fprintf( stderr, "al1: %d - %d\n", start1, end1 );
2405                 fprintf( stderr, "al2: %d - %d\n", start2, end2 );
2406 #endif
2407         }
2408
2409         fprintf( stderr, "sumscore = %f\n", sumscore );
2410         if( !divpairscore )
2411         {
2412
2413                 if( !saisho ) subnosento = subnosento->next;
2414                 for( tmppt=subnosento; tmppt; tmppt=tmppt->next )
2415                 {
2416                         tmppt->overlapaa = sumoverlap;
2417                         tmppt->opt = sumscore * 5.8 / 600 / sumoverlap;
2418                         fprintf( stderr, "tmpptr->opt = %f\n", tmppt->opt );
2419                 }
2420         }
2421 }
2422 void readlocalhomtable2( FILE*fp, int njob, LocalHom **localhomtable )
2423 {
2424         double opt;
2425         static char buff[B];
2426         int i, j, overlapaa, start1, end1, start2, end2;
2427         LocalHom *tmpptr1, *tmpptr2;
2428
2429 //      for( i=0; i<njob; i++ ) for( j=0; j<njob; j++ ) nlocalhom[i][j] = 0;
2430
2431         while ( NULL != fgets( buff, B-1, fp ) )
2432         {
2433 //              fprintf( stderr, "\n" );
2434                 sscanf( buff, "%d %d %d %lf %d %d %d %d",  &i, &j, &overlapaa, &opt, &start1, &end1, &start2, &end2 );
2435
2436 #if 0
2437                 if( start1 == end1 || start2 == end2 ) continue; //mondai ari
2438 #endif
2439
2440                 if( localhomtable[i][j].nokori++ > 0 )
2441                 {
2442                         tmpptr1 = localhomtable[i][j].last;
2443 //                      fprintf( stderr, "reallocating, localhomtable[%d][%d].nokori = %d\n", i, j, localhomtable[i][j].nokori );
2444                         tmpptr1->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2445                         tmpptr1 = tmpptr1->next;
2446                         tmpptr1->extended = -1;
2447                         tmpptr1->next = NULL;
2448                         localhomtable[i][j].last = tmpptr1;
2449 //                      fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", i, j, localhomtable[i][j].nokori );
2450                 }
2451                 else
2452                 {
2453                         tmpptr1 = localhomtable[i]+j;
2454 //                      fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", i, j, localhomtable[i][j].nokori );
2455                 }
2456
2457                 tmpptr1->start1 = start1;
2458                 tmpptr1->start2 = start2;
2459                 tmpptr1->end1 = end1;
2460                 tmpptr1->end2 = end2;
2461 //              tmpptr1->opt = ( opt / overlapaa + 0.00 ) / 5.8  * 600;
2462 //              tmpptr1->opt = opt;
2463                 tmpptr1->opt = ( opt + 0.00 ) / 5.8  * 600;
2464                 tmpptr1->overlapaa = overlapaa;
2465
2466 //              fprintf( stderr, "i=%d, j=%d, st1=%d, en1=%d, opt = %f\n", i, j, tmpptr1->start1, tmpptr1->end1, opt );
2467
2468                 if( localhomtable[j][i].nokori++ > 0 )
2469                 {
2470                         tmpptr2 = localhomtable[j][i].last;
2471                         tmpptr2->next = (LocalHom *)calloc( 1, sizeof( LocalHom ) );
2472                         tmpptr2 = tmpptr2->next;
2473                         tmpptr2->extended = -1;
2474                         tmpptr2->next = NULL;
2475                         localhomtable[j][i].last = tmpptr2;
2476 //                      fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", j, i, localhomtable[j][i].nokori );
2477                 }
2478                 else
2479                 {
2480                         tmpptr2 = localhomtable[j]+i;
2481 //                      fprintf( stderr, "### i,j = %d,%d, nokori=%d\n", j, i, localhomtable[j][i].nokori );
2482                 }
2483
2484                 tmpptr2->start2 = start1;
2485                 tmpptr2->start1 = start2;
2486                 tmpptr2->end2 = end1;
2487                 tmpptr2->end1 = end2;
2488 //              tmpptr2->opt = ( opt / overlapaa + 0.00 ) / 5.8  * 600;
2489 //              tmpptr2->opt = opt;
2490                 tmpptr2->opt = ( opt + 0.00 ) / 5.8  * 600;
2491                 tmpptr2->overlapaa = overlapaa;
2492
2493         }
2494 }