src/fastq-grep.c

   1
   2
   3 /*
   4  * fastq-grep: regular expression searches of the sequences within a fastq file
   5  *
   6  * Febuary 2011 / Daniel Jones <dcjones@cs.washington.edu>
   7  *
   8  */
   9
  10
  11 #include <stdio.h>
  12 #include <string.h>
  13 #include <getopt.h>
  14 #include <zlib.h>
  15 #include <pcre.h>
  16
  17 #include "kseq.h"
  18 KSEQ_INIT(gzFile, gzread)
  19
  20
  21 #if defined(MSDOS) || defined(OS2) || defined(WIN32) || defined(__CYGWIN__)
  22 #  include <fcntl.h>
  23 #  include <io.h>
  24 #  define SET_BINARY_MODE(file) setmode(fileno(file), O_BINARY)
  25 #else
  26 #  define SET_BINARY_MODE(file)
  27 #endif
  28
  29
  30 void print_help()
  31 {
  32     fprintf( stderr,
  33 "fastq-grep [OPTION]... PATTERN [FILE]...\n"
  34 "Search for PATTERN in the read sequences in each FILE or standard input.\n"
  35 "PATTERN, by default, is a perl compatible regular expression.\n\n"
  36 "Options:\n"
  37 "  -h, --help              print this message\n"
  38 "  -v, --invert-match      select nonmatching entries\n"
  39     );
  40 }
  41
  42 static int invert_flag;
  43 static int help_flag;
  44 static int zout_flag;
  45
  46
  47
  48 void print_fastq_entry( FILE* fout, kseq_t* seq )
  49 {
  50     fprintf(fout, "@%s\n%s\n+%s\n%s\n",
  51                   seq->name.s,
  52                   seq->seq.s,
  53                   seq->comment.s,
  54                   seq->qual.s );
  55 }
  56
  57
  58 void fastq_grep( gzFile fin, FILE* fout, pcre* re )
  59 {
  60     int rc;
  61     int ovector[3];
  62
  63     kseq_t* seq;
  64     seq = kseq_init(fin);
  65
  66     while (kseq_read(seq) >= 0) {
  67         rc = pcre_exec(re,          /* pattern */
  68                        NULL,        /* extre data */
  69                        seq->seq.s,  /* subject */
  70                        seq->seq.l,  /* subject length */
  71                        0,           /* subject offset */
  72                        0,           /* options */
  73                        ovector,     /* output vector */
  74                        3         ); /* output vector length */
  75
  76         if ((invert_flag && rc == PCRE_ERROR_NOMATCH) || rc >= 0) {
  77             print_fastq_entry( fout, seq );
  78         }
  79     }
  80
  81     kseq_destroy(seq);
  82 }
  83
  84
  85
  86 int main(int argc, char* argv[])
  87 {
  88     const char* pat;
  89     pcre* re;
  90     const char* pat_error;
  91     int pat_error_offset;
  92
  93     FILE*  fin;
  94     gzFile gzfin;
  95
  96
  97     invert_flag = 0;
  98     help_flag   = 0;
  99     zout_flag   = 0;
 100
 101     int opt;
 102     int opt_idx;
 103
 104
 105     static struct option long_options[] =
 106         {
 107           {"help", no_argument, &help_flag, 1},
 108           {"invert-match", no_argument, &invert_flag, 1},
 109           {0, 0, 0, 0}
 110         };
 111
 112     while (1) {
 113         opt = getopt_long(argc, argv, "hv", long_options, &opt_idx);
 114
 115         if( opt == -1 ) break;
 116
 117         switch (opt) {
 118             case 0:
 119                 if (long_options[opt_idx].flag != 0) break;
 120                 if (optarg) {
 121                 }
 122                 break;
 123
 124             case 'h':
 125                 help_flag = 1;
 126                 break;
 127
 128             case 'v':
 129                 invert_flag = 1;
 130                 break;
 131
 132             case '?':
 133                 return 1;
 134
 135             default:
 136                 abort();
 137         }
 138     }
 139
 140     if (help_flag) {
 141         print_help();
 142         return 0;
 143     }
 144
 145     if (optind >= argc) {
 146         fprintf(stderr, "A pattern must be specified.\n");
 147         return 1;
 148     }
 149
 150     pat = argv[optind++];
 151     re = pcre_compile( pat, PCRE_CASELESS, &pat_error, &pat_error_offset, NULL );
 152
 153
 154     if (re == NULL) {
 155         fprintf(stderr, "Syntax error in PCRE pattern at offset: %d: %s\n",
 156                 pat_error_offset, pat_error );
 157         return 1;
 158     }
 159
 160
 161     if (optind >= argc || (argc - optind == 1 && strcmp(argv[optind],"-") == 0)) {
 162         gzfin = gzdopen( fileno(stdin), "rb" );
 163         if (gzfin == NULL) {
 164             fprintf(stderr, "Malformed file 'stdin'.\n");
 165             return 1;
 166         }
 167
 168         fastq_grep(gzfin, stdout, re);
 169
 170         gzclose(gzfin);
 171     }
 172     else {
 173         for (; optind < argc; optind++) {
 174             fin = fopen(argv[optind], "rb");
 175             if (fin == NULL) {
 176                 fprintf(stderr, "No such file '%s'.\n", argv[optind]);
 177                 continue;
 178             }
 179
 180             gzfin = gzdopen(fileno(fin), "rb");
 181             if (gzfin == NULL) {
 182                 fprintf(stderr, "Malformed file '%s'.\n", argv[optind]);
 183                 continue;
 184             }
 185
 186             fastq_grep(gzfin, stdout, re);
 187
 188             gzclose(gzfin);
 189         }
 190     }
 191
 192     pcre_free(re);
 193
 194     return 0;
 195 }
 196
 197
 198