source: trunk/READSEQ/readseq.c

Last change on this file was 19452, checked in by westram, 19 months ago
  • get rid of tmpnam().
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 35.9 KB
Line 
1/* File: readseq.c
2 * main() program for ureadseq.c, ureadseq.h
3 *
4 * Reads and writes nucleic/protein sequence in various
5 * formats. Data files may have multiple sequences.
6 *
7 * Copyright 1990 by d.g.gilbert
8 * biology dept., indiana university, bloomington, in 47405
9 * e-mail: gilbertd@bio.indiana.edu
10 *
11 * This program may be freely copied and used by anyone.
12 * Developers are encourged to incorporate parts in their
13 * programs, rather than devise their own private sequence
14 * format.
15 *
16 * This should compile and run with any ANSI C compiler.
17 * Please advise me of any bugs, additions or corrections.
18 *
19 */
20
21const char *title
22    = "readSeq (1Feb93), multi-format molbio sequence reader.\n";
23
24 /*  History
25  27 Feb 90.  1st release to public.
26   4 Mar 90.  + Gary Olsen format
27              + case change
28              * minor corrections to NBRF,EMBL,others
29              * output 1 file per sequence for gcg, unknown
30              * define -DNOSTR for c-libraries w/o strstr
31              - readseq.p, pascal version, becomes out-of-date
32  24 May 90.  + Phylip 3.2 output format (no input)
33  20 Jul 90.  + Phylip 3.3 output (no input yet)
34              + interactive output re-direction
35              + verbose progress info
36              * interactive help output
37              * dropped line no.s on NBRF output
38              * patched in HyperGCG XCMD corrections,
39                - except for seq. documentation handling
40              * dropped the IG special nuc codes, as IG has
41                adopted the standard IUB codes (now if only
42                everyone would adopt a standard format !)
43  11 Oct 90.  * corrected bug in reading/writing of EMBL format
44
45  17 Oct 91.  * corrected bug in reading Olsen format
46                (serious-deletion)
47  10 Nov 91.  * corrected bug in reading some GCG format files
48                (serious-last line duplicated)
49              + add format name parsing (-fgb, -ffasta, ...)
50              + Phylip v3.4 output format (== v3.2, sequential)
51              + add checksum output to all forms that have document
52              + skip mail headers in seq file
53              + add pipe for standard input == seq file (with -p)
54              * fold in parts of MacApp Seq object
55              * strengthen format detection
56              * clarify program structure
57              * remove fixed sequence size limit (now dynamic, sizeof memory)
58              * check and fold in accumulated bug reports:
59              *   Now ANSI-C fopen(..,"w") & check open failure
60              *   Define -DFIXTOUPPER for nonANSI C libraries that mess
61                  up toupper/tolower
62              = No command-line changes; callers of readseq main() should be okay
63              - ureadseq.h functions have changed; client programs need to note.
64              + added Unix and VMS Make scripts, including validation tests
65
66   4 May 92.  + added 32 bit CRC checksum as alternative to GCG 6.5bit checksum
67                (-DBIGCHECKSUM)
68    Aug 92    = fixed Olsen format input to handle files w/ more sequences,
69                not to mess up when more than one seq has same identifier,
70                and to convert number masks to symbols.
71              = IG format fix to understand ^L
72
73  25-30 Dec 92
74              * revised command-line & interactive interface.  Suggested form is now
75                  readseq infile -format=genbank -output=outfile -item=1,3,4 ...
76                but remains compatible with prior commandlines:
77                  readseq infile -f2 -ooutfile -i3 ...
78              + added GCG MSF multi sequence file format
79              + added PIR/CODATA format
80              + added NCBI ASN.1 sequence file format
81              + added Pretty, multi sequence pretty output (only)
82              + added PAUP multi seq format
83              + added degap option
84              + added Gary Williams (GWW, G.Williams@CRC.AC.UK) reverse-complement option.
85              + added support for reading Phylip formats (interleave & sequential)
86              * string fixes, dropped need for compiler flags NOSTR, FIXTOUPPER, NEEDSTRCASECMP
87              * changed 32bit checksum to default, -DSMALLCHECKSUM for GCG version
88
89   1Feb93
90              = revert GenBank output to a fixed left number width which
91               other software depends on.
92              = fix for MSF input to handle symbols in names
93              = fix bug for possible memory overrun when truncating seqs for
94                Phylip or Paup formats (thanks Anthony Persechini)
95
96 */
97
98
99
100/*
101   Readseq has been tested with:
102      Macintosh MPW C
103      GNU gcc
104      SGI cc
105      VAX-VMS cc
106   Any ANSI C compiler should be able to handle this.
107   Old-style C compilers barf all over the source.
108
109
110How do I build the readseq program if I have an Ansi C compiler?
111#--------------------
112# Unix ANSI C
113# Use the supplied Makefile this way:
114%  make CC=name-of-c-compiler
115# OR do this...
116% gcc readseq.c ureadseq.c -o readseq
117
118#--------------------
119$!VAX-VMS cc
120$! Use the supplied Make.Com this way:
121$  @make
122$! OR, do this:
123$ cc readseq, ureadseq
124$ link readseq, ureadseq, sys$library:vaxcrtl/lib
125$ readseq :== $ MyDisk:[myacct]readseq
126
127#--------------------
128# Macintosh Simple Input/Output Window application
129# requires MPW-C and SIOW library (from APDA)
130# also uses files macinit.c, macinit.r, readseqSIOW.make
131#
132Buildprogram readseqSIOW
133
134#--------------------
135#MPW-C v3 tool
136C  ureadseq.c
137C  readseq.c
138link -w -o readseq -t MPST -c 'MPS ' ¶
139   readseq.c.o Ureadseq.c.o ¶
140    "{Libraries}"Interface.o ¶
141    "{Libraries}"ToolLibs.o ¶
142    "{Libraries}"Runtime.o ¶
143    "{CLibraries}"StdClib.o
144readseq -i1 ig.seq
145
146# MPW-C with NCBI tools
147
148set NCBI "{Boot}@molbio:ncbi:"; EXPORT NCBI
149set NCBILIB1  "{NCBI}"lib:libncbi.o; export NCBILIB1
150set NCBILIB2  "{NCBI}"lib:libncbiobj.o; export NCBILIB2
151set NCBILIB3  "{NCBI}"lib:libncbicdr.o; export NCBILIB3
152set NCBILIB4  "{NCBI}"lib:libvibrant.o; export NCBILIB4
153
154C  ureadseq.c
155C  -d NCBI -i "{NCBI}"include: ureadasn.c
156C  -d NCBI -i "{NCBI}"include: readseq.c
157link -w -o readseq -t MPST -c 'MPS ' ¶
158   ureadseq.c.o ureadasn.c.o readseq.c.o  ¶
159    {NCBILIB4} {NCBILIB2} {NCBILIB1} ¶
160    "{Libraries}"Interface.o ¶
161    "{Libraries}"ToolLibs.o ¶
162    "{Libraries}"Runtime.o ¶
163    "{CLibraries}"CSANELib.o ¶
164    "{CLibraries}"Math.o ¶
165    "{CLibraries}"StdClib.o
166
167===========================================================*/
168
169
170
171#include <stdio.h>
172#include <stdlib.h>
173#include <string.h>
174#define __NO_CTYPE
175#include <ctype.h>
176
177#include "ureadseq.h"
178#include <gets_noOverflow.h>
179
180
181static char inputfilestore[256], *inputfile = inputfilestore;
182
183const char *formats[kMaxFormat+1] = {
184    " 1. IG/Stanford",
185    " 2. GenBank/GB",
186    " 3. NBRF",
187    " 4. EMBL",
188    " 5. GCG",
189    " 6. DNAStrider",
190    " 7. Fitch",
191    " 8. Pearson/Fasta",
192    " 9. Zuker (in-only)",
193    "10. Olsen (in-only)",
194    "11. Phylip3.2",
195    "12. Phylip",
196    "13. Plain/Raw",
197    "14. PIR/CODATA",
198    "15. MSF",
199    "16. ASN.1",
200    "17. PAUP/NEXUS",
201    "18. Pretty (out-only)",
202    "" };
203
204#define kFormCount  30
205#define kMaxFormName 15
206
207const  struct formatTable {
208    const char *name;
209    short num;
210} formname[] = {
211    {"ig",  kIG},
212    {"stanford", kIG},
213    {"genbank", kGenBank},
214    {"gb", kGenBank},
215    {"nbrf", kNBRF},
216    {"embl", kEMBL},
217    {"gcg", kGCG},
218    {"uwgcg", kGCG},
219    {"dnastrider", kStrider},
220    {"strider", kStrider},
221    {"fitch", kFitch},
222    {"pearson", kPearson},
223    {"fasta", kPearson},
224    {"zuker", kZuker},
225    {"olsen", kOlsen},
226    {"phylip", kPhylip},
227    {"phylip3.2", kPhylip2},
228    {"phylip3.3", kPhylip3},
229    {"phylip3.4", kPhylip4},
230    {"phylip-interleaved", kPhylip4},
231    {"phylip-sequential", kPhylip2},
232    {"plain", kPlain},
233    {"raw", kPlain},
234    {"pir", kPIR},
235    {"codata", kPIR},
236    {"asn.1", kASN1},
237    {"msf", kMSF},
238    {"paup", kPAUP},
239    {"nexus", kPAUP},
240    {"pretty", kPretty},
241  };
242
243const char *kASN1headline = "Bioseq-set ::= {\nseq-set {\n";
244
245/* GWW table for getting the complement of a nucleotide (IUB codes) */
246/*                     ! "#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[ \]^_`abcdefghijklmnopqrstuvwxyz{|}~ */
247const char compl[] = " !\"#$%&'()*+,-./0123456789:;<=>?@TVGHNNCDNNMNKNNYRYSAABWNRN[\\]^_`tvghnncdnnmnknnyrysaabwnrn{|}~";
248
249
250
251const char *formatstr( short format)
252{
253  if (format < 1 || format > kMaxFormat) {
254    switch (format) {
255      case kASNseqentry :
256      case kASNseqset   : return formats[kASN1-1];
257      case kPhylipInterleave:
258      case kPhylipSequential: return formats[kPhylip-1];
259      default: return "(unknown)";
260      }
261    }
262  else return formats[format-1];
263}
264
265int rs_isdigit(int c){
266        return isdigit(c);
267}
268
269int parseformat( char *name2)
270{
271#define kDupmatch  -2
272  int   maxlen, i, match, matchat;
273  size_t namelen;
274  char  lname[kMaxFormName+1];
275
276  skipwhitespace(name2);
277  namelen = strlen(name2);
278  if (namelen == 0)
279    return kNoformat;
280  else if (rs_isdigit(*name2)) {
281    i = atol( name2);
282    if (i < kMinFormat || i > kMaxFormat) return kNoformat;
283    else return i;
284    }
285
286  /* else match character name */
287  maxlen = min( kMaxFormName, namelen);
288  for (i=0; i<maxlen; i++) lname[i] = to_lower(name2[i]);
289  lname[maxlen]=0;
290  matchat = kNoformat;
291
292  for (i=0; i<kFormCount; i++) {
293    match = strncmp( lname, formname[i].name, maxlen);
294    if (match == 0) {
295      if (strlen(formname[i].name) == namelen) return (formname[i].num);
296      else if (matchat == kNoformat) matchat = i;
297      else matchat = kDupmatch; /* 2 or more partial matches */
298      }
299    }
300  if (matchat == kNoformat || matchat == kDupmatch)
301    return kNoformat;
302  else
303    return formname[matchat].num;
304}
305
306
307
308static void dumpSeqList(char *list, short format)
309{
310  long i, l, listlen;
311  char s[256];
312
313  listlen = strlen(list);
314  printf("Sequences in %s  (format is %s)\n", inputfile, formatstr(format));
315  for (i=0, l=0; i < listlen; i++) {
316    if (list[i] == (char)NEWLINE) {
317      s[l] = '\0'; l = 0;
318      puts(s);
319      }
320    else if (l < 255)
321      s[l++] = list[i];
322    }
323  putchar('\n');
324}
325
326
327
328void usage()
329{
330  short   i, midi;
331
332  fputs(title, stderr);
333  fprintf(stderr,
334  "usage: readseq [-options] in.seq > out.seq\n");
335  fprintf(stderr," options\n");
336/* ? add -d[igits] to allow digits in sequence data, &/or option to specify seq charset !? */
337  fprintf(stderr, "    -a[ll]         select All sequences\n");
338  fprintf(stderr, "    -c[aselower]   change to lower case\n");
339  fprintf(stderr, "    -C[ASEUPPER]   change to UPPER CASE\n");
340  fprintf(stderr, "    -degap[=-]     remove gap symbols\n");
341  fprintf(stderr, "    -i[tem=2,3,4]  select Item number(s) from several\n");
342  fprintf(stderr, "    -l[ist]        List sequences only\n");
343  fprintf(stderr, "    -o[utput=]out.seq  redirect Output\n");
344  fprintf(stderr, "    -p[ipe]        Pipe (command line, <stdin, >stdout)\n");
345  fprintf(stderr, "    -r[everse]     change to Reverse-complement\n");
346  fprintf(stderr, "    -v[erbose]     Verbose progress\n");
347  fprintf(stderr, "    -f[ormat=]#    Format number for output,  or\n");
348  fprintf(stderr, "    -f[ormat=]Name Format name for output:\n");
349  midi = (kMaxFormat+1) / 2;
350  for (i = kMinFormat-1; i < midi; i++)
351   fprintf( stderr, "        %-20s      %-20s\n",
352    formats[i], formats[midi+i]);
353
354  /* new output format options, esp. for pretty format: */
355  fprintf(stderr, "     \n");
356  fprintf(stderr, "   Pretty format options: \n");
357  fprintf(stderr, "    -wid[th]=#            sequence line width\n");
358  fprintf(stderr, "    -tab=#                left indent\n");
359  fprintf(stderr, "    -col[space]=#         column space within sequence line on output\n");
360  fprintf(stderr, "    -gap[count]           count gap chars in sequence numbers\n");
361  fprintf(stderr, "    -nameleft, -nameright[=#]   name on left/right side [=max width]\n");
362  fprintf(stderr, "    -nametop              name at top/bottom\n");
363  fprintf(stderr, "    -numleft, -numright   seq index on left/right side\n");
364  fprintf(stderr, "    -numtop, -numbot      index on top/bottom\n");
365  fprintf(stderr, "    -match[=.]            use match base for 2..n species\n");
366  fprintf(stderr, "    -inter[line=#]        blank line(s) between sequence blocks\n");
367
368  /******  not ready yet
369  fprintf(stderr, "    -code=none,rtf,postscript,ps   code syntax\n");
370  fprintf(stderr, "    -namefont=, -numfont=, -seqfont=font   font choice\n");
371  fprintf(stderr, "       font suggestions include times,courier,helvetica\n");
372  fprintf(stderr, "    -namefontsize=, -numfontsize=, -seqfontsize=#\n");
373  fprintf(stderr, "       fontsize suggestions include 9,10,12,14\n");
374  fprintf(stderr, "    -namefontstyle=, -numfontstyle=, -seqfontstyle= style  fontstyle for names\n");
375  fprintf(stderr, "       fontstyle options are plain,italic,bold,bold-italic\n");
376  ******/
377}
378
379void erralert(short err)
380{
381  switch (err) {
382    case 0  :
383      break;
384    case eFileNotFound: fprintf(stderr, "arb_readseq: File not found: %s\n", inputfile);
385      break;
386    case eFileCreate: fprintf(stderr, "arb_readseq: Can't open output file.\n");
387      break;
388    case eASNerr: fprintf(stderr, "arb_readseq: Error in ASN.1 sequence routines.\n");
389      break;
390    case eNoData: fprintf(stderr, "arb_readseq: No data in file.\n");
391      break;
392    case eItemNotFound: fprintf(stderr, "arb_readseq: Specified item not in file.\n");
393      break;
394    case eUnequalSize:  fprintf(stderr, "arb_readseq: This format requires equal length sequences.\nSequence truncated or padded to fit.\n");
395      break;
396    case eUnknownFormat: fprintf(stderr, "arb_readseq: Error: this format is unknown to me.\n");
397      break;
398    case eOneFormat: fprintf(stderr, "arb_readseq: Warning: This format permits only 1 sequence per file.\n");
399      break;
400    case eMemFull: fprintf(stderr, "arb_readseq: Out of storage memory. Sequence truncated.\n");
401      break;
402    case ePipeStdin: fprintf(stderr, "arb_readseq: piping from stdin is prohibited.\n");
403      break;
404    default: fprintf(stderr, "arb_readseq: errorcode = %d\n", err);
405      break;
406  }
407} /* erralert */
408
409#define USERINPUT_BUFFERSIZE     128
410#define USERINPUT_BUFFERSIZE_SEQ 256
411
412int chooseFormat( boolean quietly)
413{
414  char  sform[USERINPUT_BUFFERSIZE];
415  int   midi, i, outform;
416
417    if (quietly)
418      return kPearson;  /* default */
419    else {
420      midi = (kMaxFormat+1) / 2;
421      for (i = kMinFormat-1; i < midi; i++)
422        fprintf( stderr, "        %-20s      %-20s\n",
423                        formats[i], formats[midi+i]);
424      fprintf(stderr,"\nChoose an output format (name or #): \n");
425      gets_noOverflow(sform, USERINPUT_BUFFERSIZE);
426      outform = parseformat(sform);
427      if (outform == kNoformat) outform = kPearson;
428      return outform;
429      }
430}
431
432
433
434/* read parameter(s) */
435
436boolean checkopt( boolean casesense, char *sopt, const char *smatch, short minword)
437{
438  long  lenopt, lenmatch;
439  boolean result;
440  short minmaxw;
441
442  lenopt = strlen(sopt);
443  lenmatch= strlen(smatch);
444  minmaxw= max(minword, min(lenopt, lenmatch));
445
446  if (casesense)
447    result= (!strncmp( sopt, smatch, minmaxw));
448  else
449    result= (!Strncasecmp( sopt, smatch, minmaxw ));
450  /* if (result) { */
451    /* fprintf(stderr,"true checkopt(opt=%s,match=%s,param=%s)\n", sopt, smatch, *sparam); */
452  /*  } */
453  return result;
454}
455
456
457#define   kMaxwhichlist  50
458
459/* global for readopt(), main() */
460boolean   chooseall = false, quietly = false, gotinputfile = false,
461          listonly = false, closeout = false, verbose = false,
462          manyout = false, dolower = false, doupper = false, doreverse= false,
463          askout  = true, dopipe= false, interleaved = false;
464short     nfile = 0, iwhichlist=0, nwhichlist = 0;
465short     whichlist[kMaxwhichlist+1];
466long      whichSeq = 0, outform = kNoformat;
467char      onamestore[USERINPUT_BUFFERSIZE], *oname = onamestore;
468FILE      *foo = NULL;
469
470void resetGlobals()
471/* need this when used from SIOW, as these globals are not reinited automatically
472between calls to local main() */
473{
474  chooseall = false; quietly = false; gotinputfile = false;
475  listonly = false; closeout = false; verbose = false;
476  manyout = false; dolower = false; doupper = false; doreverse= false;
477  askout  = true; dopipe= false; interleaved = false;
478  nfile = 0; iwhichlist=0; nwhichlist = 0;
479  whichSeq = 0; outform = kNoformat;
480  oname = onamestore;
481  foo = NULL;
482
483  gPrettyInit(gPretty);
484}
485
486
487#define kOptOkay  1
488#define kOptNone  0
489
490int readopt( char *sopt)
491{
492  char    sparamstore[256], *sparam= sparamstore;
493  short   n;
494
495  /* fprintf(stderr,"readopt( %s) == ", sopt); */
496
497  if (*sopt == '?') {
498    usage();
499    return kOptNone;   /*? eOptionBad or kOptNone */
500    }
501
502  else if (*sopt == '-') {
503
504      {
505          char *cp= strchr(sopt,'=');
506          *sparam= '\0';
507          if (cp) {
508              strcpy(sparam, cp+1);
509              *cp= 0;
510          }
511      }
512
513    if (checkopt( false, sopt, "-help", 2)) {
514      usage();
515      return kOptNone;
516      }
517
518    if (checkopt( false, sopt, "-all", 2)) {
519      whichSeq= 1; chooseall= true;
520      return kOptOkay;
521      }
522
523    if (checkopt( false, sopt, "-colspace", 4)) { /* test before -c[ase] */
524      n= atoi( sparam);
525      gPretty.spacer = n;
526      return kOptOkay;
527      }
528
529    if (checkopt( true, sopt, "-caselower", 2)) {
530      dolower= true;
531      return kOptOkay;
532      }
533    if (checkopt( true, sopt, "-CASEUPPER", 2)) {
534      doupper= true;
535      return kOptOkay;
536      }
537
538    if (checkopt( false, sopt, "-pipe", 2)) {
539      dopipe= true; askout= false;
540      return kOptOkay;
541      }
542
543    if (checkopt( false, sopt, "-list", 2)) {
544      listonly = true; askout = false;
545      return kOptOkay;
546      }
547
548    if (checkopt( false, sopt, "-reverse", 2)) {
549      doreverse = true;
550      return kOptOkay;
551      }
552
553    if (checkopt( false, sopt, "-verbose", 2)) {
554      verbose = true;
555      return kOptOkay;
556      }
557
558    if (checkopt( false, sopt, "-match", 5)) {
559      gPretty.domatch= true;
560      if (*sparam >= ' ') gPretty.matchchar= *sparam;
561      return kOptOkay;
562      }
563    if (checkopt( false, sopt, "-degap", 4)) {
564      gPretty.degap= true;
565      if (*sparam >= ' ') gPretty.gapchar= *sparam;
566      return kOptOkay;
567      }
568
569    if (checkopt( false, sopt, "-interline", 4)) {
570      gPretty.interline= atoi( sparam);
571      return kOptOkay;
572      }
573
574    if (checkopt( false, sopt, "-item", 2)) {
575      char  *cp = sparam;
576      nwhichlist= 0;
577      whichlist[0]= 0;
578      if (*cp == 0) cp= sopt+2; /* compatible w/ old way */
579      do {
580        while (*cp!=0 && !rs_isdigit(*cp)) cp++;
581        if (*cp!=0) {
582          n = atoi( cp);
583          whichlist[nwhichlist++]= n;
584          while (*cp!=0 && rs_isdigit(*cp)) cp++;
585          }
586      } while (*cp!=0 && n>0 && nwhichlist<kMaxwhichlist);
587      whichlist[nwhichlist++]= 0; /* 0 == stopsign for loop */
588      whichSeq= max(1,whichlist[0]); iwhichlist= 1;
589      return kOptOkay;
590      }
591
592    if (checkopt( false, sopt, "-format", 5)) {/* -format=phylip, -f2, -form=phylip */
593      if (*sparam==0) { for (sparam= sopt+2; isalpha(*sparam); sparam++) ; }
594      outform = parseformat( sparam);
595      return kOptOkay;
596      }
597    if (checkopt( false, sopt, "-f", 2)) { /* compatible w/ -fphylip prior version */
598      if (*sparam==0) sparam= sopt+2;
599      outform = parseformat( sparam);
600      return kOptOkay;
601      }
602
603    if (checkopt( false, sopt, "-output", 3)) {/* -output=myseq */
604      if (*sparam==0) { for (sparam= sopt+3; isalpha(*sparam); sparam++) ; }
605      strcpy( oname, sparam);
606      foo = fopen( oname, "w");
607      if (!foo) { erralert(eFileCreate); return eFileCreate; }
608      closeout = true;
609      askout = false;
610      return kOptOkay;
611      }
612    if (checkopt( false, sopt, "-o", 2)) {  /* compatible w/ -omyseq prior version */
613      if (*sparam==0) sparam= sopt+2;
614      strcpy( oname, sparam);
615      foo = fopen( oname, "w");
616      if (!foo) { erralert(eFileCreate); return eFileCreate; }
617      closeout = true;
618      askout = false;
619      return kOptOkay;
620      }
621
622    if (checkopt( false, sopt, "-width", 2)) {
623      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
624      n= atoi( sparam);
625      if (n>0) gPretty.seqwidth = n;
626      return kOptOkay;
627      }
628
629    if (checkopt( false, sopt, "-tab", 4)) {
630      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
631      n= atoi( sparam);
632      gPretty.tab = n;
633      return kOptOkay;
634      }
635
636    if (checkopt( false, sopt, "-gapcount", 4)) {
637      gPretty.baseonlynum = false;
638      /* if (*sparam >= ' ') gPretty.gapchar= *sparam; */
639      return kOptOkay;
640      }
641    if (checkopt( false, sopt, "-nointerleave", 8)) {
642      gPretty.noleaves = true;
643      return kOptOkay;
644      }
645
646    if (checkopt( false, sopt, "-nameleft", 7)) {
647      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
648      n= atoi( sparam);
649      if (n>0 && n<50) gPretty.namewidth =  n;
650      gPretty.nameleft= true;
651      return kOptOkay;
652      }
653    if (checkopt( false, sopt, "-nameright", 7)) {
654      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
655      n= atoi( sparam);
656      if (n>0 && n<50) gPretty.namewidth =  n;
657      gPretty.nameright= true;
658      return kOptOkay;
659      }
660    if (checkopt( false, sopt, "-nametop", 6)) {
661      gPretty.nametop= true;
662      return kOptOkay;
663      }
664
665    if (checkopt( false, sopt, "-numleft", 6)) {
666      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
667      n= atoi( sparam);
668      if (n>0 && n<50) gPretty.numwidth =  n;
669      gPretty.numleft= true;
670      return kOptOkay;
671      }
672    if (checkopt( false, sopt, "-numright", 6)) {
673      if (*sparam==0) { for (sparam= sopt+2; !rs_isdigit(*sparam) && *sparam!=0; sparam++) ; }
674      n= atoi( sparam);
675      if (n>0 && n<50) gPretty.numwidth =  n;
676      gPretty.numright= true;
677      return kOptOkay;
678      }
679
680    if (checkopt( false, sopt, "-numtop", 6)) {
681      gPretty.numtop= true;
682      return kOptOkay;
683      }
684    if (checkopt( false, sopt, "-numbottom", 6)) {
685      gPretty.numbot= true;
686      return kOptOkay;
687      }
688
689    else {
690      usage();
691      return eOptionBad;
692      }
693    }
694
695  else {
696    strcpy( inputfile, sopt);
697    gotinputfile = (*inputfile != 0);
698    nfile++;
699    return kOptOkay;
700    }
701
702 /* return kOptNone; -- never here */
703}
704
705
706
707
708/* this program suffers some as it tries to be a quiet translator pipe
709   _and_ a noisy user interactor
710*/
711
712/* return is best for SIOW, okay for others */
713#ifdef SIOW
714#define Exit(a)   return(a)
715siow_main( int argc, char *argv[])
716
717#else
718#define Exit(a)   exit(a)
719
720int main( int argc, char *argv[])
721#endif
722{
723boolean   closein = false;
724short     ifile, nseq, atseq, format, err = 0, seqtype = kDNA,
725          nlines, seqout = 0, phylvers = 2;
726long      i, skiplines, seqlen, seqlen0;
727unsigned long  checksum= 0, checkall= 0;
728char      *seq, *firstseq = NULL, *seqlist, tempname[256];
729const char *cp;
730char      seqid[256], *seqidptr = seqid;
731char      stempstore[USERINPUT_BUFFERSIZE_SEQ], *stemp = stempstore;
732FILE     *ftmp, *fin = NULL, *fout = NULL;
733long      outindexmax= 0, noutindex= 0, *outindex = NULL;
734
735#define exit_main(err) {        \
736  if (closeout) fclose(fout);   \
737  if (closein) fclose(fin);   \
738  if (*tempname!=0) remove(tempname);\
739  Exit(err); }
740
741#define indexout()  if (interleaved) {\
742  if (noutindex>=outindexmax) {\
743    outindexmax= noutindex + 20;\
744    outindex= (long*) realloc(outindex, sizeof(long)*outindexmax);\
745    if (outindex==NULL) { err= eMemFull; erralert(err); exit_main(err); }\
746    }\
747  outindex[noutindex++]= ftell(fout);\
748  }
749
750
751  resetGlobals();
752  foo = stdout;
753  *oname = 0;
754  *tempname = 0;
755  /* initialize gPretty ?? -- done in header */
756
757  for (i=1; i < argc; i++) {
758    err= readopt( argv[i]);
759    if (err <= 0) exit_main(err);
760    }
761
762                            /* pipe input from stdin !? */
763  if (dopipe && !gotinputfile) {
764      erralert(ePipeStdin);
765      exit_main(ePipeStdin);
766#if 0
767      // Disable to get rid of tmpnam warning.
768      // Case is not used by ARB.
769      int c;
770      tmpnam(tempname);
771      inputfile = tempname;
772      ftmp = fopen( inputfile, "w");
773      if (!ftmp) { erralert(eFileCreate); exit_main(eFileCreate); }
774      while ((c = getc(stdin)) != EOF) fputc(c, ftmp);
775      fclose(ftmp);
776      gotinputfile= true;
777#endif
778  }
779
780  quietly = (dopipe || (gotinputfile && (listonly || whichSeq != 0)));
781
782  if (verbose || (!quietly && !gotinputfile)) fputs(title, stderr);
783  ifile = 1;
784
785                            /* UI: Choose output */
786  if (askout && !closeout && !quietly) {
787    askout = false;
788    fprintf(stderr,"\nName of output file (?=help, defaults to display): \n");
789    gets_noOverflow(oname= onamestore, USERINPUT_BUFFERSIZE);
790    skipwhitespace(oname);
791    if (*oname == '?') { usage(); exit_main(0); }
792    else if (*oname != 0) {
793      closeout = true;
794      foo = fopen( oname, "w");
795      if (!foo) { erralert(eFileCreate); exit_main(eFileCreate); }
796      }
797    }
798
799  fout = foo;
800  if (outform == kNoformat) outform = chooseFormat(quietly);
801
802                          /* set up formats ... */
803  switch (outform) {
804    case kPhylip2:
805      interleaved= false;
806      phylvers = 2;
807      outform = kPhylip;
808      break;
809
810    case kPhylip4:
811      interleaved= true;
812      phylvers = 4;
813      outform = kPhylip;
814      break;
815
816    case kMSF:
817    case kPAUP:
818      interleaved= true;
819      break;
820
821    case kPretty:
822      gPretty.isactive= true;
823      interleaved= true;
824      break;
825
826    }
827
828  if (gPretty.isactive && gPretty.noleaves) interleaved= false;
829  if (interleaved) {
830    fout = ftmp = tmpfile();
831    outindexmax= 30; noutindex= 0;
832    outindex = (long*) malloc(outindexmax*sizeof(long));
833    if (outindex==NULL) { err= eMemFull; erralert(err); exit_main(err); }
834    }
835
836                        /* big loop over all input files */
837  do {
838                        /* select next input file */
839    gotinputfile = (*tempname != 0);
840    while ((ifile < argc) && (!gotinputfile)) {
841      if (*argv[ifile] != '-') {
842        strcpy( inputfile, argv[ifile]);
843        gotinputfile = (*inputfile != 0);
844        --nfile;
845        }
846      ifile++;
847      }
848
849    while (!gotinputfile) {
850      fprintf(stderr,"\nName an input sequence or -option: \n");
851      inputfile= inputfilestore;
852
853      gets_noOverflow(stemp= stempstore, USERINPUT_BUFFERSIZE_SEQ);
854      if (*stemp==0) goto fini;  /* !! need this to finish work during interactive use */
855      stemp= strtok(stempstore, " \n\r\t");
856      while (stemp) {
857        err= readopt( stemp); /* will read inputfile if it exists */
858        if (err<0) exit_main(err);
859        stemp= strtok( NULL, " \n\r\t");
860        }
861      }
862              /* thanks to AJB@UK.AC.DARESBURY.DLVH for this PHYLIP3 fix: */
863              /* head for end (interleave if needed) */
864    if (*inputfile == 0) break;
865
866    format = seqFileFormat( inputfile, &skiplines, &err);
867
868    if (err == 0)  {
869#ifdef NCBI
870      if (format == kASNseqentry || format == kASNseqset)
871        seqlist = listASNSeqs( inputfile, skiplines, format, &nseq, &err);
872      else
873#endif
874        seqlist = listSeqs( inputfile, skiplines, format, &nseq, &err);
875      }
876
877    if (err != 0)
878      erralert(err);
879
880    else if (listonly) {
881      dumpSeqList(seqlist,format);
882      free( seqlist);
883      }
884
885    else {
886                                /* choose whichSeq if needed */
887      if (nseq == 1 || chooseall || (quietly && whichSeq == 0)) {
888        chooseall= true;
889        whichSeq = 1;
890        quietly = true; /* no loop */
891        }
892      else if (whichSeq > nseq && quietly) {
893        erralert(eItemNotFound);
894        err= eItemNotFound;
895        }
896      else if (whichSeq > nseq || !quietly) {
897        dumpSeqList(seqlist, format);
898        fprintf(stderr,"\nChoose a sequence (# or All): \n");
899        gets_noOverflow(stemp= stempstore, USERINPUT_BUFFERSIZE_SEQ);
900        skipwhitespace(stemp);
901        if (to_lower(*stemp) == 'a') {
902          chooseall= true;
903          whichSeq = 1;
904          quietly = true; /* !? this means we don't ask for another file
905                            as well as no more whichSeqs... */
906          }
907        else if (rs_isdigit(*stemp)) whichSeq= atol(stemp);
908        else whichSeq= 1; /* default */
909        }
910      free( seqlist);
911
912      if (false /*chooseall*/) {  /* this isn't debugged yet...*/
913        fin = fopen(inputfile, "r");
914        closein= true;
915        }
916
917      while (whichSeq > 0 && whichSeq <= nseq) {
918                                /* need to open multiple output files ? */
919        manyout = ((chooseall || nwhichlist>1) && nseq > 1
920                  && (outform == kPlain || outform == kGCG));
921        if (manyout) {
922          if ( whichSeq == 1 ) erralert(eOneFormat);
923          else if (closeout) {
924            sprintf( stemp,"%s_%ld", oname, whichSeq);
925            freopen( stemp, "w", fout);
926            fprintf( stderr,"Writing sequence %ld to file %s\n", whichSeq, stemp);
927            }
928          }
929
930        if (closein) {
931          /* !! this fails... skips most seqs... */
932          /* !! in sequential read, must count seqs already read from whichSeq ... */
933          /* need major revision of ureadseq before we can do this */
934          atseq= whichSeq-1;
935          seqidptr= seqid;
936          seq = readSeqFp( whichSeq, fin, skiplines, format,
937                          &seqlen, &atseq, &err, seqidptr);
938          skiplines= 0;
939          }
940        else {
941          atseq= 0;
942          seqidptr= seqid;
943#ifdef NCBI
944          if (format == kASNseqentry || format == kASNseqset) {
945            seqidptr= NULL;
946            seq = readASNSeq( whichSeq, inputfile, skiplines, format,
947                     &seqlen, &atseq, &err, &seqidptr);
948            }
949          else
950#endif
951          seq = readSeq( whichSeq, inputfile, skiplines, format,
952                          &seqlen, &atseq, &err, seqidptr);
953          }
954
955
956        if (gPretty.degap) {
957          char *newseq;
958          long newlen;
959          newseq= compressSeq( gPretty.gapchar, seq, seqlen, &newlen);
960          if (newseq) {
961            free(seq); seq= newseq; seqlen= newlen;
962            }
963          }
964
965        if (outform == kMSF) checksum= GCGchecksum(seq, seqlen, &checkall);
966        else if (verbose) checksum= seqchecksum(seq, seqlen, &checkall);
967        if (verbose)
968          fprintf( stderr, "Sequence %ld, length= %ld, checksum= %lX, format= %s, id= %s\n",
969                whichSeq, seqlen, checksum, formatstr(format), seqidptr);
970
971        if (err != 0) erralert(err);
972        else {
973                                  /* format fixes that writeseq doesn't do */
974          switch (outform) {
975            case kPIR:
976              if (seqout == 0) fprintf( foo,"\\\\\\\n");
977              break;
978            case kASN1:
979              if (seqout == 0) fputs(kASN1headline, foo);
980              break;
981
982            case kPhylip:
983              if (seqout == 0) {
984                if (!interleaved) {  /*  bug, nseq is for 1st infile only */
985                  if (chooseall) i= nseq; else i=1;
986                  if (phylvers >= 4) fprintf(foo," %ld %ld\n", i, seqlen);
987                  else fprintf(foo," %ld %ld YF\n", i, seqlen);
988                  }
989                seqlen0 = seqlen;
990                }
991              else if (seqlen != seqlen0) {
992                erralert(eUnequalSize);
993                if (seqlen < seqlen0) seq = (char *)realloc(seq, seqlen0);
994                for (i=seqlen; i<seqlen0; i++) seq[i]= gPretty.gapchar;
995                seqlen = seqlen0;
996                seq[seqlen] = 0;
997                }
998              break;
999
1000            case kPAUP:
1001              if (seqout == 0) {
1002                seqtype= getseqtype(seq, seqlen);
1003                seqlen0 = seqlen;
1004                }
1005              else if (seqlen != seqlen0) {
1006                erralert(eUnequalSize);
1007                if (seqlen < seqlen0) seq = (char *)realloc(seq, seqlen0);
1008                for (i=seqlen; i<seqlen0; i++) seq[i]= gPretty.gapchar;
1009                seqlen = seqlen0;
1010                seq[seqlen] = 0;
1011                }
1012              break;
1013
1014            }
1015
1016          if (doupper)
1017            for (i = 0; i<seqlen; i++) seq[i] = to_upper(seq[i]);
1018          else if (dolower)
1019            for (i = 0; i<seqlen; i++) seq[i] = to_lower(seq[i]);
1020
1021          if (outform==kPhylip){
1022            for (i = 0; i<seqlen; i++) if (seq[i] == '.') seq[i] = '?';
1023          }
1024
1025          if (doreverse) {
1026            long  j, k;
1027            char  ctemp;
1028            for (j=0, k=seqlen-1; j <= k; j++, k--) {
1029              ctemp = compl[seq[j] - ' '];
1030              seq[j] = compl[seq[k] - ' '];
1031              seq[k] = ctemp;
1032              }
1033            }
1034
1035          if ((gPretty.isactive || outform==kPAUP) && gPretty.domatch && firstseq != NULL) {
1036            for (i=0; i<seqlen; i++){
1037                if (seq[i] == gPretty.matchchar) seq[i] = 'o';
1038              if (seq[i]==firstseq[i]) seq[i]= gPretty.matchchar;
1039            }
1040          }
1041
1042
1043          if (gPretty.isactive && gPretty.numtop && seqout == 0) {
1044            gPretty.numline = 1;
1045            indexout();
1046            (void) writeSeq( fout, seq, seqlen, outform, seqidptr);
1047            gPretty.numline = 2;
1048            indexout();
1049            (void) writeSeq( fout, seq, seqlen, outform, seqidptr);
1050            gPretty.numline = 0;
1051            }
1052
1053          indexout();
1054          nlines = writeSeq( fout, seq, seqlen, outform, seqidptr);
1055          seqout++;
1056          }
1057
1058        if ((gPretty.isactive || outform==kPAUP) && gPretty.domatch && firstseq == NULL) {
1059          firstseq= seq;
1060          seq = NULL;
1061          }
1062        else if (seq!=NULL) { free(seq); seq = NULL; }
1063
1064#ifdef NCBI
1065       if ( (format == kASNseqentry || format == kASNseqset)
1066          && seqidptr && seqidptr!= seqid)
1067            free(seqidptr);
1068#endif
1069        if (chooseall) whichSeq++;
1070        else if (iwhichlist<nwhichlist) whichSeq= whichlist[iwhichlist++];
1071        else whichSeq= 0;
1072        }
1073      if (closein) { fclose(fin); closein= false; }
1074      }
1075    whichSeq  = 0;
1076  } while (nfile > 0 || !quietly);
1077
1078
1079fini:
1080  if (firstseq) { free(firstseq); firstseq= NULL; }
1081  if (err || listonly) exit_main(err);
1082
1083  if (gPretty.isactive && gPretty.numbot) {
1084    gPretty.numline = 2;
1085    indexout();
1086    (void) writeSeq( fout, seq, seqlen, outform, seqidptr);
1087    gPretty.numline = 1;
1088    indexout();
1089    (void) writeSeq( fout, seq, seqlen, outform, seqidptr);
1090    gPretty.numline = 0;
1091    }
1092
1093  if (outform == kMSF) {
1094    if (*oname) cp= oname; else cp= inputfile;
1095    fprintf(foo,"\n %s  MSF: %ld  Type: N  January 01, 1776  12:00  Check: %lu ..\n\n",
1096                  cp, seqlen, checkall);
1097    }
1098
1099  if (outform == kPAUP) {
1100    fprintf(foo,"#NEXUS\n");
1101    if (*oname) cp= oname; else cp= inputfile;
1102    fprintf(foo,"[%s -- data title]\n\n", cp);
1103    /* ! now have header lines for each sequence... put them before "begin data;... */
1104    }
1105
1106  if (outform==kPhylip && interleaved) {
1107    if (phylvers >= 4) fprintf(foo," %d %ld\n", seqout, seqlen);
1108    else fprintf(foo," %d %ld YF\n", seqout, seqlen);
1109    }
1110
1111  if (interleaved) {
1112    /* interleave species lines in true output */
1113    /* nlines is # lines / sequence */
1114    short iline, j, leaf, iseq;
1115    char  *s = stempstore;
1116
1117    indexout();  noutindex--; /* mark eof */
1118
1119    for (leaf=0; leaf<nlines; leaf++) {
1120      if (outform == kMSF && leaf == 1) {
1121        fputs("//\n\n", foo);
1122        }
1123      if (outform == kPAUP && leaf==1) {
1124        switch (seqtype) {
1125          case kDNA     : cp= "dna"; break;
1126          case kRNA     : cp= "rna"; break;
1127          case kNucleic : cp= "dna"; break;
1128          case kAmino   : cp= "protein"; break;
1129          case kOtherSeq: cp= "dna"; break;
1130          }
1131        fprintf(foo,"\nbegin data;\n");
1132        fprintf(foo," dimensions ntax=%d nchar=%ld;\n", seqout, seqlen);
1133        /* fix by Ralf Westram (ARB): '-' means 'gap' ('.' means 'missing') */
1134        fprintf(foo," format datatype=%s interleave=yes missing=. gap=%c", cp, gPretty.gapchar);
1135        if (gPretty.domatch) fprintf(foo," matchchar=%c", gPretty.matchchar);
1136        fprintf(foo,";\n  matrix\n");
1137        }
1138
1139      for (iseq=0; iseq<noutindex; iseq++) {
1140        fseek(ftmp, outindex[iseq], 0);
1141        for (iline=0; iline<=leaf; iline++)
1142          if (!fgets(s, 256, ftmp)) *s= 0;
1143        if (ftell(ftmp) <= outindex[iseq+1])
1144          fputs( s, foo);
1145        }
1146
1147      for (j=0; j<gPretty.interline; j++)
1148        fputs( "\n", foo);  /* some want spacer line */
1149      }
1150    fclose(ftmp); /* tmp disappears */
1151    fout= foo;
1152    }
1153
1154  if (outform == kASN1)  fprintf( foo, "} }\n");
1155  if (outform == kPAUP)  fprintf( foo,";\n  end;\n");
1156
1157  if (outindex != NULL) free(outindex);
1158  exit_main(0);
1159}
1160
1161
Note: See TracBrowser for help on using the repository browser.