this works for a single-end test case

aln_sink.cpp

@@ -657,6 +657,7 @@ void AlnSinkWrap::finishRead(
 	ReportingMetrics&  met,         // reporting metrics
 	const PerReadMetrics& prm,      // per-read metrics
 	const Scoring& sc,              // scoring scheme
+	const bool selectBest,         // only report the best-scoring alignments
 	bool suppressSeedSummary,       // = true
 	bool suppressAlignments,        // = false
 	bool scUnMapped,                // = false
@@ -727,6 +728,7 @@ void AlnSinkWrap::finishRead(
 				&rs1_, &rs2_,
 				nconcord, select1_,
 				&rs1u_, &rs2u_,
+				selectBest,
 				bestUScore,
 				bestUDist,
 				bestP1Score,
@@ -903,6 +905,7 @@ void AlnSinkWrap::finishRead(
 				&rs1_, &rs2_,
 				ndiscord, select1_,
 				&rs1u_, &rs2u_,
+				selectBest,
 				bestUScore,
 				bestUDist,
 				bestP1Score,
@@ -1072,6 +1075,7 @@ void AlnSinkWrap::finishRead(
 			AlnScore bestUnchosenUDist, bestUnchosenP1Dist, bestUnchosenP2Dist, bestUnchosenCDist;
 			size_t off = selectByScore(
 				&rs1u_, NULL, nunpair1, select1_, NULL, NULL,
+				selectBest,
 				bestUScore,
 				bestUDist,
 				bestP1Score,
@@ -1125,6 +1129,7 @@ void AlnSinkWrap::finishRead(
 			AlnScore bestUnchosenUDist, bestUnchosenP1Dist, bestUnchosenP2Dist, bestUnchosenCDist;
 			size_t off = selectByScore(
 				&rs2u_, NULL, nunpair2, select2_, NULL, NULL,
+				selectBest,
 				bestUScore,
 				bestUDist,
 				bestP1Score,
@@ -1481,6 +1486,7 @@ size_t AlnSinkWrap::selectByScore(
 	EList<size_t>&       select, // prioritized list to put results in
 	const EList<AlnRes>* rs1u,   // alignments to select from (mate 1)
 	const EList<AlnRes>* rs2u,   // alignments to select from (mate 2, or NULL)
+	const bool           selectBest, // only select the top-scoring alignments
 	AlnScore&            bestUScore,
 	AlnScore&            bestUDist,
 	AlnScore&            bestP1Score,
@@ -1534,7 +1540,7 @@ size_t AlnSinkWrap::selectByScore(
 	if(sz == 0) {
 		return 0;
 	}
-	select.resize((size_t)num);
+
 	// Use 'selectBuf_' as a temporary list for sorting purposes
 	EList<std::pair<AlnScore, size_t> >& buf =
 		const_cast<EList<std::pair<AlnScore, size_t> >& >(selectBuf_);
@@ -1566,7 +1572,19 @@ size_t AlnSinkWrap::selectByScore(
 	if(streak > 1) {
 		buf.shufflePortion(buf.size() - streak, streak, rnd);
 	}
-
+
+	if (selectBest) {
+		// Counts the number of top-scoring alignments
+		// The final count should be <= num (# selected alignments)
+		uint64_t num_best = 1;
+		for(size_t i = 1; i < buf.size(); i++) {
+			if (buf[i].first == buf[0].first && num_best < num){
+				num_best += 1;
+			}
+		}
+		num = num_best;
+	}
+	select.resize((size_t)num);
 	// Copy the permutation into the 'select' list
 	for(size_t i = 0; i < num; i++) { select[i] = buf[i].second; }
 

aln_sink.h

@@ -219,9 +219,10 @@ struct ReportingParams {
 		THitInt pengap_,
 		bool msample_,
 		bool discord_,
-		bool mixed_)
+		bool mixed_,
+		bool kbest_)
 	{
-		init(khits_, mhits_, pengap_, msample_, discord_, mixed_);
+		init(khits_, mhits_, pengap_, msample_, discord_, mixed_, kbest_);
 	}
 
 	void init(
@@ -230,14 +231,16 @@ struct ReportingParams {
 		THitInt pengap_,
 		bool msample_,
 		bool discord_,
-		bool mixed_)
+		bool mixed_,
+		bool kbest_)
 	{
 		khits   = khits_;     // -k (or high if -a)
 		mhits   = ((mhits_ == 0) ? std::numeric_limits<THitInt>::max() : mhits_);
 		pengap  = pengap_;
 		msample = msample_;
 		discord = discord_;
 		mixed   = mixed_;
+		kbest   = kbest_;
 	}
 
 #ifndef NDEBUG
@@ -310,6 +313,10 @@ struct ReportingParams {
 	// are paired-end alignments for a paired-end read, or if the number of
 	// paired-end alignments exceeds the -m ceiling.
 	bool mixed;
+
+	// true iff both -k and --kbest are specified and we only report the best
+	// scoring alignments among the top-k hits
+	bool kbest;
 };
 
 /**
@@ -1027,7 +1034,8 @@ class AlnSinkWrap {
 		RandomSource&      rnd,         // pseudo-random generator
 		ReportingMetrics&  met,         // reporting metrics
 		const PerReadMetrics& prm,      // per-read metrics
-		const Scoring& sc,              // scoring scheme
+		const Scoring&     sc,          // scoring scheme
+		bool               selectBest,  // only report the top-scoring alignments
 		bool suppressSeedSummary = true,
 		bool suppressAlignments = false,
 		bool scUnMapped = false,
@@ -1234,6 +1242,7 @@ class AlnSinkWrap {
 		EList<size_t>&       select, // prioritized list to put results in
 		const EList<AlnRes>* rs1u,   // alignments to select from (mate 1)
 		const EList<AlnRes>* rs2u,   // alignments to select from (mate 2, or NULL)
+		const bool           selectBest, // only select the top-scoring alignments
 		AlnScore&            bestUScore,
 		AlnScore&            bestUDist,
 		AlnScore&            bestP1Score,

bt2_search.cpp

@@ -1291,7 +1291,6 @@ static void parseOption(int next_option, const char *arg) {
 	}
 	case ARG_K_BEST: {
 		kbest = true;
-		cerr << "Set -kbest to true" << endl;
 		break;
 	}
 	case ARG_VERBOSE: gVerbose = 1; break;
@@ -1802,6 +1801,11 @@ static void parseOptions(int argc, const char **argv) {
 		     << " instead" << endl;
 		multiseedMms = multiseedLen-1;
 	}
+	if (kbest && khits == 1) {
+		cerr << "Warning: --kbest is set while -k is not set or set to 1 "
+		     << "(the default value). The kbest mode would not make a "
+			 << "difference when -k is <= 1." << endl;
+	}
 	sam_print_zm = sam_print_zm && bowtie2p5;
 #ifndef NDEBUG
 	if(!gQuiet) {
@@ -4126,6 +4130,7 @@ static void multiseedSearchWorker(void *vp) {
 					rpm,                  // reporting metrics
 					prm,                  // per-read metrics
 					sc,                   // scoring scheme
+					kbest,                // only reports the top-scoring alignments
 					!seedSumm,            // suppress seed summaries?
 					seedSumm,             // suppress alignments?
 					scUnMapped,           // Consider soft-clipped bases unmapped when calculating TLEN
@@ -4483,6 +4488,7 @@ static void multiseedSearchWorker_2p5(void *vp) {
 				rpm,                  // reporting metrics
 				prm,                  // per-read metrics
 				sc,                   // scoring scheme
+				kbest,                // only reports the top-scoring alignments
 				!seedSumm,            // suppress seed summaries?
 				seedSumm,             // suppress alignments?
 				scUnMapped,           // Consider soft-clipped bases unmapped when calculating TLEN