X-Git-Url: http://woldlab.caltech.edu/gitweb/?p=samtools.git;a=blobdiff_plain;f=bam_maqcns.c;h=f36b0ee2ab443affe0635866a8d593c5cb54fdf7;hp=464288ae7da0ac74f1443bbb55c03684927cb391;hb=4a17fa7e1f91b2fe04ad334a63fc2b0d5e859d8a;hpb=b27e00385f41769d03a8cca4dbd71275fc9fa906

diff --git a/bam_maqcns.c b/bam_maqcns.c
index 464288a..f36b0ee 100644
--- a/bam_maqcns.c
+++ b/bam_maqcns.c
@@ -439,7 +439,7 @@ bam_maqindel_ret_t *bam_maqindel(int n, int pos, const bam_maqindel_opt_t *mi, c
 				}
 				if (score[i*n+j] < s) score[i*n+j] = s; // choose the higher of the two scores
 				if (pscore[i*n+j] > ps) pscore[i*n+j] = ps;
-				if (types[i] != 0) score[i*n+j] -= mi->indel_err;
+				//if (types[i] != 0) score[i*n+j] -= mi->indel_err;
 				//printf("%d, %d, %d, %d, %d, %d, %d\n", p->b->core.pos + 1, seg.qbeg, i, types[i], j,
 				//	   score[i*n+j], pscore[i*n+j]);
 			}
@@ -499,6 +499,15 @@ bam_maqindel_ret_t *bam_maqindel(int n, int pos, const bam_maqindel_opt_t *mi, c
 				if (s1 > s2) ret->gl[0] += s1 - s2 < mi->q_indel? s1 - s2 : mi->q_indel;
 				else ret->gl[1] += s2 - s1 < mi->q_indel? s2 - s1 : mi->q_indel;
 			}
+			// write cnt_ref and cnt_ambi
+			if (max1_i != 0 && max2_i != 0) {
+				for (j = 0; j < n; ++j) {
+					int diff1 = score[j] - score[max1_i * n + j];
+					int diff2 = score[j] - score[max2_i * n + j];
+					if (diff1 > 0 && diff2 > 0) ++ret->cnt_ref;
+					else if (diff1 == 0 || diff2 == 0) ++ret->cnt_ambi;
+				}
+			}
 		}
 		free(score); free(pscore); free(ref2); free(inscns);
 	}