#include <stdio.h>
#include <math.h>   /** Defines HUGE_VAL. **/
#include "sc.h"
#include "scconst.h"
#include "scio.h"
#include "scsystem.h"




void add_seq2(long pfreqs[AANUM][POSNUM], 
	long ptotals[POSNUM],
	long pfreqp[AANUM][AANUM][POSNUM][POSNUM], 
	long ptotalp[POSNUM][POSNUM],
        char reg[MAXSEQLEN],
	char seq[MAXSEQLEN], int seqlen, int bias)
/* Temporarily add or remove the sequence under consideration from the 
   positive data? (bias=-1 removes, bias=1 adds) */ 
{
  int i, posi, j, posj;

  for (i = 0; i < seqlen; i++) 
    if (seq[i] < AANUM && reg[i]!='.') {
        posi = reg[i];
	for (j = i+1; j <= i+POSNUM && j < seqlen && (reg[j] == (posj = (posi+j-i)%POSNUM)); j++) 
	  if (seq[j] < AANUM) {
	    pfreqp[seq[i]][seq[j]][posi][posj] += bias;
	    ptotalp[posi][posj] += bias;
	  }
	pfreqs[seq[i]][posi] += bias;
	ptotals[posi] += bias;
      }
}




remove_all_input_seq_from_table(long pfreqs[AANUM][POSNUM], 
	   long ptotals[POSNUM],long pfreqp[AANUM][AANUM][POSNUM][POSNUM], 
	   long ptotalp[POSNUM][POSNUM], char *testfile_name, int bias)
{
  Sequence sequence;
  char title[MAXLINE],code[MAXLINE];
  char seq[MAXSEQLEN], reg[MAXSEQLEN];
  int i, posi, j, posj;
  FILE *ftest;

 
  sequence.title = title;
  sequence.code = code;
  sequence.seq = seq;
  sequence.reg= reg;

  ftest = sopen(testfile_name,"r");
  while (getseq2(ftest,&sequence)) {
    for (i = 0; i < sequence.seqlen; i++) 
      if (seq[i] < AANUM && reg[i]!='.') {
        posi = reg[i];
	for (j = i+1; j <= i+POSNUM && j < sequence.seqlen && (reg[j] == (posj = (posi+j-i)%POSNUM)); j++) 
	  if (seq[j] < AANUM) {
	    pfreqp[seq[i]][seq[j]][posi][posj] += bias;
	    ptotalp[posi][posj] += bias;
	  }
	pfreqs[seq[i]][posi] += bias;
	ptotals[posi] += bias;
      }
  }
  sclose(ftest);
}







int getpos (FILE *fl, Sequence *sequence) 
/* Read in the next sequence from a .pos file */
{
 int c,count=0;
 char *code = sequence->code;
 char *title = sequence->title;
 char *seq = sequence->seq;
 char *reg = sequence->reg; 
 

 c=fgetc(fl);
 while (c!='>' && c!=EOF)
  c=fgetc(fl);
 if (c==EOF)
  return(0);
 /* Get the code */
 c=fgetc(fl);
 while (c!='\n')
 {
  code[count++]=c;
  c=fgetc(fl);
 }
 code[count]=0;
 count=0;
 /* Get the title info */
 c=fgetc(fl);
 while (c!='\n') {
   title[count++]=c;
   c=fgetc(fl);
 }
 title[count]=0;
 count=0;
 /* Get the sequence */
 c=fgetc(fl);
 while (c!='*' && c!=EOF)
 {
  if (aanum(c)!=(char) -1)
   seq[count++]=aanum(c);
  c=fgetc(fl);
 }
 sequence->seqlen=count;
 count=0;
 /* Get the register info */
 c=fgetc(fl);
 while (c!='*' && c!=EOF)
 {
  if (c>='a' && c<='z')
   reg[count++]=c-'a';
  else if (c=='.')
   reg[count++]='.';
  c=fgetc(fl);
 }
 if (sequence->seqlen != count) {
   fprintf(stderr, "erroneous .pos file!\n");
   /*exit(-1);*/
 }
 if (c==EOF)
  return(0);
 else
  return(1);
}