diff --git a/src/dwgsim.c b/src/dwgsim.c index 34f2077..6d2bbd5 100644 --- a/src/dwgsim.c +++ b/src/dwgsim.c @@ -498,7 +498,7 @@ void dwgsim_core(dwgsim_opt_t * opt) contigs = NULL; } - if(0 == opt->muts_only) { + if(opt->output_type != OUTPUT_TYPE_MUTS) { fprintf(stderr, "[dwgsim_core] Currently on: \n0"); } else { @@ -508,8 +508,8 @@ void dwgsim_core(dwgsim_opt_t * opt) while ((l = seq_read_fasta(opt->fp_fa, &seq, name, 0)) >= 0) { int64_t n_pairs = 0; n_ref--; - - if(1 == opt->muts_only) { + + if(opt->output_type == OUTPUT_TYPE_MUTS) { fprintf(stderr, "\r[dwgsim_core] Currently on: %s", name); if(name_len_max < strlen(name)) { name_len_max = strlen(name); @@ -616,9 +616,11 @@ void dwgsim_core(dwgsim_opt_t * opt) // generate mutations and print them out mutseq[0] = mutseq_init(); mutseq[1] = mutseq_init(); mut_diref(opt, &seq, mutseq[0], mutseq[1], contig_i, muts_input); - mut_print(name, &seq, mutseq[0], mutseq[1], opt->fp_mut, opt->fp_vcf); + if(opt->output_type != OUTPUT_TYPE_READS) { + mut_print(name, &seq, mutseq[0], mutseq[1], opt->fp_mut, opt->fp_vcf); + } - if(0 == opt->muts_only) { + if(opt->output_type != OUTPUT_TYPE_MUTS) { int num_failed = 0; for (ii = 0; ii != n_pairs; ++ii, ++ctr) { // the core loop if(0 == (ctr % 10000)) { @@ -1113,16 +1115,18 @@ int main(int argc, char *argv[]) opt->fp_fa = xopen(argv[optind+0], "r"); snprintf(fn_fai, sizeof(fn_fai), "%s.fai", argv[optind+0]); opt->fp_fai = fopen(fn_fai, "r"); // NB: depends on returning NULL; - snprintf(fn_tmp, sizeof(fn_tmp), "%s.mutations.txt", argv[optind+1]); - opt->fp_mut = xopen(fn_tmp, "w"); - snprintf(fn_tmp, sizeof(fn_tmp), "%s.mutations.vcf", argv[optind+1]); - opt->fp_vcf = xopen(fn_tmp, "w"); - if(0 == opt->muts_only) { - if (opt->output_type != OUTPUT_TYPE_BWA) { + if(opt->output_type != OUTPUT_TYPE_READS) { + snprintf(fn_tmp, sizeof(fn_tmp), "%s.mutations.txt", argv[optind+1]); + opt->fp_mut = xopen(fn_tmp, "w"); + snprintf(fn_tmp, sizeof(fn_tmp), "%s.mutations.vcf", argv[optind+1]); + opt->fp_vcf = xopen(fn_tmp, "w"); + } + if(opt->output_type != OUTPUT_TYPE_MUTS) { + if (opt->reads_output_type != READS_OUTPUT_TYPE_BWA) { snprintf(fn_tmp, sizeof(fn_tmp), "%s.bfast.fastq.gz", argv[optind+1]); opt->fp_bfast = gzopen(fn_tmp, "w"); } - if (opt->output_type != OUTPUT_TYPE_BFAST) { + if (opt->reads_output_type != READS_OUTPUT_TYPE_BFAST) { snprintf(fn_tmp, sizeof(fn_tmp), "%s.bwa.read1.fastq.gz", argv[optind+1]); opt->fp_bwa1 = gzopen(fn_tmp, "w"); snprintf(fn_tmp, sizeof(fn_tmp), "%s.bwa.read2.fastq.gz", argv[optind+1]); @@ -1134,18 +1138,20 @@ int main(int argc, char *argv[]) dwgsim_core(opt); // Close files - if(0 == opt->muts_only) { - if (opt->output_type != OUTPUT_TYPE_BWA) { - gzclose(opt->fp_bfast); + if(NULL != opt->fp_fai) fclose(opt->fp_fai); + if(opt->output_type != OUTPUT_TYPE_READS) { + fclose(opt->fp_mut); + fclose(opt->fp_vcf); + } + if(opt->output_type != OUTPUT_TYPE_MUTS) { + if (opt->reads_output_type != READS_OUTPUT_TYPE_BWA) { + gzclose(opt->fp_bfast); } - if (opt->output_type != OUTPUT_TYPE_BFAST) { - gzclose(opt->fp_bwa1); gzclose(opt->fp_bwa2); + if (opt->reads_output_type != READS_OUTPUT_TYPE_BFAST) { + gzclose(opt->fp_bwa1); gzclose(opt->fp_bwa2); } - fclose(opt->fp_fa); + fclose(opt->fp_fa); } - if(NULL != opt->fp_fai) fclose(opt->fp_fai); - fclose(opt->fp_mut); - fclose(opt->fp_vcf); dwgsim_opt_destroy(opt); diff --git a/src/dwgsim_opt.c b/src/dwgsim_opt.c index e2ff087..cf67928 100644 --- a/src/dwgsim_opt.c +++ b/src/dwgsim_opt.c @@ -63,7 +63,6 @@ dwgsim_opt_t* dwgsim_opt_init() opt->flow_order_len = 0; opt->use_base_error = 0; opt->seed = -1; - opt->muts_only = 0; opt->fixed_quality = NULL; opt->quality_std = 2.0; opt->fn_muts_input = NULL; @@ -73,6 +72,8 @@ dwgsim_opt_t* dwgsim_opt_init() opt->fp_bfast = opt->fp_bwa1 = opt->fp_bwa2 = NULL; opt->fp_fa = opt->fp_fai = NULL; opt->read_prefix = NULL; + opt->reads_output_type = READS_OUTPUT_TYPE_ALL; + opt->output_type = OUTPUT_TYPE_ALL; opt->amplicons = 0; return opt; @@ -132,7 +133,10 @@ int dwgsim_opt_usage(dwgsim_opt_t *opt) fprintf(stderr, " -B use a per-base error rate for Ion Torrent data [%s]\n", __IS_TRUE(opt->use_base_error)); fprintf(stderr, " -H haploid mode [%s]\n", __IS_TRUE(opt->is_hap)); fprintf(stderr, " -z INT random seed (-1 uses the current time) [%d]\n", opt->seed); - fprintf(stderr, " -M generate a mutations file only [%s]\n", __IS_TRUE(opt->muts_only)); + fprintf(stderr, " -M output files to generate [%d]:\n", opt->output_type); + fprintf(stderr, " 0: both reads and mutation files\n"); + fprintf(stderr, " 1: reads only\n"); + fprintf(stderr, " 2: mutations only\n"); fprintf(stderr, " -m FILE the mutations txt file to re-create [%s]\n", (MUT_INPUT_TXT != opt->fn_muts_input_type) ? "not using" : opt->fn_muts_input); fprintf(stderr, " -b FILE the bed-like file set of candidate mutations [%s]\n", (MUT_INPUT_BED == opt->fn_muts_input_type) ? "not using" : opt->fn_muts_input); fprintf(stderr, " -v FILE the vcf file set of candidate mutations (use pl tag for strand) [%s]\n", (MUT_INPUT_VCF == opt->fn_muts_input_type) ? "not using" : opt->fn_muts_input); @@ -140,7 +144,7 @@ int dwgsim_opt_usage(dwgsim_opt_t *opt) fprintf(stderr, " -P STRING a read prefix to prepend to each read name [%s]\n", (NULL == opt->read_prefix) ? "not using" : opt->read_prefix); fprintf(stderr, " -q STRING a fixed base quality to apply (single character) [%s]\n", (NULL == opt->fixed_quality) ? "not using" : opt->fixed_quality); fprintf(stderr, " -Q FLOAT standard deviation of the base quality scores [%.2lf]\n", (NULL == opt->fixed_quality) ? opt->quality_std : 0.0); - fprintf(stderr, " -o INT output type for the FASTQ files [%d]:\n", opt->output_type); + fprintf(stderr, " -o INT output type for the FASTQ files [%d]:\n", opt->reads_output_type); fprintf(stderr, " 0: interleaved (bfast) and per-read-end (bwa)\n"); fprintf(stderr, " 1: per-read-end (bwa) only\n"); fprintf(stderr, " 2: interleaved (bfast) only\n"); @@ -204,7 +208,7 @@ dwgsim_opt_parse(dwgsim_opt_t *opt, int argc, char *argv[]) int c; int muts_input_type = 0; - while ((c = getopt(argc, argv, "id:s:N:C:1:2:e:E:r:F:R:X:I:c:S:A:n:y:BHf:z:Mm:b:v:x:P:q:Q:o:ah")) >= 0) { + while ((c = getopt(argc, argv, "id:s:N:C:1:2:e:E:r:F:R:X:I:c:S:A:n:y:BHf:z:M:m:b:v:x:P:q:Q:o:ah")) >= 0) { switch (c) { case 'i': opt->is_inner = 1; break; case 'd': opt->dist = dwgsim_atoi(optarg, 'd', 0); break; @@ -237,7 +241,7 @@ dwgsim_opt_parse(dwgsim_opt_t *opt, int argc, char *argv[]) case 'H': opt->is_hap = 1; break; case 'h': return 0; case 'z': opt->seed = dwgsim_atoi(optarg, 'z', 1); break; - case 'M': opt->muts_only = 1; break; + case 'M': opt->output_type = dwgsim_atoi(optarg, 'M', 0); break; case 'm': free(opt->fn_muts_input); opt->fn_muts_input = strdup(optarg); @@ -293,7 +297,7 @@ dwgsim_opt_parse(dwgsim_opt_t *opt, int argc, char *argv[]) } break; case 'Q': opt->quality_std = atof(optarg); break; - case 'o': opt->output_type = atoi(optarg); break; + case 'o': opt->reads_output_type = atoi(optarg); break; case 'a': opt->amplicons = 1; break; default: fprintf(stderr, "Unrecognized option: -%c\n", c); return 0; } @@ -364,7 +368,7 @@ dwgsim_opt_parse(dwgsim_opt_t *opt, int argc, char *argv[]) fprintf(stderr, "Warning: remember to use the -P option with dwgsim_eval\n"); } - __check_option(opt->output_type, 0, 2, "-o"); + __check_option(opt->reads_output_type, 0, 2, "-o"); switch(muts_input_type) { case 0x0: @@ -448,7 +452,7 @@ dwgsim_opt_parse(dwgsim_opt_t *opt, int argc, char *argv[]) opt->e[1].by = (opt->e[1].end - opt->e[1].start) / opt->length[1]; } - __check_option(opt->muts_only, 0, 1, "-M"); + __check_option(opt->output_type, OUTPUT_TYPE_ALL, OUTPUT_TYPE_MUTS, "-M"); __check_option(opt->amplicons, 0, 1, "-a"); if (opt->amplicons == 1 && opt->fn_regions_bed != NULL) { diff --git a/src/dwgsim_opt.h b/src/dwgsim_opt.h index d8725b4..d9242d3 100644 --- a/src/dwgsim_opt.h +++ b/src/dwgsim_opt.h @@ -5,9 +5,13 @@ #define ERROR_RATE_NUM_RANDOM_READS 1000000 +#define READS_OUTPUT_TYPE_ALL 0 +#define READS_OUTPUT_TYPE_BWA 1 +#define READS_OUTPUT_TYPE_BFAST 2 + #define OUTPUT_TYPE_ALL 0 -#define OUTPUT_TYPE_BWA 1 -#define OUTPUT_TYPE_BFAST 2 +#define OUTPUT_TYPE_READS 1 +#define OUTPUT_TYPE_MUTS 2 typedef struct { @@ -37,7 +41,6 @@ typedef struct { int32_t use_base_error; int32_t is_hap; int32_t seed; - int32_t muts_only; char *fixed_quality; double quality_std; char *fn_muts_input; @@ -51,6 +54,7 @@ typedef struct { FILE *fp_fa; FILE *fp_fai; char *read_prefix; + int32_t reads_output_type; int32_t output_type; int32_t amplicons; } dwgsim_opt_t;