[med-svn] [Git][med-team/kma][upstream] New upstream version 1.3.11

Nilesh Patra gitlab at salsa.debian.org
Sun Feb 14 15:24:30 GMT 2021



Nilesh Patra pushed to branch upstream at Debian Med / kma


Commits:
4ac3c304 by Nilesh Patra at 2021-02-14T20:50:41+05:30
New upstream version 1.3.11
- - - - -


10 changed files:

- align.c
- alnfrags.c
- assembly.c
- nw.c
- nw.h
- runkma.c
- sam.c
- sparse.c
- version.h
- xml.c


Changes:

=====================================
align.c
=====================================
@@ -201,7 +201,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 		Stat.score = 0;
 		Stat.len = 1;
 		Stat.match = 0;
-		Stat.gaps = 0;
+		Stat.tGaps = 0;
+		Stat.qGaps = 0;
 		Stat.pos = 0;
 		aligned->s[0] = 0;
 		aligned->len = 0;
@@ -216,7 +217,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 		Stat.score = 0;
 		Stat.len = 1;
 		Stat.match = 0;
-		Stat.gaps = 0;
+		Stat.tGaps = 0;
+		Stat.qGaps = 0;
 		Stat.pos = 0;
 		aligned->s[0] = 0;
 		aligned->len = 0;
@@ -228,7 +230,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 	Stat.len = 0;
 	Stat.score = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	value = points->tStart[start] - 1;
 	Stat.pos = value;
 	i = points->qStart[start];
@@ -263,8 +266,10 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 			if(t_s == 0) {
 				while(bias < NWstat.len && (Frag_align->t[bias] == 5 || Frag_align->q[bias] == 5)) {
 					if(Frag_align->t[bias] == 5) {
-						--NWstat.gaps;
+						--NWstat.tGaps;
 						++(Frag_align->start);
+					} else {
+						--NWstat.qGaps;
 					}
 					++bias;
 				}
@@ -278,11 +283,12 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 			memcpy(aligned->s, Frag_align->s + bias, NWstat.len);
 			memcpy(aligned->q, Frag_align->q + bias, NWstat.len);
 			aligned->start = q_s + Frag_align->start;
-			Stat.pos -= (NWstat.len - NWstat.gaps);
+			Stat.pos -= (NWstat.len - NWstat.tGaps);
 			Stat.score = NWstat.score;
 			Stat.len = NWstat.len;
 			Stat.match = NWstat.match;
-			Stat.gaps = NWstat.gaps;
+			Stat.tGaps = NWstat.tGaps;
+			Stat.qGaps = NWstat.qGaps;
 		} else {
 			aligned->start = q_s;
 		}
@@ -340,7 +346,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 				Stat.score = 0;
 				Stat.len = 1;
 				Stat.match = 0;
-				Stat.gaps = 0;
+				Stat.tGaps = 0;
+				Stat.qGaps = 0;
 				aligned->s[0] = 0;
 				aligned->len = 0;
 				points->len = 0;
@@ -361,7 +368,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 				Stat.score += NWstat.score;
 				Stat.len += NWstat.len;
 				Stat.match += NWstat.match;
-				Stat.gaps += NWstat.gaps;
+				Stat.tGaps += NWstat.tGaps;
+				Stat.qGaps += NWstat.qGaps;
 			}
 		} else {
 			stop = 0;
@@ -398,8 +406,10 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 			bias = NWstat.len - 1;
 			while(bias && (Frag_align->t[bias] == 5 || Frag_align->q[bias] == 5)) {
 				if(Frag_align->t[bias] == 5) {
-					--NWstat.gaps;
+					--NWstat.tGaps;
 					++(Frag_align->end);
+				} else {
+					--NWstat.qGaps;
 				}
 				--bias;
 			}
@@ -418,7 +428,8 @@ AlnScore KMA(const HashMapCCI *template_index, const unsigned char *qseq, int q_
 		Stat.score += NWstat.score;
 		Stat.len += NWstat.len;
 		Stat.match += NWstat.match;
-		Stat.gaps += NWstat.gaps;
+		Stat.tGaps += NWstat.tGaps;
+		Stat.qGaps += NWstat.qGaps;
 	} else {
 		Frag_align->end = 0;
 	}
@@ -569,7 +580,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 		Stat.score = 0;
 		Stat.len = 1;
 		Stat.match = 0;
-		Stat.gaps = 0;
+		Stat.tGaps = 0;
+		Stat.qGaps = 0;
 		Stat.pos = 0;
 		points->len = 0;
 		return Stat;
@@ -584,7 +596,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 		Stat.score = 0;
 		Stat.len = 1;
 		Stat.match = 0;
-		Stat.gaps = 0;
+		Stat.tGaps = 0;
+		Stat.qGaps = 0;
 		Stat.pos = 0;
 		points->len = 0;
 		return Stat;
@@ -594,7 +607,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 	Stat.len = 0;
 	Stat.score = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	value = points->tStart[start] - 1;
 	Stat.pos = value;
 	i = points->qStart[start];
@@ -623,11 +637,12 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 				NWstat = NW_band_score(template_index->seq, qseq, -1 - (t_s == 0), t_s, t_e, q_s, q_e, band, matrices, t_len);
 				//NWstat = NW_score(template_index->seq, qseq, -1 - (t_s == 0), t_s, t_e, q_s, q_e, matrices, t_len);
 			}
-			Stat.pos -= (NWstat.len - NWstat.gaps);
+			Stat.pos -= (NWstat.len - NWstat.tGaps);
 			Stat.score = NWstat.score;
 			Stat.len = NWstat.len;
 			Stat.match = NWstat.match;
-			Stat.gaps = NWstat.gaps;
+			Stat.tGaps = NWstat.tGaps;
+			Stat.qGaps = NWstat.qGaps;
 		}
 	}
 	
@@ -679,7 +694,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 				Stat.score = 0;
 				Stat.len = 1;
 				Stat.match = 0;
-				Stat.gaps = 0;
+				Stat.tGaps = 0;
+				Stat.qGaps = 0;
 				points->len = 0;
 				return Stat;
 			}
@@ -693,7 +709,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 				Stat.score += NWstat.score;
 				Stat.len += NWstat.len;
 				Stat.match += NWstat.match;
-				Stat.gaps += NWstat.gaps;
+				Stat.tGaps += NWstat.tGaps;
+				Stat.qGaps += NWstat.qGaps;
 			}
 		} else {
 			stop = 0;
@@ -727,7 +744,8 @@ AlnScore KMA_score(const HashMapCCI *template_index, const unsigned char *qseq,
 		Stat.score += NWstat.score;
 		Stat.len += NWstat.len;
 		Stat.match += NWstat.match;
-		Stat.gaps += NWstat.gaps;
+		Stat.tGaps += NWstat.tGaps;
+		Stat.qGaps += NWstat.qGaps;
 	}
 	points->len = 0;
 	


=====================================
alnfrags.c
=====================================
@@ -34,6 +34,7 @@
 #include "stdstat.h"
 #include "threader.h"
 #include "updatescores.h"
+#define mrcheck(mrc, Stat, q_len, t_len) ((mrc * q_len <= Stat.len - Stat.qGaps) || (mrc * t_len <= Stat.len - Stat.tGaps))
 
 int (*alnFragsPE)(HashMapCCI**, int*, int*, int, double, double, int, CompDNA*, CompDNA*, CompDNA*, CompDNA*, unsigned char*, unsigned char*, unsigned char*, unsigned char*, Qseqs*, Qseqs*, int, int*, int*, long unsigned*, long unsigned*, int*, int*, int*, int*, int*, int*, int, long*, FILE*, AlnPoints*, NWmat*, volatile int*, volatile int*) = alnFragsUnionPE;
 
@@ -103,7 +104,8 @@ int alnFragsSE(HashMapCCI **templates_index, int *matched_templates, int *templa
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 				points->len = 0;
 			}
 		} else {
@@ -117,9 +119,9 @@ int alnFragsSE(HashMapCCI **templates_index, int *matched_templates, int *templa
 		/* get read score */
 		aln_len = alnStat.len;
 		start = alnStat.pos;
-		end = start + aln_len - alnStat.gaps;
+		end = start + aln_len - alnStat.tGaps;
 		t_len = template_lengths[abs(template)];
-		if(template_lengths[abs(template)] < end && ((t_len < q_len) ? mrc * t_len : mrc * q_len) <= alnStat.match) {
+		if(template_lengths[abs(template)] < end && mrcheck(mrc, alnStat, q_len, t_len)) {
 			end -= template_lengths[abs(template)];
 		}
 		
@@ -247,7 +249,8 @@ int alnFragsUnionPE(HashMapCCI **templates_index, int *matched_templates, int *t
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 				points->len = 0;
 			}
 		} else {
@@ -258,9 +261,9 @@ int alnFragsUnionPE(HashMapCCI **templates_index, int *matched_templates, int *t
 		aln_len = alnStat.len;
 		read_score = alnStat.score;
 		t_len = template_lengths[abs(template)];
-		if(minlen <= aln_len && 0 < read_score && ((t_len < qseq_comp->seqlen) ? mrc * t_len : mrc * qseq_comp->seqlen) <= alnStat.match) {
+		if(minlen <= aln_len && 0 < read_score && mrcheck(mrc, alnStat, qseq_comp->seqlen, t_len)) {
 			start = alnStat.pos;
-			end = alnStat.pos + alnStat.len - alnStat.gaps;
+			end = alnStat.pos + alnStat.len - alnStat.tGaps;
 			if(start == 0 && end == t_len) {
 				read_score += abs(W1);
 			}
@@ -294,7 +297,8 @@ int alnFragsUnionPE(HashMapCCI **templates_index, int *matched_templates, int *t
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 			}
 			rc = 1;
 		} else {
@@ -304,9 +308,9 @@ int alnFragsUnionPE(HashMapCCI **templates_index, int *matched_templates, int *t
 		/* get read score */
 		aln_len = alnStat.len;
 		read_score = alnStat.score;
-		if(minlen <= aln_len && 0 < read_score && ((t_len < qseq_r_comp->seqlen) ? mrc * t_len : mrc * qseq_r_comp->seqlen) <= alnStat.match) {
+		if(minlen <= aln_len && 0 < read_score && mrcheck(mrc, alnStat, qseq_r_comp->seqlen, t_len)) {
 			start = alnStat.pos;
-			end = alnStat.pos + alnStat.len - alnStat.gaps;
+			end = alnStat.pos + alnStat.len - alnStat.tGaps;
 			
 			if(start == 0 && end == t_len) {
 				read_score += abs(W1);
@@ -576,7 +580,8 @@ int alnFragsPenaltyPE(HashMapCCI **templates_index, int *matched_templates, int
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 				points->len = 0;
 			}
 		} else {
@@ -587,9 +592,9 @@ int alnFragsPenaltyPE(HashMapCCI **templates_index, int *matched_templates, int
 		aln_len = alnStat.len;
 		read_score = alnStat.score;
 		t_len = template_lengths[abs(template)];
-		if(minlen <= aln_len && 0 < read_score && ((t_len < qseq_comp->seqlen) ? mrc * t_len : mrc * qseq_comp->seqlen) <= alnStat.match) {
+		if(minlen <= aln_len && 0 < read_score && mrcheck(mrc, alnStat, qseq_comp->seqlen, t_len)) {
 			start = alnStat.pos;
-			end = alnStat.pos + alnStat.len - alnStat.gaps;
+			end = alnStat.pos + alnStat.len - alnStat.tGaps;
 			
 			if(start == 0 && end == t_len) {
 				read_score += abs(W1);
@@ -622,7 +627,8 @@ int alnFragsPenaltyPE(HashMapCCI **templates_index, int *matched_templates, int
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 			}
 			rc = 1;
 		} else {
@@ -632,9 +638,9 @@ int alnFragsPenaltyPE(HashMapCCI **templates_index, int *matched_templates, int
 		/* get read score */
 		aln_len = alnStat.len;
 		read_score = alnStat.score;
-		if(minlen <= aln_len && 0 < read_score && ((t_len < qseq_r_comp->seqlen) ? mrc * t_len : mrc * qseq_r_comp->seqlen) <= alnStat.match) {
+		if(minlen <= aln_len && 0 < read_score && mrcheck(mrc, alnStat, qseq_r_comp->seqlen, t_len)) {
 			start = alnStat.pos;
-			end = alnStat.pos + alnStat.len - alnStat.gaps;
+			end = alnStat.pos + alnStat.len - alnStat.tGaps;
 			
 			if(start == 0 && end == t_len) {
 				read_score += abs(W1);
@@ -899,7 +905,8 @@ int alnFragsForcePE(HashMapCCI **templates_index, int *matched_templates, int *t
 				alnStat.score = 0;
 				alnStat.pos = 0;
 				alnStat.len = 0;
-				alnStat.gaps = 0;
+				alnStat.tGaps = 0;
+				alnStat.qGaps = 0;
 				points->len = 0;
 			}
 		} else {
@@ -907,7 +914,7 @@ int alnFragsForcePE(HashMapCCI **templates_index, int *matched_templates, int *t
 		}
 		
 		t_len = template_lengths[abs(template)];
-		if(0 < alnStat.score && minlen <= alnStat.len && ((t_len < qseq_comp->seqlen) ? mrc * t_len : mrc * qseq_comp->seqlen) <= alnStat.match) {
+		if(0 < alnStat.score && minlen <= alnStat.len && mrcheck(mrc, alnStat, qseq_comp->seqlen, t_len)) {
 			if(arc) {
 				if(rc < 0) {
 					/* rc */
@@ -919,7 +926,8 @@ int alnFragsForcePE(HashMapCCI **templates_index, int *matched_templates, int *t
 					alnStat_r.score = 0;
 					alnStat_r.pos = 0;
 					alnStat_r.len = 0;
-					alnStat_r.gaps = 0;
+					alnStat_r.tGaps = 0;
+					alnStat_r.qGaps = 0;
 				}
 				rc = 1;
 			} else {
@@ -927,17 +935,17 @@ int alnFragsForcePE(HashMapCCI **templates_index, int *matched_templates, int *t
 			}
 			
 			/* get read score */
-			if(0 < alnStat_r.score && minlen <= alnStat_r.len) {
+			if(0 < alnStat_r.score && minlen <= alnStat_r.len && mrcheck(mrc, alnStat_r, qseq_r_comp->seqlen, t_len)) {
 				aln_len = alnStat.len + alnStat_r.len;
 				
 				/* Handle negative insertsizes caused by trimming,
 				user stupidity or sample error. */
 				if(alnStat.pos < alnStat_r.pos) {
 					start = alnStat.pos;
-					end = alnStat_r.pos + alnStat_r.len - alnStat_r.gaps;
+					end = alnStat_r.pos + alnStat_r.len - alnStat_r.tGaps;
 				} else {
 					start = alnStat_r.pos;
-					end = alnStat.pos + alnStat.len - alnStat.gaps;
+					end = alnStat.pos + alnStat.len - alnStat.tGaps;
 				}
 				
 				read_score = alnStat.score + alnStat_r.score;
@@ -953,7 +961,7 @@ int alnFragsForcePE(HashMapCCI **templates_index, int *matched_templates, int *t
 		}
 		
 		/* save best match(es) */
-		if(read_score > kmersize && score >= scoreT && ((t_len < qseq_r_comp->seqlen) ? mrc * t_len : mrc * qseq_r_comp->seqlen) <= alnStat.match) {
+		if(read_score > kmersize && score >= scoreT) {
 			if(score > bestScore) { // save as best match
 				bestScore = score;
 				*best_read_score = read_score;
@@ -1117,7 +1125,7 @@ void * alnFrags_threaded(void * arg) {
 			unmapped = 0;
 		}
 		
-		if(sam && unmapped) {
+		if(sam && !(sam & 2096) && unmapped) {
 			if(unmapped & 1) {
 				stats[1] = flag;
 				nibble2base(qseq->seq, qseq->len);


=====================================
assembly.c
=====================================
@@ -39,6 +39,7 @@
 #include "stdstat.h"
 #include "threader.h"
 #include "xml.h"
+#define mrcheck(mrc, Stat, q_len, t_len) ((mrc * q_len <= Stat.len - Stat.qGaps) || (mrc * t_len <= Stat.len - Stat.tGaps))
 
 void * (*assembly_KMA_Ptr)(void *) = &assemble_KMA;
 int (*significantBase)(int, int, double) = &significantNuc;
@@ -432,11 +433,11 @@ void * assemble_KMA_threaded(void *arg) {
 						/* get read score */
 						aln_len = alnStat.len;
 						start = alnStat.pos;
-						end = start + aln_len - alnStat.gaps;
+						end = start + aln_len - alnStat.tGaps;
 						
 						/* Get normed score check read coverage */
 						read_score = alnStat.score;
-						if(minlen <= aln_len && ((t_len < qseq->len) ? mrc * t_len : mrc * qseq->len) <= alnStat.match) {
+						if(minlen <= aln_len && mrcheck(mrc, alnStat, qseq->len, t_len)) {
 							score = 1.0 * read_score / aln_len;
 						} else {
 							read_score = 0;
@@ -986,11 +987,11 @@ void * assemble_KMA_dense_threaded(void *arg) {
 						/* get read score */
 						aln_len = alnStat.len;
 						start = alnStat.pos;
-						end = start + aln_len - alnStat.gaps;
+						end = start + aln_len - alnStat.tGaps;
 						
 						/* Get normed score check read coverage */
 						read_score = alnStat.score;
-						if(minlen <= aln_len && ((t_len < qseq->len) ? mrc * t_len : mrc * qseq->len) <= alnStat.match) {
+						if(minlen <= aln_len && mrcheck(mrc, alnStat, qseq->len, t_len)) {
 							score = 1.0 * read_score / aln_len;
 						} else {
 							read_score = 0;
@@ -1293,7 +1294,7 @@ void * skip_assemble_KMA(void *arg) {
 	}
 	
 	aligned_assem->cover = 0; 
-	aligned_assem->aln_len = (1 - exp((-1.0) * aligned_assem->depth / t_len)) * t_len; // expected coverage from depth
+	aligned_assem->aln_len = 0;//(1 - exp((-1.0) * aligned_assem->depth / t_len)) * t_len; // expected coverage from depth
 	
 	return NULL;
 }
@@ -1858,11 +1859,11 @@ void * assemble_KMA(void *arg) {
 						/* get read score */
 						aln_len = alnStat.len;
 						start = alnStat.pos;
-						end = start + aln_len - alnStat.gaps;
+						end = start + aln_len - alnStat.tGaps;
 						
 						/* Get normed score check read coverage */
 						read_score = alnStat.score;
-						if(minlen <= aln_len && ((t_len < qseq->len) ? mrc * t_len : mrc * qseq->len) <= alnStat.match) {
+						if(minlen <= aln_len && mrcheck(mrc, alnStat, qseq->len, t_len)) {
 							score = 1.0 * read_score / aln_len;
 						} else {
 							read_score = 0;


=====================================
nw.c
=====================================
@@ -50,13 +50,15 @@ AlnScore NW(const long unsigned *template, const unsigned char *queryOrg, int k,
 		if(t_len == q_len) {
 			Stat.len = 0;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = 0;
 			Stat.score = 0;
 			aligned->s[0] = 0;
 		} else if(t_len == 0) {
 			Stat.len = q_len;
 			Stat.match = 0;
-			Stat.gaps = q_len;
+			Stat.tGaps = q_len;
+			Stat.qGaps = 0;
 			Stat.score = W1 + (q_len - 1) * U;
 			memset(aligned->s, '_', q_len);
 			aligned->s[q_len] = 0;
@@ -65,7 +67,8 @@ AlnScore NW(const long unsigned *template, const unsigned char *queryOrg, int k,
 		} else {
 			Stat.len = t_len;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = t_len;
 			Stat.score = W1 + (t_len - 1) * U;
 			memset(aligned->s, '_', t_len);
 			aligned->s[t_len] = 0;
@@ -252,7 +255,8 @@ AlnScore NW(const long unsigned *template, const unsigned char *queryOrg, int k,
 	nuc_pos = m + t_s;
 	Stat.len = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	while(E_ptr[n] != 0) {
 		if(nuc_pos == template_length) {
 			nuc_pos = 0;
@@ -273,26 +277,28 @@ AlnScore NW(const long unsigned *template, const unsigned char *queryOrg, int k,
 				++nuc_pos;
 				E_ptr += (q_len + 1);
 				++Stat.len;
+				++Stat.qGaps;
 			}
 			aligned->t[Stat.len] = getNuc(template, nuc_pos);
 			aligned->q[Stat.len] = 5;
 			aligned->s[Stat.len] = '_';		
 			++nuc_pos;
 			E_ptr += (q_len + 1);
+			++Stat.qGaps;
 		} else {
 			while(!(E_ptr[n] >> 3)) {
 				aligned->t[Stat.len] = 5;
 				aligned->q[Stat.len] = query[n];
 				aligned->s[Stat.len] = '_';
-				++Stat.gaps;
 				++n;
 				++Stat.len;
+				++Stat.tGaps;
 			}
 			aligned->t[Stat.len] = 5;
 			aligned->q[Stat.len] = query[n];
 			aligned->s[Stat.len] = '_';
-			++Stat.gaps;
 			++n;
+			++Stat.tGaps;
 		}
 		++Stat.len;
 	}
@@ -328,15 +334,17 @@ AlnScore NW_band(const long unsigned *template, const unsigned char *queryOrg, i
 	
 	if(t_len == 0 || q_len == 0) {
 		if(t_len == q_len) {
+			Stat.score = 0;
 			Stat.len = 0;
 			Stat.match = 0;
-			Stat.gaps = 0;
-			Stat.score = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = 0;
 			aligned->s[0] = 0;
 		} else if(t_len == 0) {
 			Stat.len = q_len;
 			Stat.match = 0;
-			Stat.gaps = q_len;
+			Stat.tGaps = q_len;
+			Stat.qGaps = 0;
 			Stat.score = W1 + (q_len - 1) * U;
 			memset(aligned->s, '_', q_len);
 			aligned->s[q_len] = 0;
@@ -345,7 +353,8 @@ AlnScore NW_band(const long unsigned *template, const unsigned char *queryOrg, i
 		} else {
 			Stat.len = t_len;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = t_len;
 			Stat.score = W1 + (t_len - 1) * U;
 			memset(aligned->s, '_', t_len);
 			aligned->s[t_len] = 0;
@@ -576,7 +585,8 @@ AlnScore NW_band(const long unsigned *template, const unsigned char *queryOrg, i
 	nuc_pos = m + t_s;
 	Stat.len = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	while(E_ptr[n] != 0) {
 		if(nuc_pos == template_length) {
 			nuc_pos = 0;
@@ -598,6 +608,7 @@ AlnScore NW_band(const long unsigned *template, const unsigned char *queryOrg, i
 				E_ptr += (bq_len + 1);
 				--n;
 				++Stat.len;
+				++Stat.qGaps;
 			}
 			aligned->t[Stat.len] = getNuc(template, nuc_pos);
 			aligned->q[Stat.len] = 5;
@@ -605,22 +616,23 @@ AlnScore NW_band(const long unsigned *template, const unsigned char *queryOrg, i
 			++nuc_pos;
 			E_ptr += (bq_len + 1);
 			--n;
+			++Stat.qGaps;
 		} else {
 			while(!(E_ptr[n] >> 3)) {
 				aligned->t[Stat.len] = 5;
 				aligned->q[Stat.len] = query[q_pos];
 				aligned->s[Stat.len] = '_';
-				++Stat.gaps;
 				++n;
 				++q_pos;
 				++Stat.len;
+				++Stat.tGaps;
 			}
 			aligned->t[Stat.len] = 5;
 			aligned->q[Stat.len] = query[q_pos];
 			aligned->s[Stat.len] = '_';
-			++Stat.gaps;
 			++n;
 			++q_pos;
+			++Stat.tGaps;
 		}
 		++Stat.len;
 	}
@@ -654,17 +666,20 @@ AlnScore NW_score(const long unsigned *template, const unsigned char *queryOrg,
 		if(t_len == q_len) {
 			Stat.len = 0;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = 0;
 			Stat.score = 0;
 		} else if(t_len == 0) {
 			Stat.len = q_len;
 			Stat.match = 0;
-			Stat.gaps = q_len;
+			Stat.tGaps = q_len;
+			Stat.qGaps = 0;
 			Stat.score = W1 + (q_len - 1) * U;
 		} else {
 			Stat.len = t_len;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = t_len;
 			Stat.score = W1 + (t_len - 1) * U;
 		}
 		return Stat;
@@ -841,7 +856,8 @@ AlnScore NW_score(const long unsigned *template, const unsigned char *queryOrg,
 	nuc_pos = m + t_s;
 	Stat.len = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	while(E_ptr[n] != 0) {
 		if(nuc_pos == template_length) {
 			nuc_pos = 0;
@@ -856,16 +872,18 @@ AlnScore NW_score(const long unsigned *template, const unsigned char *queryOrg,
 				++nuc_pos;
 				E_ptr += (q_len + 1);
 				++Stat.len;
+				++Stat.qGaps;
 			}
+			++Stat.qGaps;
 			++nuc_pos;
 			E_ptr += (q_len + 1);
 		} else {
 			while(!(E_ptr[n] >> 3)) {
-				++Stat.gaps;
 				++n;
 				++Stat.len;
+				++Stat.tGaps;
 			}
-			++Stat.gaps;
+			++Stat.tGaps;
 			++n;
 		}
 		++Stat.len;
@@ -899,17 +917,20 @@ AlnScore NW_band_score(const long unsigned *template, const unsigned char *query
 		if(t_len == q_len) {
 			Stat.len = 0;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = 0;
 			Stat.score = 0;
 		} else if(t_len == 0) {
 			Stat.len = q_len;
 			Stat.match = 0;
-			Stat.gaps = q_len;
+			Stat.tGaps = q_len;
+			Stat.qGaps = 0;
 			Stat.score = W1 + (q_len - 1) * U;
 		} else {
 			Stat.len = t_len;
 			Stat.match = 0;
-			Stat.gaps = 0;
+			Stat.tGaps = 0;
+			Stat.qGaps = t_len;
 			Stat.score = W1 + (t_len - 1) * U;
 		}
 		return Stat;
@@ -1129,7 +1150,8 @@ AlnScore NW_band_score(const long unsigned *template, const unsigned char *query
 	nuc_pos = m + t_s;
 	Stat.len = 0;
 	Stat.match = 0;
-	Stat.gaps = 0;
+	Stat.tGaps = 0;
+	Stat.qGaps = 0;
 	while(E_ptr[n] != 0) {
 		if(nuc_pos == template_length) {
 			nuc_pos = 0;
@@ -1145,18 +1167,20 @@ AlnScore NW_band_score(const long unsigned *template, const unsigned char *query
 				E_ptr += (bq_len + 1);
 				--n;
 				++Stat.len;
+				++Stat.qGaps;
 			}
+			++Stat.qGaps;
 			++nuc_pos;
 			E_ptr += (bq_len + 1);
 			--n;
 		} else {
 			while(!(E_ptr[n] >> 3)) {
-				++Stat.gaps;
 				++n;
 				++q_pos;
 				++Stat.len;
+				++Stat.tGaps;
 			}
-			++Stat.gaps;
+			++Stat.tGaps;
 			++n;
 			++q_pos;
 		}


=====================================
nw.h
=====================================
@@ -38,7 +38,8 @@ struct alnScore {
 	int len;
 	int pos;
 	int match;
-	int gaps;
+	int tGaps;
+	int qGaps;
 };
 
 struct aln {


=====================================
runkma.c
=====================================
@@ -335,7 +335,7 @@ int runKMA(char *templatefilename, char *outputfilename, char *exePrev, int ConC
 		alnThread->template_lengths = template_lengths;
 		alnThread->templates_index = templates_index;
 		alnThread->mq = mq;
-		alnThread->sam = (sam == 1) ? 1 : 0;
+		alnThread->sam = sam;//(sam == 1) ? 1 : 0;
 		alnThread->scoreT = scoreT;
 		alnThread->next = alnThreads;
 		alnThreads = alnThread;
@@ -405,7 +405,7 @@ int runKMA(char *templatefilename, char *outputfilename, char *exePrev, int ConC
 	alnThread->kmersize = kmersize;
 	alnThread->minlen = minlen;
 	alnThread->mq = mq;
-	alnThread->sam = (sam == 1) ? 1 : 0;
+	alnThread->sam = sam;//(sam == 1) ? 1 : 0;
 	alnThread->scoreT = scoreT;
 	alnThread->template_lengths = template_lengths;
 	alnThread->templates_index = templates_index;
@@ -1258,11 +1258,11 @@ int runKMA(char *templatefilename, char *outputfilename, char *exePrev, int ConC
 					}
 				}
 			} else {
-				if(sam || ID_t == 0.0) {
+				if((sam && !(sam & 2096)) || ID_t == 0.0) {
 					thread->template_index = templates_index[template];
 					thread->template_name = nameLoad(template_name, name_file);
 					thread->template = template;
-					assembly_KMA_Ptr(thread);
+					skip_assemble_KMA(thread);
 					//skip_assemble_KMA(template, sam, t_len, thread->template_name, fileCount, template_fragments, aligned_assem, qseq, header);
 					if(ID_t == 0.0) {
 						depth = aligned_assem->depth;
@@ -2424,7 +2424,7 @@ int runKMA_MEM(char *templatefilename, char *outputfilename, char *exePrev, int
 					}
 				}
 			} else {
-				if(sam || ID_t == 0.0) {
+				if((sam && !(sam & 2096)) || ID_t == 0.0) {
 					/* load DB */
 					seq_seeker *= sizeof(long unsigned);
 					lseek(seq_in_no, seq_seeker, SEEK_CUR);


=====================================
sam.c
=====================================
@@ -142,7 +142,6 @@ int samwrite(const Qseqs *qseq, const Qseqs *header, const Qseqs *Qual, char *rn
 	2048	supplementary alignment
 	*/
 	
-	
 	qname = (char *) header->seq;
 	seq = qseq->seq;
 	if(Qual) {
@@ -150,6 +149,7 @@ int samwrite(const Qseqs *qseq, const Qseqs *header, const Qseqs *Qual, char *rn
 	} else {
 		qual = "*";
 	}
+	
 	if(aligned) {
 		mapQ = 254 < aligned->mapQ ? 254 : aligned->mapQ;
 		et = *stats;
@@ -172,11 +172,11 @@ int samwrite(const Qseqs *qseq, const Qseqs *header, const Qseqs *Qual, char *rn
 	}
 	rnext = "*";
 	pnext = 0;
-	
 	tab = 0;
 	if(qname) {
 		while(*qname) {
 			if(*qname == '\t') {
+				tab = -tab;
 				*qname = 0;
 			} else {
 				++qname;
@@ -192,15 +192,16 @@ int samwrite(const Qseqs *qseq, const Qseqs *header, const Qseqs *Qual, char *rn
 	}
 	if(aligned) {
 		cigar = makeCigar(Cigar, aligned);
-	} else if(2 * sizeof(int) + 1 < header->len && header->seq[header->len - 2 * sizeof(int) - 1] == 0) {
-		cigar = (char *) Cigar->seq;
-		sprintf(cigar, "%dS", *((int*) (header->seq + (header->len - 2 * sizeof(int)))));
+		if(2 * sizeof(int) + 1 < header->len && header->seq[header->len - 2 * sizeof(int) - 1] == 0) {
+			cigar = (char *) Cigar->seq;
+			sprintf(cigar, "%dS", *((int*) (header->seq + (header->len - 2 * sizeof(int)))));
+		}
 	}
 	size = fprintf(stdout, "%s\t%d\t%s\t%d\t%d\t%s\t%s\t%d\t%d\t%s\t%s\tET:i:%d\tAS:i:%d\n", qname, flag, rname, pos, mapQ, cigar, rnext, pnext, tlen, (char *) seq, qual, et, score);
 	unlock(lock);
 	
-	if(tab) {
-		qname[tab] = '\t';
+	if(tab < 0) {
+		qname[-tab] = '\t';
 	}
 	
 	return size;


=====================================
sparse.c
=====================================
@@ -153,8 +153,12 @@ char ** load_DBs_Sparse(char *templatefilename, unsigned **template_lengths, uns
 		
 		/* load lengths */
 		fseek(DB_file, DB_size * sizeof(int), SEEK_CUR);
-		fread(*template_lengths, sizeof(int), DB_size, DB_file);
-		fread(*template_ulengths, sizeof(int), DB_size, DB_file);	
+		i = fread(*template_lengths, sizeof(int), DB_size, DB_file);
+		i += fread(*template_ulengths, sizeof(int), DB_size, DB_file);	
+		if(!i) {
+			fprintf(stderr, "DB needs to sparse indexed, to run a sparse mapping.\n");
+			exit(1);
+		}
 	}
 	templatefilename[file_len] = 0;
 	fclose(DB_file);


=====================================
version.h
=====================================
@@ -17,4 +17,4 @@
  * limitations under the License.
 */
 
-#define KMA_VERSION "1.3.10"
+#define KMA_VERSION "1.3.11"


=====================================
xml.c
=====================================
@@ -212,7 +212,7 @@ void hitXML(FILE *out, const int template, const char *template_name, const Aln
 	fprintf(out, "\t\t\t<Hsp_query-from>%d</Hsp_query-from>\n", ((flag & 16) ? (aligned->end) : (aligned->start)) + 1);
 	fprintf(out, "\t\t\t<Hsp_query-to>%d</Hsp_query-to>\n", ((flag & 16) ? (aligned->start) : (aligned->end)) + 1);
 	fprintf(out, "\t\t\t<Hsp_hit-from>%d</Hsp_hit-from>\n", alnStat->pos + 1);
-	fprintf(out, "\t\t\t<Hsp_hit-to>%d</Hsp_hit-to>\n", alnStat->pos + alnStat->len - alnStat->gaps + 1);
+	fprintf(out, "\t\t\t<Hsp_hit-to>%d</Hsp_hit-to>\n", alnStat->pos + alnStat->len - alnStat->tGaps + 1);
 	fprintf(out, "\t\t\t<Hsp_query-frame>%d</Hsp_query-frame>\n", aligned->start % 3);
 	fprintf(out, "\t\t\t<Hsp_hit-frame>%d</Hsp_hit-frame>\n", alnStat->pos % 3);
 	fprintf(out, "\t\t\t<Hsp_identity>%d</Hsp_identity>\n", Ms);



View it on GitLab: https://salsa.debian.org/med-team/kma/-/commit/4ac3c30413e380eb61bbe03a1ecd5d8310fee40a

-- 
View it on GitLab: https://salsa.debian.org/med-team/kma/-/commit/4ac3c30413e380eb61bbe03a1ecd5d8310fee40a
You're receiving this email because of your account on salsa.debian.org.


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/debian-med-commit/attachments/20210214/a2c03997/attachment-0001.html>


More information about the debian-med-commit mailing list