From dda9466bb1eb65630a6ce7e4be956f767525070e Mon Sep 17 00:00:00 2001 From: Archieyoung <750826700@qq.com> Date: Mon, 3 Dec 2018 12:36:18 +0800 Subject: [PATCH 1/4] SA full lenght CIGAR --- format.c | 33 ++++++++++++++++++++++----------- 1 file changed, 22 insertions(+), 11 deletions(-) diff --git a/format.c b/format.c index e81ef32a..53fc2e98 100644 --- a/format.c +++ b/format.c @@ -506,18 +506,29 @@ void mm_write_sam2(kstring_t *s, const mm_idx_t *mi, const mm_bseq1_t *t, int se mm_sprintf_lite(s, "\tSA:Z:"); for (i = 0; i < n_regs; ++i) { const mm_reg1_t *q = ®s[i]; - int l_M, l_I = 0, l_D = 0, clip5 = 0, clip3 = 0; if (r == q || q->parent != q->id || q->p == 0) continue; - if (q->qe - q->qs < q->re - q->rs) l_M = q->qe - q->qs, l_D = (q->re - q->rs) - l_M; - else l_M = q->re - q->rs, l_I = (q->qe - q->qs) - l_M; - clip5 = q->rev? t->l_seq - q->qe : q->qs; - clip3 = q->rev? q->qs : t->l_seq - q->qe; - mm_sprintf_lite(s, "%s,%d,%c,", mi->seq[q->rid].name, q->rs+1, "+-"[q->rev]); - if (clip5) mm_sprintf_lite(s, "%dS", clip5); - if (l_M) mm_sprintf_lite(s, "%dM", l_M); - if (l_I) mm_sprintf_lite(s, "%dI", l_I); - if (l_D) mm_sprintf_lite(s, "%dD", l_D); - if (clip3) mm_sprintf_lite(s, "%dS", clip3); + + uint32_t k, clip_len[2]; + + clip_len[0] = q->rev? t->l_seq - q->qe : q->qs; + clip_len[1] = q->rev? q->qs : t->l_seq - q->qe; + + int clip_char = !(opt_flag&MM_F_SOFTCLIP)? 'H' : 'S'; + if (clip_len[0]) mm_sprintf_lite(s, "%d%c", clip_len[0], clip_char); + for (k = 0; k < r->p->n_cigar; ++k) + mm_sprintf_lite(s, "%d%c", r->p->cigar[k]>>4, "MIDNSHP=XB"[r->p->cigar[k]&0xf]); + if (clip_len[1]) mm_sprintf_lite(s, "%d%c", clip_len[1], clip_char); + + // if (q->qe - q->qs < q->re - q->rs) l_M = q->qe - q->qs, l_D = (q->re - q->rs) - l_M; + // else l_M = q->re - q->rs, l_I = (q->qe - q->qs) - l_M; + // clip5 = q->rev? t->l_seq - q->qe : q->qs; + // clip3 = q->rev? q->qs : t->l_seq - q->qe; + // mm_sprintf_lite(s, "%s,%d,%c,", mi->seq[q->rid].name, q->rs+1, "+-"[q->rev]); + // if (clip5) mm_sprintf_lite(s, "%dS", clip5); + // if (l_M) mm_sprintf_lite(s, "%dM", l_M); + // if (l_I) mm_sprintf_lite(s, "%dI", l_I); + // if (l_D) mm_sprintf_lite(s, "%dD", l_D); + // if (clip3) mm_sprintf_lite(s, "%dS", clip3); mm_sprintf_lite(s, ",%d,%d;", q->mapq, q->blen - q->mlen + q->p->n_ambi); } } From 70bb6ede5abbebd643a00ae847b18c8bbd216a5c Mon Sep 17 00:00:00 2001 From: Archieyoung <750826700@qq.com> Date: Mon, 3 Dec 2018 13:11:14 +0800 Subject: [PATCH 2/4] fix SA full lenght CIGAR by repace r with q --- format.c | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/format.c b/format.c index 53fc2e98..e7bd8be2 100644 --- a/format.c +++ b/format.c @@ -515,8 +515,8 @@ void mm_write_sam2(kstring_t *s, const mm_idx_t *mi, const mm_bseq1_t *t, int se int clip_char = !(opt_flag&MM_F_SOFTCLIP)? 'H' : 'S'; if (clip_len[0]) mm_sprintf_lite(s, "%d%c", clip_len[0], clip_char); - for (k = 0; k < r->p->n_cigar; ++k) - mm_sprintf_lite(s, "%d%c", r->p->cigar[k]>>4, "MIDNSHP=XB"[r->p->cigar[k]&0xf]); + for (k = 0; k < q->p->n_cigar; ++k) + mm_sprintf_lite(s, "%d%c", q->p->cigar[k]>>4, "MIDNSHP=XB"[q->p->cigar[k]&0xf]); if (clip_len[1]) mm_sprintf_lite(s, "%d%c", clip_len[1], clip_char); // if (q->qe - q->qs < q->re - q->rs) l_M = q->qe - q->qs, l_D = (q->re - q->rs) - l_M; From 4600430770c80e6c93fe0a94110a1d32034a1ee0 Mon Sep 17 00:00:00 2001 From: Archieyoung <750826700@qq.com> Date: Mon, 3 Dec 2018 13:52:41 +0800 Subject: [PATCH 3/4] fix SA full lenght CIGAR by restore the missing print of chr and pos --- format.c | 2 ++ 1 file changed, 2 insertions(+) diff --git a/format.c b/format.c index e7bd8be2..1fb514c7 100644 --- a/format.c +++ b/format.c @@ -513,6 +513,8 @@ void mm_write_sam2(kstring_t *s, const mm_idx_t *mi, const mm_bseq1_t *t, int se clip_len[0] = q->rev? t->l_seq - q->qe : q->qs; clip_len[1] = q->rev? q->qs : t->l_seq - q->qe; + mm_sprintf_lite(s, "%s,%d,%c,", mi->seq[q->rid].name, q->rs+1, "+-"[q->rev]); + int clip_char = !(opt_flag&MM_F_SOFTCLIP)? 'H' : 'S'; if (clip_len[0]) mm_sprintf_lite(s, "%d%c", clip_len[0], clip_char); for (k = 0; k < q->p->n_cigar; ++k) From 0b954ce6cba543eb45ea015398349b67991ae46e Mon Sep 17 00:00:00 2001 From: Archieyoung <750826700@qq.com> Date: Mon, 3 Dec 2018 14:34:03 +0800 Subject: [PATCH 4/4] remove comment --- format.c | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/format.c b/format.c index 1fb514c7..51cc14c8 100644 --- a/format.c +++ b/format.c @@ -521,16 +521,6 @@ void mm_write_sam2(kstring_t *s, const mm_idx_t *mi, const mm_bseq1_t *t, int se mm_sprintf_lite(s, "%d%c", q->p->cigar[k]>>4, "MIDNSHP=XB"[q->p->cigar[k]&0xf]); if (clip_len[1]) mm_sprintf_lite(s, "%d%c", clip_len[1], clip_char); - // if (q->qe - q->qs < q->re - q->rs) l_M = q->qe - q->qs, l_D = (q->re - q->rs) - l_M; - // else l_M = q->re - q->rs, l_I = (q->qe - q->qs) - l_M; - // clip5 = q->rev? t->l_seq - q->qe : q->qs; - // clip3 = q->rev? q->qs : t->l_seq - q->qe; - // mm_sprintf_lite(s, "%s,%d,%c,", mi->seq[q->rid].name, q->rs+1, "+-"[q->rev]); - // if (clip5) mm_sprintf_lite(s, "%dS", clip5); - // if (l_M) mm_sprintf_lite(s, "%dM", l_M); - // if (l_I) mm_sprintf_lite(s, "%dI", l_I); - // if (l_D) mm_sprintf_lite(s, "%dD", l_D); - // if (clip3) mm_sprintf_lite(s, "%dS", clip3); mm_sprintf_lite(s, ",%d,%d;", q->mapq, q->blen - q->mlen + q->p->n_ambi); } }