diff --git a/bwamem.c b/bwamem.c index 4fffe38..ae1886a 100644 --- a/bwamem.c +++ b/bwamem.c @@ -554,7 +554,9 @@ void bwa_hit2sam(kstring_t *str, const int8_t mat[25], int q, int r, int w, cons } else kputw(0, str); kputc('\t', str); } else kputsn("\t*\t0\t0\t", 7, str); - if (!(p->flag&0x10)) { // print SEQ and QUAL, the forward strand + if (p->flag&0x100) { // for secondary alignments, don't write SEQ and QUAL + kputsn("*\t*", 3, str); + } else if (!(p->flag&0x10)) { // print SEQ and QUAL, the forward strand int i, qb = 0, qe = s->l_seq; if (!(p->flag&4) && is_hard) qb = p->qb, qe = p->qe; ks_resize(str, str->l + (qe - qb) + 1); @@ -610,7 +612,7 @@ void mem_alnreg2hit(const mem_alnreg_t *a, bwahit_t *h) { h->rb = a->rb; h->re = a->re; h->qb = a->qb; h->qe = a->qe; h->score = a->score; - h->sub = a->sub > a->csub? a->sub : a->csub; + h->sub = a->secondary >= 0? -1 : a->sub > a->csub? a->sub : a->csub; h->qual = 0; // quality unset h->flag = a->secondary >= 0? 0x100 : 0; // only the "secondary" bit is set } @@ -623,10 +625,10 @@ void mem_sam_se(const mem_opt_t *opt, const bntseq_t *bns, const uint8_t *pac, b if (a->n > 0) { for (k = 0; k < a->n; ++k) { bwahit_t h; - if (a->a[k].secondary >= 0) continue; + if (a->a[k].secondary >= 0 && !(opt->flag&MEM_F_ALL)) continue; mem_alnreg2hit(&a->a[k], &h); h.flag |= extra_flag; - h.qual = mem_approx_mapq_se(opt, &a->a[k]); + h.qual = a->a[k].secondary >= 0? 0 : mem_approx_mapq_se(opt, &a->a[k]); bwa_hit2sam(&str, opt->mat, opt->q, opt->r, opt->w, bns, pac, s, &h, opt->flag&MEM_F_HARDCLIP, m); } } else bwa_hit2sam(&str, opt->mat, opt->q, opt->r, opt->w, bns, pac, s, 0, opt->flag&MEM_F_HARDCLIP, m); diff --git a/bwamem.h b/bwamem.h index fa55b44..5cf3ac5 100644 --- a/bwamem.h +++ b/bwamem.h @@ -14,6 +14,7 @@ typedef struct __smem_i smem_i; #define MEM_F_HARDCLIP 0x1 #define MEM_F_PE 0x2 #define MEM_F_NOPAIRING 0x4 +#define MEM_F_ALL 0x8 typedef struct { int a, b, q, r, w; diff --git a/fastmap.c b/fastmap.c index 819e301..49c46fd 100644 --- a/fastmap.c +++ b/fastmap.c @@ -26,11 +26,12 @@ int main_mem(int argc, char *argv[]) void *ko = 0, *ko2 = 0; opt = mem_opt_init(); - while ((c = getopt(argc, argv, "CPHk:c:v:s:r:t:R:")) >= 0) { + while ((c = getopt(argc, argv, "aCPHk:c:v:s:r:t:R:")) >= 0) { if (c == 'k') opt->min_seed_len = atoi(optarg); else if (c == 't') opt->n_threads = atoi(optarg), opt->n_threads = opt->n_threads > 1? opt->n_threads : 1; else if (c == 'P') opt->flag |= MEM_F_NOPAIRING; else if (c == 'H') opt->flag |= MEM_F_HARDCLIP; + else if (c == 'a') opt->flag |= MEM_F_ALL; else if (c == 'c') opt->max_occ = atoi(optarg); else if (c == 'v') bwa_verbose = atoi(optarg); else if (c == 'r') opt->split_factor = atof(optarg); @@ -49,6 +50,9 @@ int main_mem(int argc, char *argv[]) fprintf(stderr, " -r FLOAT look for internal seeds inside a seed longer than {-k} * FLOAT [%g]\n", opt->split_factor); fprintf(stderr, " -R STR read group header line such as '@RG\tID:foo\tSM:bar' [null]\n"); fprintf(stderr, " -v INT verbose level [%d]\n", bwa_verbose); + fprintf(stderr, " -a output all alignments for SE or unpaired PE\n"); + fprintf(stderr, " -P perform mate SW only but skip pairing\n"); + fprintf(stderr, " -H hard clipping\n"); fprintf(stderr, " -C append FASTA/FASTQ comment to SAM output\n"); fprintf(stderr, "\n"); free(opt);