diff --git a/.vscode/extensions.json b/.vscode/extensions.json
new file mode 100644
index 0000000..8b3d6e9
--- /dev/null
+++ b/.vscode/extensions.json
@@ -0,0 +1,5 @@
+{
+    "recommendations": [
+        "mhutchie.git-graph"
+    ]
+}
\ No newline at end of file
diff --git a/.vscode/launch.json b/.vscode/launch.json
index dca4429..22af17e 100644
--- a/.vscode/launch.json
+++ b/.vscode/launch.json
@@ -10,12 +10,16 @@
             "type": "cppdbg",
             "request": "launch",
             "program": "${workspaceRoot}/minimap2",
-            "args": ["-x", 
-                    "ava-ont", 
-                    "-t", 
-                    "1", 
-                    "/public/home/zzh/work/3gseq/TGM-2021YFF/Acinetobacter_pittii.fastq", 
-                    "/public/home/zzh/work/3gseq/TGM-2021YFF/Acinetobacter_pittii.fastq"],
+            "args": [
+                "-x",
+                "ava-ont",
+                "-t",
+                "1",
+                "/public/home/zzh/work/3gseq/TGM-2021YFF/Acinetobacter_pittii.fastq",
+                "/public/home/zzh/work/3gseq/TGM-2021YFF/Acinetobacter_pittii.fastq",
+                "-o",
+                "reads.paf"
+            ],
             "cwd": "${workspaceFolder}", // 当前工作路径：当前文件所在的工作空间
         }
     ]
diff --git a/.vscode/settings.json b/.vscode/settings.json
new file mode 100644
index 0000000..8909a6f
--- /dev/null
+++ b/.vscode/settings.json
@@ -0,0 +1,6 @@
+{
+    "files.associations": {
+        "minimap.h": "c",
+        "time.h": "c"
+    }
+}
\ No newline at end of file
diff --git a/Makefile b/Makefile
index 17b13b6..fa26a31 100644
--- a/Makefile
+++ b/Makefile
@@ -1,4 +1,5 @@
-CFLAGS=		-g -Wall -O2 -Wc++-compat #-Wextra
+CFLAGS=		-g -Wall -Wc++-compat #-Wextra
+#CFLAGS=		-g -Wall -O2 -Wc++-compat #-Wextra
 CPPFLAGS=	-DHAVE_KALLOC
 INCLUDES=
 OBJS=		kthread.o kalloc.o misc.o bseq.o sketch.o sdust.o options.o index.o \
diff --git a/index.c b/index.c
index 0d8a2ae..8a3e309 100644
--- a/index.c
+++ b/index.c
@@ -15,6 +15,11 @@
 #include "kvec.h"
 #include "khash.h"
 
+#ifdef ANALYSIS_PERF
+extern int64_t get_mseconds();
+extern int64_t time_mm_idx_reader_read;
+#endif
+
 #define idx_hash(a) ((a)>>1)
 #define idx_eq(a, b) ((a)>>1 == (b)>>1)
 KHASH_INIT(idx, uint64_t, uint64_t, 1, idx_hash, idx_eq)
@@ -608,6 +613,9 @@ void mm_idx_reader_close(mm_idx_reader_t *r)
 mm_idx_t *mm_idx_reader_read(mm_idx_reader_t *r, int n_threads)
 {
 	mm_idx_t *mi;
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds();
+#endif
 	if (r->is_idx) {
 		mi = mm_idx_load(r->fp.idx);
 		if (mi && mm_verbose >= 2 && (mi->k != r->opt.k || mi->w != r->opt.w || (mi->flag&MM_I_HPC) != (r->opt.flag&MM_I_HPC)))
@@ -618,6 +626,9 @@ mm_idx_t *mm_idx_reader_read(mm_idx_reader_t *r, int n_threads)
 		if (r->fp_out) mm_idx_dump(r->fp_out, mi);
 		mi->index = r->n_parts++;
 	}
+#ifdef ANALYSIS_PERF
+	time_mm_idx_reader_read += get_mseconds() - tmp_cur_time;
+#endif
 	return mi;
 }
 
diff --git a/lchain.c b/lchain.c
index 7df5cab..ee72d68 100644
--- a/lchain.c
+++ b/lchain.c
@@ -5,7 +5,11 @@
 #include "mmpriv.h"
 #include "kalloc.h"
 #include "krmq.h"
-
+#ifdef ANALYSIS_PERF
+extern int64_t get_mseconds();
+extern int64_t time_mg_lchain_dp,
+	time_mg_chain_backtrack;
+#endif
 static int64_t mg_chain_bk_end(int32_t max_drop, const mm128_t *z, const int32_t *f, const int64_t *p, int32_t *t, int64_t k)
 {
 	int64_t i = z[k].y, end_i = -1, max_i = i;
@@ -151,7 +155,9 @@ mm128_t *mg_lchain_dp(int max_dist_x, int max_dist_y, int bw, int max_skip, int
 	int32_t *f, *t, *v, n_u, n_v, mmax_f = 0, max_drop = bw;
 	int64_t *p, i, j, max_ii, st = 0, n_iter = 0;
 	uint64_t *u;
-
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds(), tmp_diff = 0;
+#endif
 	if (_u) *_u = 0, *n_u_ = 0;
 	if (n == 0 || a == 0) {
 		kfree(km, a);
@@ -205,10 +211,20 @@ mm128_t *mg_lchain_dp(int max_dist_x, int max_dist_y, int bw, int max_skip, int
 			max_ii = i;
 		if (mmax_f < max_f) mmax_f = max_f;
 	}
-
+#ifdef ANALYSIS_PERF
+	int64_t tmp_inner_time = get_mseconds();
+#endif
 	u = mg_chain_backtrack(km, n, f, p, v, t, min_cnt, min_sc, max_drop, &n_u, &n_v);
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_inner_time;
+	__sync_fetch_and_add(&time_mg_chain_backtrack, tmp_diff);
+#endif
 	*n_u_ = n_u, *_u = u; // NB: note that u[] may not be sorted by score here
 	kfree(km, p); kfree(km, f); kfree(km, t);
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mg_lchain_dp, tmp_diff);
+#endif
 	if (n_u == 0) {
 		kfree(km, a); kfree(km, v);
 		return 0;
diff --git a/main.c b/main.c
index d50b565..6f899b6 100644
--- a/main.c
+++ b/main.c
@@ -7,6 +7,38 @@
 #include "mmpriv.h"
 #include "ketopt.h"
 
+// 用来调试，计算感兴趣部分的运行时间
+#include "sys/time.h"
+// 获取当前毫秒数
+int64_t get_mseconds()
+{
+	struct timeval tv;
+	gettimeofday(&tv, NULL);
+	return (int64_t)1000 * (tv.tv_sec + ((1e-6) * tv.tv_usec));
+}
+
+// 记录运行时间的变量
+#ifdef ANALYSIS_PERF
+
+int64_t time_mm_idx_reader_read,
+	time_mm_map_file_frag,
+	time_map_work_for_block_1,
+	time_map_work_for_block_2,
+	time_map_work_for_block_3,
+	time_mm_map_frag_b1,
+	time_mm_map_frag_b2,
+	time_mm_map_frag_b3,
+	time_mm_map_frag_b4,
+	time_mm_map_frag_b5,
+	time_mm_map_frag_b6,
+	time_mg_lchain_dp = 0,
+	time_collect_seed_hits_heap = 0,
+	time_collect_seed_hits = 0,
+	time_mg_chain_backtrack = 0;
+
+#endif
+//////////////////////////////////
+
 #ifdef __linux__
 #include <sys/resource.h>
 #include <sys/time.h>
@@ -18,85 +50,90 @@ void liftrlimit()
 	setrlimit(RLIMIT_AS, &r);
 }
 #else
-void liftrlimit() {}
+void liftrlimit()
+{
+}
 #endif
 
 static ko_longopt_t long_options[] = {
-	{ "bucket-bits",    ko_required_argument, 300 },
-	{ "mb-size",        ko_required_argument, 'K' },
-	{ "seed",           ko_required_argument, 302 },
-	{ "no-kalloc",      ko_no_argument,       303 },
-	{ "print-qname",    ko_no_argument,       304 },
-	{ "no-self",        ko_no_argument,       'D' },
-	{ "print-seeds",    ko_no_argument,       306 },
-	{ "max-chain-skip", ko_required_argument, 307 },
-	{ "min-dp-len",     ko_required_argument, 308 },
-	{ "print-aln-seq",  ko_no_argument,       309 },
-	{ "splice",         ko_no_argument,       310 },
-	{ "cost-non-gt-ag", ko_required_argument, 'C' },
-	{ "no-long-join",   ko_no_argument,       312 },
-	{ "sr",             ko_no_argument,       313 },
-	{ "frag",           ko_required_argument, 314 },
-	{ "secondary",      ko_required_argument, 315 },
-	{ "cs",             ko_optional_argument, 316 },
-	{ "end-bonus",      ko_required_argument, 317 },
-	{ "no-pairing",     ko_no_argument,       318 },
-	{ "splice-flank",   ko_required_argument, 319 },
-	{ "idx-no-seq",     ko_no_argument,       320 },
-	{ "end-seed-pen",   ko_required_argument, 321 },
-	{ "for-only",       ko_no_argument,       322 },
-	{ "rev-only",       ko_no_argument,       323 },
-	{ "heap-sort",      ko_required_argument, 324 },
-	{ "all-chain",      ko_no_argument,       'P' },
-	{ "dual",           ko_required_argument, 326 },
-	{ "max-clip-ratio", ko_required_argument, 327 },
-	{ "min-occ-floor",  ko_required_argument, 328 },
-	{ "MD",             ko_no_argument,       329 },
-	{ "lj-min-ratio",   ko_required_argument, 330 },
-	{ "score-N",        ko_required_argument, 331 },
-	{ "eqx",            ko_no_argument,       332 },
-	{ "paf-no-hit",     ko_no_argument,       333 },
-	{ "split-prefix",   ko_required_argument, 334 },
-	{ "no-end-flt",     ko_no_argument,       335 },
-	{ "hard-mask-level",ko_no_argument,       336 },
-	{ "cap-sw-mem",     ko_required_argument, 337 },
-	{ "max-qlen",       ko_required_argument, 338 },
-	{ "max-chain-iter", ko_required_argument, 339 },
-	{ "junc-bed",       ko_required_argument, 340 },
-	{ "junc-bonus",     ko_required_argument, 341 },
-	{ "sam-hit-only",   ko_no_argument,       342 },
-	{ "chain-gap-scale",ko_required_argument, 343 },
-	{ "alt",            ko_required_argument, 344 },
-	{ "alt-drop",       ko_required_argument, 345 },
-	{ "mask-len",       ko_required_argument, 346 },
-	{ "rmq",            ko_optional_argument, 347 },
-	{ "qstrand",        ko_no_argument,       348 },
-	{ "cap-kalloc",     ko_required_argument, 349 },
-	{ "q-occ-frac",     ko_required_argument, 350 },
-	{ "chain-skip-scale",ko_required_argument,351 },
-	{ "print-chains",   ko_no_argument,       352 },
-	{ "no-hash-name",   ko_no_argument,       353 },
-	{ "secondary-seq",  ko_no_argument,       354 },
-	{ "help",           ko_no_argument,       'h' },
-	{ "max-intron-len", ko_required_argument, 'G' },
-	{ "version",        ko_no_argument,       'V' },
-	{ "min-count",      ko_required_argument, 'n' },
-	{ "min-chain-score",ko_required_argument, 'm' },
-	{ "mask-level",     ko_required_argument, 'M' },
-	{ "min-dp-score",   ko_required_argument, 's' },
-	{ "sam",            ko_no_argument,       'a' },
-	{ 0, 0, 0 }
-};
+	{"bucket-bits", ko_required_argument, 300},
+	{"mb-size", ko_required_argument, 'K'},
+	{"seed", ko_required_argument, 302},
+	{"no-kalloc", ko_no_argument, 303},
+	{"print-qname", ko_no_argument, 304},
+	{"no-self", ko_no_argument, 'D'},
+	{"print-seeds", ko_no_argument, 306},
+	{"max-chain-skip", ko_required_argument, 307},
+	{"min-dp-len", ko_required_argument, 308},
+	{"print-aln-seq", ko_no_argument, 309},
+	{"splice", ko_no_argument, 310},
+	{"cost-non-gt-ag", ko_required_argument, 'C'},
+	{"no-long-join", ko_no_argument, 312},
+	{"sr", ko_no_argument, 313},
+	{"frag", ko_required_argument, 314},
+	{"secondary", ko_required_argument, 315},
+	{"cs", ko_optional_argument, 316},
+	{"end-bonus", ko_required_argument, 317},
+	{"no-pairing", ko_no_argument, 318},
+	{"splice-flank", ko_required_argument, 319},
+	{"idx-no-seq", ko_no_argument, 320},
+	{"end-seed-pen", ko_required_argument, 321},
+	{"for-only", ko_no_argument, 322},
+	{"rev-only", ko_no_argument, 323},
+	{"heap-sort", ko_required_argument, 324},
+	{"all-chain", ko_no_argument, 'P'},
+	{"dual", ko_required_argument, 326},
+	{"max-clip-ratio", ko_required_argument, 327},
+	{"min-occ-floor", ko_required_argument, 328},
+	{"MD", ko_no_argument, 329},
+	{"lj-min-ratio", ko_required_argument, 330},
+	{"score-N", ko_required_argument, 331},
+	{"eqx", ko_no_argument, 332},
+	{"paf-no-hit", ko_no_argument, 333},
+	{"split-prefix", ko_required_argument, 334},
+	{"no-end-flt", ko_no_argument, 335},
+	{"hard-mask-level", ko_no_argument, 336},
+	{"cap-sw-mem", ko_required_argument, 337},
+	{"max-qlen", ko_required_argument, 338},
+	{"max-chain-iter", ko_required_argument, 339},
+	{"junc-bed", ko_required_argument, 340},
+	{"junc-bonus", ko_required_argument, 341},
+	{"sam-hit-only", ko_no_argument, 342},
+	{"chain-gap-scale", ko_required_argument, 343},
+	{"alt", ko_required_argument, 344},
+	{"alt-drop", ko_required_argument, 345},
+	{"mask-len", ko_required_argument, 346},
+	{"rmq", ko_optional_argument, 347},
+	{"qstrand", ko_no_argument, 348},
+	{"cap-kalloc", ko_required_argument, 349},
+	{"q-occ-frac", ko_required_argument, 350},
+	{"chain-skip-scale", ko_required_argument, 351},
+	{"print-chains", ko_no_argument, 352},
+	{"no-hash-name", ko_no_argument, 353},
+	{"secondary-seq", ko_no_argument, 354},
+	{"help", ko_no_argument, 'h'},
+	{"max-intron-len", ko_required_argument, 'G'},
+	{"version", ko_no_argument, 'V'},
+	{"min-count", ko_required_argument, 'n'},
+	{"min-chain-score", ko_required_argument, 'm'},
+	{"mask-level", ko_required_argument, 'M'},
+	{"min-dp-score", ko_required_argument, 's'},
+	{"sam", ko_no_argument, 'a'},
+	{0, 0, 0}};
 
 static inline int64_t mm_parse_num2(const char *str, char **q)
 {
 	double x;
 	char *p;
 	x = strtod(str, &p);
-	if (*p == 'G' || *p == 'g') x *= 1e9, ++p;
-	else if (*p == 'M' || *p == 'm') x *= 1e6, ++p;
-	else if (*p == 'K' || *p == 'k') x *= 1e3, ++p;
-	if (q) *q = p;
+	if (*p == 'G' || *p == 'g')
+		x *= 1e9, ++p;
+	else if (*p == 'M' || *p == 'm')
+		x *= 1e6, ++p;
+	else if (*p == 'K' || *p == 'k')
+		x *= 1e3, ++p;
+	if (q)
+		*q = p;
 	return (int64_t)(x + .499);
 }
 
@@ -107,14 +144,23 @@ static inline int64_t mm_parse_num(const char *str)
 
 static inline void yes_or_no(mm_mapopt_t *opt, int64_t flag, int long_idx, const char *arg, int yes_to_set)
 {
-	if (yes_to_set) {
-		if (strcmp(arg, "yes") == 0 || strcmp(arg, "y") == 0) opt->flag |= flag;
-		else if (strcmp(arg, "no") == 0 || strcmp(arg, "n") == 0) opt->flag &= ~flag;
-		else fprintf(stderr, "[WARNING]\033[1;31m option '--%s' only accepts 'yes' or 'no'.\033[0m\n", long_options[long_idx].name);
-	} else {
-		if (strcmp(arg, "yes") == 0 || strcmp(arg, "y") == 0) opt->flag &= ~flag;
-		else if (strcmp(arg, "no") == 0 || strcmp(arg, "n") == 0) opt->flag |= flag;
-		else fprintf(stderr, "[WARNING]\033[1;31m option '--%s' only accepts 'yes' or 'no'.\033[0m\n", long_options[long_idx].name);
+	if (yes_to_set)
+	{
+		if (strcmp(arg, "yes") == 0 || strcmp(arg, "y") == 0)
+			opt->flag |= flag;
+		else if (strcmp(arg, "no") == 0 || strcmp(arg, "n") == 0)
+			opt->flag &= ~flag;
+		else
+			fprintf(stderr, "[WARNING]\033[1;31m option '--%s' only accepts 'yes' or 'no'.\033[0m\n", long_options[long_idx].name);
+	}
+	else
+	{
+		if (strcmp(arg, "yes") == 0 || strcmp(arg, "y") == 0)
+			opt->flag &= ~flag;
+		else if (strcmp(arg, "no") == 0 || strcmp(arg, "n") == 0)
+			opt->flag |= flag;
+		else
+			fprintf(stderr, "[WARNING]\033[1;31m option '--%s' only accepts 'yes' or 'no'.\033[0m\n", long_options[long_idx].name);
 	}
 }
 
@@ -135,194 +181,360 @@ int main(int argc, char *argv[])
 	mm_realtime0 = realtime();
 	mm_set_opt(0, &ipt, &opt);
 
-	while ((c = ketopt(&o, argc, argv, 1, opt_str, long_options)) >= 0) { // test command line options and apply option -x/preset first
-		if (c == 'x') {
-			if (mm_set_opt(o.arg, &ipt, &opt) < 0) {
+#ifdef ANALYSIS_PERF
+
+	time_mm_idx_reader_read = 0;
+	time_mm_map_file_frag = 0;
+	time_map_work_for_block_1 = 0;
+	time_map_work_for_block_2 = 0;
+	time_map_work_for_block_3 = 0;
+	time_mm_map_frag_b1 = 0;
+	time_mm_map_frag_b2 = 0;
+	time_mm_map_frag_b3 = 0;
+	time_mm_map_frag_b4 = 0;
+	time_mm_map_frag_b5 = 0;
+	time_mm_map_frag_b6 = 0;
+
+#endif
+	while ((c = ketopt(&o, argc, argv, 1, opt_str, long_options)) >= 0)
+	{ // test command line options and apply option -x/preset first
+		if (c == 'x')
+		{
+			if (mm_set_opt(o.arg, &ipt, &opt) < 0)
+			{
 				fprintf(stderr, "[ERROR] unknown preset '%s'\n", o.arg);
 				return 1;
 			}
-		} else if (c == ':') {
+		}
+		else if (c == ':')
+		{
 			fprintf(stderr, "[ERROR] missing option argument\n");
 			return 1;
-		} else if (c == '?') {
+		}
+		else if (c == '?')
+		{
 			fprintf(stderr, "[ERROR] unknown option in \"%s\"\n", argv[o.i - 1]);
 			return 1;
 		}
 	}
 	o = KETOPT_INIT;
 
-	while ((c = ketopt(&o, argc, argv, 1, opt_str, long_options)) >= 0) {
-		if (c == 'w') ipt.w = atoi(o.arg);
-		else if (c == 'k') ipt.k = atoi(o.arg);
-		else if (c == 'H') ipt.flag |= MM_I_HPC;
-		else if (c == 'd') fnw = o.arg; // the above are indexing related options, except -I
-		else if (c == 't') n_threads = atoi(o.arg);
-		else if (c == 'v') mm_verbose = atoi(o.arg);
-		else if (c == 'g') opt.max_gap = (int)mm_parse_num(o.arg);
-		else if (c == 'G') mm_mapopt_max_intron_len(&opt, (int)mm_parse_num(o.arg));
-		else if (c == 'F') opt.max_frag_len = (int)mm_parse_num(o.arg);
-		else if (c == 'N') old_best_n = opt.best_n, opt.best_n = atoi(o.arg);
-		else if (c == 'p') opt.pri_ratio = atof(o.arg);
-		else if (c == 'M') opt.mask_level = atof(o.arg);
-		else if (c == 'c') opt.flag |= MM_F_OUT_CG | MM_F_CIGAR;
-		else if (c == 'D') opt.flag |= MM_F_NO_DIAG;
-		else if (c == 'P') opt.flag |= MM_F_ALL_CHAINS;
-		else if (c == 'X') opt.flag |= MM_F_ALL_CHAINS | MM_F_NO_DIAG | MM_F_NO_DUAL | MM_F_NO_LJOIN; // -D -P --no-long-join --dual=no
-		else if (c == 'a') opt.flag |= MM_F_OUT_SAM | MM_F_CIGAR;
-		else if (c == 'Q') opt.flag |= MM_F_NO_QUAL;
-		else if (c == 'Y') opt.flag |= MM_F_SOFTCLIP;
-		else if (c == 'L') opt.flag |= MM_F_LONG_CIGAR;
-		else if (c == 'y') opt.flag |= MM_F_COPY_COMMENT;
-		else if (c == 'T') opt.sdust_thres = atoi(o.arg);
-		else if (c == 'n') opt.min_cnt = atoi(o.arg);
-		else if (c == 'm') opt.min_chain_score = atoi(o.arg);
-		else if (c == 'A') opt.a = atoi(o.arg);
-		else if (c == 'B') opt.b = atoi(o.arg);
-		else if (c == 'b') opt.transition = atoi(o.arg);
-		else if (c == 's') opt.min_dp_max = atoi(o.arg);
-		else if (c == 'C') opt.noncan = atoi(o.arg);
-		else if (c == 'I') ipt.batch_size = mm_parse_num(o.arg);
-		else if (c == 'K') opt.mini_batch_size = mm_parse_num(o.arg);
-		else if (c == 'e') opt.occ_dist = mm_parse_num(o.arg);
-		else if (c == 'R') rg = o.arg;
-		else if (c == 'h') fp_help = stdout;
-		else if (c == '2') opt.flag |= MM_F_2_IO_THREADS;
-		else if (c == 'J') {
+	while ((c = ketopt(&o, argc, argv, 1, opt_str, long_options)) >= 0)
+	{
+		if (c == 'w')
+			ipt.w = atoi(o.arg);
+		else if (c == 'k')
+			ipt.k = atoi(o.arg);
+		else if (c == 'H')
+			ipt.flag |= MM_I_HPC;
+		else if (c == 'd')
+			fnw = o.arg; // the above are indexing related options, except -I
+		else if (c == 't')
+			n_threads = atoi(o.arg);
+		else if (c == 'v')
+			mm_verbose = atoi(o.arg);
+		else if (c == 'g')
+			opt.max_gap = (int)mm_parse_num(o.arg);
+		else if (c == 'G')
+			mm_mapopt_max_intron_len(&opt, (int)mm_parse_num(o.arg));
+		else if (c == 'F')
+			opt.max_frag_len = (int)mm_parse_num(o.arg);
+		else if (c == 'N')
+			old_best_n = opt.best_n, opt.best_n = atoi(o.arg);
+		else if (c == 'p')
+			opt.pri_ratio = atof(o.arg);
+		else if (c == 'M')
+			opt.mask_level = atof(o.arg);
+		else if (c == 'c')
+			opt.flag |= MM_F_OUT_CG | MM_F_CIGAR;
+		else if (c == 'D')
+			opt.flag |= MM_F_NO_DIAG;
+		else if (c == 'P')
+			opt.flag |= MM_F_ALL_CHAINS;
+		else if (c == 'X')
+			opt.flag |= MM_F_ALL_CHAINS | MM_F_NO_DIAG | MM_F_NO_DUAL | MM_F_NO_LJOIN; // -D -P --no-long-join --dual=no
+		else if (c == 'a')
+			opt.flag |= MM_F_OUT_SAM | MM_F_CIGAR;
+		else if (c == 'Q')
+			opt.flag |= MM_F_NO_QUAL;
+		else if (c == 'Y')
+			opt.flag |= MM_F_SOFTCLIP;
+		else if (c == 'L')
+			opt.flag |= MM_F_LONG_CIGAR;
+		else if (c == 'y')
+			opt.flag |= MM_F_COPY_COMMENT;
+		else if (c == 'T')
+			opt.sdust_thres = atoi(o.arg);
+		else if (c == 'n')
+			opt.min_cnt = atoi(o.arg);
+		else if (c == 'm')
+			opt.min_chain_score = atoi(o.arg);
+		else if (c == 'A')
+			opt.a = atoi(o.arg);
+		else if (c == 'B')
+			opt.b = atoi(o.arg);
+		else if (c == 'b')
+			opt.transition = atoi(o.arg);
+		else if (c == 's')
+			opt.min_dp_max = atoi(o.arg);
+		else if (c == 'C')
+			opt.noncan = atoi(o.arg);
+		else if (c == 'I')
+			ipt.batch_size = mm_parse_num(o.arg);
+		else if (c == 'K')
+			opt.mini_batch_size = mm_parse_num(o.arg);
+		else if (c == 'e')
+			opt.occ_dist = mm_parse_num(o.arg);
+		else if (c == 'R')
+			rg = o.arg;
+		else if (c == 'h')
+			fp_help = stdout;
+		else if (c == '2')
+			opt.flag |= MM_F_2_IO_THREADS;
+		else if (c == 'J')
+		{
 			int t;
 			t = atoi(o.arg);
-			if (t == 0) opt.flag |= MM_F_SPLICE_OLD;
-			else if (t == 1) opt.flag &= ~MM_F_SPLICE_OLD;
-		} else if (c == 'o') {
-			if (strcmp(o.arg, "-") != 0) {
-				if (freopen(o.arg, "wb", stdout) == NULL) {
+			if (t == 0)
+				opt.flag |= MM_F_SPLICE_OLD;
+			else if (t == 1)
+				opt.flag &= ~MM_F_SPLICE_OLD;
+		}
+		else if (c == 'o')
+		{
+			if (strcmp(o.arg, "-") != 0)
+			{
+				if (freopen(o.arg, "wb", stdout) == NULL)
+				{
 					fprintf(stderr, "[ERROR]\033[1;31m failed to write the output to file '%s'\033[0m: %s\n", o.arg, strerror(errno));
 					exit(1);
 				}
 			}
 		}
-		else if (c == 300) ipt.bucket_bits = atoi(o.arg); // --bucket-bits
-		else if (c == 302) opt.seed = atoi(o.arg); // --seed
-		else if (c == 303) mm_dbg_flag |= MM_DBG_NO_KALLOC; // --no-kalloc
-		else if (c == 304) mm_dbg_flag |= MM_DBG_PRINT_QNAME; // --print-qname
-		else if (c == 306) mm_dbg_flag |= MM_DBG_PRINT_QNAME | MM_DBG_PRINT_SEED, n_threads = 1; // --print-seed
-		else if (c == 307) opt.max_chain_skip = atoi(o.arg); // --max-chain-skip
-		else if (c == 339) opt.max_chain_iter = atoi(o.arg); // --max-chain-iter
-		else if (c == 308) opt.min_ksw_len = atoi(o.arg); // --min-dp-len
-		else if (c == 309) mm_dbg_flag |= MM_DBG_PRINT_QNAME | MM_DBG_PRINT_ALN_SEQ, n_threads = 1; // --print-aln-seq
-		else if (c == 310) opt.flag |= MM_F_SPLICE; // --splice
-		else if (c == 312) opt.flag |= MM_F_NO_LJOIN; // --no-long-join
-		else if (c == 313) opt.flag |= MM_F_SR; // --sr
-		else if (c == 317) opt.end_bonus = atoi(o.arg); // --end-bonus
-		else if (c == 318) opt.flag |= MM_F_INDEPEND_SEG; // --no-pairing
-		else if (c == 320) ipt.flag |= MM_I_NO_SEQ; // --idx-no-seq
-		else if (c == 321) opt.anchor_ext_shift = atoi(o.arg); // --end-seed-pen
-		else if (c == 322) opt.flag |= MM_F_FOR_ONLY; // --for-only
-		else if (c == 323) opt.flag |= MM_F_REV_ONLY; // --rev-only
-		else if (c == 327) opt.max_clip_ratio = atof(o.arg); // --max-clip-ratio
-		else if (c == 328) opt.min_mid_occ = atoi(o.arg); // --min-occ-floor
-		else if (c == 329) opt.flag |= MM_F_OUT_MD; // --MD
-		else if (c == 331) opt.sc_ambi = atoi(o.arg); // --score-N
-		else if (c == 332) opt.flag |= MM_F_EQX; // --eqx
-		else if (c == 333) opt.flag |= MM_F_PAF_NO_HIT; // --paf-no-hit
-		else if (c == 334) opt.split_prefix = o.arg; // --split-prefix
-		else if (c == 335) opt.flag |= MM_F_NO_END_FLT; // --no-end-flt
-		else if (c == 336) opt.flag |= MM_F_HARD_MLEVEL; // --hard-mask-level
-		else if (c == 337) opt.max_sw_mat = mm_parse_num(o.arg); // --cap-sw-mat
-		else if (c == 338) opt.max_qlen = mm_parse_num(o.arg); // --max-qlen
-		else if (c == 340) junc_bed = o.arg; // --junc-bed
-		else if (c == 341) opt.junc_bonus = atoi(o.arg); // --junc-bonus
-		else if (c == 342) opt.flag |= MM_F_SAM_HIT_ONLY; // --sam-hit-only
-		else if (c == 343) opt.chain_gap_scale = atof(o.arg); // --chain-gap-scale
-		else if (c == 351) opt.chain_skip_scale = atof(o.arg); // --chain-skip-scale
-		else if (c == 344) alt_list = o.arg; // --alt
-		else if (c == 345) opt.alt_drop = atof(o.arg); // --alt-drop
-		else if (c == 346) opt.mask_len = mm_parse_num(o.arg); // --mask-len
-		else if (c == 348) opt.flag |= MM_F_QSTRAND | MM_F_NO_INV; // --qstrand
-		else if (c == 349) opt.cap_kalloc = mm_parse_num(o.arg); // --cap-kalloc
-		else if (c == 350) opt.q_occ_frac = atof(o.arg); // --q-occ-frac
-		else if (c == 352) mm_dbg_flag |= MM_DBG_PRINT_CHAIN; // --print-chains
-		else if (c == 353) opt.flag |= MM_F_NO_HASH_NAME; // --no-hash-name
-		else if (c == 354) opt.flag |= MM_F_SECONDARY_SEQ; // --secondary-seq
-		else if (c == 330) {
+		else if (c == 300)
+			ipt.bucket_bits = atoi(o.arg); // --bucket-bits
+		else if (c == 302)
+			opt.seed = atoi(o.arg); // --seed
+		else if (c == 303)
+			mm_dbg_flag |= MM_DBG_NO_KALLOC; // --no-kalloc
+		else if (c == 304)
+			mm_dbg_flag |= MM_DBG_PRINT_QNAME; // --print-qname
+		else if (c == 306)
+			mm_dbg_flag |= MM_DBG_PRINT_QNAME | MM_DBG_PRINT_SEED, n_threads = 1; // --print-seed
+		else if (c == 307)
+			opt.max_chain_skip = atoi(o.arg); // --max-chain-skip
+		else if (c == 339)
+			opt.max_chain_iter = atoi(o.arg); // --max-chain-iter
+		else if (c == 308)
+			opt.min_ksw_len = atoi(o.arg); // --min-dp-len
+		else if (c == 309)
+			mm_dbg_flag |= MM_DBG_PRINT_QNAME | MM_DBG_PRINT_ALN_SEQ, n_threads = 1; // --print-aln-seq
+		else if (c == 310)
+			opt.flag |= MM_F_SPLICE; // --splice
+		else if (c == 312)
+			opt.flag |= MM_F_NO_LJOIN; // --no-long-join
+		else if (c == 313)
+			opt.flag |= MM_F_SR; // --sr
+		else if (c == 317)
+			opt.end_bonus = atoi(o.arg); // --end-bonus
+		else if (c == 318)
+			opt.flag |= MM_F_INDEPEND_SEG; // --no-pairing
+		else if (c == 320)
+			ipt.flag |= MM_I_NO_SEQ; // --idx-no-seq
+		else if (c == 321)
+			opt.anchor_ext_shift = atoi(o.arg); // --end-seed-pen
+		else if (c == 322)
+			opt.flag |= MM_F_FOR_ONLY; // --for-only
+		else if (c == 323)
+			opt.flag |= MM_F_REV_ONLY; // --rev-only
+		else if (c == 327)
+			opt.max_clip_ratio = atof(o.arg); // --max-clip-ratio
+		else if (c == 328)
+			opt.min_mid_occ = atoi(o.arg); // --min-occ-floor
+		else if (c == 329)
+			opt.flag |= MM_F_OUT_MD; // --MD
+		else if (c == 331)
+			opt.sc_ambi = atoi(o.arg); // --score-N
+		else if (c == 332)
+			opt.flag |= MM_F_EQX; // --eqx
+		else if (c == 333)
+			opt.flag |= MM_F_PAF_NO_HIT; // --paf-no-hit
+		else if (c == 334)
+			opt.split_prefix = o.arg; // --split-prefix
+		else if (c == 335)
+			opt.flag |= MM_F_NO_END_FLT; // --no-end-flt
+		else if (c == 336)
+			opt.flag |= MM_F_HARD_MLEVEL; // --hard-mask-level
+		else if (c == 337)
+			opt.max_sw_mat = mm_parse_num(o.arg); // --cap-sw-mat
+		else if (c == 338)
+			opt.max_qlen = mm_parse_num(o.arg); // --max-qlen
+		else if (c == 340)
+			junc_bed = o.arg; // --junc-bed
+		else if (c == 341)
+			opt.junc_bonus = atoi(o.arg); // --junc-bonus
+		else if (c == 342)
+			opt.flag |= MM_F_SAM_HIT_ONLY; // --sam-hit-only
+		else if (c == 343)
+			opt.chain_gap_scale = atof(o.arg); // --chain-gap-scale
+		else if (c == 351)
+			opt.chain_skip_scale = atof(o.arg); // --chain-skip-scale
+		else if (c == 344)
+			alt_list = o.arg; // --alt
+		else if (c == 345)
+			opt.alt_drop = atof(o.arg); // --alt-drop
+		else if (c == 346)
+			opt.mask_len = mm_parse_num(o.arg); // --mask-len
+		else if (c == 348)
+			opt.flag |= MM_F_QSTRAND | MM_F_NO_INV; // --qstrand
+		else if (c == 349)
+			opt.cap_kalloc = mm_parse_num(o.arg); // --cap-kalloc
+		else if (c == 350)
+			opt.q_occ_frac = atof(o.arg); // --q-occ-frac
+		else if (c == 352)
+			mm_dbg_flag |= MM_DBG_PRINT_CHAIN; // --print-chains
+		else if (c == 353)
+			opt.flag |= MM_F_NO_HASH_NAME; // --no-hash-name
+		else if (c == 354)
+			opt.flag |= MM_F_SECONDARY_SEQ; // --secondary-seq
+		else if (c == 330)
+		{
 			fprintf(stderr, "[WARNING] \033[1;31m --lj-min-ratio has been deprecated.\033[0m\n");
-		} else if (c == 314) { // --frag
+		}
+		else if (c == 314)
+		{ // --frag
 			yes_or_no(&opt, MM_F_FRAG_MODE, o.longidx, o.arg, 1);
-		} else if (c == 315) { // --secondary
+		}
+		else if (c == 315)
+		{ // --secondary
 			yes_or_no(&opt, MM_F_NO_PRINT_2ND, o.longidx, o.arg, 0);
-		} else if (c == 316) { // --cs
+		}
+		else if (c == 316)
+		{ // --cs
 			opt.flag |= MM_F_OUT_CS | MM_F_CIGAR;
-			if (o.arg == 0 || strcmp(o.arg, "short") == 0) {
+			if (o.arg == 0 || strcmp(o.arg, "short") == 0)
+			{
 				opt.flag &= ~MM_F_OUT_CS_LONG;
-			} else if (strcmp(o.arg, "long") == 0) {
+			}
+			else if (strcmp(o.arg, "long") == 0)
+			{
 				opt.flag |= MM_F_OUT_CS_LONG;
-			} else if (strcmp(o.arg, "none") == 0) {
+			}
+			else if (strcmp(o.arg, "none") == 0)
+			{
 				opt.flag &= ~MM_F_OUT_CS;
-			} else if (mm_verbose >= 2) {
+			}
+			else if (mm_verbose >= 2)
+			{
 				fprintf(stderr, "[WARNING]\033[1;31m --cs only takes 'short' or 'long'. Invalid values are assumed to be 'short'.\033[0m\n");
 			}
-		} else if (c == 319) { // --splice-flank
+		}
+		else if (c == 319)
+		{ // --splice-flank
 			yes_or_no(&opt, MM_F_SPLICE_FLANK, o.longidx, o.arg, 1);
-		} else if (c == 324) { // --heap-sort
+		}
+		else if (c == 324)
+		{ // --heap-sort
 			yes_or_no(&opt, MM_F_HEAP_SORT, o.longidx, o.arg, 1);
-		} else if (c == 326) { // --dual
+		}
+		else if (c == 326)
+		{ // --dual
 			yes_or_no(&opt, MM_F_NO_DUAL, o.longidx, o.arg, 0);
-		} else if (c == 347) { // --rmq
-			if (o.arg) yes_or_no(&opt, MM_F_RMQ, o.longidx, o.arg, 1);
-			else opt.flag |= MM_F_RMQ;
-		} else if (c == 'S') {
+		}
+		else if (c == 347)
+		{ // --rmq
+			if (o.arg)
+				yes_or_no(&opt, MM_F_RMQ, o.longidx, o.arg, 1);
+			else
+				opt.flag |= MM_F_RMQ;
+		}
+		else if (c == 'S')
+		{
 			opt.flag |= MM_F_OUT_CS | MM_F_CIGAR | MM_F_OUT_CS_LONG;
 			if (mm_verbose >= 2)
 				fprintf(stderr, "[WARNING]\033[1;31m option -S is deprecated and may be removed in future. Please use --cs=long instead.\033[0m\n");
-		} else if (c == 'V') {
+		}
+		else if (c == 'V')
+		{
 			puts(MM_VERSION);
 			return 0;
-		} else if (c == 'r') {
+		}
+		else if (c == 'r')
+		{
 			opt.bw = (int)mm_parse_num2(o.arg, &s);
-			if (*s == ',') opt.bw_long = (int)mm_parse_num2(s + 1, &s);
-		} else if (c == 'U') {
+			if (*s == ',')
+				opt.bw_long = (int)mm_parse_num2(s + 1, &s);
+		}
+		else if (c == 'U')
+		{
 			opt.min_mid_occ = strtol(o.arg, &s, 10);
-			if (*s == ',') opt.max_mid_occ = strtol(s + 1, &s, 10);
-		} else if (c == 'f') {
+			if (*s == ',')
+				opt.max_mid_occ = strtol(s + 1, &s, 10);
+		}
+		else if (c == 'f')
+		{
 			double x;
 			char *p;
 			x = strtod(o.arg, &p);
-			if (x < 1.0) opt.mid_occ_frac = x, opt.mid_occ = 0;
-			else opt.mid_occ = (int)(x + .499);
-			if (*p == ',') opt.max_occ = (int)(strtod(p+1, &p) + .499);
-		} else if (c == 'u') {
-			if (*o.arg == 'b') opt.flag |= MM_F_SPLICE_FOR|MM_F_SPLICE_REV; // both strands
-			else if (*o.arg == 'f') opt.flag |= MM_F_SPLICE_FOR, opt.flag &= ~MM_F_SPLICE_REV; // match GT-AG
-			else if (*o.arg == 'r') opt.flag |= MM_F_SPLICE_REV, opt.flag &= ~MM_F_SPLICE_FOR; // match CT-AC (reverse complement of GT-AG)
-			else if (*o.arg == 'n') opt.flag &= ~(MM_F_SPLICE_FOR|MM_F_SPLICE_REV); // don't try to match the GT-AG signal
-			else {
+			if (x < 1.0)
+				opt.mid_occ_frac = x, opt.mid_occ = 0;
+			else
+				opt.mid_occ = (int)(x + .499);
+			if (*p == ',')
+				opt.max_occ = (int)(strtod(p + 1, &p) + .499);
+		}
+		else if (c == 'u')
+		{
+			if (*o.arg == 'b')
+				opt.flag |= MM_F_SPLICE_FOR | MM_F_SPLICE_REV; // both strands
+			else if (*o.arg == 'f')
+				opt.flag |= MM_F_SPLICE_FOR, opt.flag &= ~MM_F_SPLICE_REV; // match GT-AG
+			else if (*o.arg == 'r')
+				opt.flag |= MM_F_SPLICE_REV, opt.flag &= ~MM_F_SPLICE_FOR; // match CT-AC (reverse complement of GT-AG)
+			else if (*o.arg == 'n')
+				opt.flag &= ~(MM_F_SPLICE_FOR | MM_F_SPLICE_REV); // don't try to match the GT-AG signal
+			else
+			{
 				fprintf(stderr, "[ERROR]\033[1;31m unrecognized cDNA direction\033[0m\n");
 				return 1;
 			}
-		} else if (c == 'z') {
+		}
+		else if (c == 'z')
+		{
 			opt.zdrop = opt.zdrop_inv = strtol(o.arg, &s, 10);
-			if (*s == ',') opt.zdrop_inv = strtol(s + 1, &s, 10);
-		} else if (c == 'O') {
+			if (*s == ',')
+				opt.zdrop_inv = strtol(s + 1, &s, 10);
+		}
+		else if (c == 'O')
+		{
 			opt.q = opt.q2 = strtol(o.arg, &s, 10);
-			if (*s == ',') opt.q2 = strtol(s + 1, &s, 10);
-		} else if (c == 'E') {
+			if (*s == ',')
+				opt.q2 = strtol(s + 1, &s, 10);
+		}
+		else if (c == 'E')
+		{
 			opt.e = opt.e2 = strtol(o.arg, &s, 10);
-			if (*s == ',') opt.e2 = strtol(s + 1, &s, 10);
+			if (*s == ',')
+				opt.e2 = strtol(s + 1, &s, 10);
 		}
 	}
-	if ((opt.flag & MM_F_SPLICE) && (opt.flag & MM_F_FRAG_MODE)) {
+	if ((opt.flag & MM_F_SPLICE) && (opt.flag & MM_F_FRAG_MODE))
+	{
 		fprintf(stderr, "[ERROR]\033[1;31m --splice and --frag should not be specified at the same time.\033[0m\n");
 		return 1;
 	}
-	if (!fnw && !(opt.flag&MM_F_CIGAR))
+	if (!fnw && !(opt.flag & MM_F_CIGAR))
 		ipt.flag |= MM_I_NO_SEQ;
 	if (mm_check_opt(&ipt, &opt) < 0)
 		return 1;
-	if (opt.best_n == 0) {
+	if (opt.best_n == 0)
+	{
 		fprintf(stderr, "[WARNING]\033[1;31m changed '-N 0' to '-N %d --secondary=no'.\033[0m\n", old_best_n);
 		opt.best_n = old_best_n, opt.flag |= MM_F_NO_PRINT_2ND;
 	}
 
-	if (argc == o.ind || fp_help == stdout) {
+	if (argc == o.ind || fp_help == stdout)
+	{
 		fprintf(fp_help, "Usage: minimap2 [options] <target.fa>|<target.idx> [query.fa] [...]\n");
 		fprintf(fp_help, "Options:\n");
 		fprintf(fp_help, "  Indexing:\n");
@@ -339,7 +551,7 @@ int main(int argc, char *argv[])
 		fprintf(fp_help, "    -r NUM[,NUM] chaining/alignment bandwidth and long-join bandwidth [%d,%d]\n", opt.bw, opt.bw_long);
 		fprintf(fp_help, "    -n INT       minimal number of minimizers on a chain [%d]\n", opt.min_cnt);
 		fprintf(fp_help, "    -m INT       minimal chaining score (matching bases minus log gap penalty) [%d]\n", opt.min_chain_score);
-//		fprintf(fp_help, "    -T INT       SDUST threshold; 0 to disable SDUST [%d]\n", opt.sdust_thres); // TODO: this option is never used; might be buggy
+		//		fprintf(fp_help, "    -T INT       SDUST threshold; 0 to disable SDUST [%d]\n", opt.sdust_thres); // TODO: this option is never used; might be buggy
 		fprintf(fp_help, "    -X           skip self and dual mappings (for the all-vs-all mode)\n");
 		fprintf(fp_help, "    -p FLOAT     min secondary-to-primary score ratio [%g]\n", opt.pri_ratio);
 		fprintf(fp_help, "    -N INT       retain at most INT secondary alignments [%d]\n", opt.best_n);
@@ -364,7 +576,7 @@ int main(int argc, char *argv[])
 		fprintf(fp_help, "    -Y           use soft clipping for supplementary alignments\n");
 		fprintf(fp_help, "    -t INT       number of threads [%d]\n", n_threads);
 		fprintf(fp_help, "    -K NUM       minibatch size for mapping [500M]\n");
-//		fprintf(fp_help, "    -v INT       verbose level [%d]\n", mm_verbose);
+		//		fprintf(fp_help, "    -v INT       verbose level [%d]\n", mm_verbose);
 		fprintf(fp_help, "    --version    show version number\n");
 		fprintf(fp_help, "  Preset:\n");
 		fprintf(fp_help, "    -x STR       preset (always applied before other options; see minimap2.1 for details) []\n");
@@ -375,45 +587,55 @@ int main(int argc, char *argv[])
 		fprintf(fp_help, "                 - splice/splice:hq - long-read/Pacbio-CCS spliced alignment\n");
 		fprintf(fp_help, "                 - sr - genomic short-read mapping\n");
 		fprintf(fp_help, "\nSee `man ./minimap2.1' for detailed description of these and other advanced command-line options.\n");
-		return fp_help == stdout? 0 : 1;
+		return fp_help == stdout ? 0 : 1;
 	}
 
-	if ((opt.flag & MM_F_SR) && argc - o.ind > 3) {
+	if ((opt.flag & MM_F_SR) && argc - o.ind > 3)
+	{
 		fprintf(stderr, "[ERROR] incorrect input: in the sr mode, please specify no more than two query files.\n");
 		return 1;
 	}
 	idx_rdr = mm_idx_reader_open(argv[o.ind], &ipt, fnw);
-	if (idx_rdr == 0) {
+	if (idx_rdr == 0)
+	{
 		fprintf(stderr, "[ERROR] failed to open file '%s': %s\n", argv[o.ind], strerror(errno));
 		return 1;
 	}
-	if (!idx_rdr->is_idx && fnw == 0 && argc - o.ind < 2) {
+	if (!idx_rdr->is_idx && fnw == 0 && argc - o.ind < 2)
+	{
 		fprintf(stderr, "[ERROR] missing input: please specify a query file to map or option -d to keep the index\n");
 		mm_idx_reader_close(idx_rdr);
 		return 1;
 	}
-	if (opt.best_n == 0 && (opt.flag&MM_F_CIGAR) && mm_verbose >= 2)
+	if (opt.best_n == 0 && (opt.flag & MM_F_CIGAR) && mm_verbose >= 2)
 		fprintf(stderr, "[WARNING]\033[1;31m `-N 0' reduces alignment accuracy. Please use --secondary=no to suppress secondary alignments.\033[0m\n");
-	while ((mi = mm_idx_reader_read(idx_rdr, n_threads)) != 0) {
+	while ((mi = mm_idx_reader_read(idx_rdr, n_threads)) != 0)
+	{
 		int ret;
-		if ((opt.flag & MM_F_CIGAR) && (mi->flag & MM_I_NO_SEQ)) {
+		if ((opt.flag & MM_F_CIGAR) && (mi->flag & MM_I_NO_SEQ))
+		{
 			fprintf(stderr, "[ERROR] the prebuilt index doesn't contain sequences.\n");
 			mm_idx_destroy(mi);
 			mm_idx_reader_close(idx_rdr);
 			return 1;
 		}
-		if ((opt.flag & MM_F_OUT_SAM) && idx_rdr->n_parts == 1) {
-			if (mm_idx_reader_eof(idx_rdr)) {
+		if ((opt.flag & MM_F_OUT_SAM) && idx_rdr->n_parts == 1)
+		{
+			if (mm_idx_reader_eof(idx_rdr))
+			{
 				if (opt.split_prefix == 0)
 					ret = mm_write_sam_hdr(mi, rg, MM_VERSION, argc, argv);
 				else
 					ret = mm_write_sam_hdr(0, rg, MM_VERSION, argc, argv);
-			} else {
+			}
+			else
+			{
 				ret = mm_write_sam_hdr(0, rg, MM_VERSION, argc, argv);
 				if (opt.split_prefix == 0 && mm_verbose >= 2)
 					fprintf(stderr, "[WARNING]\033[1;31m For a multi-part index, no @SQ lines will be outputted. Please use --split-prefix.\033[0m\n");
 			}
-			if (ret != 0) {
+			if (ret != 0)
+			{
 				mm_idx_destroy(mi);
 				mm_idx_reader_close(idx_rdr);
 				return 1;
@@ -422,25 +644,37 @@ int main(int argc, char *argv[])
 		if (mm_verbose >= 3)
 			fprintf(stderr, "[M::%s::%.3f*%.2f] loaded/built the index for %d target sequence(s)\n",
 					__func__, realtime() - mm_realtime0, cputime() / (realtime() - mm_realtime0), mi->n_seq);
-		if (argc != o.ind + 1) mm_mapopt_update(&opt, mi);
-		if (mm_verbose >= 3) mm_idx_stat(mi);
-		if (junc_bed) mm_idx_bed_read(mi, junc_bed, 1);
-		if (alt_list) mm_idx_alt_read(mi, alt_list);
-		if (argc - (o.ind + 1) == 0) {
+		if (argc != o.ind + 1)
+			mm_mapopt_update(&opt, mi);
+		if (mm_verbose >= 3)
+			mm_idx_stat(mi);
+		if (junc_bed)
+			mm_idx_bed_read(mi, junc_bed, 1);
+		if (alt_list)
+			mm_idx_alt_read(mi, alt_list);
+		if (argc - (o.ind + 1) == 0)
+		{
 			mm_idx_destroy(mi);
 			continue; // no query files
 		}
 		ret = 0;
-		if (!(opt.flag & MM_F_FRAG_MODE)) {
-			for (i = o.ind + 1; i < argc; ++i) {
+		if (!(opt.flag & MM_F_FRAG_MODE))
+		{
+			// 如果是read / read overlap操作，那么把第一个输入的read.fq当做索引，后续的read比对到该索引
+			for (i = o.ind + 1; i < argc; ++i)
+			{
 				ret = mm_map_file(mi, argv[i], &opt, n_threads);
-				if (ret < 0) break;
+				if (ret < 0)
+					break;
 			}
-		} else {
-			ret = mm_map_file_frag(mi, argc - (o.ind + 1), (const char**)&argv[o.ind + 1], &opt, n_threads);
+		}
+		else
+		{
+			ret = mm_map_file_frag(mi, argc - (o.ind + 1), (const char **)&argv[o.ind + 1], &opt, n_threads);
 		}
 		mm_idx_destroy(mi);
-		if (ret < 0) {
+		if (ret < 0)
+		{
 			fprintf(stderr, "ERROR: failed to map the query file\n");
 			exit(EXIT_FAILURE);
 		}
@@ -449,19 +683,42 @@ int main(int argc, char *argv[])
 	mm_idx_reader_close(idx_rdr);
 
 	if (opt.split_prefix)
-		mm_split_merge(argc - (o.ind + 1), (const char**)&argv[o.ind + 1], &opt, n_parts);
+		mm_split_merge(argc - (o.ind + 1), (const char **)&argv[o.ind + 1], &opt, n_parts);
 
-	if (fflush(stdout) == EOF) {
+	if (fflush(stdout) == EOF)
+	{
 		perror("[ERROR] failed to write the results");
 		exit(EXIT_FAILURE);
 	}
 
-	if (mm_verbose >= 3) {
+	if (mm_verbose >= 3)
+	{
 		fprintf(stderr, "[M::%s] Version: %s\n", __func__, MM_VERSION);
 		fprintf(stderr, "[M::%s] CMD:", __func__);
 		for (i = 0; i < argc; ++i)
 			fprintf(stderr, " %s", argv[i]);
 		fprintf(stderr, "\n[M::%s] Real time: %.3f sec; CPU: %.3f sec; Peak RSS: %.3f GB\n", __func__, realtime() - mm_realtime0, cputime(), peakrss() / 1024.0 / 1024.0 / 1024.0);
 	}
+#ifdef ANALYSIS_PERF
+	fprintf(stderr, "\n");
+
+	fprintf(stderr, "time_mm_idx_reader_read: %f s\n", time_mm_idx_reader_read / 1000.0);
+	fprintf(stderr, "time_mm_map_file_frag: %f s\n", time_mm_map_file_frag / 1000.0);
+	fprintf(stderr, "time_map_work_for_block_1: %f s\n", time_map_work_for_block_1 / 1000.0 / n_threads);
+	fprintf(stderr, "time_map_work_for_block_2: %f s\n", time_map_work_for_block_2 / 1000.0 / n_threads);
+	fprintf(stderr, "time_map_work_for_block_3: %f s\n", time_map_work_for_block_3 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b1: %f s\n", time_mm_map_frag_b1 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b2: %f s\n", time_mm_map_frag_b2 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b3: %f s\n", time_mm_map_frag_b3 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b4: %f s\n", time_mm_map_frag_b4 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b5: %f s\n", time_mm_map_frag_b5 / 1000.0 / n_threads);
+	fprintf(stderr, "time_mm_map_frag_b6: %f s\n", time_mm_map_frag_b6 / 1000.0 / n_threads);
+	fprintf(stderr, "time_collect_seed_hits_heap: %f s\n", time_collect_seed_hits_heap / 1000.0 / n_threads);
+	fprintf(stderr, "time_collect_seed_hits: %f s\n", time_collect_seed_hits / 1000.0 / n_threads);
+	fprintf(stderr, "time_mg_lchain_dp: %f s\n", time_mg_lchain_dp / 1000.0 / n_threads);
+	fprintf(stderr, "time_mg_chain_backtrack: %f s\n", time_mg_chain_backtrack / 1000.0 / n_threads);
+	
+	fprintf(stderr, "\n");
+#endif
 	return 0;
 }
diff --git a/map.c b/map.c
index 038888f..e68567b 100644
--- a/map.c
+++ b/map.c
@@ -10,17 +10,37 @@
 #include "bseq.h"
 #include "khash.h"
 
+#ifdef ANALYSIS_PERF
+extern int64_t get_mseconds();
+extern int64_t time_mm_map_file_frag,
+	time_map_work_for_block_1,
+	time_map_work_for_block_2,
+	time_map_work_for_block_3,
+	time_mm_map_frag_b1,
+	time_mm_map_frag_b2,
+	time_mm_map_frag_b3,
+	time_mm_map_frag_b4,
+	time_mm_map_frag_b5,
+	time_mm_map_frag_b6,
+	time_mg_lchain_dp,
+	time_collect_seed_hits_heap,
+	time_collect_seed_hits,
+	time_mg_chain_backtrack;
+#endif
+
 mm_tbuf_t *mm_tbuf_init(void)
 {
 	mm_tbuf_t *b;
-	b = (mm_tbuf_t*)calloc(1, sizeof(mm_tbuf_t));
-	if (!(mm_dbg_flag & 1)) b->km = km_init();
+	b = (mm_tbuf_t *)calloc(1, sizeof(mm_tbuf_t));
+	if (!(mm_dbg_flag & 1))
+		b->km = km_init();
 	return b;
 }
 
 void mm_tbuf_destroy(mm_tbuf_t *b)
 {
-	if (b == 0) return;
+	if (b == 0)
+		return;
 	km_destroy(b->km);
 	free(b);
 }
@@ -35,22 +55,30 @@ static int mm_dust_minier(void *km, int n, mm128_t *a, int l_seq, const char *se
 	int n_dreg, j, k, u = 0;
 	const uint64_t *dreg;
 	sdust_buf_t *sdb;
-	if (sdust_thres <= 0) return n;
+	if (sdust_thres <= 0)
+		return n;
 	sdb = sdust_buf_init(km);
-	dreg = sdust_core((const uint8_t*)seq, l_seq, sdust_thres, 64, &n_dreg, sdb);
-	for (j = k = 0; j < n; ++j) { // squeeze out minimizers that significantly overlap with LCRs
-		int32_t qpos = (uint32_t)a[j].y>>1, span = a[j].x&0xff;
+	dreg = sdust_core((const uint8_t *)seq, l_seq, sdust_thres, 64, &n_dreg, sdb);
+	for (j = k = 0; j < n; ++j)
+	{ // squeeze out minimizers that significantly overlap with LCRs
+		int32_t qpos = (uint32_t)a[j].y >> 1, span = a[j].x & 0xff;
 		int32_t s = qpos - (span - 1), e = s + span;
-		while (u < n_dreg && (int32_t)dreg[u] <= s) ++u;
-		if (u < n_dreg && (int32_t)(dreg[u]>>32) < e) {
+		while (u < n_dreg && (int32_t)dreg[u] <= s)
+			++u;
+		if (u < n_dreg && (int32_t)(dreg[u] >> 32) < e)
+		{
 			int v, l = 0;
-			for (v = u; v < n_dreg && (int32_t)(dreg[v]>>32) < e; ++v) { // iterate over LCRs overlapping this minimizer
-				int ss = s > (int32_t)(dreg[v]>>32)? s : dreg[v]>>32;
-				int ee = e < (int32_t)dreg[v]? e : (uint32_t)dreg[v];
+			for (v = u; v < n_dreg && (int32_t)(dreg[v] >> 32) < e; ++v)
+			{ // iterate over LCRs overlapping this minimizer
+				int ss = s > (int32_t)(dreg[v] >> 32) ? s : dreg[v] >> 32;
+				int ee = e < (int32_t)dreg[v] ? e : (uint32_t)dreg[v];
 				l += ee - ss;
 			}
-			if (l <= span>>1) a[k++] = a[j]; // keep the minimizer if less than half of it falls in masked region
-		} else a[k++] = a[j];
+			if (l <= span >> 1)
+				a[k++] = a[j]; // keep the minimizer if less than half of it falls in masked region
+		}
+		else
+			a[k++] = a[j];
 	}
 	sdust_buf_destroy(sdb);
 	return k; // the new size
@@ -60,9 +88,10 @@ static void collect_minimizers(void *km, const mm_mapopt_t *opt, const mm_idx_t
 {
 	int i, n, sum = 0;
 	mv->n = 0;
-	for (i = n = 0; i < n_segs; ++i) {
+	for (i = n = 0; i < n_segs; ++i)
+	{
 		size_t j;
-		mm_sketch(km, seqs[i], qlens[i], mi->w, mi->k, i, mi->flag&MM_I_HPC, mv);
+		mm_sketch(km, seqs[i], qlens[i], mi->w, mi->k, i, mi->flag & MM_I_HPC, mv);
 		for (j = n; j < mv->n; ++j)
 			mv->a[j].y += sum << 1;
 		if (opt->sdust_thres > 0) // mask low-complexity minimizers
@@ -78,72 +107,96 @@ KSORT_INIT(heap, mm128_t, heap_lt)
 static inline int skip_seed(int flag, uint64_t r, const mm_seed_t *q, const char *qname, int qlen, const mm_idx_t *mi, int *is_self)
 {
 	*is_self = 0;
-	if (qname && (flag & (MM_F_NO_DIAG|MM_F_NO_DUAL))) {
-		const mm_idx_seq_t *s = &mi->seq[r>>32];
+	if (qname && (flag & (MM_F_NO_DIAG | MM_F_NO_DUAL)))
+	{
+		const mm_idx_seq_t *s = &mi->seq[r >> 32];
 		int cmp;
 		cmp = strcmp(qname, s->name);
-		if ((flag&MM_F_NO_DIAG) && cmp == 0 && (int)s->len == qlen) {
-			if ((uint32_t)r>>1 == (q->q_pos>>1)) return 1; // avoid the diagnonal anchors
-			if ((r&1) == (q->q_pos&1)) *is_self = 1; // this flag is used to avoid spurious extension on self chain
+		if ((flag & MM_F_NO_DIAG) && cmp == 0 && (int)s->len == qlen)
+		{
+			if ((uint32_t)r >> 1 == (q->q_pos >> 1))
+				return 1; // avoid the diagnonal anchors
+			if ((r & 1) == (q->q_pos & 1))
+				*is_self = 1; // this flag is used to avoid spurious extension on self chain
 		}
-		if ((flag&MM_F_NO_DUAL) && cmp > 0) // all-vs-all mode: map once
+		if ((flag & MM_F_NO_DUAL) && cmp > 0) // all-vs-all mode: map once
 			return 1;
 	}
-	if (flag & (MM_F_FOR_ONLY|MM_F_REV_ONLY)) {
-		if ((r&1) == (q->q_pos&1)) { // forward strand
-			if (flag & MM_F_REV_ONLY) return 1;
-		} else {
-			if (flag & MM_F_FOR_ONLY) return 1;
+	if (flag & (MM_F_FOR_ONLY | MM_F_REV_ONLY))
+	{
+		if ((r & 1) == (q->q_pos & 1))
+		{ // forward strand
+			if (flag & MM_F_REV_ONLY)
+				return 1;
+		}
+		else
+		{
+			if (flag & MM_F_FOR_ONLY)
+				return 1;
 		}
 	}
 	return 0;
 }
 
 static mm128_t *collect_seed_hits_heap(void *km, const mm_mapopt_t *opt, int max_occ, const mm_idx_t *mi, const char *qname, const mm128_v *mv, int qlen, int64_t *n_a, int *rep_len,
-								  int *n_mini_pos, uint64_t **mini_pos)
+									   int *n_mini_pos, uint64_t **mini_pos)
 {
 	int i, n_m, heap_size = 0;
 	int64_t j, n_for = 0, n_rev = 0;
 	mm_seed_t *m;
 	mm128_t *a, *heap;
-
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds(), tmp_diff = 0;
+#endif
 	m = mm_collect_matches(km, &n_m, qlen, max_occ, opt->max_max_occ, opt->occ_dist, mi, mv, n_a, rep_len, n_mini_pos, mini_pos);
 
-	heap = (mm128_t*)kmalloc(km, n_m * sizeof(mm128_t));
-	a = (mm128_t*)kmalloc(km, *n_a * sizeof(mm128_t));
+	heap = (mm128_t *)kmalloc(km, n_m * sizeof(mm128_t));
+	a = (mm128_t *)kmalloc(km, *n_a * sizeof(mm128_t));
 
-	for (i = 0, heap_size = 0; i < n_m; ++i) {
-		if (m[i].n > 0) {
+	for (i = 0, heap_size = 0; i < n_m; ++i)
+	{
+		if (m[i].n > 0)
+		{
 			heap[heap_size].x = m[i].cr[0];
-			heap[heap_size].y = (uint64_t)i<<32;
+			heap[heap_size].y = (uint64_t)i << 32;
 			++heap_size;
 		}
 	}
 	ks_heapmake_heap(heap_size, heap);
-	while (heap_size > 0) {
-		mm_seed_t *q = &m[heap->y>>32];
+	while (heap_size > 0)
+	{
+		mm_seed_t *q = &m[heap->y >> 32];
 		mm128_t *p;
 		uint64_t r = heap->x;
 		int32_t is_self, rpos = (uint32_t)r >> 1;
-		if (!skip_seed(opt->flag, r, q, qname, qlen, mi, &is_self)) {
-			if ((r&1) == (q->q_pos&1)) { // forward strand
+		if (!skip_seed(opt->flag, r, q, qname, qlen, mi, &is_self))
+		{
+			if ((r & 1) == (q->q_pos & 1))
+			{ // forward strand
 				p = &a[n_for++];
-				p->x = (r&0xffffffff00000000ULL) | rpos;
+				p->x = (r & 0xffffffff00000000ULL) | rpos;
 				p->y = (uint64_t)q->q_span << 32 | q->q_pos >> 1;
-			} else { // reverse strand
+			}
+			else
+			{ // reverse strand
 				p = &a[(*n_a) - (++n_rev)];
-				p->x = 1ULL<<63 | (r&0xffffffff00000000ULL) | rpos;
-				p->y = (uint64_t)q->q_span << 32 | (qlen - ((q->q_pos>>1) + 1 - q->q_span) - 1);
+				p->x = 1ULL << 63 | (r & 0xffffffff00000000ULL) | rpos;
+				p->y = (uint64_t)q->q_span << 32 | (qlen - ((q->q_pos >> 1) + 1 - q->q_span) - 1);
 			}
 			p->y |= (uint64_t)q->seg_id << MM_SEED_SEG_SHIFT;
-			if (q->is_tandem) p->y |= MM_SEED_TANDEM;
-			if (is_self) p->y |= MM_SEED_SELF;
+			if (q->is_tandem)
+				p->y |= MM_SEED_TANDEM;
+			if (is_self)
+				p->y |= MM_SEED_SELF;
 		}
 		// update the heap
-		if ((uint32_t)heap->y < q->n - 1) {
+		if ((uint32_t)heap->y < q->n - 1)
+		{
 			++heap[0].y;
-			heap[0].x = m[heap[0].y>>32].cr[(uint32_t)heap[0].y];
-		} else {
+			heap[0].x = m[heap[0].y >> 32].cr[(uint32_t)heap[0].y];
+		}
+		else
+		{
 			heap[0] = heap[heap_size - 1];
 			--heap_size;
 		}
@@ -153,15 +206,21 @@ static mm128_t *collect_seed_hits_heap(void *km, const mm_mapopt_t *opt, int max
 	kfree(km, heap);
 
 	// reverse anchors on the reverse strand, as they are in the descending order
-	for (j = 0; j < n_rev>>1; ++j) {
+	for (j = 0; j < n_rev >> 1; ++j)
+	{
 		mm128_t t = a[(*n_a) - 1 - j];
 		a[(*n_a) - 1 - j] = a[(*n_a) - (n_rev - j)];
 		a[(*n_a) - (n_rev - j)] = t;
 	}
-	if (*n_a > n_for + n_rev) {
+	if (*n_a > n_for + n_rev)
+	{
 		memmove(a + n_for, a + (*n_a) - n_rev, n_rev * sizeof(mm128_t));
 		*n_a = n_for + n_rev;
 	}
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_collect_seed_hits_heap, tmp_diff);
+#endif
 	return a;
 }
 
@@ -171,31 +230,48 @@ static mm128_t *collect_seed_hits(void *km, const mm_mapopt_t *opt, int max_occ,
 	int i, n_m;
 	mm_seed_t *m;
 	mm128_t *a;
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds(), tmp_diff = 0;
+#endif
 	m = mm_collect_matches(km, &n_m, qlen, max_occ, opt->max_max_occ, opt->occ_dist, mi, mv, n_a, rep_len, n_mini_pos, mini_pos);
-	a = (mm128_t*)kmalloc(km, *n_a * sizeof(mm128_t));
-	for (i = 0, *n_a = 0; i < n_m; ++i) {
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_collect_seed_hits, tmp_diff);
+#endif	
+	a = (mm128_t *)kmalloc(km, *n_a * sizeof(mm128_t));
+	for (i = 0, *n_a = 0; i < n_m; ++i)
+	{
 		mm_seed_t *q = &m[i];
 		const uint64_t *r = q->cr;
 		uint32_t k;
-		for (k = 0; k < q->n; ++k) {
+		for (k = 0; k < q->n; ++k)
+		{
 			int32_t is_self, rpos = (uint32_t)r[k] >> 1;
 			mm128_t *p;
-			if (skip_seed(opt->flag, r[k], q, qname, qlen, mi, &is_self)) continue;
+			if (skip_seed(opt->flag, r[k], q, qname, qlen, mi, &is_self))
+				continue;
 			p = &a[(*n_a)++];
-			if ((r[k]&1) == (q->q_pos&1)) { // forward strand
-				p->x = (r[k]&0xffffffff00000000ULL) | rpos;
+			if ((r[k] & 1) == (q->q_pos & 1))
+			{ // forward strand
+				p->x = (r[k] & 0xffffffff00000000ULL) | rpos;
 				p->y = (uint64_t)q->q_span << 32 | q->q_pos >> 1;
-			} else if (!(opt->flag & MM_F_QSTRAND)) { // reverse strand and not in the query-strand mode
-				p->x = 1ULL<<63 | (r[k]&0xffffffff00000000ULL) | rpos;
-				p->y = (uint64_t)q->q_span << 32 | (qlen - ((q->q_pos>>1) + 1 - q->q_span) - 1);
-			} else { // reverse strand; query-strand
-				int32_t len = mi->seq[r[k]>>32].len;
-				p->x = 1ULL<<63 | (r[k]&0xffffffff00000000ULL) | (len - (rpos + 1 - q->q_span) - 1); // coordinate only accurate for non-HPC seeds
+			}
+			else if (!(opt->flag & MM_F_QSTRAND))
+			{ // reverse strand and not in the query-strand mode
+				p->x = 1ULL << 63 | (r[k] & 0xffffffff00000000ULL) | rpos;
+				p->y = (uint64_t)q->q_span << 32 | (qlen - ((q->q_pos >> 1) + 1 - q->q_span) - 1);
+			}
+			else
+			{ // reverse strand; query-strand
+				int32_t len = mi->seq[r[k] >> 32].len;
+				p->x = 1ULL << 63 | (r[k] & 0xffffffff00000000ULL) | (len - (rpos + 1 - q->q_span) - 1); // coordinate only accurate for non-HPC seeds
 				p->y = (uint64_t)q->q_span << 32 | q->q_pos >> 1;
 			}
 			p->y |= (uint64_t)q->seg_id << MM_SEED_SEG_SHIFT;
-			if (q->is_tandem) p->y |= MM_SEED_TANDEM;
-			if (is_self) p->y |= MM_SEED_SELF;
+			if (q->is_tandem)
+				p->y |= MM_SEED_TANDEM;
+			if (is_self)
+				p->y |= MM_SEED_SELF;
 		}
 	}
 	kfree(km, m);
@@ -205,20 +281,25 @@ static mm128_t *collect_seed_hits(void *km, const mm_mapopt_t *opt, int max_occ,
 
 static void chain_post(const mm_mapopt_t *opt, int max_chain_gap_ref, const mm_idx_t *mi, void *km, int qlen, int n_segs, const int *qlens, int *n_regs, mm_reg1_t *regs, mm128_t *a)
 {
-	if (!(opt->flag & MM_F_ALL_CHAINS)) { // don't choose primary mapping(s)
-		mm_set_parent(km, opt->mask_level, opt->mask_len, *n_regs, regs, opt->a * 2 + opt->b, opt->flag&MM_F_HARD_MLEVEL, opt->alt_drop);
-		if (n_segs <= 1) mm_select_sub(km, opt->pri_ratio, mi->k*2, opt->best_n, 1, opt->max_gap * 0.8, n_regs, regs);
-		else mm_select_sub_multi(km, opt->pri_ratio, 0.2f, 0.7f, max_chain_gap_ref, mi->k*2, opt->best_n, n_segs, qlens, n_regs, regs);
+	if (!(opt->flag & MM_F_ALL_CHAINS))
+	{ // don't choose primary mapping(s)
+		mm_set_parent(km, opt->mask_level, opt->mask_len, *n_regs, regs, opt->a * 2 + opt->b, opt->flag & MM_F_HARD_MLEVEL, opt->alt_drop);
+		if (n_segs <= 1)
+			mm_select_sub(km, opt->pri_ratio, mi->k * 2, opt->best_n, 1, opt->max_gap * 0.8, n_regs, regs);
+		else
+			mm_select_sub_multi(km, opt->pri_ratio, 0.2f, 0.7f, max_chain_gap_ref, mi->k * 2, opt->best_n, n_segs, qlens, n_regs, regs);
 	}
 }
 
 static mm_reg1_t *align_regs(const mm_mapopt_t *opt, const mm_idx_t *mi, void *km, int qlen, const char *seq, int *n_regs, mm_reg1_t *regs, mm128_t *a)
 {
-	if (!(opt->flag & MM_F_CIGAR)) return regs;
+	if (!(opt->flag & MM_F_CIGAR))
+		return regs;
 	regs = mm_align_skeleton(km, opt, mi, qlen, seq, n_regs, regs, a); // this calls mm_filter_regs()
-	if (!(opt->flag & MM_F_ALL_CHAINS)) { // don't choose primary mapping(s)
-		mm_set_parent(km, opt->mask_level, opt->mask_len, *n_regs, regs, opt->a * 2 + opt->b, opt->flag&MM_F_HARD_MLEVEL, opt->alt_drop);
-		mm_select_sub(km, opt->pri_ratio, mi->k*2, opt->best_n, 0, opt->max_gap * 0.8, n_regs, regs);
+	if (!(opt->flag & MM_F_ALL_CHAINS))
+	{ // don't choose primary mapping(s)
+		mm_set_parent(km, opt->mask_level, opt->mask_len, *n_regs, regs, opt->a * 2 + opt->b, opt->flag & MM_F_HARD_MLEVEL, opt->alt_drop);
+		mm_select_sub(km, opt->pri_ratio, mi->k * 2, opt->best_n, 0, opt->max_gap * 0.8, n_regs, regs);
 		mm_set_sam_pri(*n_regs, regs);
 	}
 	return regs;
@@ -232,139 +313,203 @@ void mm_map_frag(const mm_idx_t *mi, int n_segs, const int *qlens, const char **
 	int64_t n_a;
 	uint64_t *u, *mini_pos;
 	mm128_t *a;
-	mm128_v mv = {0,0,0};
+	mm128_v mv = {0, 0, 0};
 	mm_reg1_t *regs0;
 	km_stat_t kmst;
 	float chn_pen_gap, chn_pen_skip;
-
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds(), tmp_diff = 0;
+#endif
 	for (i = 0, qlen_sum = 0; i < n_segs; ++i)
 		qlen_sum += qlens[i], n_regs[i] = 0, regs[i] = 0;
 
-	if (qlen_sum == 0 || n_segs <= 0 || n_segs > MM_MAX_SEG) return;
-	if (opt->max_qlen > 0 && qlen_sum > opt->max_qlen) return;
+	if (qlen_sum == 0 || n_segs <= 0 || n_segs > MM_MAX_SEG)
+		return;
+	if (opt->max_qlen > 0 && qlen_sum > opt->max_qlen)
+		return;
 
-	hash  = qname && !(opt->flag & MM_F_NO_HASH_NAME)? __ac_X31_hash_string(qname) : 0;
+	hash = qname && !(opt->flag & MM_F_NO_HASH_NAME) ? __ac_X31_hash_string(qname) : 0;
 	hash ^= __ac_Wang_hash(qlen_sum) + __ac_Wang_hash(opt->seed);
-	hash  = __ac_Wang_hash(hash);
+	hash = __ac_Wang_hash(hash);
 
 	collect_minimizers(b->km, opt, mi, n_segs, qlens, seqs, &mv);
-	if (opt->q_occ_frac > 0.0f) mm_seed_mz_flt(b->km, &mv, opt->mid_occ, opt->q_occ_frac);
-	if (opt->flag & MM_F_HEAP_SORT) a = collect_seed_hits_heap(b->km, opt, opt->mid_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
-	else a = collect_seed_hits(b->km, opt, opt->mid_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b1, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
+	if (opt->q_occ_frac > 0.0f)
+		mm_seed_mz_flt(b->km, &mv, opt->mid_occ, opt->q_occ_frac);
+	if (opt->flag & MM_F_HEAP_SORT)
+		a = collect_seed_hits_heap(b->km, opt, opt->mid_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
+	else
+		a = collect_seed_hits(b->km, opt, opt->mid_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
 
-	if (mm_dbg_flag & MM_DBG_PRINT_SEED) {
+	if (mm_dbg_flag & MM_DBG_PRINT_SEED)
+	{
 		fprintf(stderr, "RS\t%d\n", rep_len);
 		for (i = 0; i < n_a; ++i)
-			fprintf(stderr, "SD\t%s\t%d\t%c\t%d\t%d\t%d\n", mi->seq[a[i].x<<1>>33].name, (int32_t)a[i].x, "+-"[a[i].x>>63], (int32_t)a[i].y, (int32_t)(a[i].y>>32&0xff),
-					i == 0? 0 : ((int32_t)a[i].y - (int32_t)a[i-1].y) - ((int32_t)a[i].x - (int32_t)a[i-1].x));
+			fprintf(stderr, "SD\t%s\t%d\t%c\t%d\t%d\t%d\n", mi->seq[a[i].x << 1 >> 33].name, (int32_t)a[i].x, "+-"[a[i].x >> 63], (int32_t)a[i].y, (int32_t)(a[i].y >> 32 & 0xff),
+					i == 0 ? 0 : ((int32_t)a[i].y - (int32_t)a[i - 1].y) - ((int32_t)a[i].x - (int32_t)a[i - 1].x));
 	}
-
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b2, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
 	// set max chaining gap on the query and the reference sequence
 	if (is_sr)
-		max_chain_gap_qry = qlen_sum > opt->max_gap? qlen_sum : opt->max_gap;
-	else max_chain_gap_qry = opt->max_gap;
-	if (opt->max_gap_ref > 0) {
+		max_chain_gap_qry = qlen_sum > opt->max_gap ? qlen_sum : opt->max_gap;
+	else
+		max_chain_gap_qry = opt->max_gap;
+	if (opt->max_gap_ref > 0)
+	{
 		max_chain_gap_ref = opt->max_gap_ref; // always honor mm_mapopt_t::max_gap_ref if set
-	} else if (opt->max_frag_len > 0) {
+	}
+	else if (opt->max_frag_len > 0)
+	{
 		max_chain_gap_ref = opt->max_frag_len - qlen_sum;
-		if (max_chain_gap_ref < opt->max_gap) max_chain_gap_ref = opt->max_gap;
-	} else max_chain_gap_ref = opt->max_gap;
+		if (max_chain_gap_ref < opt->max_gap)
+			max_chain_gap_ref = opt->max_gap;
+	}
+	else
+		max_chain_gap_ref = opt->max_gap;
 
-	chn_pen_gap  = opt->chain_gap_scale * 0.01 * mi->k;
+	chn_pen_gap = opt->chain_gap_scale * 0.01 * mi->k;
 	chn_pen_skip = opt->chain_skip_scale * 0.01 * mi->k;
-	if (opt->flag & MM_F_RMQ) {
+	if (opt->flag & MM_F_RMQ)
+	{
 		a = mg_lchain_rmq(opt->max_gap, opt->rmq_inner_dist, opt->bw, opt->max_chain_skip, opt->rmq_size_cap, opt->min_cnt, opt->min_chain_score,
 						  chn_pen_gap, chn_pen_skip, n_a, a, &n_regs0, &u, b->km);
-	} else {
+	}
+	else
+	{
 		a = mg_lchain_dp(max_chain_gap_ref, max_chain_gap_qry, opt->bw, opt->max_chain_skip, opt->max_chain_iter, opt->min_cnt, opt->min_chain_score,
 						 chn_pen_gap, chn_pen_skip, is_splice, n_segs, n_a, a, &n_regs0, &u, b->km);
 	}
-
-	if (opt->bw_long > opt->bw && (opt->flag & (MM_F_SPLICE|MM_F_SR|MM_F_NO_LJOIN)) == 0 && n_segs == 1 && n_regs0 > 1) { // re-chain/long-join for long sequences
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b3, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
+	if (opt->bw_long > opt->bw && (opt->flag & (MM_F_SPLICE | MM_F_SR | MM_F_NO_LJOIN)) == 0 && n_segs == 1 && n_regs0 > 1)
+	{ // re-chain/long-join for long sequences
 		int32_t st = (int32_t)a[0].y, en = (int32_t)a[(int32_t)u[0] - 1].y;
-		if (qlen_sum - (en - st) > opt->rmq_rescue_size || en - st > qlen_sum * opt->rmq_rescue_ratio) {
+		if (qlen_sum - (en - st) > opt->rmq_rescue_size || en - st > qlen_sum * opt->rmq_rescue_ratio)
+		{
 			int32_t i;
-			for (i = 0, n_a = 0; i < n_regs0; ++i) n_a += (int32_t)u[i];
+			for (i = 0, n_a = 0; i < n_regs0; ++i)
+				n_a += (int32_t)u[i];
 			kfree(b->km, u);
 			radix_sort_128x(a, a + n_a);
 			a = mg_lchain_rmq(opt->max_gap, opt->rmq_inner_dist, opt->bw_long, opt->max_chain_skip, opt->rmq_size_cap, opt->min_cnt, opt->min_chain_score,
 							  chn_pen_gap, chn_pen_skip, n_a, a, &n_regs0, &u, b->km);
 		}
-	} else if (opt->max_occ > opt->mid_occ && rep_len > 0 && !(opt->flag & MM_F_RMQ)) { // re-chain, mostly for short reads
+	}
+	else if (opt->max_occ > opt->mid_occ && rep_len > 0 && !(opt->flag & MM_F_RMQ))
+	{ // re-chain, mostly for short reads
 		int rechain = 0;
-		if (n_regs0 > 0) { // test if the best chain has all the segments
+		if (n_regs0 > 0)
+		{ // test if the best chain has all the segments
 			int n_chained_segs = 1, max = 0, max_i = -1, max_off = -1, off = 0;
-			for (i = 0; i < n_regs0; ++i) { // find the best chain
-				if (max < (int)(u[i]>>32)) max = u[i]>>32, max_i = i, max_off = off;
+			for (i = 0; i < n_regs0; ++i)
+			{ // find the best chain
+				if (max < (int)(u[i] >> 32))
+					max = u[i] >> 32, max_i = i, max_off = off;
 				off += (uint32_t)u[i];
 			}
 			for (i = 1; i < (int32_t)u[max_i]; ++i) // count the number of segments in the best chain
-				if ((a[max_off+i].y&MM_SEED_SEG_MASK) != (a[max_off+i-1].y&MM_SEED_SEG_MASK))
+				if ((a[max_off + i].y & MM_SEED_SEG_MASK) != (a[max_off + i - 1].y & MM_SEED_SEG_MASK))
 					++n_chained_segs;
 			if (n_chained_segs < n_segs)
 				rechain = 1;
-		} else rechain = 1;
-		if (rechain) { // redo chaining with a higher max_occ threshold
+		}
+		else
+			rechain = 1;
+		if (rechain)
+		{ // redo chaining with a higher max_occ threshold
 			kfree(b->km, a);
 			kfree(b->km, u);
 			kfree(b->km, mini_pos);
-			if (opt->flag & MM_F_HEAP_SORT) a = collect_seed_hits_heap(b->km, opt, opt->max_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
-			else a = collect_seed_hits(b->km, opt, opt->max_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
+			if (opt->flag & MM_F_HEAP_SORT)
+				a = collect_seed_hits_heap(b->km, opt, opt->max_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
+			else
+				a = collect_seed_hits(b->km, opt, opt->max_occ, mi, qname, &mv, qlen_sum, &n_a, &rep_len, &n_mini_pos, &mini_pos);
 			a = mg_lchain_dp(max_chain_gap_ref, max_chain_gap_qry, opt->bw, opt->max_chain_skip, opt->max_chain_iter, opt->min_cnt, opt->min_chain_score,
 							 chn_pen_gap, chn_pen_skip, is_splice, n_segs, n_a, a, &n_regs0, &u, b->km);
 		}
 	}
 	b->frag_gap = max_chain_gap_ref;
 	b->rep_len = rep_len;
-
-	regs0 = mm_gen_regs(b->km, hash, qlen_sum, n_regs0, u, a, !!(opt->flag&MM_F_QSTRAND));
-	if (mi->n_alt) {
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b4, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
+	regs0 = mm_gen_regs(b->km, hash, qlen_sum, n_regs0, u, a, !!(opt->flag & MM_F_QSTRAND));
+	if (mi->n_alt)
+	{
 		mm_mark_alt(mi, n_regs0, regs0);
 		mm_hit_sort(b->km, &n_regs0, regs0, opt->alt_drop); // this step can be merged into mm_gen_regs(); will do if this shows up in profile
 	}
 
-	if (mm_dbg_flag & (MM_DBG_PRINT_SEED|MM_DBG_PRINT_CHAIN))
+	if (mm_dbg_flag & (MM_DBG_PRINT_SEED | MM_DBG_PRINT_CHAIN))
 		for (j = 0; j < n_regs0; ++j)
 			for (i = regs0[j].as; i < regs0[j].as + regs0[j].cnt; ++i)
-				fprintf(stderr, "CN\t%d\t%s\t%d\t%c\t%d\t%d\t%d\n", j, mi->seq[a[i].x<<1>>33].name, (int32_t)a[i].x, "+-"[a[i].x>>63], (int32_t)a[i].y, (int32_t)(a[i].y>>32&0xff),
-						i == regs0[j].as? 0 : ((int32_t)a[i].y - (int32_t)a[i-1].y) - ((int32_t)a[i].x - (int32_t)a[i-1].x));
+				fprintf(stderr, "CN\t%d\t%s\t%d\t%c\t%d\t%d\t%d\n", j, mi->seq[a[i].x << 1 >> 33].name, (int32_t)a[i].x, "+-"[a[i].x >> 63], (int32_t)a[i].y, (int32_t)(a[i].y >> 32 & 0xff),
+						i == regs0[j].as ? 0 : ((int32_t)a[i].y - (int32_t)a[i - 1].y) - ((int32_t)a[i].x - (int32_t)a[i - 1].x));
 
 	chain_post(opt, max_chain_gap_ref, mi, b->km, qlen_sum, n_segs, qlens, &n_regs0, regs0, a);
-	if (!is_sr && !(opt->flag&MM_F_QSTRAND)) {
+	if (!is_sr && !(opt->flag & MM_F_QSTRAND))
+	{
 		mm_est_err(mi, qlen_sum, n_regs0, regs0, a, n_mini_pos, mini_pos);
 		n_regs0 = mm_filter_strand_retained(n_regs0, regs0);
 	}
-
-	if (n_segs == 1) { // uni-segment
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b5, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
+	if (n_segs == 1)
+	{ // uni-segment
 		regs0 = align_regs(opt, mi, b->km, qlens[0], seqs[0], &n_regs0, regs0, a);
-		regs0 = (mm_reg1_t*)realloc(regs0, sizeof(*regs0) * n_regs0);
+		regs0 = (mm_reg1_t *)realloc(regs0, sizeof(*regs0) * n_regs0);
 		mm_set_mapq(b->km, n_regs0, regs0, opt->min_chain_score, opt->a, rep_len, is_sr);
 		n_regs[0] = n_regs0, regs[0] = regs0;
-	} else { // multi-segment
+	}
+	else
+	{ // multi-segment
 		mm_seg_t *seg;
 		seg = mm_seg_gen(b->km, hash, n_segs, qlens, n_regs0, regs0, n_regs, regs, a); // split fragment chain to separate segment chains
 		free(regs0);
-		for (i = 0; i < n_segs; ++i) {
-			mm_set_parent(b->km, opt->mask_level, opt->mask_len, n_regs[i], regs[i], opt->a * 2 + opt->b, opt->flag&MM_F_HARD_MLEVEL, opt->alt_drop); // update mm_reg1_t::parent
+		for (i = 0; i < n_segs; ++i)
+		{
+			mm_set_parent(b->km, opt->mask_level, opt->mask_len, n_regs[i], regs[i], opt->a * 2 + opt->b, opt->flag & MM_F_HARD_MLEVEL, opt->alt_drop); // update mm_reg1_t::parent
 			regs[i] = align_regs(opt, mi, b->km, qlens[i], seqs[i], &n_regs[i], regs[i], seg[i].a);
 			mm_set_mapq(b->km, n_regs[i], regs[i], opt->min_chain_score, opt->a, rep_len, is_sr);
 		}
 		mm_seg_free(b->km, n_segs, seg);
-		if (n_segs == 2 && opt->pe_ori >= 0 && (opt->flag&MM_F_CIGAR))
+		if (n_segs == 2 && opt->pe_ori >= 0 && (opt->flag & MM_F_CIGAR))
 			mm_pair(b->km, max_chain_gap_ref, opt->pe_bonus, opt->a * 2 + opt->b, opt->a, qlens, n_regs, regs); // pairing
 	}
-
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_mm_map_frag_b6, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
 	kfree(b->km, mv.a);
 	kfree(b->km, a);
 	kfree(b->km, u);
 	kfree(b->km, mini_pos);
 
-	if (b->km) {
+	if (b->km)
+	{
 		km_stat(b->km, &kmst);
 		if (mm_dbg_flag & MM_DBG_PRINT_QNAME)
 			fprintf(stderr, "QM\t%s\t%d\tcap=%ld,nCore=%ld,largest=%ld\n", qname, qlen_sum, kmst.capacity, kmst.n_cores, kmst.largest);
 		assert(kmst.n_blocks == kmst.n_cores); // otherwise, there is a memory leak
-		if (kmst.largest > 1U<<28 || (opt->cap_kalloc > 0 && kmst.capacity > opt->cap_kalloc)) {
+		if (kmst.largest > 1U << 28 || (opt->cap_kalloc > 0 && kmst.capacity > opt->cap_kalloc))
+		{
 			if (mm_dbg_flag & MM_DBG_PRINT_QNAME)
 				fprintf(stderr, "[W::%s] reset thread-local memory after read %s\n", __func__, qname);
 			km_destroy(b->km);
@@ -384,7 +529,8 @@ mm_reg1_t *mm_map(const mm_idx_t *mi, int qlen, const char *seq, int *n_regs, mm
  * Multi-threaded mapping *
  **************************/
 
-typedef struct {
+typedef struct
+{
 	int n_processed, n_threads, n_fp;
 	int64_t mini_batch_size;
 	const mm_mapopt_t *opt;
@@ -397,9 +543,10 @@ typedef struct {
 	FILE *fp_split, **fp_parts;
 } pipeline_t;
 
-typedef struct {
+typedef struct
+{
 	const pipeline_t *p;
-    int n_seq, n_frag;
+	int n_seq, n_frag;
 	mm_bseq1_t *seq;
 	int *n_reg, *seg_off, *n_seg, *rep_len, *frag_gap;
 	mm_reg1_t **reg;
@@ -408,40 +555,62 @@ typedef struct {
 
 static void worker_for(void *_data, long i, int tid) // kt_for() callback
 {
-    step_t *s = (step_t*)_data;
+	step_t *s = (step_t *)_data;
 	int qlens[MM_MAX_SEG], j, off = s->seg_off[i], pe_ori = s->p->opt->pe_ori;
 	const char *qseqs[MM_MAX_SEG];
 	double t = 0.0;
 	mm_tbuf_t *b = s->buf[tid];
 	assert(s->n_seg[i] <= MM_MAX_SEG);
-	if (mm_dbg_flag & MM_DBG_PRINT_QNAME) {
+	if (mm_dbg_flag & MM_DBG_PRINT_QNAME)
+	{
 		fprintf(stderr, "QR\t%s\t%d\t%d\n", s->seq[off].name, tid, s->seq[off].l_seq);
 		t = realtime();
 	}
-	for (j = 0; j < s->n_seg[i]; ++j) {
-		if (s->n_seg[i] == 2 && ((j == 0 && (pe_ori>>1&1)) || (j == 1 && (pe_ori&1))))
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds(), tmp_diff = 0;
+#endif
+	for (j = 0; j < s->n_seg[i]; ++j)
+	{
+		if (s->n_seg[i] == 2 && ((j == 0 && (pe_ori >> 1 & 1)) || (j == 1 && (pe_ori & 1))))
 			mm_revcomp_bseq(&s->seq[off + j]);
 		qlens[j] = s->seq[off + j].l_seq;
 		qseqs[j] = s->seq[off + j].seq;
 	}
-	if (s->p->opt->flag & MM_F_INDEPEND_SEG) {
-		for (j = 0; j < s->n_seg[i]; ++j) {
-			mm_map_frag(s->p->mi, 1, &qlens[j], &qseqs[j], &s->n_reg[off+j], &s->reg[off+j], b, s->p->opt, s->seq[off+j].name);
-			s->rep_len[off + j] = b->rep_len;
-			s->frag_gap[off + j] = b->frag_gap;
-		}
-	} else {
-		mm_map_frag(s->p->mi, s->n_seg[i], qlens, qseqs, &s->n_reg[off], &s->reg[off], b, s->p->opt, s->seq[off].name);
-		for (j = 0; j < s->n_seg[i]; ++j) {
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_map_work_for_block_1, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
+	if (s->p->opt->flag & MM_F_INDEPEND_SEG)
+	{
+		for (j = 0; j < s->n_seg[i]; ++j)
+		{
+			mm_map_frag(s->p->mi, 1, &qlens[j], &qseqs[j], &s->n_reg[off + j], &s->reg[off + j], b, s->p->opt, s->seq[off + j].name);
 			s->rep_len[off + j] = b->rep_len;
 			s->frag_gap[off + j] = b->frag_gap;
 		}
 	}
+	else
+	{
+		mm_map_frag(s->p->mi, s->n_seg[i], qlens, qseqs, &s->n_reg[off], &s->reg[off], b, s->p->opt, s->seq[off].name);
+		for (j = 0; j < s->n_seg[i]; ++j)
+		{
+			s->rep_len[off + j] = b->rep_len;
+			s->frag_gap[off + j] = b->frag_gap;
+		}
+	}
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_map_work_for_block_2, tmp_diff);
+	tmp_cur_time = get_mseconds();
+#endif
 	for (j = 0; j < s->n_seg[i]; ++j) // flip the query strand and coordinate to the original read strand
-		if (s->n_seg[i] == 2 && ((j == 0 && (pe_ori>>1&1)) || (j == 1 && (pe_ori&1)))) {
+		if (s->n_seg[i] == 2 && ((j == 0 && (pe_ori >> 1 & 1)) || (j == 1 && (pe_ori & 1))))
+		{
 			int k, t;
 			mm_revcomp_bseq(&s->seq[off + j]);
-			for (k = 0; k < s->n_reg[off + j]; ++k) {
+			for (k = 0; k < s->n_reg[off + j]; ++k)
+			{
 				mm_reg1_t *r = &s->reg[off + j][k];
 				t = r->qs;
 				r->qs = qlens[j] - r->qe;
@@ -449,6 +618,10 @@ static void worker_for(void *_data, long i, int tid) // kt_for() callback
 				r->rev = !r->rev;
 			}
 		}
+#ifdef ANALYSIS_PERF
+	tmp_diff = get_mseconds() - tmp_cur_time;
+	__sync_fetch_and_add(&time_map_work_for_block_3, tmp_diff);
+#endif
 	if (mm_dbg_flag & MM_DBG_PRINT_QNAME)
 		fprintf(stderr, "QT\t%s\t%d\t%.6f\n", s->seq[off].name, tid, realtime() - t);
 }
@@ -462,56 +635,65 @@ static void merge_hits(step_t *s)
 
 	km = km_init();
 	for (f = 0; f < s->n_frag; ++f)
-		max_seg = max_seg > s->n_seg[f]? max_seg : s->n_seg[f];
+		max_seg = max_seg > s->n_seg[f] ? max_seg : s->n_seg[f];
 	qlens = CALLOC(int, max_seg + s->p->n_parts * 3);
 	n_reg_part = qlens + max_seg;
 	rep_len_part = n_reg_part + s->p->n_parts;
 	frag_gap_part = rep_len_part + s->p->n_parts;
-	for (f = 0, k = k0 = 0; f < s->n_frag; ++f) {
+	for (f = 0, k = k0 = 0; f < s->n_frag; ++f)
+	{
 		k0 = k;
-		for (i = 0; i < s->n_seg[f]; ++i, ++k) {
+		for (i = 0; i < s->n_seg[f]; ++i, ++k)
+		{
 			int j, l, t, rep_len = 0;
 			qlens[i] = s->seq[k].l_seq;
-			for (j = 0, s->n_reg[k] = 0; j < s->p->n_parts; ++j) {
-				mm_err_fread(&n_reg_part[j],    sizeof(int), 1, fp[j]);
-				mm_err_fread(&rep_len_part[j],  sizeof(int), 1, fp[j]);
+			for (j = 0, s->n_reg[k] = 0; j < s->p->n_parts; ++j)
+			{
+				mm_err_fread(&n_reg_part[j], sizeof(int), 1, fp[j]);
+				mm_err_fread(&rep_len_part[j], sizeof(int), 1, fp[j]);
 				mm_err_fread(&frag_gap_part[j], sizeof(int), 1, fp[j]);
 				s->n_reg[k] += n_reg_part[j];
 				if (rep_len < rep_len_part[j])
 					rep_len = rep_len_part[j];
 			}
 			s->reg[k] = CALLOC(mm_reg1_t, s->n_reg[k]);
-			for (j = 0, l = 0; j < s->p->n_parts; ++j) {
-				for (t = 0; t < n_reg_part[j]; ++t, ++l) {
+			for (j = 0, l = 0; j < s->p->n_parts; ++j)
+			{
+				for (t = 0; t < n_reg_part[j]; ++t, ++l)
+				{
 					mm_reg1_t *r = &s->reg[k][l];
 					uint32_t capacity;
 					mm_err_fread(r, sizeof(mm_reg1_t), 1, fp[j]);
 					r->rid += s->p->rid_shift[j];
-					if (opt->flag & MM_F_CIGAR) {
+					if (opt->flag & MM_F_CIGAR)
+					{
 						mm_err_fread(&capacity, 4, 1, fp[j]);
-						r->p = (mm_extra_t*)calloc(capacity, 4);
+						r->p = (mm_extra_t *)calloc(capacity, 4);
 						r->p->capacity = capacity;
 						mm_err_fread(r->p, r->p->capacity, 4, fp[j]);
 					}
 				}
 			}
-			if (!(opt->flag&MM_F_SR) && s->seq[k].l_seq >= opt->rank_min_len)
+			if (!(opt->flag & MM_F_SR) && s->seq[k].l_seq >= opt->rank_min_len)
 				mm_update_dp_max(s->seq[k].l_seq, s->n_reg[k], s->reg[k], opt->rank_frac, opt->a, opt->b);
-			for (j = 0; j < s->n_reg[k]; ++j) {
+			for (j = 0; j < s->n_reg[k]; ++j)
+			{
 				mm_reg1_t *r = &s->reg[k][j];
-				if (r->p) r->p->dp_max2 = 0; // reset ->dp_max2 as mm_set_parent() doesn't clear it; necessary with mm_update_dp_max()
-				r->subsc = 0; // this may not be necessary
-				r->n_sub = 0; // n_sub will be an underestimate as we don't see all the chains now, but it can't be accurate anyway
+				if (r->p)
+					r->p->dp_max2 = 0; // reset ->dp_max2 as mm_set_parent() doesn't clear it; necessary with mm_update_dp_max()
+				r->subsc = 0;		   // this may not be necessary
+				r->n_sub = 0;		   // n_sub will be an underestimate as we don't see all the chains now, but it can't be accurate anyway
 			}
 			mm_hit_sort(km, &s->n_reg[k], s->reg[k], opt->alt_drop);
-			mm_set_parent(km, opt->mask_level, opt->mask_len, s->n_reg[k], s->reg[k], opt->a * 2 + opt->b, opt->flag&MM_F_HARD_MLEVEL, opt->alt_drop);
-			if (!(opt->flag & MM_F_ALL_CHAINS)) {
-				mm_select_sub(km, opt->pri_ratio, s->p->mi->k*2, opt->best_n, 0, opt->max_gap * 0.8, &s->n_reg[k], s->reg[k]);
+			mm_set_parent(km, opt->mask_level, opt->mask_len, s->n_reg[k], s->reg[k], opt->a * 2 + opt->b, opt->flag & MM_F_HARD_MLEVEL, opt->alt_drop);
+			if (!(opt->flag & MM_F_ALL_CHAINS))
+			{
+				mm_select_sub(km, opt->pri_ratio, s->p->mi->k * 2, opt->best_n, 0, opt->max_gap * 0.8, &s->n_reg[k], s->reg[k]);
 				mm_set_sam_pri(s->n_reg[k], s->reg[k]);
 			}
 			mm_set_mapq(km, s->n_reg[k], s->reg[k], opt->min_chain_score, opt->a, rep_len, !!(opt->flag & MM_F_SR));
 		}
-		if (s->n_seg[f] == 2 && opt->pe_ori >= 0 && (opt->flag&MM_F_CIGAR))
+		if (s->n_seg[f] == 2 && opt->pe_ori >= 0 && (opt->flag & MM_F_CIGAR))
 			mm_pair(km, frag_gap_part[0], opt->pe_bonus, opt->a * 2 + opt->b, opt->a, qlens, &s->n_reg[k0], &s->reg[k0]);
 	}
 	free(qlens);
@@ -521,107 +703,141 @@ static void merge_hits(step_t *s)
 static void *worker_pipeline(void *shared, int step, void *in)
 {
 	int i, j, k;
-    pipeline_t *p = (pipeline_t*)shared;
-    if (step == 0) { // step 0: read sequences
+	pipeline_t *p = (pipeline_t *)shared;
+	if (step == 0)
+	{ // step 0: read sequences
 		int with_qual = (!!(p->opt->flag & MM_F_OUT_SAM) && !(p->opt->flag & MM_F_NO_QUAL));
 		int with_comment = !!(p->opt->flag & MM_F_COPY_COMMENT);
 		int frag_mode = (p->n_fp > 1 || !!(p->opt->flag & MM_F_FRAG_MODE));
-        step_t *s;
-        s = (step_t*)calloc(1, sizeof(step_t));
-		if (p->n_fp > 1) s->seq = mm_bseq_read_frag2(p->n_fp, p->fp, p->mini_batch_size, with_qual, with_comment, &s->n_seq);
-		else s->seq = mm_bseq_read3(p->fp[0], p->mini_batch_size, with_qual, with_comment, frag_mode, &s->n_seq);
-		if (s->seq) {
+		step_t *s;
+		s = (step_t *)calloc(1, sizeof(step_t));
+		if (p->n_fp > 1)
+			s->seq = mm_bseq_read_frag2(p->n_fp, p->fp, p->mini_batch_size, with_qual, with_comment, &s->n_seq);
+		else
+			s->seq = mm_bseq_read3(p->fp[0], p->mini_batch_size, with_qual, with_comment, frag_mode, &s->n_seq);
+		if (s->seq)
+		{
 			s->p = p;
 			for (i = 0; i < s->n_seq; ++i)
 				s->seq[i].rid = p->n_processed++;
-			s->buf = (mm_tbuf_t**)calloc(p->n_threads, sizeof(mm_tbuf_t*));
+			s->buf = (mm_tbuf_t **)calloc(p->n_threads, sizeof(mm_tbuf_t *));
 			for (i = 0; i < p->n_threads; ++i)
 				s->buf[i] = mm_tbuf_init();
-			s->n_reg = (int*)calloc(5 * s->n_seq, sizeof(int));
+			s->n_reg = (int *)calloc(5 * s->n_seq, sizeof(int));
 			s->seg_off = s->n_reg + s->n_seq; // seg_off, n_seg, rep_len and frag_gap are allocated together with n_reg
 			s->n_seg = s->seg_off + s->n_seq;
 			s->rep_len = s->n_seg + s->n_seq;
 			s->frag_gap = s->rep_len + s->n_seq;
-			s->reg = (mm_reg1_t**)calloc(s->n_seq, sizeof(mm_reg1_t*));
+			s->reg = (mm_reg1_t **)calloc(s->n_seq, sizeof(mm_reg1_t *));
 			for (i = 1, j = 0; i <= s->n_seq; ++i)
-				if (i == s->n_seq || !frag_mode || !mm_qname_same(s->seq[i-1].name, s->seq[i].name)) {
+				if (i == s->n_seq || !frag_mode || !mm_qname_same(s->seq[i - 1].name, s->seq[i].name))
+				{
 					s->n_seg[s->n_frag] = i - j;
 					s->seg_off[s->n_frag++] = j;
 					j = i;
 				}
 			return s;
-		} else free(s);
-    } else if (step == 1) { // step 1: map
-		if (p->n_parts > 0) merge_hits((step_t*)in);
-		else kt_for(p->n_threads, worker_for, in, ((step_t*)in)->n_frag);
+		}
+		else
+			free(s);
+	}
+	else if (step == 1)
+	{ // step 1: map
+		if (p->n_parts > 0)
+			merge_hits((step_t *)in);
+		else
+			kt_for(p->n_threads, worker_for, in, ((step_t *)in)->n_frag);
 		return in;
-    } else if (step == 2) { // step 2: output
+	}
+	else if (step == 2)
+	{ // step 2: output
 		void *km = 0;
-        step_t *s = (step_t*)in;
+		step_t *s = (step_t *)in;
 		const mm_idx_t *mi = p->mi;
-		for (i = 0; i < p->n_threads; ++i) mm_tbuf_destroy(s->buf[i]);
+		for (i = 0; i < p->n_threads; ++i)
+			mm_tbuf_destroy(s->buf[i]);
 		free(s->buf);
-		if ((p->opt->flag & MM_F_OUT_CS) && !(mm_dbg_flag & MM_DBG_NO_KALLOC)) km = km_init();
-		for (k = 0; k < s->n_frag; ++k) {
+		if ((p->opt->flag & MM_F_OUT_CS) && !(mm_dbg_flag & MM_DBG_NO_KALLOC))
+			km = km_init();
+		for (k = 0; k < s->n_frag; ++k)
+		{
 			int seg_st = s->seg_off[k], seg_en = s->seg_off[k] + s->n_seg[k];
-			for (i = seg_st; i < seg_en; ++i) {
+			for (i = seg_st; i < seg_en; ++i)
+			{
 				mm_bseq1_t *t = &s->seq[i];
-				if (p->opt->split_prefix && p->n_parts == 0) { // then write to temporary files
-					mm_err_fwrite(&s->n_reg[i],    sizeof(int), 1, p->fp_split);
-					mm_err_fwrite(&s->rep_len[i],  sizeof(int), 1, p->fp_split);
+				if (p->opt->split_prefix && p->n_parts == 0)
+				{ // then write to temporary files
+					mm_err_fwrite(&s->n_reg[i], sizeof(int), 1, p->fp_split);
+					mm_err_fwrite(&s->rep_len[i], sizeof(int), 1, p->fp_split);
 					mm_err_fwrite(&s->frag_gap[i], sizeof(int), 1, p->fp_split);
-					for (j = 0; j < s->n_reg[i]; ++j) {
+					for (j = 0; j < s->n_reg[i]; ++j)
+					{
 						mm_reg1_t *r = &s->reg[i][j];
 						mm_err_fwrite(r, sizeof(mm_reg1_t), 1, p->fp_split);
-						if (p->opt->flag & MM_F_CIGAR) {
+						if (p->opt->flag & MM_F_CIGAR)
+						{
 							mm_err_fwrite(&r->p->capacity, 4, 1, p->fp_split);
 							mm_err_fwrite(r->p, r->p->capacity, 4, p->fp_split);
 						}
 					}
-				} else if (s->n_reg[i] > 0) { // the query has at least one hit
-					for (j = 0; j < s->n_reg[i]; ++j) {
+				}
+				else if (s->n_reg[i] > 0)
+				{ // the query has at least one hit
+					for (j = 0; j < s->n_reg[i]; ++j)
+					{
 						mm_reg1_t *r = &s->reg[i][j];
 						assert(!r->sam_pri || r->id == r->parent);
 						if ((p->opt->flag & MM_F_NO_PRINT_2ND) && r->id != r->parent)
 							continue;
 						if (p->opt->flag & MM_F_OUT_SAM)
-							mm_write_sam3(&p->str, mi, t, i - seg_st, j, s->n_seg[k], &s->n_reg[seg_st], (const mm_reg1_t*const*)&s->reg[seg_st], km, p->opt->flag, s->rep_len[i]);
+							mm_write_sam3(&p->str, mi, t, i - seg_st, j, s->n_seg[k], &s->n_reg[seg_st], (const mm_reg1_t *const *)&s->reg[seg_st], km, p->opt->flag, s->rep_len[i]);
 						else
 							mm_write_paf3(&p->str, mi, t, r, km, p->opt->flag, s->rep_len[i]);
 						mm_err_puts(p->str.s);
 					}
-				} else if ((p->opt->flag & MM_F_PAF_NO_HIT) || ((p->opt->flag & MM_F_OUT_SAM) && !(p->opt->flag & MM_F_SAM_HIT_ONLY))) { // output an empty hit, if requested
+				}
+				else if ((p->opt->flag & MM_F_PAF_NO_HIT) || ((p->opt->flag & MM_F_OUT_SAM) && !(p->opt->flag & MM_F_SAM_HIT_ONLY)))
+				{ // output an empty hit, if requested
 					if (p->opt->flag & MM_F_OUT_SAM)
-						mm_write_sam3(&p->str, mi, t, i - seg_st, -1, s->n_seg[k], &s->n_reg[seg_st], (const mm_reg1_t*const*)&s->reg[seg_st], km, p->opt->flag, s->rep_len[i]);
+						mm_write_sam3(&p->str, mi, t, i - seg_st, -1, s->n_seg[k], &s->n_reg[seg_st], (const mm_reg1_t *const *)&s->reg[seg_st], km, p->opt->flag, s->rep_len[i]);
 					else
 						mm_write_paf3(&p->str, mi, t, 0, 0, p->opt->flag, s->rep_len[i]);
 					mm_err_puts(p->str.s);
 				}
 			}
-			for (i = seg_st; i < seg_en; ++i) {
-				for (j = 0; j < s->n_reg[i]; ++j) free(s->reg[i][j].p);
+			for (i = seg_st; i < seg_en; ++i)
+			{
+				for (j = 0; j < s->n_reg[i]; ++j)
+					free(s->reg[i][j].p);
 				free(s->reg[i]);
-				free(s->seq[i].seq); free(s->seq[i].name);
-				if (s->seq[i].qual) free(s->seq[i].qual);
-				if (s->seq[i].comment) free(s->seq[i].comment);
+				free(s->seq[i].seq);
+				free(s->seq[i].name);
+				if (s->seq[i].qual)
+					free(s->seq[i].qual);
+				if (s->seq[i].comment)
+					free(s->seq[i].comment);
 			}
 		}
-		free(s->reg); free(s->n_reg); free(s->seq); // seg_off, n_seg, rep_len and frag_gap were allocated with reg; no memory leak here
+		free(s->reg);
+		free(s->n_reg);
+		free(s->seq); // seg_off, n_seg, rep_len and frag_gap were allocated with reg; no memory leak here
 		km_destroy(km);
 		if (mm_verbose >= 3)
 			fprintf(stderr, "[M::%s::%.3f*%.2f] mapped %d sequences\n", __func__, realtime() - mm_realtime0, cputime() / (realtime() - mm_realtime0), s->n_seq);
 		free(s);
 	}
-    return 0;
+	return 0;
 }
 
 static mm_bseq_file_t **open_bseqs(int n, const char **fn)
 {
 	mm_bseq_file_t **fp;
 	int i, j;
-	fp = (mm_bseq_file_t**)calloc(n, sizeof(mm_bseq_file_t*));
-	for (i = 0; i < n; ++i) {
-		if ((fp[i] = mm_bseq_open(fn[i])) == 0) {
+	fp = (mm_bseq_file_t **)calloc(n, sizeof(mm_bseq_file_t *));
+	for (i = 0; i < n; ++i)
+	{
+		if ((fp[i] = mm_bseq_open(fn[i])) == 0)
+		{
 			if (mm_verbose >= 1)
 				fprintf(stderr, "ERROR: failed to open file '%s': %s\n", fn[i], strerror(errno));
 			for (j = 0; j < i; ++j)
@@ -635,26 +851,36 @@ static mm_bseq_file_t **open_bseqs(int n, const char **fn)
 
 int mm_map_file_frag(const mm_idx_t *idx, int n_segs, const char **fn, const mm_mapopt_t *opt, int n_threads)
 {
+#ifdef ANALYSIS_PERF
+	int64_t tmp_cur_time = get_mseconds();
+#endif
 	int i, pl_threads;
 	pipeline_t pl;
-	if (n_segs < 1) return -1;
+	if (n_segs < 1)
+		return -1;
 	memset(&pl, 0, sizeof(pipeline_t));
 	pl.n_fp = n_segs;
 	pl.fp = open_bseqs(pl.n_fp, fn);
-	if (pl.fp == 0) return -1;
+	if (pl.fp == 0)
+		return -1;
 	pl.opt = opt, pl.mi = idx;
-	pl.n_threads = n_threads > 1? n_threads : 1;
+	pl.n_threads = n_threads > 1 ? n_threads : 1;
 	pl.mini_batch_size = opt->mini_batch_size;
 	if (opt->split_prefix)
 		pl.fp_split = mm_split_init(opt->split_prefix, idx);
-	pl_threads = n_threads == 1? 1 : (opt->flag&MM_F_2_IO_THREADS)? 3 : 2;
+	pl_threads = n_threads == 1 ? 1 : (opt->flag & MM_F_2_IO_THREADS) ? 3
+																	  : 2;
 	kt_pipeline(pl_threads, worker_pipeline, &pl, 3);
 
 	free(pl.str.s);
-	if (pl.fp_split) fclose(pl.fp_split);
+	if (pl.fp_split)
+		fclose(pl.fp_split);
 	for (i = 0; i < pl.n_fp; ++i)
 		mm_bseq_close(pl.fp[i]);
 	free(pl.fp);
+#ifdef ANALYSIS_PERF
+	time_mm_map_file_frag += get_mseconds() - tmp_cur_time;
+#endif
 	return 0;
 }
 
@@ -668,19 +894,22 @@ int mm_split_merge(int n_segs, const char **fn, const mm_mapopt_t *opt, int n_sp
 	int i;
 	pipeline_t pl;
 	mm_idx_t *mi;
-	if (n_segs < 1 || n_split_idx < 1) return -1;
+	if (n_segs < 1 || n_split_idx < 1)
+		return -1;
 	memset(&pl, 0, sizeof(pipeline_t));
 	pl.n_fp = n_segs;
 	pl.fp = open_bseqs(pl.n_fp, fn);
-	if (pl.fp == 0) return -1;
+	if (pl.fp == 0)
+		return -1;
 	pl.opt = opt;
 	pl.mini_batch_size = opt->mini_batch_size;
 
 	pl.n_parts = n_split_idx;
-	pl.fp_parts  = CALLOC(FILE*, pl.n_parts);
+	pl.fp_parts = CALLOC(FILE *, pl.n_parts);
 	pl.rid_shift = CALLOC(uint32_t, pl.n_parts);
 	pl.mi = mi = mm_split_merge_prep(opt->split_prefix, n_split_idx, pl.fp_parts, pl.rid_shift);
-	if (pl.mi == 0) {
+	if (pl.mi == 0)
+	{
 		free(pl.fp_parts);
 		free(pl.rid_shift);
 		return -1;
diff --git a/minimap.h b/minimap.h
index 59f6fa0..18f5b8f 100644
--- a/minimap.h
+++ b/minimap.h
@@ -7,6 +7,9 @@
 
 #define MM_VERSION "2.26-r1175"
 
+// 用来开关调试性能分析，运行时间等信息
+#define ANALYSIS_PERF 1
+
 #define MM_F_NO_DIAG       (0x001LL) // no exact diagonal hit
 #define MM_F_NO_DUAL       (0x002LL) // skip pairs where query name is lexicographically larger than target name
 #define MM_F_CIGAR         (0x004LL)
diff --git a/sketch.c b/sketch.c
index f830693..30230a2 100644
--- a/sketch.c
+++ b/sketch.c
@@ -106,7 +106,7 @@ void mm_sketch(void *km, const char *str, int len, int w, int k, uint32_t rid, i
 			kmer[0] = (kmer[0] << 2 | c) & mask;           // forward k-mer
 			kmer[1] = (kmer[1] >> 2) | (3ULL^c) << shift1; // reverse k-mer
 			if (kmer[0] == kmer[1]) continue; // skip "symmetric k-mers" as we don't know it strand
-			z = kmer[0] < kmer[1]? 0 : 1; // strand
+			z = kmer[0] < kmer[1]? 0 : 1; // strand // kmer的strand到底是什么意思？为什么通过比较就能确定正反？
 			++l;
 			if (l >= k && kmer_span < 256) {
 				info.x = hash64(kmer[z], mask) << 8 | kmer_span;