// SRS-FastA Notes
// d.gilbert, oct. 1995
// gilbertd@bio.indiana.edu


SRS-FastA modification notes:

This note describes modifications to Thure Etzold's Sequence Retrieval System
(SRS) and William Pearson's FastA software that allow one to perform
biosequence similarity searches on subsets of biosequence data libraries.

*) Perl and HTML Glue
*) FastA
*) SRS



*) Perl and HTML Glue

 Perl
  sfsubsearch.pl  -- do general SRS query, or select section query, then 
    Fasta similarity search.  Normally you would have this available
    from a HTTP server's "cgi-bin" folder, where it is called by 
    'sfsubsearch.html'
    
  sfsectionlib.pl -- generate index libraries for fixed section queries
   (e.g., species subsets).  Run this when SRS indices are updated.
   It will update data files and update the 'sfsectionsearch.html'
   list of data.
  
 HTML
  sfmain.html  -- main document with description of functions, and links to
    two forms below.  Publish this on an HTTP server.
  sfsubsearch.html  -- form for general SRS-Fasta query, calls sfsubsearch.pl
  sfsectionsearch.html -- form for fixed section queries, calls sfsubsearch.pl
  
  
-- in sfsubsearch.pl, variables defined at the top will need to be localized
   for each server.  These include

###  Server specific vars
$SRSROOT = "/b4/srs";  ## path to SRS software
$seqsearchpath  = 'nice /b4/srsfasta/bin/';  # [t]fasta app is appended to path
$SRSserver  = "http://iubio.bio.indiana.edu:81/srs/srsc";


-- in sfsectionlib.pl, the variables defined at the top will need to be localized 
   for each server.  These include 
   
###  Server specific vars
$SRSROOT = "/b4/srs";  ## path to SRS software

# SRS Databanks to use
@nadatabanks=('genbank','gbnew', 'gbest');  
@aadatabanks=('pir','swissprot', 'swissnew'); 

# SRS queries to generate predefined subset libraries
# format is 'data-field=srs-query=index-file-name=description'
push(@queries, 'org=arabidopsis=arabidopsis');
push(@queries, 'org=caenorhabditis=caenorhabditis');
...

  
*) FastA

 The main change to FastA software is the addition of sequence library
 type code that lets FastA use the index files (described in SRS section)
 to read only specified sequence records from a set of data library files.
 The files changed include small changes to FFFASTA.C, NXGETAA.C, ALTLIB.H,
 and addition of a new source file SUBGETAA.C with the bulk of code to
 read sequences from the subset indices.  This file format may be generated
 by software other than SRS fairly easily, so isn't dependent on SRS software.
 
 For using with Fasta, this subset library type is defined as #12, so when
 Fasta/Tfasta ask for a library, you provide the name of the subset indices file
 (without suffix), followed by the library type number (12), as in
 
    test sequence file name: my.seq
     library file name: []srslibs  12
 
 where the files srslibs.names and srslibs.indices were generated with 
 an SRS command like 'getz [databank-section:query] -fse srslibs' 
 
 
*) SRS

 The main change to SRS software is to add an option to write a file of
 indices for the sequences that match a given query.  This amounts to
 adding a new command argument (ListFSE) to GETZ &/or WGETZ (I use it from Getz),
 and a new function in the shared SEQLIB.C file called "SlbListFSE".
 The ListFSE function will produce two files, "queryname.indices" and
 "queryname.names".  
 
 The version of SRS that I am working from, srs4_02,
 is not the current one, so you may need to do some judicious cut & paste
 to get these patches into the current release.  I've made various other
 changes to my working version of SRS software, so the difference files
 supplied show many changes not relevant to this modification.
 
 
 
-- File structure
The two files .names and .indices are created by calling the new
ListFSE function in getz or wgetz (using the "-fse dataname" command
argument).  This is their current structure.

The .names file is a list of data file names, in this form
full-pathname        file-format   file-number
devc:srs:data:pir:pir1.dat 2 0
devc:srs:data:pir:pir2.dat 2 1

The .indices file has an entry for each sequence matching the SRS query,
in this format:
file-number   start-of-document  start-of-sequence  end-of-sequence
0 33946 34851 35946
0 35946 36851 1104
1 12972 14168 15074



-- The changes to GETZ.C or WGETZ.C are small, adding merely a new
 function call:

src/getz.c: 
 else if (ParGetNum ("printText") || ParGetNum ("printData") ||
     *ParGetStr ("fieldList")
     || *ParGetStr("listFSE") /* dgg <<< */
     ) {
        /*  ParDefNum ("listFOSN", 0); nicht mehr notwendig, neues seqlib */
     EntryPrintSet (setName, 1, 0);
  }

src/wgetz.c:

  /* dgg --- add for subset indices */
  else if ((tmp = ParGetStr ("listFSE")) && *tmp) { 
    if (QryDo (query, "Q")) { 
   ENTRYo   *entry;
   SETo     *set;
   char  *setName;
   IDoENTRY id;
   INT4     setEntryN, k;
   
   setName= "Q";
   set = SetGet (setName);
   setEntryN = SetSize (setName);
   for (k=1;  k <= setEntryN;  k++) {
    SetGetID (set, k, &id); 
    entry = EntryOpen (&id);
    if (entry) SlbListFSE( entry,  setName,  setEntryN,  k);  
    EntryClose (&entry);
    }
   }
  }

 
-- Add items to the odd/parameters.sdl and odd/arglist.sdl to
 allow this new function to be called in getz and wgetz.  Also
 the SRSTEMP folder variable is defined to specify where the
 new files *.names and *.indices will be created and stored.
 

odd/parameter.sdl:
#partable 
... 
! dgg addition
    #parameter /name="listFSE" /type=string /str=""
         /comment="produce a File of Sequence Entries"
    #parameter /name=tempDirName /type=string /str="SRSTEMP:"  
         /comment="name of temporary files directory"  ! dgg

odd/arglist.sdl:
 
#command /name=getz
...
    #arg /name="-fse"  /parameter="listFSE"     ! dgg added
    #arg /name="-tdir" /parameter=tempDirName     ! dgg added
 
 
#command /name=wgetz
...
    #arg /name="-fse"  /parameter="listFSE"     ! dgg added
    #arg /name="-tdir" /parameter=tempDirName     ! dgg added
 
 
 
-- each sequence data ODD file that defines a #filetype must add
 an entry "/seqtype=genbank" for instance.  I did this for
 EMBL/SwissProt which share a #filetype, GenBank and PIR in this way.
 Please note that I use the CODATA format of PIR, while others may use
 the NBRF formatted PIR dataset.

odd/embl.sdl:
#filetype /ID=%EMBLSEQ_FILE 
          /begstr="SQ   " /advance=1 
          /seqtype=embl ! SLBxEMBL, dgg addition
          /typename="dat" /maxline=90 /exitstr="//"

odd/genbank.sdl:
#filetype /ID=%GBSEQ_FILE 
          /begstr="ORIGIN" /advance=1
          /seqtype=genbank ! SLBxGENBANK, dgg addition
          /typename="seq" /maxline=200 /exitstr="//"

odd/pir.sdl:
#filetype /ID=%PIRASEQ_FILE 
          /typename="dat" /maxline=200 
          /seqtype=codata ! SLBxCODATA, dgg addition
          /begstr="SEQUENCE" /advance=1 /exitstr="///"
 
 
-- in srsgeneral-def.sdl, the allowable seqtypes are defined in the
  %%O_FILTYP file type definition.  The seqtypes must
 match in with the values (+10) used by FastA for filetypes.  This file
 type value is written into the 'sequence.names' file in the ListFSE function.

odd/srsgeneral_def.sdl:

#object /id=%%O_FILTYP /ppname=SRS
...
    #attribute /name=seqtype /declname="type" /type=uchar
        /valtype=name /defaultval=0
        /rem="Type of sequence file.  Distinguish base file format for various
              uses." !dgg modified
        !/rem="Type of sequence file. This was only important on VMS and with
        !      GCG. GCG version 8.0 does not have a record separator as with 
        !      version 7.0."
    {
       #value /name=normal /val=0
        #value /name=gcg7 /val=1
        #value /name=gcg8 /val=2
        ! dgg additions for fasta lib reader format values, val-10 for pearson value
        #value /name=pearson /val=10
        #value /name=genbank /val=11
        #value /name=codata  /val=12
        #value /name=embl    /val=13
        #value /name=ig      /val=14
        #value /name=nbrf    /val=15
    }
 

-- the ListFSE routine added to srs/src/seqlib.c 

 /* --- 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];
  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);
   _FilLN(outFile);
      file = fopen( outFile, "r+");
      if (!file) file = fopen( outFile, "w");
      else fseek(file, 0, 2);

   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 *********************************************************

: :

*/

static INT4 gSeqFormat = -1;

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

 if ( *ParGetStr("listFSE")) {
  /* a self-contained function - can be called from other places - dgg */
  SlbListFSE( entry,  setName,  entryN,  entryCurrN);
  }
  
.
