From ef5dd318ca295243aaaf147600f7cebf650a3b0e Mon Sep 17 00:00:00 2001
From: Heng Li <lh3@me.com>
Date: Sat, 24 Jun 2017 22:57:43 -0400
Subject: [PATCH] implemented chain splitting; NOT tested!!!

---
 align.c  | 16 ++++++++++++----
 map.c    |  2 +-
 mmpriv.h |  2 +-
 3 files changed, 14 insertions(+), 6 deletions(-)

diff --git a/align.c b/align.c
index 05af553..4d55871 100644
--- a/align.c
+++ b/align.c
@@ -224,10 +224,10 @@ static void mm_align1(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int
 	kfree(km, tseq);
 }
 
-void mm_align_skeleton(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int qlen, const char *qstr, int n_regs, mm_reg1_t *regs, mm128_t *a)
+mm_reg1_t *mm_align_skeleton(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int qlen, const char *qstr, int *n_regs_, mm_reg1_t *regs, mm128_t *a)
 {
 	extern unsigned char seq_nt4_table[256];
-	int i, reg;
+	int i, r, n_regs = *n_regs_;
 	uint8_t *qseq0[2];
 	ksw_extz_t ez;
 
@@ -239,10 +239,18 @@ void mm_align_skeleton(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int
 		qseq0[1][qlen - 1 - i] = qseq0[0][i] < 4? 3 - qseq0[0][i] : 4;
 	}
 
-	for (reg = 0; reg < n_regs; ++reg) {
+	for (r = 0; r < n_regs; ++r) {
 		mm_reg1_t r2;
-		mm_align1(km, opt, mi, qlen, qseq0, &regs[reg], &r2, a, &ez);
+		mm_align1(km, opt, mi, qlen, qseq0, &regs[r], &r2, a, &ez);
+		if (r2.cnt > 0) {
+			regs = (mm_reg1_t*)realloc(regs, n_regs + 1); // this should be very rare
+			if (r + 1 != n_regs)
+				memmove(&regs[r + 2], &regs[r + 1], sizeof(mm_reg1_t) * (n_regs - r - 1));
+			regs[r + 1] = r2;
+			++n_regs;
+		}
 	}
 
 	kfree(km, qseq0[0]); kfree(km, qseq0[1]); kfree(km, ez.cigar);
+	return regs;
 }
diff --git a/map.c b/map.c
index b01d4e7..f068c6d 100644
--- a/map.c
+++ b/map.c
@@ -291,7 +291,7 @@ mm_reg1_t *mm_map_frag(const mm_mapopt_t *opt, const mm_idx_t *mi, mm_tbuf_t *b,
 	regs = mm_gen_reg(qlen, n_u, u, a);
 	*n_regs = n_u;
 	mm_select_sub(opt->mask_level, opt->pri_ratio, n_regs, regs, b->km);
-	mm_align_skeleton(b->km, opt, mi, qlen, seq, *n_regs, regs, a);
+	regs = mm_align_skeleton(b->km, opt, mi, qlen, seq, n_regs, regs, a);
 
 	// free
 	kfree(b->km, a);
diff --git a/mmpriv.h b/mmpriv.h
index 1314b99..c5456cd 100644
--- a/mmpriv.h
+++ b/mmpriv.h
@@ -16,7 +16,7 @@ uint32_t ks_ksmall_uint32_t(size_t n, uint32_t arr[], size_t kk);
 
 int mm_chain_dp(int max_dist, int bw, int max_skip, int min_sc, int n, mm128_t *a, uint64_t **_u, void *km);
 
-void mm_align_skeleton(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int qlen, const char *qstr, int n_regs, mm_reg1_t *regs, mm128_t *a);
+mm_reg1_t *mm_align_skeleton(void *km, const mm_mapopt_t *opt, const mm_idx_t *mi, int qlen, const char *qstr, int *n_regs_, mm_reg1_t *regs, mm128_t *a);
 
 #ifdef __cplusplus
 }