/*****************************************************************************/
char seqlib_srs_ID[] = "$Id: seqlib.c,v 1.60 1995/01/31 22:01:33 etzold Exp $";
/*
**
**    $Source: /kappa1/srs4_0/src/RCS/seqlib.c,v $
**    $Revision: 1.60 $
**    $Date: 1995/01/31 22:01:33 $
**    $Author: etzold $
**
**    $Locker:  $
**    $State: Exp $
**
**
**                  SRS V3.x Copyright by Thure Etzold
**
**             Author: Gerald Schaefer
**                     EMBL, Meyerhofstrasse 1
**                     D-69012 Heidelberg, Germany
**                     Tel: 06221 387372  
**                     Email: schaefer@embl-heidelberg.de
**
**
**    Requires: list,expr,stack
**
** * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 
**
**    Description    sequence computation using generic list
**    =========== 
**
**    Global Parameters:
**        "printf"
**        "writeFOSN"
**        "shiftbeginpos"
**        "shiftendpos"
**        "isshiftbeginrelend"
**        "isshiftendrelbegin"
**        "takeuncompletefeature"
**        "takeuncompleteshiftfeature"
**
** * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 
*/


/*
#define TRACE
#define DUMMY
#define TRACE1
*/

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <stdarg.h>
#include <ctype.h>
#include "srsos.h"

#include "message.h"
#include "sm.h"
#include "lst.h"
#include "futil.h"
#include "par.h"
#include "parser.h"
#include "id.h"
#include "library.h"
#include "entry.h"
#include "seq.h"
#include "query.h"
#include "seqlib.h"
#include "stack.h"
#include "expr.h"

#define _SRS
#define _SLB
#define _CONSTANTS
#include SRSINCLUDE

#define vabs(x)  ((x)<0 ? -(x) : (x))
#define NAMELEN  40
enum slbfile_type {SLBxNORMAL=0, SLBxGCG7, SLBxGCG8,
	/* dgg additions */
	SLBxPEARSON=10, SLBxGENBANK, SLBxCODATA, SLBxEMBL,
	SLBxIG, SLBxNBRF
	};


#ifdef sun
extern INT4 printf ();
#endif


/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
**  external, global and module wide variables
*/

extern FILE        *protokoll;
ENTRYo             *global_entry; 
SEQo               *globalSeq = NULL;  /* remembers, if entry sequence was 
					  already read */
SEQo               *seq = NULL;
SEQo               *param_list = NULL;
static SLBoFEATURE feature = { 0, 0, 1, NULL, NULL, NULL, 0, NULL, NULL, NULL};
static FILE        *file=NULL;
static FILE				* fofn= NULL; /* dgg */
static SMoBUFF     *featureBuff=NULL;
INT4               global_before_f = 0;
INT4               global_after_f = 0;

/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
**  static functions
*/

static INT4 join(STKo *parameter_stack, STKo *eval_stack);
static INT4 complement(STKo *parameter_stack, STKo *eval_stack);
static INT4 replace(STKo *parameter_stack, STKo *eval_stack);
static INT4 one_of(STKo *parameter_stack, STKo *eval_stack);
static void SlbInitFeature (SLBoFEATURE *feature);

static void SlbJoinList();
static void SlbComplementList();
static void SlbReplaceList(char *replace);
static INT4 SlbInsert(char *acc, INT4 bpos, INT4 epos);
static INT4 SlbModifyAccno(char *acc);
static INT4 SlbShiftAndCheck();
static INT4 SlbEvaluate();
static void SlbWriteFosn (SEQo *seq);
static SEQo *SlbRetrieveAccNo (char *accno);


typedef struct FUNC_DESC_ {
           char func_name[20];
           INT4 (*func) ();
        } FUNC_DESC;
#define  MAXFUNC  4
FUNC_DESC func_table[MAXFUNC] = {
             { "join", join },
             { "complement", complement },
             { "replace", replace },
             { "one-of", one_of }
          };



/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
**  default function for printing lines
*/

static INT4 SlbPrintF (char *formatStr, ...)
{
    va_list ap;

    va_start (ap, formatStr);
    vfprintf (file, formatStr, ap);
    va_end (ap);

    if (formatStr[strlen (formatStr) -1] != '\n')
    fprintf (file, "\n");

    return 1;
}



	/* --- ListFSE parts --------------*/

/* dgg - from hash.c, for hashing file names... */
static  long HashString(void *s)
{
  register unsigned long h=0, g;
  register char *p = (char *) s;
  while (*p) {
    h= (h << 4) + *p++;
    if ((g= (h & 0xf0000000)))
      h= (h ^ (g >> 24)) ^ g;
  	}
  return (long) h;
}


static long   lastdocstart = -1;
static short  lastfileid = -1;
static long   nFofn = 0, maxFofn = 0;
static long * hashFofn = NULL;
char * gOutentrySuffix = ".indices";   
char * gOutnamesSuffix = ".names";

void SlbListFSE(ENTRYo *entry, char *setName, int entryN, int entryCurrN)
{
	/* dgg addition 
	- list sequence position in file 
	*/
	
	if ( EntryOpenData (entry)) {
		short   entrynuminfile, seqformat, fileid;
		long	  docstart, datastart, dataend;
		char	* filename;
		IDoENTRY  nextid;
	  ENTRYo   *nextentry;
	  
		filename= entry->file[1]->nam;
		fileid= 0; /* this must be set to filename entry in file of filenames ! */
		entrynuminfile= entry->id->fip; /* need to look at next entry for end index !*/
		docstart= entry->fip[0];
		datastart= entry->fip[1];
		/* !! NOTE: for GENBANK (? and others) the entry->fip
		currently points at the "ORIGIN" line above data. 
		In SwissProt,EMBL it points to "SQ   SEQUENCE   544 AA..." above data
		In PIR, it points to 'SEQUENCE' line above data.
		Must bypass this line here or in index users progs... 
		*/
 		datastart= ftell(entry->file[1]->fil); /* this is after 1st line read of EntryOpenData() */

		dataend= 0; /* flag we don't know end */
		/* we want real file format, in fasta/pearson lib file values */
		seqformat= entry->lib->form->fil_t[1]->type - SLBxPEARSON;
			    
    if (!entryCurrN || !file) {
			/* force output to files -- file of entry indices & file of filenames */
    	char  *tmp, *outentry;
    	char  outDirName[FILxXNAM+1], outFile[FILxXNAM+1];
	    
	    lastdocstart= 0; lastfileid= 0;
			nFofn = 0; 
			maxFofn= 20;
			if (hashFofn) free(hashFofn);
			hashFofn = (long*) malloc( maxFofn * sizeof(long));
	    
	    if ((tmp = ParGetStr ("tempDirName"))) sprintf (outDirName, "%s", tmp);    
	    else  *outDirName = '\0';
    
    	outentry= ParGetStr("listFSE");
			sprintf (outFile, "%s%s%s", outDirName, outentry, gOutentrySuffix);
#if 1
			_FilLN(outFile);
      file = fopen( outFile, "r+");
      if (!file) file = fopen( outFile, "w");
      else fseek(file, 0, 2);
#else
      file = FilOpenU(outFile, &errCode);
      _ErrMsg2 (errCode, outFile);
#endif

			sprintf (outFile, "%s%s%s", outDirName, outentry, gOutnamesSuffix);
			_FilLN(outFile);
      fofn = fopen( outFile, "r+");
      if (!fofn) 
      	fofn = fopen( outFile, "w");
      else {
      		/* read & store hashFofn values */
      	char *ep, aline[FILxXNAM+20];
      	*aline= '\0';
     		while ( fgets( aline, sizeof(aline), fofn) ) {
     			ep= aline;
     			while (*ep && *ep != '\t' && *ep != ' ' && *ep != '\n') ep++;
      	  if (*ep) *ep= '\0';
		   		if (nFofn >= maxFofn) {
		   			maxFofn = nFofn + 20;
		   			hashFofn= ( long*) realloc(hashFofn, maxFofn * sizeof( long));
		   			}
      		if (hashFofn) hashFofn[nFofn]= HashString( aline);
      		nFofn++;
      		}
      	}
      }
      
  	if (fofn) {
			/* !! save time by storing filenames in local array !? */
	    short linenum, more;
	    long  namelen;
			long filehash;
			
			more= 1;
			filehash= HashString( filename);
			if (hashFofn) for (linenum=0; linenum<nFofn; linenum++) {
				if (filehash == hashFofn[linenum]) { 
					fileid= linenum;
					more= 0;
					break;
					}
				}
				
      if (more) {
	char * libname;
      	libname= LibGetName(entry->lib,"full"); 
      	/*libname= entry->lib->lnam[0]; /* short name !? */

	fileid= nFofn;
     	fseek( fofn, 0, 2);
      	fprintf( fofn,"%s\t%d\t%d\t%s\n", filename, seqformat, fileid,libname);
	   		if (nFofn >= maxFofn) {
	   			maxFofn = nFofn + 20;
	   			hashFofn= ( long*) realloc(hashFofn, maxFofn * sizeof( long));
	   			}
	   		if (hashFofn) hashFofn[nFofn]= filehash;
	   		nFofn++;
   			linenum++;
    		}
    	}

	  if (fileid == lastfileid && docstart == lastdocstart)
	  	goto skipListFSE;
    
		/* find data end -- shouldn't this be in entry record somewhere ! */
		IdCopy( &nextid, entry->id);
		nextid.fip++; /* skip to next entry in library */
		nextentry = EntryOpen (&nextid);
   	if (nextentry) { 
		  if ( EntryOpenText ( nextentry) ) dataend= nextentry->fip[0]; 
    	EntryClose (&nextentry);
			}
			
		/* now write it all out -- use binary output to file for compactness? */			
		fprintf( file, "%d\t%ld\t%ld\t%ld\n", fileid, docstart, datastart, dataend);

skipListFSE: 
		 /* save to screen out dups */
		/* strncpy(lastename, EntryGetFullName(entry), sizeof(lastename)); */
		lastfileid= fileid;
		lastdocstart= docstart;
		}

  if ( entryCurrN == entryN - 1 ) {
		if (file) { fclose (file); file = NULL; }  
		if (fofn) { fclose(fofn); fofn= NULL; } /* dgg */
		if (hashFofn) { free(hashFofn); hashFofn= NULL; }
		}
}	

	/*--------end ListFSE----------- */	




/****** SlbEntryPrint *********************************************************
**
**      prints a single sequence entry using various options; function needs
**      to know how many entries have to be printed and are currently printed
**      in order to decide if the output file (if any) has to be opened or
**      closed.
**
**      INPUT:      entry object [R]
**                  set Name (for constructing output file name) [R]
**                  number of entries to print in whole list [R]
**                  number of current entry [R]
**              IMPLICIT:
**
**      RETURNS:    
*/

INT4 gSeqFormat = -1;

void SlbEntryPrint (ENTRYo *entry, char *setName, int entryN, int entryCurrN)
{
  SLBoFEATURE *feature;
  SEQo        *seq = NULL;
  char        *fileExt;
  char	      * forms, * tmp;
  INT4        errCode, seqFormat, (*print)(char*,...);

#if 1
     tmp= ParGetStr("listFSE");
	if ( tmp && *tmp) {
		/* a self-contained function - can be called from other places - dgg */
		SlbListFSE( entry,  setName,  entryN,  entryCurrN);
		}
#endif
		
	forms= ParGetStr("seqFormat");
        switch (tolower (forms[0])) {
	  case 'a':  /* dgg */
	  	seqFormat= SEQxAsIs;
	  	fileExt= ".seq";
	  	break;
	  case 'c':  /* dgg */
	  	seqFormat= SEQxCODATA;
	  	fileExt= ".cdt";
	  	break;
	  case 'g':  /* dgg */
		  if (tolower(forms[1]) == 'c') {
		    seqFormat = SEQxGCG;
		    fileExt = ".gcg";
		  	}
			else { 
		  	seqFormat= SEQxGENBANK;
		  	fileExt= ".gb";
		  	}
	    break;
	    
	  case 'p':
	    seqFormat = SEQxPIR;
	    fileExt = ".pir";
	    break;
	  case 'e':
	    seqFormat = SEQxEMBL;
	    fileExt = ".dat";
	    break;
	  case 's':
	    seqFormat = SEQxSWISS;
	    fileExt = ".dat";
	    break;
	  case 'f':
	    seqFormat = SEQxFASTA;
	    fileExt = ".fas";
	    break;
	  }
  
  gSeqFormat = seqFormat; /* pass to SlbGetSequence w/o parameter change - dgg */

  /* 
  ** the next deals with the file name... part of which may be specified by 
  ** user and open the file ...if not already open
  */

  if (ParGetNum ("printToFile")) {
    char *tmp, *outFileName, outDirName[FILxXNAM+1], fileName[FILxXNAM+1];

    ParDefFunction ("printf", (INT4 (*)()) SlbPrintF);
    outFileName = ParGetStr ("outFileName");	  
    
    if ((tmp = ParGetStr ("outDirName")))
	    sprintf (outDirName, "%s%c", tmp, DIRDELIMCHR);    /* '/' was unix-ism */
    else
      *outDirName = '\0';
    
    if (!entryCurrN && seqFormat != SEQxGCG) {
      if (outFileName && *outFileName)
				sprintf (fileName, "%s%s", outDirName, outFileName);
      else
				sprintf (fileName, "%s%s%s", outDirName, setName, fileExt);
    	}
    else if (seqFormat == SEQxGCG) {
      if (entryN == 1 && outFileName && *outFileName)
				sprintf (fileName, "%s%s", outDirName, outFileName);
      else
				sprintf (fileName, "%s%s%s", outDirName, entry->entry_nm, fileExt);
    	}
    if (*fileName) {
      file = FilOpenW (fileName, &errCode);
      _ErrMsg2 (errCode, fileName);
    	}
 	 }    

  /*  
  ** print text (annotation)
  */
  
  if (!(print = (INT4 (*)(char*,...)) ParGetFunction ("printf")))
    print = (INT4 (*)(char*,...)) printf;

  if (ParGetNum ("printText")) 
    EntryPrint (entry, seqFormat == SEQxGCG? 1 : 0); 
  else if (ParGetNum ("parSetFields")) {
    printf ("entry: %s \n", entry->full_nm);
    EntryPrintFields (entry); 
  }
  else if (seqFormat != SEQxFASTA)
    print ("entry: %s \n", entry->full_nm);

	
  /* 
  ** print data (sequence)
  */

  if (ParGetNum ("printData")) {
    if (IdIsSub (entry->id)) 
      errCode = SlbGetFeature (entry, &seq, &feature);
    else
      errCode = SlbGetSequence (entry, &seq);
    
    _ErrIf (errCode) 
      _ErrMsg2 (e__failureinslbgetx, entry->full_nm);

    else {
      switch (seqFormat) {
      	
      case SEQxAsIs:
      	{
      	char *cp, *maxp, *ep;
      	char line[SEQxXLN+1];
      	char hasnl;
      	long nline;
      	
      	cp= seq->seq;
      	maxp= seq->seq + seq->len;
      	while (cp<maxp) {
      		ep= strchr(cp, '\n');
      		if (!ep) { 
      			 hasnl= 0;
      			 if (maxp-cp < 60) ep= maxp; 
      			 else ep= cp+60; 
      			 }
      		else hasnl= 1;
      		nline= ep - cp; 
      		if (nline > SEQxXLN) nline= SEQxXLN;
      		strncpy( line, cp, nline);
      		line[nline]= 0;
     			print( "%s\n", line);
       		cp = ep;
      		if (hasnl) cp++;
      		}
  			}
      	break;
      	
      case SEQxGENBANK:
 				SeqWriteGenbank(seq, entry->entry_nm, print);
     		break;
      case SEQxCODATA:
				SeqWriteCODATA(seq, entry->entry_nm, print);
      	break;

      case SEQxPIR:
	SeqWritePIR (seq, entry->entry_nm, 0, print);
	break;
      case SEQxEMBL:
	SeqWriteEMBL (seq, entry->entry_nm, print);
	break;
      case SEQxSWISS:
	SeqWriteSwiss (seq, entry->entry_nm, print);
	break;
      case SEQxGCG:
	SeqWriteGCG (seq, entry->entry_nm, print);
	break;
      case SEQxFASTA:
	SeqWriteFasta (seq, entry->entry_nm, print);
	break;
      }
    }
  }

  /*
  ** close file 
  */

  if (file && (entryCurrN == entryN - 1 || seqFormat == SEQxGCG)) {
    fclose (file);
    file = NULL;  /* set back to default value /stdout/ */
 	 }
}


/****** SlbGetNamedFunction ***************************************************
**
**      looks for the named function in the local function table
**
**      INPUT:      address of function name [R]
**              IMPLICIT:
**
**      RETURNS:    address of function descriptor, if function exists,
**                  otherwise NULL
*/

static FUNC_DESC* SlbGetNamedFunction(char* name)
{
  INT4 i;
  for (i=0; i<MAXFUNC; i++)
      if (!strcmp(func_table[i].func_name, name)) return &func_table[i];

  return NULL; /* Function undefined */
}


/****** join ******************************************************************
**
**      implements the join function of a sequence expression.
**      As the parameter stack and evaluation stack is not used any more
**      SlbJoinList could be called directly.
**
**      INPUT:      address of parameter stack [R]
**                  address of evaluation stack [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
*/

static INT4 join(STKo *parameter_stack, STKo *eval_stack)
{
#ifdef TRACE
  printf("Function *join* called\n");   /* only test */
  printf("parameterstack should be empty (new list version) !!\n");

  StkPrintAll(parameter_stack);   /* only test */
#endif
    if (!StkIsEmpty(parameter_stack)) 
        _ErrRet2(e__joinillegalpar, "in stack");

    /*
    ** call function that organizes the list structures
    */
    SlbJoinList();

    return 1;
}

/****** complement ************************************************************
**
**      implements the complement function of a sequence expression.
**      As the parameter stack and evaluation stack is not used any more
**      SlbComplementList could be called directly.
**
**      INPUT:      address of parameter stack [R]
**                  address of evaluation stack [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
*/

static INT4 complement(STKo *parameter_stack, STKo *eval_stack)
{
#ifdef TRACE
  printf("Function *complement* called\n");   /* only test */
  printf("parameterstack should be empty (new list version) !!\n");
  StkPrintAll(parameter_stack);   /* only test */
#endif

    /*
    ** call function that organizes the list structures
    */
    SlbComplementList();

    return 1;
}

/****** replace ***************************************************************
**
**      implements the replace function of a sequence expression.
**      The string "replace" has to be set in the last inserted
**      sequence object.
**
**      INPUT:      address of parameter stack [R]
**                  address of evaluation stack [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
*/

static INT4 replace(STKo *parameter_stack, STKo *eval_stack)
{
  char repstring[100];

  strcpy (repstring, "_ERROR_");

#if defined(TRACE)
  printf("Function *replace* called\n");
  printf("parameterstack should contain replace string !!\n");
  StkPrintAll(parameter_stack); 
#endif

    /* get replace string from parameter_stack */
    ExpPopVarname(parameter_stack, repstring);
#if defined(TRACE)
    printf("** repstring popped from StacK: *%s*\n", repstring);
#endif

    /*
    ** call function that organizes the list structures
    */
    SlbReplaceList(repstring);

    return 1;
}

/****** one_of ****************************************************************
**
**      implements the one_of function of a sequence expression.
**
**      INPUT:      address of parameter stack [R]
**                  address of evaluation stack [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
**                  e__foundoneof <entry name>   
*/

static INT4 one_of(STKo *parameter_stack, STKo *eval_stack)
{
    STKoSYM op;
    STKoSYM first;
    INT4 param_count = 0;
    INT4 rv=1;

#if defined(TRACE)
    printf("Function *one_of* called\n");
#endif

    if (ParGetNum("discardoneof")) {
        _ErrRet(e__foundoneof);
    }

    /*
    ** auf dem Stack koennen Zahlen (pos) liegen oder es kann eine
    ** Parameterliste mit Ranges existieren !!
    */

    while (StkPop(parameter_stack, &op) != e__eostack)
    {
        param_count++;
        /* printf("Parameter: %d\n", param_count); */
        if (param_count==1) first = op;
    }

    if (param_count > 0) {
        rv = StkPush(eval_stack, first);
        _ErrRet2(rv, "when evaluating one-of");
    }
    else {   /* param_list with ranges */
        LstFirst ((void**)&param_list);
        if (LstNext ((void**)&param_list)) 
            LstDeleteToEndOfList((void**)&param_list);
    }

    return 1;
}

/****** SlbExe ****************************************************************
**
**      this function implements the operator EXE.
**      Symbols are popped from the evaluation stack until the next
**      function name is popped. If there are values on the stack,
**      the values are pushed into a private parameter stack. If the
**      function name is popped, the corresponding function is called
**      and supplied with the parameter stack.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  reads the content of the evaluation stack
**                  writes the content of the parameter stack
**
**      RETURNS:    1 on success
**                  e__doesnotexistfunction
*/

INT4 SlbExe (STKo *eval_stack)
{
  /* pop symbols from eval_stack,
     if operator: execute them (not possible here!)
     if value: put it into private parameter_stack
     if variable: put it's value into private parameter_stack
     if functionname: execute function giving it a pointer
                      to the private parameter_stack
     if OP_SEQ: nothing to do
  */

  STKo private_stack;
  STKoSYM sym;
  FUNC_DESC *func_desc_ptr;
  float float_val;
  INT4 int_flag=1, rv=1;

  if (!StkIniSym(&private_stack, 20)) {
      _ErrRet2(e__allocfail, 
               "in initialization of private parameter stack");
  }
   
  while (1) {
    rv = StkPop(eval_stack, &sym);

    if (sym.any.code == STKxFNAME) {
      /* get function_description pointer ... */
      func_desc_ptr = SlbGetNamedFunction(sym.var_name.name);
 /* ??? Should there be a number that decides what parameter list to use ? */
      if (func_desc_ptr)
	rv = (func_desc_ptr->func) (&private_stack, eval_stack);
      else rv = e__doesnotexistfunction;
      _ErrIf(rv) {
	StkDelSym(&private_stack);
	_ErrRet2(rv, sym.var_name.name);
      }
      break;
    }
    
    switch (sym.any.code) {
    case 0:
      break;
    case STKxNAME:
#ifdef TRACE
      printf(".. trying to put *%s* in parameter_stack\n", sym.var_name.name);
#endif
      rv = StkPush(&private_stack, sym);
      _ErrIf(rv) {
	StkDelSym(&private_stack);
	_ErrMsg(rv); /* because _ErrRet does not print message */
	_ErrRet(rv);
      }
      break;
      
    case STKxNUMBER:  /* put value into private parameter stack */
    case STKxFLOAT:
    case STKxINUM:
    case STKxFNUM:
      int_flag = 1;
      ExpGetVal(&sym,&float_val,&int_flag);
#ifdef TRACE
if (int_flag)
  printf(".. trying to put %d in parameter_stack\n", (INT4)float_val);
else
  printf(".. trying to put %f in parameter_stack\n", float_val);
#endif

      ExpGetSym(&sym,float_val,int_flag);
      rv = StkPush(&private_stack, sym);
      _ErrIf(rv) {
	StkDelSym(&private_stack);
	_ErrMsg(rv); /* because _ErrRet does not print message */
	_ErrRet(rv);
      }
      break;
    }

  }

  StkDelSym(&private_stack);
  return 1;
}


/****** SlbOpAfter ************************************************************
**
**      this function implements the operator AFTER ( > ).
**      An integer value is popped from the evaluation stack. 
**      Negative positions are rejected.
**      Legal positions are made negative afterwards to mark the
**      existence of the 'after' operator.
**      The value is pushed onto the evaluation stack again. 
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**
**      RETURNS:    1
*/

INT4 SlbOpAfter (STKo *eval_stack)
{
  float op1val, ergval;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpAfter* called\n");
  if (protokoll)
      fprintf(protokoll, "***** Operator AFTER *****\n");
#endif

  if (!ParGetNum("takeuncompletefeature"))
      _ErrRet(e__operatorafternotallowed);
  else global_after_f = 1;

  rv = ExpPopVal(eval_stack, &op1val, &int_flag);
  if (!int_flag) 
      _ErrRet2(e__illegalvalueonstack, op1val);

  /* negative positions must not exist, parser should signal error before */
  if (op1val<0)
      _ErrRet(e__illegalposition);

  ergval = op1val>0 ? op1val : 1;  /* default */

  /*
  **  make ergval negative to use it later as a sign
  **  for the existence of after_f in range
  */

  ergval = -ergval;

#ifdef TRACE
  printf(" > %d = %d\n", (INT4)op1val, (INT4)ergval);
#endif
  ExpPushVal(eval_stack, ergval, int_flag);
  return rv;
}


/****** SlbOpBefore ***********************************************************
**
**      this function implements the operator BEFORE ( < ).
**      An integer value is popped from the evaluation stack. 
**      Negative positions are rejected.
**      Legal positions are made negative afterwards to mark the
**      existence of the 'before' operator.
**      The value is pushed onto the evaluation stack again. 
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**
**      RETURNS:    1
*/

INT4 SlbOpBefore (STKo *eval_stack)
{
  float op1val, ergval;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpBefore* called\n");
#endif

  if (!ParGetNum("takeuncompletefeature"))
      _ErrRet(e__operatorbeforenotallowed);
  else global_before_f = 1;

  rv = ExpPopVal(eval_stack, &op1val, &int_flag);
  if (!int_flag)
      _ErrRet2(e__illegalvalueonstack, op1val);

/* negative positions must not exist, parser should signal error before */
  if (op1val<0)
      _ErrRet(e__illegalposition);

  ergval = op1val>1 ? op1val : 1;  /* default */

/*
**  make ergval negative to use it later as a sign
**  for the existence of before_f in range
*/
  ergval = -ergval;

#ifdef TRACE
  printf(" < %d = %d\n", (INT4)op1val, (INT4)ergval);
#endif
  ExpPushVal(eval_stack, ergval, int_flag);
  return rv;
}

/****** SlbOpRange ************************************************************
**
**      this function implements the operator RANGE ( .. ).
**      Two integer values are popped from the evaluation stack.
**      A new sequence object is inserted into the param_list, using
**      the popped value as begin and end position. The accession number
**      of the new sequence object is initially an empty string.
**      If a position is negative it means that an operator (before, after)
**      existed. Therefore a flag (before_f, after_f) is set TRUE in SlbInsert.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**                  inserts a new sequence object into param_list
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbOpRange (STKo *eval_stack)
{
  float op1val, op2val;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpRange* called\n");
#endif
  rv = ExpPopVals(eval_stack, &op1val, &op2val, &int_flag);
  if (!int_flag)
      _ErrRet3(e__illegalvaluesonstack, op1val, op2val);

#ifdef TRACE1
  printf("computing range %d .. %d \n", (INT4)op1val, (INT4)op2val);
#endif

  if (vabs(op1val) > vabs(op2val))  /* remember the existence of operators */
      _ErrRet(e__illegaldirectioninrange);

  /* Generate new range and insert it into list */
  rv = SlbInsert("", (INT4)op1val, (INT4)op2val);
  return rv;
}

/****** SlbOpMakeRange ********************************************************
**
**      this function implements the operator MAKERANGE ( .. ).
**      One integer value gets popped from the evaluation stack.
**      A new sequence object is inserted into the param_list, using
**      the popped value as begin and end position. The accession number
**      of the new sequence object is initially an empty string.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**                  inserts a new sequence object into param_list
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbOpMakeRange (STKo *eval_stack)
{
  float op1val;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpMakeRange* called\n");
#endif
  rv = ExpPopVal(eval_stack, &op1val, &int_flag);
  if (!int_flag)
      _ErrRet2(e__illegalvalueonstack, op1val);

  if (op1val<=0)
      _ErrRet(e__illegalposition);

#ifdef TRACE
  printf("computing range %d .. %d \n", (INT4)op1val, (INT4)op1val);
#endif

  /* Generate new range and insert it into list */
  rv = SlbInsert("", (INT4)op1val, (INT4)op1val);
  return rv;
}


/****** SlbOpSelect ***********************************************************
**
**      this function implements the operator SELECT ( . ).
**      Two integer values are popped from the evaluation stack.
**      The user is asked to specify a position between that positions.
**      The user-specified or default position is pushed onto stack again.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**
**      RETURNS:    1
*/

INT4 SlbOpSelect (STKo *eval_stack)
{
  char answer[20];
  float op1val, op2val, ergval, userval;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpSelect* called\n");
  if (protokoll)
      fprintf(protokoll, "***** Operator SELECT *****\n");
#endif

  rv = ExpPopVals(eval_stack, &op1val, &op2val, &int_flag);
  if (!int_flag)
      _ErrRet3(e__illegalvaluesonstack, op1val, op2val);

  ergval = op1val;  /* default */

  if (ParGetNum("interactiveselect")) {
      printf("what position do you want BETWEEN %d and %d [%d]: ",
              (INT4)op1val, (INT4)op2val, (INT4)ergval);
      gets(answer);
      if (strlen(answer) > 0) {
         if (sscanf(answer, "%f", &userval) == 1) {
            if (userval < op1val || userval > op2val)
               printf("** position too low or too high: default used !!\n");
            else ergval = userval;
         }
         else printf("** illegal input: default used !!\n");
      }
  }

#ifdef TRACE
  printf(" %d.%d = %d\n", (INT4)op1val, (INT4)op2val, (INT4)ergval);
#endif
  ExpPushVal(eval_stack, ergval, int_flag);
  return rv;
}

/****** SlbOpCut **************************************************************
**
**      this function implements the operator CUT ( ^ ).
**      Two integer values are popped from the evaluation stack.
**      A new sequence object is inserted into the param_list, using
**      the popped values as begin and end position. Because there should
**      be only a site specified, begin and end position must be adjacent.
**      At the moment the end position is set according to begin position.
**      The accession number of the new sequence object is initially an
**      empty string.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**                  inserts a new sequence object into param_list
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbOpCut (STKo *eval_stack)
{
  char answer[20];
  float op1val, op2val, ergval, userval;
  INT4 int_flag=1;
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpCut* called\n");
  if (protokoll)
      fprintf(protokoll, "***** Operator CUT *****\n");
#endif
  rv = ExpPopVals(eval_stack, &op1val, &op2val, &int_flag);
  if (!int_flag)
      _ErrRet3(e__illegalvaluesonstack, op1val, op2val);
  ergval = op1val;  /* default */

  if (op2val-op1val > 1 && ParGetNum("interactivecut")) {
      printf("after what position do you want to CUT %d..%d [%d]: ",
              (INT4)op1val, (INT4)op2val, (INT4)ergval);
      gets(answer);
      if (strlen(answer) > 0) {
         if (sscanf(answer, "%f", &userval) == 1) {
            if (userval < op1val || userval > op2val-1)
               printf("** position too low or too high: default used !!\n");
            else ergval = userval;
         }
         else printf("** illegal input: default used !!\n");
      }
  }

#ifdef TRACE
  printf(" %d^%d = %d^%d\n",(INT4)op1val,(INT4)op2val,(INT4)ergval,(INT4)ergval+1);
#endif
  op1val = ergval;
  op2val = op1val + 1;  /* the two positions must be adjacent */
  /* Generate new sequence object and insert it into list */
  rv = SlbInsert("", (INT4)op1val, (INT4)op2val);
  if (rv) 
      param_list->cut_f = 1; /* set flag (site_flag) */
  return rv;
}


/****** SlbOpSetAccno *********************************************************
**
**      this function implements the operator SET_ACCNO ( : ).
**      The accession number is popped from the evaluation stack.
**      Modifies (sets) the accession number of the sequence object,
**      that was last inserted into the param_list.
**
**      INPUT:      address of evaluation stack [R]
**              IMPLICIT:
**                  modifies the content of the evaluation stack
**                  modifies last inserted sequence object in param_list
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbOpSetAccno (STKo *eval_stack)
{
  char accno[NAMELEN];
  INT4 rv;

#ifdef TRACE1
  printf("Function *SlbOpSetAccno* called ");
#endif

  /*
  ** get accession number from evaluation stack
  */
  rv = ExpPopVarname(eval_stack, accno);
#ifdef TRACE
  printf("(%s)\n", accno);
#endif

  /*
  ** set accession number of last inserted range in list
  */
  rv = SlbModifyAccno(accno);

  return rv;
}


/****** SlbPrintFeature *******************************************************
**
**      this function is used to concatenate all lines belonging to
**      the same feature description in a database. After concatenating
**      all lines into one string, the resulting string will be parsed
**      and evaluated.
**
**      INPUT:      address of line in feature table [R]
**              IMPLICIT:
**                  modifies the global featureline
**
**      RETURNS:
*/

static void SlbPrintFeature(char *ln)
{
  INT4 len;

#if defined(TRACE)
  printf ("%s", ln);
#endif
  memcpy (ln, "  ", 2);   /* blank FT */
  len = strlen (ln);
  if (ln[ len-1 ] == '\n')
    ln[len-1] = '\0';
  BuffCopyString (featureBuff, ln);
}


/**api* SlbMakeFeature *****************************************************
**
**      accesses the feature table of an entry from the sequence library
**      and returns the evaluated sequence object; also returns a pointer
**      to the global feature object, wherein translation exceptions and
**      special codons may be stored.
**
**      INPUT:      buffer with the feature description [R]
**                  address of entry object [R]
**                  pointer to address of sequence object [W]
**                  pointer to address of feature description object [W]
**              IMPLICIT:
**                  featureBuff (SMoBUFF *) [W]
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbMakeFeature (SMoBUFF *featureBuff, ENTRYo *entry, SEQo *entry_seq,
		     SEQo **res_seq, SLBoFEATURE **featurePtr)
{
  static SMoBUFF *buff=NULL;

#if defined(TRACE1)
  static INT4 feature_count = 1;
#endif
  EXPo  *expression;
  INT4   rv;
  char  *tmp;

  if (!buff) 
    buff = BuffNew (100);
  else
    BuffReset (buff);

  global_entry = entry;   /* needed in SlbSeqFromAccNo to get library name */ 
  globalSeq = entry_seq;  /* entry sequence already read */
  SlbInitFeature (&feature); 
  if (featurePtr)
    *featurePtr = &feature;

  expression = (EXPo*) LibObjByName("expression", "embl_feature");
  if (!expression)
    _ErrRet4(e__objnotfound, 5, "seqlib.c", "embl_feature");

  ExpInit(expression);
 
  SlbListReset();  /* maybe at another place ??? */

  /*
  **  ExpEval builds only a list structure (param_list), that has
  **  to be evaluated later. The sequence of the entry was passed
  **  to SlbMakeFeature via the attribute entry_seq.
  **  After ExpEval parsing is already completely done and 
  **  feature.transl_except points to one of maybe many translation
  **  exception objects in a list, if there exists one.
  */

  BuffCopyString (buff, BuffGetPtr (featureBuff)+5);
  while ((tmp = BuffSearch (buff,"\nFT")))
    strncpy(tmp,"   ",3);

  rv = ExpEval (BuffGetPtr (buff), "location", EXPxEVALUATE);
  if (rv == e__parsefail)
    rv = e__entryparsefail;
  _ErrRet2(rv, entry->entry_nm);
  rv = SlbShiftAndCheck();
  _ErrMsg2(rv, entry->entry_nm);
  _ErrRet(rv);

  rv = SlbEvaluate();
  _ErrRet2(rv, entry->entry_nm);
    
  /*
  **  the last inserted one is the RESULT
  **  move it to the result sequence list
  */

  LstMove((void**) res_seq, (void**)&seq);
  globalSeq = NULL;  /* reset to initial state */

  return 1;
}


/**api* SlbGetFeature *********************************************************
**
**      accesses the feature table of an entry from the sequence library
**      and returns the evaluated sequence object; also returns a pointer
**      to the global feature object, wherein translation exceptions and
**      special codons may be stored.
**
**      INPUT:      address of entry object [R]
**                  pointer to address of sequence object [W]
**                  pointer to address of feature description object [W]
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbGetFeature (ENTRYo *entry, SEQo **res_seq, SLBoFEATURE **featurePtr)
{
  EXPo  *expression;
  INT4   rv, (*printSave)();
  
  SlbInitFeature (&feature);
  if (featurePtr)
    *featurePtr = &feature;
  expression = entry->lib->form->featureExpression;
  ExpInit(expression);
   SlbListReset(); 
  global_entry = entry; 
  
  /*
  ** goto begin of feature in entry
  */
  rv=EntryOpenText(entry);
  EntryBeginSub(entry);
  while ( (rv = EntryNextSub (entry)) < entry->id->subentry_n) 
    if (_ErrIs (rv) || rv == 0) 
      _ErrRet2 (e__entryparsefail, entry->entry_nm);

  /*
  ** print feature to the feature buffer
  */
  if (!featureBuff)
    featureBuff = BuffNew (1000);
  else
    BuffReset (featureBuff);
  printSave = ParGetFunction ("printf");
  ParDefFunction ("printf", (INT4(*)()) SlbPrintFeature);
  EntryPrintSub (entry);  /* copy feature to featureBuff */
  ParDefFunction ("printf", printSave);  /* reset print function */

  /*
  **  ExpEval builds only a list structure (param_list), that has
  **  to be evaluated later. The sequence of the entry has not
  **  been read yet. After ExpEval parsing is already completely
  **  done and feature.transl_except points to one of maybe many
  **  translation exception objects in a list, if there exists one.
  */

  rv = ExpEval(BuffGetPtr (featureBuff), "location", EXPxEVALUATE);
  if (rv == e__parsefail)
    rv = e__entryparsefail;
  _ErrRet2(rv, entry->entry_nm);
  rv = SlbShiftAndCheck();
  _ErrRet(rv); 

  if (ParGetNum ("listFOSN") && !ParGetNum ("printData")) 
    SlbWriteFosn (param_list);
  else {    
    rv = SlbEvaluate();
    _ErrRet(rv); 
    
    /*
    **  the last inserted one is the RESULT
    **  move it to the result sequence list
    */
    LstMove ((void**) res_seq, (void**)&seq);
  }
  return 1;
}


/**api* SlbWriteFosn **********************************************************
**
**      writes list of features as extended FOSN file with begin and end
**      positions within parent sequence, complement...join
**
**      INPUT:      address of sequence object [R]
**              IMPLICIT:
**
**      RETURNS:    nothing
*/

static void SlbWriteFosn (SEQo *seq)
{
  SEQo *joinSeq;
  char hexNumStr[20], l[133];
  INT4  (*print)();

  IdToStr (global_entry->id, hexNumStr);
  
  if (!(print = ParGetFunction ("printf")))
    print = (INT4 (*)()) printf;
  
  /*
  ** check first if there is a sublist (-> join) and process that
  */ 

  if (seq->sublist && ParGetNum ("fosnWithPos")) {  
    joinSeq = seq->sublist;
    
    /*
    ** start with the last sequence if the WHOLE join must be complemented
    ** in that case all individual joined sequences must be complemented 
    */
    
    if (seq->compl_f) 
      LstLast ((void **) &joinSeq);
    else
      LstFirst ((void **) &joinSeq);
    
    do {
      if (seq->compl_f)
	joinSeq->compl_f = !joinSeq->compl_f;
      
      if (!*joinSeq->name)
	sprintf (l, "%s:%s", global_entry->lib->nam, 
		 global_entry->entry_nm);
      else
	sprintf (l, "%s", joinSeq->name);
      sprintf (l, "%s  Begin: %d End: %d", l, joinSeq->bpos, joinSeq->epos);
      sprintf (l, "%s Strand: %c", l, joinSeq->compl_f ? '-' : '+');
      sprintf (l, "%s Join: %s_%d", l,  
	       global_entry->entry_nm, global_entry->id->subentry_n);
      sprintf (l, "%s !Id: %s\n", l, hexNumStr);
      print (l);
    } while (seq->compl_f ? LstPrev ((void **) &joinSeq) :
	     LstNext ((void **) &joinSeq));
  }
  
  /*
  ** normal sequence
  */

  else {
    sprintf (l, "%s:%s", global_entry->lib->nam, global_entry->entry_nm);
    if (ParGetNum ("fosnWithPos")) {
      sprintf (l, "%s  Begin: %d End: %d", l, seq->bpos, seq->epos);
      sprintf (l, "%s Strand: %c", l, seq->compl_f ? '-' : '+');
    }
    sprintf (l, "%s !Id: %s\n", l, hexNumStr);
    print (l);
  }
}



/**api* SlbNextSequence *******************************************************
**
**      'Walks' to the next sequence in the file;
**
**      INPUT:      address of entry object [W]
**              IMPLICIT:
**
**      RETURNS:    1 on success
**                  0 if at the end of the file
*/

INT4 SlbNextSequence (ENTRYo *entry)
{
  SRSoFILTYP *fileType=entry->lib->form->fil_t[1];
  FILo       *file=entry->file[1];
  INT4       k;

  while (1) {				     
    if (SmXLoc (file->l, fileType->fsts, fileType->fstl, fileType->fstr) != -1)
      break;
    if (FilURead (file) == e__eof)
      return 0;
  }
  
  /* 
  ** advance lines 
  */

  for (k=1; k < fileType->adv; k++) 
    FilURead (file);

  entry->isDataBegin = 1;
  entry->fip[1] = file->crfdd;
  return 1;
}


/**api* SlbGetSequence ********************************************************
**
**      accesses the sequence part of an entry from the sequence library
**      and returns a sequence object;
**
**      INPUT:      address of entry object [R]
**                  pointer to address of sequence object [W]
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
*/

INT4 SlbGetSequence (ENTRYo *entry, SEQo **seqPtr)
{
  SEQo       *seq=NULL;
  SRSoFILTYP *fileType;
  FIP        fipSave;
  INT4       isGetSeqLenOnly = ParGetNum ("getSeqLenOnly");  

#if defined(TRACE)
  printf("SlbGetSequence: %d\n", sequence_count++);
#endif

  if (entry->isDataBegin)
    entry->isDataBegin = 0;  /* stream access - file pointer is positioned */
  else
    if (!EntryOpenData (entry))
      return 0;

  /*
  ** get the file type object and assume that the first line after the
  ** sequence title is a comment line only if a production "find" for
  ** the entry name is specified
  */

  fileType = entry->lib->form->fil_t[ENTxDATA];

  SeqSetRecordSeparator (fileType->type == SLBxGCG7 ? 1 : 0); /* gcg format 
								 for vms */
  SeqFromLib (seqPtr, entry->file[ENTxDATA],
	      (*fileType->find) ? fileType->find : NULL, /* scanf gets name */
	      fileType->exit, 
	      (*fileType->find ? 1 : 0) | /* no title line -> no comment ... */
	      (isGetSeqLenOnly ? 2 : 0)
	      | ((gSeqFormat == SEQxAsIs) ? 4 : 0) /* dgg */
	      );   

  /*
  ** if the sequence is split into two files (GCG format) then get the next 
  ** sequence which is the continuation - a file is considered to be split
  ** if it has exactly the split size.
  */

  if ((*seqPtr)->len && (*seqPtr)->len == fileType->gcgSplitSize) {
    fipSave = entry->fip[1];
    SlbNextSequence (entry);
    SlbGetSequence (entry, &seq);
    if (isGetSeqLenOnly)
      (*seqPtr)->len += seq->len;
    else
      SeqCat (seqPtr, seq);
    SeqClr (&seq);
    entry->fip[1] = fipSave;
  }

  if (!*fileType->find) /* no title line -> no name either */
    strcpy ((*seqPtr)->name, entry->entry_nm);
  return 1;
}


/****** SlbJoinList ***********************************************************
**
**      implements the join function of a sequence expression.
**      Creates a new function object as first member of a local list.
**      The current parameter list is assigned to the sublist pointer
**      of the new function object. As there is no parameter list any more,
**      the local list can be assigned to the (new) parameter list.
**
**      INPUT:
**              IMPLICIT:
**                  param_list changes
**
**      RETURNS:    nothing
*/

static void SlbJoinList()
{
    SEQo *local_list = NULL;
    
    if (SeqNew(&local_list)) {
      /*sprintf (local_list->name, "%s", global_entry->entry_nm);*/
      strcpy(local_list->name, "JOIN");
      
      /* 
       **  no hashing  if LstHashSearch or LstDelete... is never used
       */
      LstAssign((void**)&local_list->sublist, (void**)&param_list);
      LstAssign((void**)&param_list, (void**)&local_list);
    }
    else printf("ERROR: JOIN object could not be created in SlbJoinList !\n");
}

/****** SlbComplementList *****************************************************
**
**      implements the complement function of a sequence expression.
**      Simply inverts the complement flag of the current sequence object
**      in the param_list.
**
**      INPUT:
**              IMPLICIT:
**                  inverts the complement flag of the current sequence object
**
**      RETURNS:    nothing
*/

static void SlbComplementList()
{
    param_list->compl_f = !param_list->compl_f;
}

/****** SlbReplaceList ********************************************************
**
**      implements the replace function of a sequence expression.
**      Simply sets the replace string of the current sequence object
**      in the param_list.
**
**      INPUT:
**                  address of replace string [R]
**              IMPLICIT:
**                  sets the replace string of the current sequence object
**
**      RETURNS:    nothing
*/

static void SlbReplaceList(char *replace)
{
    param_list->replace = (char*) malloc(strlen(replace) + 1);
    strcpy(param_list->replace, replace);
}

/****** SlbInsert *************************************************************
**
**      inserts a new (current) subsequence object (range) into the param_list.
**      The accession number of a sequence object is up to now initially
**      set to an empty string, i.e. the accession number equals the
**      accession number of the current entry.
**      In case of the CUT operator the cut_f flag has to be set afterwards.
**
**      INPUT:      address of accession number [R]
**                  start position of subsequence [R]
**                  end position of subsequence [R]
**              IMPLICIT:
**                  new sequence object in param_list
**
**      RETURNS:    1 on success, 0 otherwise
*/

static INT4 SlbInsert(char *acc, INT4 bpos, INT4 epos)
{
    if (!SeqNew(&param_list))  return 0;
    strcpy(param_list->accno, acc);

    if (bpos<0) {
        param_list->before_f = 1;
        param_list->bpos = -bpos;
    }
    else {
        param_list->before_f = 0;
        param_list->bpos = bpos;
    }

    if (epos<0) {
        param_list->after_f = 1;
        param_list->epos = -epos;
    }
    else {
        param_list->after_f = 0;
        param_list->epos = epos;
    }

    return 1;
}

/****** SlbModifyAccno ********************************************************
**
**      sets accession number of last (current) sequence object in the
**      param_list. The accession number of a sequence object is initially
**      set to an empty string, i.e. the accession number equals the
**      accession number of the current entry.
**
**      INPUT:      address of accession number [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
*/

static INT4 SlbModifyAccno(char *acc)
{
  if (param_list)
    strcpy(param_list->accno, acc);
  else {
    _ErrMsg(e__noseqinlist);
    return 0;
  }

  return 1;
}

/****** SlbModifyBegin ********************************************************
**
**      modifies the begin positiom of the specified sequence object.
**      If an illegal position (<1) would result, an error is
**      reported.
**
**      INPUT:      address of a sequence object [R]
**                  delta for shift [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
*/

static INT4 SlbModifyBegin(SEQo *range, INT4 delta)
{
  if (range)
    if (range->bpos + delta > 0)
      range->bpos += delta;
    else {
      if (!ParGetNum("takeuncompleteshiftfeature")) {
	_ErrMsg(e__illegalrangemod);
	return 0;
      }
      else range->bpos = 1;
    }
  else {
    _ErrMsg(e__noseqinlist);
    return 0;
  }
  
  return 1;
}

/****** SlbModifyEnd **********************************************************
**
**      modifies the end positiom of the specified sequence object.
**      If an illegal position (<1) would result, an error is
**      reported.
**
**      INPUT:      address of a sequence object [R]
**                  delta for shift [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
*/

static INT4 SlbModifyEnd(SEQo *range, INT4 delta)
{
  if (range)
    if (range->epos + delta > 0)
      range->epos += delta;
    else {
      if (!ParGetNum("takeuncompleteshiftfeature")) {
	_ErrMsg(e__illegalrangemod);
	return 0;
      }
      else range->epos = 1;
    }
  else {
    _ErrMsg(e__noseqinlist);
    return 0;
  }
  
  return 1;
}

void SlbDoTranslExcept(SEQo *seq, SLBoFEATURE *feature)
{
  char *except;
  if (LstFirst((void**)&feature->transl_except)) {
    do  {
      if (feature->transl_except->beginPos+2 > strlen(seq->seq))
	fprintf(stderr, "FATAL_ERROR in SlbDoTranslExcept\n");
      else {
	except = seq->seq + feature->transl_except->beginPos;
	*(except-1) = '@';
	*(except)   = *feature->transl_except->aa;
	*(except+1) = '@';
      }
    }
    while (LstNext((void**)&feature->transl_except));
    LstDeleteAll((void**)&feature->transl_except);
  }
}
	  
/****** SlbSeqFromAccNo *******************************************************
**
**      gets the complete sequence string belonging to the entry with the
**      specified accession number. To retrieve the sequence string the
**      unique accession number is transformed into a unique id_name.
**      If the sequence string is not already in the global list (seq)
**      it will be retrieved by a query and the sequence object will be
**      inserted into the global list.
**
**      INPUT:      address of accession number [R]
**              IMPLICIT:
**
**      RETURNS:    address of the complete sequence string
**                  0, in case of an error (message already printed)
*/

static SEQo *SlbSeqFromAccNo (char *accno)
{
    ENTRYo *entry = global_entry;
    SEQo  *local_list = NULL;
    char id_name[40], *tmp;

    /*
    **  get id_name
    */

    if (strlen(accno) > 0) {
      IDoENTRY *id;

      if ((tmp = ParGetStr ("findAccNoCommand")) && *tmp)
	return SlbRetrieveAccNo (accno);

      if (!(id = QryFindID ("embl", "acc", accno))/* &&
	  !(id = QryFindID ("emnew", "acc", accno) ) */) {
	_ErrMsg3(e__couldnotgetsequence,accno, global_entry->entry_nm);
	return NULL;
      }
      entry = EntryOpen(id);
    }
    else if (globalSeq)
      return globalSeq;  /* see SlbMakeFeature */

    sprintf (id_name, "%s:%s", entry->lib->nam, entry->entry_nm);
        
    /*
    ** search for the sequence with the name id_name in sequence list;
    ** if sequence not in list, get the new sequence object
    ** by a query, and store it into the sequence list.
    **
    ** If the sequence was fetched the first time, it may be possible that
    ** a couple of translation exceptions (feature.transl_except) exist.
    ** They are only important if accno is empty string, because otherwise 
    ** you would modify the sequence of another entry. The modification has
    ** to be done only once per entry, therefore the list of translation
    ** exceptions has to be deleted afterwards.
    */
    
#ifdef TRACE1
    printf("\n*** Looking for *%s* in list of sequences\n", id_name);
    LstPrintAllNames((void*)seq);
#endif
    
    if (!LstHashSearch((void**)&seq, id_name))  {
        if (!SlbGetSequence(entry, &local_list))  {
            if (entry != global_entry)
	      EntryClose(&entry);
            LstDeleteAll((void**)&local_list);
            _ErrMsg2(e__couldnotgetsequence, accno);
            return NULL;
        }
        if (entry != global_entry)
	  EntryClose(&entry);

        if (!LstMove((void**)&seq, (void**)&local_list))  {
            LstDeleteAll((void**)&local_list);
            _ErrExit(f__failedtomovelist);
        }
	strcpy(seq->name, id_name);

	/* 
	** check translation exceptions (see comment above)
	*/
	if (feature.transl_except && ! *accno) {
	  SlbDoTranslExcept(seq, &feature); 
	}   /* what to do else ??? */

    }
#ifdef TRACE1
    else printf("***  *%s* in list of sequences found !!!\n", id_name);
#endif

    if (entry != global_entry)
      EntryClose(&entry);  /*PURIFY*/
    return seq;
}


/****** SlbRetrieveAccNo *******************************************************
**
**      Uses the command specified in "findAccNoCommand" to get entry by 
**      accession number.
**
**      INPUT:      address of accession number [R]
**              IMPLICIT:
**
**      RETURNS:    address of the complete sequence string
**                  0, in case of an error (message already printed)
*/

static SEQo *SlbRetrieveAccNo (char *accno)
{
  SEQo *seq=NULL;
  char command[256], fileName[256];

  sprintf (command, ParGetStr ("findAccNoCommand"), accno, accno, accno);
  sprintf (fileName, ParGetStr ("inputSeqFile"), accno);
  system (command);
  return SeqReadAll (&seq, fileName) ? seq : NULL;
}


/****** SeqReplace ************************************************************
**
**
**      INPUT:      address of a sequence object [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
**                  e__incorrectlength
*/

INT4 SeqReplace (SEQo *new, SEQo *seq, char *replace, 
		 INT4 begPos, INT4 endPos, INT4 doInsert)
{
  INT4 length;

  /* calculate length of output sequence */
  length = strlen (seq->seq) + strlen (replace);
  if (!doInsert)  
    length -= endPos - begPos + 1;

  new->seq = (char*) malloc (length + 1);
  new->len = length;
  
  if (doInsert) {
    strncpy (new->seq, seq->seq, begPos);
    new->seq[begPos] = '\0';
    strcat (new->seq, replace);
    strcat (new->seq, &seq->seq[begPos]); /* endPos = begPos +1 */
  }
  else {
    strncpy (new->seq, seq->seq, begPos - 1);
    new->seq[begPos - 1] = '\0';
    strcat (new->seq, replace);
    strcat (new->seq, &seq->seq[endPos]);
  }
  
  if (new->len != strlen (new->seq)) 
    _ErrRet1 (e__incorrectlength);
  return 1;
}


/****** SeqSubSeq *************************************************************
**
**
**      INPUT:      address of a sequence object [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
**                  e__incorrectlength
*/

INT4 SeqGetSubSeq (SEQo *new, SEQo *seq, INT4 begPos, INT4 endPos)
{
  INT4 length;

  length = endPos - begPos + 1;
  if (length <= 0)  
    _ErrRet(e__negativelength);

  new->seq = (char*) malloc (length + 1);
  strncpy (new->seq, &(seq->seq[ begPos-1 ]), length);
  new->seq[length] = '\0';
  new->len = length;
  return 1;
}


/****** SlbEvalSequence *******************************************************
**
**      gets the sequence string belonging to the specified sequence object.
**      If the sequence string is not already calculated the calculation
**      will be done here using the complete sequence of an entry with
**      an unique accession number and the begin and end of a partial
**      sequence. At last the sequence string may be complemented.
**
**      INPUT:      address of a sequence object [R]
**              IMPLICIT:
**
**      RETURNS:    1 on success
**                  e__error, in case of an already printed error message
**                  e__postoohigh, position is greater than length of sequence
*/

static INT4 SlbEvalSequence (SEQo *new)
{
  SEQo *seq;
  INT4 rv;
  
  if (new->seq)
    return 1;
  
  /* fetch pointer (tmp) to complete sequence */
  if ((seq = SlbSeqFromAccNo (new->accno)) == NULL) 
    _ErrRet (e__error);
  
  /* 
  ** check translation exceptions (see comment above)
  */

  if (feature.transl_except && ! *new->accno)
    SlbDoTranslExcept(seq, &feature);
  
  strcpy (new->name, seq->name);
  
  /*
  ** report error if begin and end positions are greater than sequence length  
  */

  if (new->bpos > seq->len)  {
    if (!ParGetNum ("takeuncompletefeature"))
      _ErrRet1 (e__postoohigh);
    else new->bpos = seq->len;
  }
  if (new->epos > seq->len)  {
    if (!ParGetNum ("takeuncompletefeature")) 
      _ErrRet1 (e__postoohigh);
    else new->epos = seq->len;
  }
  

  if (new->replace)
    rv = SeqReplace (new, seq, new->replace, new->epos, new->bpos, new->cut_f);
  else 
    rv = SeqGetSubSeq (new, seq, new->bpos, new->epos);
  
  _ErrRet (rv);

  if (new->compl_f)
    SeqComplement (new);
  return 1;
}




/****** SlbShiftAndCheck *****************************************************
**
**      looks for the current values of the modification parameters (shift) and
**      checks the global param_list. In the param_list there could be
**      one sequence object generated by an expression like "135 .. 177" or
**      one function object that has a sublist of sequence objects generated
**      by an expression like "join(135 .. 177, complement(23 .. 97))".
**
**      INPUT:
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
**
*/

static INT4 SlbShiftAndCheck()
{
  SEQo *range;
  INT4 begPos=0, endPos=0, isBegRelEnd=0, isEndRelBeg=0;
  INT4 head_len=0, tail_len=0;
  INT4 rv;

  if (!param_list)  {
    _ErrMsg(e__nothingtoevaluate);
    return 0;
  }

  begPos = ParGetNum ("shiftbeginpos");
  endPos = ParGetNum ("shiftendPos");
  isBegRelEnd = ParGetNum ("isshiftbeginrelend");
  isEndRelBeg = ParGetNum ("isshiftendrelbegin");

  /*
  ** return error message for all combinations of
  ** modification parameters that should not occur
  */

  if ((isBegRelEnd && isEndRelBeg) || (isBegRelEnd && begPos < endPos) ||
      (isEndRelBeg && begPos < endPos)) {
    _ErrMsg(e__illegalmodparam);
    return 0;
  }

  /*
  ** compute length of head or tail subsequence
  */

  if (isEndRelBeg)
    head_len = endPos - begPos + 1;
  if (isBegRelEnd)
    tail_len = endPos - begPos + 1;

  param_list->part_len = 0;  /* whole sequence */
  if (head_len > 0) param_list->part_len = -head_len;
  if (tail_len > 0) param_list->part_len = tail_len;
  
  range = param_list->sublist;  /* for shorter source code */
  
  if (!range) {
    if (begPos && !isBegRelEnd)
      rv = SlbModifyBegin (param_list, begPos);
    if (endPos && !isEndRelBeg)
      rv = SlbModifyEnd (param_list, endPos);
    /* in case of a 0 returned there was not made a modification */
    
    rv = SlbEvalSequence(param_list);
    _ErrRet (rv);  /* already printed */
  }
  else  {

    /*
    **  range is the sublist of a (Join-) function
    */

    if (begPos && !isBegRelEnd)
      if (LstFirst((void**)&range))
	rv = SlbModifyBegin (range, begPos);
    
    if (endPos && !isEndRelBeg)
      if (LstLast((void**)&range))
	rv = SlbModifyEnd (range, endPos);
    
    if (LstFirst((void**)&range)) {
      do  {
	rv = SlbEvalSequence(range);
	_ErrRet (rv);  /* already printed */
      }
      while (LstNext((void**)&range));
    }
    else {
      _ErrMsg(e__noparamlist);
      return 0;
    }
  } 

  return 1;
}


/****** SlbEvaluate ***********************************************************
**
**      evaluates the global param_list. In the param_list there could be
**      one sequence object generated by an expression like "135 .. 177" or
**      one function object that has a sublist of sequence objects generated
**      by an expression like "join(135 .. 177, complement(23 .. 97))".
**
**      INPUT:
**              IMPLICIT:
**
**      RETURNS:    1 on success, 0 otherwise
**
*/

static INT4 SlbEvaluate()
{
  SEQo *range;
  
  if (!param_list)  {
    _ErrMsg(e__nothingtoevaluate);
    return 0;
  }
  
  range = param_list->sublist;  /* for shorter source code */
  
  if (!range) {
    /* the single range (param_list) has been evaluated already */
    SeqJoin (&seq, param_list);
  }
  else  {

    /*
    **  range is the sublist of a (Join-) function.
    **  Make sure that all ranges in this list
    **  are already evaluated (SlbEvalSequence)
    **  before joining them (see SlbShiftAndCheck)
    */
    
    if (LstFirst((void**)&range)) {
      SeqJoin (&seq, range);
    }
    else {
      _ErrMsg(e__noparamlist);
      return 0;
    }
  }
  
  /*
  ** give sequence a name, and complement it if necessary
  */
  if (param_list->compl_f && param_list->sublist)
    SeqComplement(seq);
  strcpy(seq->name, "RESULT");
  strcpy(seq->accno, "RESULT");
  
  if (param_list->part_len < 0)
    SeqHead(seq, -param_list->part_len);
  if (param_list->part_len > 0)
    SeqTail(seq, param_list->part_len);
  
  return 1;
}



/**api* SlbListReset **********************************************************
**
**      deletes all local lists (seq, param_list) created during evaluation
**
**      INPUT:
**              IMPLICIT:
**
**      RETURNS:    1
*/

INT4 SlbListReset()
{
  LstDeleteAll((void**)&param_list);
  if (seq)
    LstDeleteAll((void**)&seq);
  return 1;
}


/****** SlbCodonDump ****************************************
**
**      dumps all information belonging to the specified object.
**      This is the print function set through LstManageClass
**      for codontype objects (see SlbCodonNew).
**
**      INPUT:      address of an object [R]
**                  file pointer [W] or NULL
**              IMPLICIT:
**
**      RETURNS:    nothing
*/

static void SlbCodonDump (FILE *file, void *obj)
{
  CODON *tmp = (CODON*)obj;
  
  if (!file)
    file = stderr;
  
  fprintf (file, "/codon:  %s  ->  %s\n", tmp->seq, tmp->aa);
}


/****** SlbTranslExceptDump ****************************************
**
**      dumps all information belonging to the specified object.
**      This is the print function set through LstManageClass
**      for trextype objects (see SlbTranslExceptNew).
**
**      INPUT:      address of an object [R]
**                  file pointer [W] or NULL
**              IMPLICIT:
**
**      RETURNS:    nothing
*/

static void SlbTranslExceptDump (FILE *file, void *obj)
{
  TRANSL_EXCEPT *tmp = (TRANSL_EXCEPT*)obj;
  
  if (!file)
    file = stderr;
  
  fprintf (file, "/transl_except:  %d..%d  ->  %s\n",
	   tmp->beginPos, tmp->endPos, tmp->aa);
  
  if ((tmp->endPos - tmp->beginPos) != 2)
    fprintf (file, "ERROR: /transl_except:  %d..%d  ->  %s\n",
	     tmp->beginPos, tmp->endPos, tmp->aa);
  
}


/****** SlbCitationDump ****************************************
**
**      dumps all information belonging to the specified object.
**      This is the print function set through LstManageClass
**      for citationtype objects (see SlbCitationNew).
**
**      INPUT:      address of an object [R]
**                  file pointer [W] or NULL
**              IMPLICIT:
**
**      RETURNS:    nothing
*/

static void SlbCitationDump (FILE *file, void *obj)
{
  CITATION *tmp = (CITATION*)obj;

  if (!file)
    file = stderr;
  fprintf (file, "*** /CITATION %d\n", tmp->citationNum);
}


/**api* SlbCodonNew *************************************************
**
**      gets a new initialized codon object.
**      The codontype_id tells, if the class CODON is already
**      registered by the list module (LstManageClass).
**
**      INPUT:      pointer to address of codon object [W]
**              IMPLICIT:
**
**      RETURNS:
*/

INT4 SlbCodonNew (CODON **codon)
{
  static INT4 codontype_id=0;

  if (!codontype_id)
    LstManageClass (&codontype_id, sizeof (CODON),
		    NULL, NULL, SlbCodonDump);

  return LstNew ((void **) codon, codontype_id);
}

/**api* SlbTranslExceptNew *************************************************
**
**      gets a new initialized translation exception object.
**      The trextype_id tells, if the class TRANSL_EXCEPT is already
**      registered by the list module (LstManageClass).
**
**      INPUT:      pointer to address of translation exception object [W]
**              IMPLICIT:
**
**      RETURNS:
*/

INT4 SlbTranslExceptNew (TRANSL_EXCEPT **transl_except)
{
  static INT4 trextype_id=0;

  if (!trextype_id)
    LstManageClass (&trextype_id, sizeof (TRANSL_EXCEPT),
		    NULL, NULL, SlbTranslExceptDump);

  return LstNew ((void **) transl_except, trextype_id);
}

/**api* SlbCitationNew *************************************************
**
**      gets a new initialized citation object.
**      The citationtype_id tells, if the class CITATION is already
**      registered by the list module (LstManageClass).
**
**      INPUT:      pointer to address of citation object [W]
**              IMPLICIT:
**
**      RETURNS:
*/

INT4 SlbCitationNew (CITATION **citation)
{
  static INT4 citationtype_id=0;

  if (!citationtype_id)
    LstManageClass (&citationtype_id, sizeof (CITATION),
		    NULL, NULL, SlbCitationDump);

  return LstNew ((void **) citation, citationtype_id);
}


/****** SlbInitFeature ********************************************************
**
**      initializes global feature object ...functions that modify it can
**      be called directly form the parser.     
**      For managing translation exceptions the list module has been used.
**      Therefore, when resetting the feature object, an already existing
**      list has to be destroyed.
**
**      INPUT:      feature object [W]
**              IMPLICIT:
**
**      RETURNS:   
*/

static void SlbInitFeature (SLBoFEATURE *feature)
{
  feature->codeTableNum = 0;       /* default translation table */
  feature->codonStartShifted = 0;  /* normally there should be no shifting */
  feature->codonStart = 1;         /* default start position */
  feature->isPseudo = 0;

  /*
  **  if lists of special treated codons, translation exceptions or 
  **  citations exist, free the lists
  */
  if (feature->codon)
    LstDeleteAll((void**)&feature->codon);
  feature->codon = NULL;           /* no list */
  if (feature->transl_except)
    LstDeleteAll((void**)&feature->transl_except);
  feature->transl_except = NULL;   /* no list */
  if (feature->citation)
    LstDeleteAll((void**)&feature->citation);
  feature->citation = NULL;        /* no list */

  if (!feature->gene) {
    feature->gene = BuffNew (20);
    feature->product = BuffNew (20);
    feature->note = BuffNew (20);
  }
  else {
    BuffReset (feature->gene);
    BuffReset (feature->product);
    BuffReset (feature->note);
  }
}


/**api* SlbSetCodonStart *****************************************************
**
**            
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetCodonStart (PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  tokstr[20];

  PrsUnSym (tokList, tokstr, &tokval);
  feature.codonStart = atoi (tokstr);
  return 1;
}


/**api* SlbSetCodon **********************************************************
**
**            
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetCodon (PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  tokstr[20];

  SlbCodonNew (&feature.codon);  /* get new object */
  PrsUnSym (tokList, tokstr, &tokval);            /*  )  */
  PrsUnSym (tokList, feature.codon->aa, &tokval);   /* aa  */
  PrsUnSym (tokList, feature.codon->seq, &tokval);  /* seq */

  return 1;
}


/**api* SlbSetTranslTable  ***************************************************
**
**            
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetTranslTable(PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  tokstr[20];

  PrsUnSym (tokList, tokstr, &tokval);
  feature.codeTableNum = atoi (tokstr);

  return 1;
}


/**api* SlbSetTranslExcept ****************************************************
**
**            
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetTranslExcept (PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  tokstr[20];

#if defined(IGNORE_TRANSL_EXCEPT)
  PrsUnSym (tokList, tokstr, &tokval);  /*  )  */
  PrsUnSym (tokList, tokstr, &tokval);  /* aa  */
  PrsUnSym (tokList, tokstr, &tokval);  /* posend   */
  PrsUnSym (tokList, tokstr, &tokval);  /* posbegin */
  return 1; 
#endif

  SlbTranslExceptNew (&feature.transl_except);  /* get new object */
  PrsUnSym (tokList, tokstr, &tokval);  /*  )  */
  PrsUnSym (tokList, feature.transl_except->aa, &tokval);  /* aa  */
  PrsUnSym (tokList, tokstr, &tokval);  /* posend   */
  feature.transl_except->endPos = atoi (tokstr);
  PrsUnSym (tokList, tokstr, &tokval);  /* posbegin */
  feature.transl_except->beginPos = atoi (tokstr);

  return 1;
}


/**api* SlbSetCitation *****************************************************
**
**            
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetCitation (PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  tokstr[20];

  SlbCitationNew (&feature.citation);  /* get new object */
  PrsUnSym (tokList, tokstr, &tokval);
  feature.citation->citationNum = atoi (tokstr);
  return 1;
}


/**api* SlbSetPseudo **********************************************************
**
**       Called during feature parsing. Sets the "pseudo gene" flag.           
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSetPseudo (PRSoST *tokList,  FILo *file)
{
  feature.isPseudo = 1;
  return 1;
}

/**api* SlbSaveQualText *******************************************************
**
**      Called during feature parsing. Saves the text for certain qualifiers.
**
**      INPUT:
**              IMPLICIT:
**                  feature (SLBoFEATURE *) [W]
**
**      RETURNS:    1
*/

INT4 SlbSaveQualText (PRSoST *tokList,  FILo *file)
{
  UINT4 tokval;
  char  text[500], qualName[30];

  PrsUnSym (tokList, text, &tokval);
  PrsUnSym (tokList, qualName, &tokval);

  switch (tolower (qualName[0])) {
  case 'g':  /* gene */
    BuffCopyString (feature.gene, text);
    break;
  case 'p':  /* product */
    BuffCopyString (feature.product, text);
    break;
  case 'n':  /* note */
    BuffCopyString (feature.note, text);
    break;
  default:
    _ErrExit2 (e__unknownoption, qualName);
  }
  return 1;
}


