src/parse.c

   1 /*
   2  * This file is part of fastq-tools.
   3  *
   4  * Copyright (c) 2011 by Daniel C. Jones <dcjones@cs.washington.edu>
   5  *
   6  */
   7
   8 #include "parse.h"
   9 #include "common.h"
  10 #include <stdlib.h>
  11 #include <ctype.h>
  12
  13
  14 static const size_t init_str_size  = 128;
  15 static const size_t fastq_buf_size = 4096;
  16
  17 static void fastq_alloc_str(str_t* s)
  18 {
  19     s->s = malloc_or_die(init_str_size);
  20     s->s[0] = '\0';
  21     s->n = 0;
  22     s->size = init_str_size;
  23 }
  24
  25
  26 static void fastq_expand_str(str_t* s)
  27 {
  28     s->size *= 2;
  29     realloc_or_die(s->s, s->size);
  30 }
  31
  32
  33 seq_t* fastq_alloc_seq()
  34 {
  35     seq_t* seq = malloc_or_die(sizeof(seq_t));
  36     fastq_alloc_str(&seq->id1);
  37     fastq_alloc_str(&seq->seq);
  38     fastq_alloc_str(&seq->id2);
  39     fastq_alloc_str(&seq->qual);
  40
  41     return seq;
  42 }
  43
  44
  45 void fastq_free_seq(seq_t* seq)
  46 {
  47     free(seq->id1.s);
  48     free(seq->seq.s);
  49     free(seq->id2.s);
  50     free(seq->qual.s);
  51     free(seq);
  52 }
  53
  54
  55 typedef enum
  56 {
  57     STATE_EOF,
  58     STATE_ERROR,
  59     STATE_ID1,
  60     STATE_SEQ,
  61     STATE_ID2,
  62     STATE_QUAL
  63
  64 } fastq_state;
  65
  66
  67 fastq_t* fastq_open(FILE* f)
  68 {
  69     fastq_t* fqf = malloc_or_die(sizeof(fastq_t));
  70     or_die((int)(fqf->file = gzdopen(fileno(f), "rb")),
  71            "Can not open gzip file.");
  72
  73     fqf->state = STATE_ID1;
  74     fqf->buf = malloc_or_die(fastq_buf_size);
  75     fqf->buf[0] = '\0';
  76     fqf->c = fqf->buf;
  77
  78     return fqf;
  79 }
  80
  81
  82 void fastq_close(fastq_t* fqf)
  83 {
  84     gzclose(fqf->file);
  85     free(fqf->buf);
  86     free(fqf);
  87 }
  88
  89
  90 void fastq_refill(fastq_t* f)
  91 {
  92     int errnum;
  93     const char* errmsg;
  94
  95     int n = gzread(f->file, f->buf, fastq_buf_size - 1);
  96
  97     if (n <= 0) {
  98         if (gzeof(f->file)) {
  99             f->state = STATE_EOF;
 100             n = 0;
 101         }
 102         else {
 103             errmsg = gzerror(f->file, &errnum);
 104             fprintf(stderr, "I/O error: %s\n", errmsg);
 105             exit(1);
 106         }
 107     }
 108
 109     f->buf[n] = '\0';
 110     f->c = f->buf;
 111 }
 112
 113
 114 void fastq_get_line(fastq_t* f, str_t* s)
 115 {
 116     int i = 0;
 117
 118     if (f->state == STATE_EOF) goto fastq_get_line_done;
 119
 120     while (1) {
 121         switch (*f->c) {
 122             case '\0':
 123                 fastq_refill(f);
 124                 if (f->state == STATE_EOF) goto fastq_get_line_done;
 125                 break;
 126
 127             case '\r':
 128                 f->c++;
 129                 break;
 130
 131             case '\n':
 132                 goto fastq_get_line_done;
 133
 134             default:
 135                 while (s->size < i + 2) {
 136                     fastq_expand_str(s);
 137                 }
 138                 if (s) s->s[i++] = *f->c;
 139                 f->c++;
 140         }
 141
 142     }
 143
 144 fastq_get_line_done:
 145     if (s) {
 146         s->s[i] = '\0';
 147         s->n = i;
 148     }
 149 }
 150
 151
 152
 153 int fastq_next(fastq_t* f, seq_t* seq)
 154 {
 155     if (f->state == STATE_EOF) return 0;
 156
 157     while (1) {
 158
 159         /* read more, if needed */
 160         if (*f->c == '\0' ) {
 161             fastq_refill(f);
 162             if (f->state == STATE_EOF) return 0;
 163             continue;
 164         }
 165
 166         /* skip over leading whitespace */
 167         else if (isspace(*f->c)) {
 168             /* do nothing */
 169         }
 170
 171         /* skip comments */
 172         else if (*f->c == ';') {
 173             fastq_get_line(f, NULL);
 174             if (f->state == STATE_EOF) return 0;
 175         }
 176
 177         /* read id1 */
 178         else if (f->state == STATE_ID1) {
 179             if (*f->c == '@' || *f->c == '>') {
 180                 f->c++;
 181                 fastq_get_line(f, &seq->id1);
 182                 if (f->state == STATE_EOF) return 0;
 183
 184                 f->state = STATE_SEQ;
 185             }
 186             else {
 187                 fprintf(stderr,
 188                         "Malformed FASTQ file: expecting an '@' or '>', saw a '%c'\n",
 189                         *f->c);
 190                 exit(1);
 191             }
 192         }
 193
 194         /* read sequence */
 195         else if (f->state == STATE_SEQ) {
 196             fastq_get_line(f, &seq->seq);
 197             if (f->state == STATE_EOF) return 0;
 198
 199             f->state = STATE_ID2;
 200         }
 201
 202         /* read id2 */
 203         else if (f->state == STATE_ID2) {
 204             if (*f->c == '+') {
 205                 f->c++;
 206                 fastq_get_line(f, &seq->id2);
 207                 if (f->state == STATE_EOF) return 0;
 208
 209                 f->state = STATE_QUAL;
 210             }
 211             else {
 212                 /* fasta style entry */
 213                 seq->id2.s[0]  = '\0';
 214                 seq->qual.s[0] = '\0';
 215
 216                 f->state = STATE_ID1;
 217                 break;
 218             }
 219         }
 220
 221         /* read quality string */
 222         else if (f->state == STATE_QUAL) {
 223             fastq_get_line(f, &seq->qual);
 224             if (f->state == STATE_EOF) return 1;
 225
 226             f->state = STATE_ID1;
 227             break;
 228         }
 229
 230         else {
 231             fputs("Inexplicable error in fastq parser.\n", stderr);
 232             exit(1);
 233         }
 234
 235         f->c++;
 236     }
 237
 238     return 1;
 239 }
 240
 241
 242 void fastq_print(FILE* fout, seq_t* seq)
 243 {
 244     /* FASTQ */
 245     if (seq->qual.n > 0) {
 246         fprintf(fout, "@%s\n%s\n+%s\n%s\n",
 247                       seq->id1.s,
 248                       seq->seq.s,
 249                       seq->id2.s,
 250                       seq->qual.s );
 251     }
 252
 253     /* FASTA */
 254     else {
 255         fprintf(fout, ">%s\n%s\n",
 256                       seq->id1.s,
 257                       seq->seq.s );
 258     }
 259 }
 260
 261