|
9
|
1 #include <assert.h>
|
|
|
2 #include <ctype.h>
|
|
|
3 #include <stdlib.h>
|
|
|
4 #include <zlib.h>
|
|
|
5 #include <stdio.h>
|
|
|
6 #include <getopt.h>
|
|
|
7 #include <unistd.h>
|
|
|
8 #include "sickle.h"
|
|
|
9 #include "kseq.h"
|
|
|
10 #include "print_record.h"
|
|
|
11
|
|
|
12 __KS_GETC(gzread, BUFFER_SIZE)
|
|
|
13 __KS_GETUNTIL(gzread, BUFFER_SIZE)
|
|
|
14 __KSEQ_READ
|
|
|
15
|
|
|
16 int paired_qual_threshold = 20;
|
|
|
17 int paired_length_threshold = 20;
|
|
|
18
|
|
|
19 static struct option paired_long_options[] = {
|
|
|
20 {"qual-type", required_argument, 0, 't'},
|
|
|
21 {"pe-file1", required_argument, 0, 'f'},
|
|
|
22 {"pe-file2", required_argument, 0, 'r'},
|
|
|
23 {"pe-combo", required_argument, 0, 'c'},
|
|
|
24 {"output-pe1", required_argument, 0, 'o'},
|
|
|
25 {"output-pe2", required_argument, 0, 'p'},
|
|
|
26 {"output-single", required_argument, 0, 's'},
|
|
|
27 {"output-combo", required_argument, 0, 'm'},
|
|
|
28 {"qual-threshold", required_argument, 0, 'q'},
|
|
|
29 {"length-threshold", required_argument, 0, 'l'},
|
|
|
30 {"no-fiveprime", no_argument, 0, 'x'},
|
|
|
31 {"truncate-n", no_argument, 0, 'n'},
|
|
|
32 {"gzip-output", no_argument, 0, 'g'},
|
|
|
33 {"output-combo-all", required_argument, 0, 'M'},
|
|
|
34 {"quiet", no_argument, 0, 'z'},
|
|
|
35 {GETOPT_HELP_OPTION_DECL},
|
|
|
36 {GETOPT_VERSION_OPTION_DECL},
|
|
|
37 {NULL, 0, NULL, 0}
|
|
|
38 };
|
|
|
39
|
|
|
40 void paired_usage (int status, char *msg) {
|
|
|
41
|
|
|
42 fprintf(stderr, "\nIf you have separate files for forward and reverse reads:\n");
|
|
|
43 fprintf(stderr, "Usage: %s pe [options] -f <paired-end forward fastq file> -r <paired-end reverse fastq file> -t <quality type> -o <trimmed PE forward file> -p <trimmed PE reverse file> -s <trimmed singles file>\n\n", PROGRAM_NAME);
|
|
|
44 fprintf(stderr, "If you have one file with interleaved forward and reverse reads:\n");
|
|
|
45 fprintf(stderr, "Usage: %s pe [options] -c <interleaved input file> -t <quality type> -m <interleaved trimmed paired-end output> -s <trimmed singles file>\n\n\
|
|
|
46 If you have one file with interleaved reads as input and you want ONLY one interleaved file as output:\n\
|
|
|
47 Usage: %s pe [options] -c <interleaved input file> -t <quality type> -M <interleaved trimmed output>\n\n", PROGRAM_NAME, PROGRAM_NAME);
|
|
|
48 fprintf(stderr, "Options:\n\
|
|
|
49 Paired-end separated reads\n\
|
|
|
50 --------------------------\n\
|
|
|
51 -f, --pe-file1, Input paired-end forward fastq file (Input files must have same number of records)\n\
|
|
|
52 -r, --pe-file2, Input paired-end reverse fastq file\n\
|
|
|
53 -o, --output-pe1, Output trimmed forward fastq file\n\
|
|
|
54 -p, --output-pe2, Output trimmed reverse fastq file. Must use -s option.\n\n\
|
|
|
55 Paired-end interleaved reads\n\
|
|
|
56 ----------------------------\n");
|
|
|
57 fprintf(stderr,"-c, --pe-combo, Combined (interleaved) input paired-end fastq\n\
|
|
|
58 -m, --output-combo, Output combined (interleaved) paired-end fastq file. Must use -s option.\n\
|
|
|
59 -M, --output-combo-all, Output combined (interleaved) paired-end fastq file with any discarded read written to output file as a single N. Cannot be used with the -s option.\n\n\
|
|
|
60 Global options\n\
|
|
|
61 --------------\n\
|
|
|
62 -t, --qual-type, Type of quality values (solexa (CASAVA < 1.3), illumina (CASAVA 1.3 to 1.7), sanger (which is CASAVA >= 1.8)) (required)\n");
|
|
|
63 fprintf(stderr, "-s, --output-single, Output trimmed singles fastq file\n\
|
|
|
64 -q, --qual-threshold, Threshold for trimming based on average quality in a window. Default 20.\n\
|
|
|
65 -l, --length-threshold, Threshold to keep a read based on length after trimming. Default 20.\n\
|
|
|
66 -x, --no-fiveprime, Don't do five prime trimming.\n\
|
|
|
67 -n, --truncate-n, Truncate sequences at position of first N.\n");
|
|
|
68
|
|
|
69
|
|
|
70 fprintf(stderr, "-g, --gzip-output, Output gzipped files.\n--quiet, do not output trimming info\n\
|
|
|
71 --help, display this help and exit\n\
|
|
|
72 --version, output version information and exit\n\n");
|
|
|
73
|
|
|
74 if (msg) fprintf(stderr, "%s\n\n", msg);
|
|
|
75 exit(status);
|
|
|
76 }
|
|
|
77
|
|
|
78
|
|
|
79 int paired_main(int argc, char *argv[]) {
|
|
|
80
|
|
|
81 gzFile pe1 = NULL; /* forward input file handle */
|
|
|
82 gzFile pe2 = NULL; /* reverse input file handle */
|
|
|
83 gzFile pec = NULL; /* combined input file handle */
|
|
|
84 kseq_t *fqrec1 = NULL;
|
|
|
85 kseq_t *fqrec2 = NULL;
|
|
|
86 int l1, l2;
|
|
|
87 FILE *outfile1 = NULL; /* forward output file handle */
|
|
|
88 FILE *outfile2 = NULL; /* reverse output file handle */
|
|
|
89 FILE *combo = NULL; /* combined output file handle */
|
|
|
90 FILE *single = NULL; /* single output file handle */
|
|
|
91 gzFile outfile1_gzip = NULL;
|
|
|
92 gzFile outfile2_gzip = NULL;
|
|
|
93 gzFile combo_gzip = NULL;
|
|
|
94 gzFile single_gzip = NULL;
|
|
|
95 int debug = 0;
|
|
|
96 int optc;
|
|
|
97 extern char *optarg;
|
|
|
98 int qualtype = -1;
|
|
|
99 cutsites *p1cut;
|
|
|
100 cutsites *p2cut;
|
|
|
101 char *outfn1 = NULL; /* forward file out name */
|
|
|
102 char *outfn2 = NULL; /* reverse file out name */
|
|
|
103 char *outfnc = NULL; /* combined file out name */
|
|
|
104 char *sfn = NULL; /* single/combined file out name */
|
|
|
105 char *infn1 = NULL; /* forward input filename */
|
|
|
106 char *infn2 = NULL; /* reverse input filename */
|
|
|
107 char *infnc = NULL; /* combined input filename */
|
|
|
108 int kept_p = 0;
|
|
|
109 int discard_p = 0;
|
|
|
110 int kept_s1 = 0;
|
|
|
111 int kept_s2 = 0;
|
|
|
112 int discard_s1 = 0;
|
|
|
113 int discard_s2 = 0;
|
|
|
114 int quiet = 0;
|
|
|
115 int no_fiveprime = 0;
|
|
|
116 int trunc_n = 0;
|
|
|
117 int gzip_output = 0;
|
|
|
118 int combo_all=0;
|
|
|
119 int combo_s=0;
|
|
|
120 int total=0;
|
|
|
121
|
|
|
122 while (1) {
|
|
|
123 int option_index = 0;
|
|
|
124 optc = getopt_long(argc, argv, "df:r:c:t:o:p:m:M:s:q:l:xng", paired_long_options, &option_index);
|
|
|
125
|
|
|
126 if (optc == -1)
|
|
|
127 break;
|
|
|
128
|
|
|
129 switch (optc) {
|
|
|
130 if (paired_long_options[option_index].flag != 0)
|
|
|
131 break;
|
|
|
132
|
|
|
133 case 'f':
|
|
|
134 infn1 = (char *) malloc(strlen(optarg) + 1);
|
|
|
135 strcpy(infn1, optarg);
|
|
|
136 break;
|
|
|
137
|
|
|
138 case 'r':
|
|
|
139 infn2 = (char *) malloc(strlen(optarg) + 1);
|
|
|
140 strcpy(infn2, optarg);
|
|
|
141 break;
|
|
|
142
|
|
|
143 case 'c':
|
|
|
144 infnc = (char *) malloc(strlen(optarg) + 1);
|
|
|
145 strcpy(infnc, optarg);
|
|
|
146 break;
|
|
|
147
|
|
|
148 case 't':
|
|
|
149 if (!strcmp(optarg, "illumina")) qualtype = ILLUMINA;
|
|
|
150 else if (!strcmp(optarg, "solexa")) qualtype = SOLEXA;
|
|
|
151 else if (!strcmp(optarg, "sanger")) qualtype = SANGER;
|
|
|
152 else {
|
|
|
153 fprintf(stderr, "Error: Quality type '%s' is not a valid type.\n", optarg);
|
|
|
154 return EXIT_FAILURE;
|
|
|
155 }
|
|
|
156 break;
|
|
|
157
|
|
|
158 case 'o':
|
|
|
159 outfn1 = (char *) malloc(strlen(optarg) + 1);
|
|
|
160 strcpy(outfn1, optarg);
|
|
|
161 break;
|
|
|
162
|
|
|
163 case 'p':
|
|
|
164 outfn2 = (char *) malloc(strlen(optarg) + 1);
|
|
|
165 strcpy(outfn2, optarg);
|
|
|
166 break;
|
|
|
167
|
|
|
168 case 'm':
|
|
|
169 outfnc = (char *) malloc(strlen(optarg) + 1);
|
|
|
170 strcpy(outfnc, optarg);
|
|
|
171 combo_s = 1;
|
|
|
172 break;
|
|
|
173
|
|
|
174 case 'M':
|
|
|
175 outfnc = (char *) malloc(strlen(optarg) + 1);
|
|
|
176 strcpy(outfnc, optarg);
|
|
|
177 combo_all = 1;
|
|
|
178 break;
|
|
|
179
|
|
|
180 case 's':
|
|
|
181 sfn = (char *) malloc(strlen(optarg) + 1);
|
|
|
182 strcpy(sfn, optarg);
|
|
|
183 break;
|
|
|
184
|
|
|
185 case 'q':
|
|
|
186 paired_qual_threshold = atoi(optarg);
|
|
|
187 if (paired_qual_threshold < 0) {
|
|
|
188 fprintf(stderr, "Quality threshold must be >= 0\n");
|
|
|
189 return EXIT_FAILURE;
|
|
|
190 }
|
|
|
191 break;
|
|
|
192
|
|
|
193 case 'l':
|
|
|
194 paired_length_threshold = atoi(optarg);
|
|
|
195 if (paired_length_threshold < 0) {
|
|
|
196 fprintf(stderr, "Length threshold must be >= 0\n");
|
|
|
197 return EXIT_FAILURE;
|
|
|
198 }
|
|
|
199 break;
|
|
|
200
|
|
|
201 case 'x':
|
|
|
202 no_fiveprime = 1;
|
|
|
203 break;
|
|
|
204
|
|
|
205 case 'n':
|
|
|
206 trunc_n = 1;
|
|
|
207 break;
|
|
|
208
|
|
|
209 case 'g':
|
|
|
210 gzip_output = 1;
|
|
|
211 break;
|
|
|
212
|
|
|
213 case 'z':
|
|
|
214 quiet = 1;
|
|
|
215 break;
|
|
|
216
|
|
|
217 case 'd':
|
|
|
218 debug = 1;
|
|
|
219 break;
|
|
|
220
|
|
|
221 case_GETOPT_HELP_CHAR(paired_usage);
|
|
|
222 case_GETOPT_VERSION_CHAR(PROGRAM_NAME, VERSION, AUTHORS);
|
|
|
223
|
|
|
224 case '?':
|
|
|
225 paired_usage(EXIT_FAILURE, NULL);
|
|
|
226 break;
|
|
|
227
|
|
|
228 default:
|
|
|
229 paired_usage(EXIT_FAILURE, NULL);
|
|
|
230 break;
|
|
|
231 }
|
|
|
232 }
|
|
|
233
|
|
|
234 /* required: qualtype */
|
|
|
235 if (qualtype == -1) {
|
|
|
236 paired_usage(EXIT_FAILURE, "****Error: Quality type is required.");
|
|
|
237 }
|
|
|
238
|
|
|
239 /* make sure minimum input filenames are specified */
|
|
|
240 if (!infn1 && !infnc) {
|
|
|
241 paired_usage(EXIT_FAILURE, "****Error: Must have either -f OR -c argument.");
|
|
|
242 }
|
|
|
243
|
|
|
244 if (infnc) { /* using combined input file */
|
|
|
245
|
|
|
246 if (infn1 || infn2 || outfn1 || outfn2) {
|
|
|
247 paired_usage(EXIT_FAILURE, "****Error: Cannot have -f, -r, -o, or -p options with -c.");
|
|
|
248 }
|
|
|
249
|
|
|
250 if ((combo_all && combo_s) || (!combo_all && !combo_s)) {
|
|
|
251 paired_usage(EXIT_FAILURE, "****Error: Must have only one of either -m or -M options with -c.");
|
|
|
252 }
|
|
|
253
|
|
|
254 if ((combo_s && !sfn) || (combo_all && sfn)) {
|
|
|
255 paired_usage(EXIT_FAILURE, "****Error: -m option must have -s option, and -M option cannot have -s option.");
|
|
|
256 }
|
|
|
257
|
|
|
258 /* check for duplicate file names */
|
|
|
259 if (!strcmp(infnc, outfnc) || (combo_s && (!strcmp(infnc, sfn) || !strcmp(outfnc, sfn)))) {
|
|
|
260 fprintf(stderr, "****Error: Duplicate filename between combo input, combo output, and/or single output file names.\n\n");
|
|
|
261 return EXIT_FAILURE;
|
|
|
262 }
|
|
|
263
|
|
|
264 /* get combined output file */
|
|
|
265 if (!gzip_output) {
|
|
|
266 combo = fopen(outfnc, "w");
|
|
|
267 if (!combo) {
|
|
|
268 fprintf(stderr, "****Error: Could not open combo output file '%s'.\n\n", outfnc);
|
|
|
269 return EXIT_FAILURE;
|
|
|
270 }
|
|
|
271 } else {
|
|
|
272 combo_gzip = gzopen(outfnc, "w");
|
|
|
273 if (!combo_gzip) {
|
|
|
274 fprintf(stderr, "****Error: Could not open combo output file '%s'.\n\n", outfnc);
|
|
|
275 return EXIT_FAILURE;
|
|
|
276 }
|
|
|
277 }
|
|
|
278
|
|
|
279 pec = gzopen(infnc, "r");
|
|
|
280 if (!pec) {
|
|
|
281 fprintf(stderr, "****Error: Could not open combined input file '%s'.\n\n", infnc);
|
|
|
282 return EXIT_FAILURE;
|
|
|
283 }
|
|
|
284
|
|
|
285 } else { /* using forward and reverse input files */
|
|
|
286
|
|
|
287 if (infn1 && (!infn2 || !outfn1 || !outfn2 || !sfn)) {
|
|
|
288 paired_usage(EXIT_FAILURE, "****Error: Using the -f option means you must have the -r, -o, -p, and -s options.");
|
|
|
289 }
|
|
|
290
|
|
|
291 if (infn1 && (infnc || combo_all || combo_s)) {
|
|
|
292 paired_usage(EXIT_FAILURE, "****Error: The -f option cannot be used in combination with -c, -m, or -M.");
|
|
|
293 }
|
|
|
294
|
|
|
295 if (!strcmp(infn1, infn2) || !strcmp(infn1, outfn1) || !strcmp(infn1, outfn2) ||
|
|
|
296 !strcmp(infn1, sfn) || !strcmp(infn2, outfn1) || !strcmp(infn2, outfn2) ||
|
|
|
297 !strcmp(infn2, sfn) || !strcmp(outfn1, outfn2) || !strcmp(outfn1, sfn) || !strcmp(outfn2, sfn)) {
|
|
|
298
|
|
|
299 fprintf(stderr, "****Error: Duplicate input and/or output file names.\n\n");
|
|
|
300 return EXIT_FAILURE;
|
|
|
301 }
|
|
|
302
|
|
|
303 pe1 = gzopen(infn1, "r");
|
|
|
304 if (!pe1) {
|
|
|
305 fprintf(stderr, "****Error: Could not open input file '%s'.\n\n", infn1);
|
|
|
306 return EXIT_FAILURE;
|
|
|
307 }
|
|
|
308
|
|
|
309 pe2 = gzopen(infn2, "r");
|
|
|
310 if (!pe2) {
|
|
|
311 fprintf(stderr, "****Error: Could not open input file '%s'.\n\n", infn2);
|
|
|
312 return EXIT_FAILURE;
|
|
|
313 }
|
|
|
314
|
|
|
315 if (!gzip_output) {
|
|
|
316 outfile1 = fopen(outfn1, "w");
|
|
|
317 if (!outfile1) {
|
|
|
318 fprintf(stderr, "****Error: Could not open output file '%s'.\n\n", outfn1);
|
|
|
319 return EXIT_FAILURE;
|
|
|
320 }
|
|
|
321
|
|
|
322 outfile2 = fopen(outfn2, "w");
|
|
|
323 if (!outfile2) {
|
|
|
324 fprintf(stderr, "****Error: Could not open output file '%s'.\n\n", outfn2);
|
|
|
325 return EXIT_FAILURE;
|
|
|
326 }
|
|
|
327 } else {
|
|
|
328 outfile1_gzip = gzopen(outfn1, "w");
|
|
|
329 if (!outfile1_gzip) {
|
|
|
330 fprintf(stderr, "****Error: Could not open output file '%s'.\n\n", outfn1);
|
|
|
331 return EXIT_FAILURE;
|
|
|
332 }
|
|
|
333
|
|
|
334 outfile2_gzip = gzopen(outfn2, "w");
|
|
|
335 if (!outfile2_gzip) {
|
|
|
336 fprintf(stderr, "****Error: Could not open output file '%s'.\n\n", outfn2);
|
|
|
337 return EXIT_FAILURE;
|
|
|
338 }
|
|
|
339
|
|
|
340 }
|
|
|
341 }
|
|
|
342
|
|
|
343 /* get singles output file handle */
|
|
|
344 if (sfn && !combo_all) {
|
|
|
345 if (!gzip_output) {
|
|
|
346 single = fopen(sfn, "w");
|
|
|
347 if (!single) {
|
|
|
348 fprintf(stderr, "****Error: Could not open single output file '%s'.\n\n", sfn);
|
|
|
349 return EXIT_FAILURE;
|
|
|
350 }
|
|
|
351 } else {
|
|
|
352 single_gzip = gzopen(sfn, "w");
|
|
|
353 if (!single_gzip) {
|
|
|
354 fprintf(stderr, "****Error: Could not open single output file '%s'.\n\n", sfn);
|
|
|
355 return EXIT_FAILURE;
|
|
|
356 }
|
|
|
357 }
|
|
|
358 }
|
|
|
359
|
|
|
360 if (pec) {
|
|
|
361 fqrec1 = kseq_init(pec);
|
|
|
362 fqrec2 = (kseq_t *) malloc(sizeof(kseq_t));
|
|
|
363 fqrec2->f = fqrec1->f;
|
|
|
364 } else {
|
|
|
365 fqrec1 = kseq_init(pe1);
|
|
|
366 fqrec2 = kseq_init(pe2);
|
|
|
367 }
|
|
|
368
|
|
|
369 while ((l1 = kseq_read(fqrec1)) >= 0) {
|
|
|
370
|
|
|
371 l2 = kseq_read(fqrec2);
|
|
|
372 if (l2 < 0) {
|
|
|
373 fprintf(stderr, "Warning: PE file 2 is shorter than PE file 1. Disregarding rest of PE file 1.\n");
|
|
|
374 break;
|
|
|
375 }
|
|
|
376
|
|
|
377 p1cut = sliding_window(fqrec1, qualtype, paired_length_threshold, paired_qual_threshold, no_fiveprime, trunc_n, debug);
|
|
|
378 p2cut = sliding_window(fqrec2, qualtype, paired_length_threshold, paired_qual_threshold, no_fiveprime, trunc_n, debug);
|
|
|
379 total += 2;
|
|
|
380
|
|
|
381 if (debug) printf("p1cut: %d,%d\n", p1cut->five_prime_cut, p1cut->three_prime_cut);
|
|
|
382 if (debug) printf("p2cut: %d,%d\n", p2cut->five_prime_cut, p2cut->three_prime_cut);
|
|
|
383
|
|
|
384 /* The sequence and quality print statements below print out the sequence string starting from the 5' cut */
|
|
|
385 /* and then only print out to the 3' cut, however, we need to adjust the 3' cut */
|
|
|
386 /* by subtracting the 5' cut because the 3' cut was calculated on the original sequence */
|
|
|
387
|
|
|
388 /* if both sequences passed quality and length filters, then output both records */
|
|
|
389 if (p1cut->three_prime_cut >= 0 && p2cut->three_prime_cut >= 0) {
|
|
|
390 if (!gzip_output) {
|
|
|
391 if (pec) {
|
|
|
392 print_record (combo, fqrec1, p1cut);
|
|
|
393 print_record (combo, fqrec2, p2cut);
|
|
|
394 } else {
|
|
|
395 print_record (outfile1, fqrec1, p1cut);
|
|
|
396 print_record (outfile2, fqrec2, p2cut);
|
|
|
397 }
|
|
|
398 } else {
|
|
|
399 if (pec) {
|
|
|
400 print_record_gzip (combo_gzip, fqrec1, p1cut);
|
|
|
401 print_record_gzip (combo_gzip, fqrec2, p2cut);
|
|
|
402 } else {
|
|
|
403 print_record_gzip (outfile1_gzip, fqrec1, p1cut);
|
|
|
404 print_record_gzip (outfile2_gzip, fqrec2, p2cut);
|
|
|
405 }
|
|
|
406 }
|
|
|
407
|
|
|
408 kept_p += 2;
|
|
|
409 }
|
|
|
410
|
|
|
411 /* if only one sequence passed filter, then put its record in singles and discard the other */
|
|
|
412 /* or put an "N" record in if that option was chosen. */
|
|
|
413 else if (p1cut->three_prime_cut >= 0 && p2cut->three_prime_cut < 0) {
|
|
|
414 if (!gzip_output) {
|
|
|
415 if (combo_all) {
|
|
|
416 print_record (combo, fqrec1, p1cut);
|
|
|
417 print_record_N (combo, fqrec2, qualtype);
|
|
|
418 } else {
|
|
|
419 print_record (single, fqrec1, p1cut);
|
|
|
420 }
|
|
|
421 } else {
|
|
|
422 if (combo_all) {
|
|
|
423 print_record_gzip (combo_gzip, fqrec1, p1cut);
|
|
|
424 print_record_N_gzip (combo_gzip, fqrec2, qualtype);
|
|
|
425 } else {
|
|
|
426 print_record_gzip (single_gzip, fqrec1, p1cut);
|
|
|
427 }
|
|
|
428 }
|
|
|
429
|
|
|
430 kept_s1++;
|
|
|
431 discard_s2++;
|
|
|
432 }
|
|
|
433
|
|
|
434 else if (p1cut->three_prime_cut < 0 && p2cut->three_prime_cut >= 0) {
|
|
|
435 if (!gzip_output) {
|
|
|
436 if (combo_all) {
|
|
|
437 print_record_N (combo, fqrec1, qualtype);
|
|
|
438 print_record (combo, fqrec2, p2cut);
|
|
|
439 } else {
|
|
|
440 print_record (single, fqrec2, p2cut);
|
|
|
441 }
|
|
|
442 } else {
|
|
|
443 if (combo_all) {
|
|
|
444 print_record_N_gzip (combo_gzip, fqrec1, qualtype);
|
|
|
445 print_record_gzip (combo_gzip, fqrec2, p2cut);
|
|
|
446 } else {
|
|
|
447 print_record_gzip (single_gzip, fqrec2, p2cut);
|
|
|
448 }
|
|
|
449 }
|
|
|
450
|
|
|
451 kept_s2++;
|
|
|
452 discard_s1++;
|
|
|
453
|
|
|
454 } else {
|
|
|
455
|
|
|
456 /* If both records are to be discarded, but the -M option */
|
|
|
457 /* is being used, then output two "N" records */
|
|
|
458 if (combo_all) {
|
|
|
459 if (!gzip_output) {
|
|
|
460 print_record_N (combo, fqrec1, qualtype);
|
|
|
461 print_record_N (combo, fqrec2, qualtype);
|
|
|
462 } else {
|
|
|
463 print_record_N_gzip (combo_gzip, fqrec1, qualtype);
|
|
|
464 print_record_N_gzip (combo_gzip, fqrec2, qualtype);
|
|
|
465 }
|
|
|
466 }
|
|
|
467
|
|
|
468 discard_p += 2;
|
|
|
469 }
|
|
|
470
|
|
|
471 free(p1cut);
|
|
|
472 free(p2cut);
|
|
|
473 } /* end of while ((l1 = kseq_read (fqrec1)) >= 0) */
|
|
|
474
|
|
|
475 if (l1 < 0) {
|
|
|
476 l2 = kseq_read(fqrec2);
|
|
|
477 if (l2 >= 0) {
|
|
|
478 fprintf(stderr, "Warning: PE file 1 is shorter than PE file 2. Disregarding rest of PE file 2.\n");
|
|
|
479 }
|
|
|
480 }
|
|
|
481
|
|
|
482 if (!quiet) {
|
|
|
483 if (infn1 && infn2) fprintf(stdout, "\nPE forward file: %s\nPE reverse file: %s\n", infn1, infn2);
|
|
|
484 if (infnc) fprintf(stdout, "\nPE interleaved file: %s\n", infnc);
|
|
|
485 fprintf(stdout, "\nTotal input FastQ records: %d (%d pairs)\n", total, (total / 2));
|
|
|
486 fprintf(stdout, "\nFastQ paired records kept: %d (%d pairs)\n", kept_p, (kept_p / 2));
|
|
|
487 if (pec) fprintf(stdout, "FastQ single records kept: %d\n", (kept_s1 + kept_s2));
|
|
|
488 else fprintf(stdout, "FastQ single records kept: %d (from PE1: %d, from PE2: %d)\n", (kept_s1 + kept_s2), kept_s1, kept_s2);
|
|
|
489
|
|
|
490 fprintf(stdout, "FastQ paired records discarded: %d (%d pairs)\n", discard_p, (discard_p / 2));
|
|
|
491
|
|
|
492 if (pec) fprintf(stdout, "FastQ single records discarded: %d\n\n", (discard_s1 + discard_s2));
|
|
|
493 else fprintf(stdout, "FastQ single records discarded: %d (from PE1: %d, from PE2: %d)\n\n", (discard_s1 + discard_s2), discard_s1, discard_s2);
|
|
|
494 }
|
|
|
495
|
|
|
496 kseq_destroy(fqrec1);
|
|
|
497 if (pec) free(fqrec2);
|
|
|
498 else kseq_destroy(fqrec2);
|
|
|
499
|
|
|
500 if (sfn && !combo_all) {
|
|
|
501 if (!gzip_output) fclose(single);
|
|
|
502 else gzclose(single_gzip);
|
|
|
503 }
|
|
|
504
|
|
|
505 if (pec) {
|
|
|
506 gzclose(pec);
|
|
|
507 if (!gzip_output) fclose(combo);
|
|
|
508 else gzclose(combo_gzip);
|
|
|
509 } else {
|
|
|
510 gzclose(pe1);
|
|
|
511 gzclose(pe2);
|
|
|
512 if (!gzip_output) {
|
|
|
513 fclose(outfile1);
|
|
|
514 fclose(outfile2);
|
|
|
515 } else {
|
|
|
516 gzclose(outfile1_gzip);
|
|
|
517 gzclose(outfile2_gzip);
|
|
|
518 }
|
|
|
519 }
|
|
|
520
|
|
|
521 return EXIT_SUCCESS;
|
|
|
522 } /* end of paired_main() */
|
