 +			kputsn("NA", 2, str);
 +		}
 +		kputc(',', str);
-+		kputw(p->pair_score, str); // ztz5
++		kputw(p->osub, str); // ztz5
++		kputc(',', str);
++		kputw(p->pair_score, str); // ztz6
 +		kputc(',', str);
 +		if(p->pair_sub >= 0) {
-+			kputw(p->pair_score - p->pair_sub, str); // ztz6
++			kputw(p->pair_score - p->pair_sub, str); // ztz7
 +		} else {
 +			kputsn("NA", 2, str);
 +		}
 +		kputc(',', str);
-+		kputw(p->n_comp, str); // ztz7
++		kputw(p->n_comp, str); // ztz8
 +		kputc(',', str);
-+		kputw(p->seedlen0, str); // ztz8
++		kputw(p->seedlen0, str); // ztz9
 +		kputc(',', str);
-+		kputw(p->sub_n, str); // ztz9
++		kputw(p->sub_n, str); // ztz10
 +		kputc(',', str);
-+		kputw(p->pair_nsub, str); // ztz10
++		kputw(p->osub_n, str); // ztz11
 +		kputc(',', str);
-+		kputw((int)(p->frac_rep * 1000.0), str); // ztz11
++		kputw(p->pair_nsub, str); // ztz12
 +		kputc(',', str);
-+		kputw((int)(p->pair_frac_rep * 1000.0), str); // ztz12
++		kputw((int)(p->frac_rep * 1000.0), str); // ztz13
 +		kputc(',', str);
-+		kputw((int)(p->seedcov * 1000.0), str); // ztz13
++		kputw((int)(p->pair_frac_rep * 1000.0), str); // ztz14
 +		kputc(',', str);
-+		kputw((int)(p->pair_seedcov * 1000.0), str); // ztz14
++		kputw((int)(p->seedcov * 1000.0), str); // ztz15
++		kputc(',', str);
++		kputw((int)(p->pair_seedcov * 1000.0), str); // ztz16
 +	}
  	kputc('\n', str);
  }
  
-@@ -1131,6 +1180,13 @@ mem_aln_t mem_reg2aln(const mem_opt_t *o
+@@ -1129,8 +1184,17 @@ mem_aln_t mem_reg2aln(const mem_opt_t *o
+ 	a.rid = bns_pos2rid(bns, pos);
+ 	assert(a.rid == ar->rid);
  	a.pos = pos - bns->anns[a.rid].offset;
- 	a.score = ar->score; a.sub = ar->sub > ar->csub? ar->sub : ar->csub;
+-	a.score = ar->score; a.sub = ar->sub > ar->csub? ar->sub : ar->csub;
++	a.score = ar->score;
++	a.sub = ar->sub > ar->csub? ar->sub : ar->csub;
++	a.osub = ar->osub;
  	a.is_alt = ar->is_alt; a.alt_sc = ar->alt_sc;
 +	a.seedlen0 = ar->seedlen0;
 +	a.n_comp = ar->n_comp;
 +	a.seedcov = ((float)ar->seedcov) / l_query;
 +	a.frac_rep = ar->frac_rep;
-+	a.sub_n = ar->sub_n;
++	a.sub_n = ar->sub_n; a.osub_n = ar->osub_n;
 +	a.pair_score = a.pair_sub = a.pair_nsub = 0;
 +	a.pair_frac_rep = a.pair_seedcov = 0.;
  	free(query);

diff -rupN -x Xcode -x .git bwa/bwamem.c bwa-mod/bwamem.c
--- bwa/bwamem.c	2015-07-08 13:14:28.000000000 -0400
+++ bwa-mod/bwamem.c	2016-07-29 22:21:24.000000000 -0400
@@ -838,6 +838,7 @@ void mem_aln2sam(const mem_opt_t *opt, c
 	ks_resize(str, str->l + s->l_seq + l_name + (s->qual? s->l_seq : 0) + 20);
 	kputsn(s->name, l_name, str); kputc('\t', str); // QNAME
 	kputw((p->flag&0xffff) | (p->flag&0x10000? 0x100 : 0), str); kputc('\t', str); // FLAG

 }
diff -rupN -x Xcode -x .git bwa/bwamem.h bwa-mod/bwamem.h
--- bwa/bwamem.h	2015-07-08 13:14:28.000000000 -0400
+++ bwa-mod/bwamem.h	2016-07-30 14:46:21.000000000 -0400
@@ -63,9 +63,11 @@ typedef struct {
 	int score;      // best local SW score
 	int truesc;     // actual score corresponding to the aligned region; possibly smaller than $score
 	int sub;        // 2nd best SW score
+	int osub;       // more liberal 2nd best SW score
 	int alt_sc;
 	int csub;       // SW score of a tandem hit
 	int sub_n;      // approximate number of suboptimal hits
+	int osub_n;     // more liberal approximate number of suboptimal hits
 	int w;          // actual band width used in extension
 	int seedcov;    // length of regions coverged by seeds
 	int secondary;  // index of the parent hit shadowing the current hit; <0 if primary
@@ -93,7 +95,11 @@ typedef struct { // This struct is only 
 	uint32_t *cigar; // CIGAR in the BAM encoding: opLen<<4|op; op to integer mapping: MIDSH=>01234
 	char *XA;        // alternative mappings
 
-	int score, sub, alt_sc;
+	int score, sub, osub, alt_sc;
+	float seedcov, pair_seedcov;
+	int sub_n, osub_n, pair_score, pair_sub, pair_nsub;
+	float frac_rep, pair_frac_rep;
+	int seedlen0, n_comp;
 } mem_aln_t;
