sources/multicoil/scpos.c

   1 #include <stdio.h>
   2 #include <math.h>   /** Defines HUGE_VAL. **/
   3 #include "sc.h"
   4 #include "scconst.h"
   5 #include "scio.h"
   6 #include "scsystem.h"
   7
   8
   9
  10
  11 void add_seq2(long pfreqs[AANUM][POSNUM],
  12         long ptotals[POSNUM],
  13         long pfreqp[AANUM][AANUM][POSNUM][POSNUM],
  14         long ptotalp[POSNUM][POSNUM],
  15         char reg[MAXSEQLEN],
  16         char seq[MAXSEQLEN], int seqlen, int bias)
  17 /* Temporarily add or remove the sequence under consideration from the
  18    positive data? (bias=-1 removes, bias=1 adds) */
  19 {
  20   int i, posi, j, posj;
  21
  22   for (i = 0; i < seqlen; i++)
  23     if (seq[i] < AANUM && reg[i]!='.') {
  24         posi = reg[i];
  25         for (j = i+1; j <= i+POSNUM && j < seqlen && (reg[j] == (posj = (posi+j-i)%POSNUM)); j++)
  26           if (seq[j] < AANUM) {
  27             pfreqp[seq[i]][seq[j]][posi][posj] += bias;
  28             ptotalp[posi][posj] += bias;
  29           }
  30         pfreqs[seq[i]][posi] += bias;
  31         ptotals[posi] += bias;
  32       }
  33 }
  34
  35
  36
  37
  38 remove_all_input_seq_from_table(long pfreqs[AANUM][POSNUM],
  39            long ptotals[POSNUM],long pfreqp[AANUM][AANUM][POSNUM][POSNUM],
  40            long ptotalp[POSNUM][POSNUM], char *testfile_name, int bias)
  41 {
  42   Sequence sequence;
  43   char title[MAXLINE],code[MAXLINE];
  44   char seq[MAXSEQLEN], reg[MAXSEQLEN];
  45   int i, posi, j, posj;
  46   FILE *ftest;
  47
  48
  49   sequence.title = title;
  50   sequence.code = code;
  51   sequence.seq = seq;
  52   sequence.reg= reg;
  53
  54   ftest = sopen(testfile_name,"r");
  55   while (getseq2(ftest,&sequence)) {
  56     for (i = 0; i < sequence.seqlen; i++)
  57       if (seq[i] < AANUM && reg[i]!='.') {
  58         posi = reg[i];
  59         for (j = i+1; j <= i+POSNUM && j < sequence.seqlen && (reg[j] == (posj = (posi+j-i)%POSNUM)); j++)
  60           if (seq[j] < AANUM) {
  61             pfreqp[seq[i]][seq[j]][posi][posj] += bias;
  62             ptotalp[posi][posj] += bias;
  63           }
  64         pfreqs[seq[i]][posi] += bias;
  65         ptotals[posi] += bias;
  66       }
  67   }
  68   sclose(ftest);
  69 }
  70
  71
  72
  73
  74
  75
  76
  77 int getpos (FILE *fl, Sequence *sequence)
  78 /* Read in the next sequence from a .pos file */
  79 {
  80  int c,count=0;
  81  char *code = sequence->code;
  82  char *title = sequence->title;
  83  char *seq = sequence->seq;
  84  char *reg = sequence->reg;
  85
  86
  87  c=fgetc(fl);
  88  while (c!='>' && c!=EOF)
  89   c=fgetc(fl);
  90  if (c==EOF)
  91   return(0);
  92  /* Get the code */
  93  c=fgetc(fl);
  94  while (c!='\n')
  95  {
  96   code[count++]=c;
  97   c=fgetc(fl);
  98  }
  99  code[count]=0;
 100  count=0;
 101  /* Get the title info */
 102  c=fgetc(fl);
 103  while (c!='\n') {
 104    title[count++]=c;
 105    c=fgetc(fl);
 106  }
 107  title[count]=0;
 108  count=0;
 109  /* Get the sequence */
 110  c=fgetc(fl);
 111  while (c!='*' && c!=EOF)
 112  {
 113   if (aanum(c)!=(char) -1)
 114    seq[count++]=aanum(c);
 115   c=fgetc(fl);
 116  }
 117  sequence->seqlen=count;
 118  count=0;
 119  /* Get the register info */
 120  c=fgetc(fl);
 121  while (c!='*' && c!=EOF)
 122  {
 123   if (c>='a' && c<='z')
 124    reg[count++]=c-'a';
 125   else if (c=='.')
 126    reg[count++]='.';
 127   c=fgetc(fl);
 128  }
 129  if (sequence->seqlen != count) {
 130    fprintf(stderr, "erroneous .pos file!\n");
 131    /*exit(-1);*/
 132  }
 133  if (c==EOF)
 134   return(0);
 135  else
 136   return(1);
 137 }
 138