Thread scaling improvements

.gitignore

@@ -24,3 +24,7 @@ bowtie2-align-l-debug
 bowtie2-inspect-s-debug
 bowtie2-inspect-l-debug
 Xcode
+reads_*.fq
+reads_*.sam
+.DS_Store
+.tmp.simple_tests.*.pl

aligner_cache.cpp

@@ -51,16 +51,14 @@ bool SAVal::repOk(const AlignmentCache& ac) const {
  * Add a new association between a read sequnce ('seq') and a
  * reference sequence ('')
  */
-bool AlignmentCache::addOnTheFly(
+bool AlignmentCache::addOnTheFlyImpl(
 	QVal& qv,         // qval that points to the range of reference substrings
 	const SAKey& sak, // the key holding the reference substring
 	TIndexOffU topf,    // top range elt in BWT index
 	TIndexOffU botf,    // bottom range elt in BWT index
 	TIndexOffU topb,    // top range elt in BWT' index
-	TIndexOffU botb,    // bottom range elt in BWT' index
-	bool getLock)
+	TIndexOffU botb)    // bottom range elt in BWT' index
 {
-    ThreadSafe ts(lockPtr(), shared_ && getLock);
 	bool added = true;
 	// If this is the first reference sequence we're associating with
 	// the query sequence, initialize the QVal.
@@ -104,6 +102,23 @@ bool AlignmentCache::addOnTheFly(
 	return true; 
 }
 
+bool AlignmentCache::addOnTheFly(
+	QVal& qv,         // qval that points to the range of reference substrings
+	const SAKey& sak, // the key holding the reference substring
+	TIndexOffU topf,    // top range elt in BWT index
+	TIndexOffU botf,    // bottom range elt in BWT index
+	TIndexOffU topb,    // top range elt in BWT' index
+	TIndexOffU botb,    // bottom range elt in BWT' index
+	bool getLock)
+{
+	if(shared_ && getLock) {
+		ThreadSafe ts(mutex_m);
+		return addOnTheFlyImpl(qv, sak, topf, botf, topb, botb);
+	} else {
+		return addOnTheFlyImpl(qv, sak, topf, botf, topb, botb);
+	}
+}
+
 #ifdef ALIGNER_CACHE_MAIN
 
 #include <iostream>

aligner_cache.h

@@ -444,10 +444,8 @@ class AlignmentCache {
 		samap_(CACHE_PAGE_SZ, CA_CAT),
 		salist_(CA_CAT),
 		shared_(shared),
-        mutex_m(),
-		version_(0)
-	{
-	}
+		mutex_m(),
+		version_(0) { }
 
 	/**
 	 * Given a QVal, populate the given EList of SATuples with records
@@ -462,37 +460,11 @@ class AlignmentCache {
 		size_t& nelt,
 		bool getLock = true)
 	{
-        ThreadSafe ts(lockPtr(), shared_ && getLock);
-		assert(qv.repOk(*this));
-		const size_t refi = qv.offset();
-		const size_t reff = refi + qv.numRanges();
-		// For each reference sequence sufficiently similar to the
-		// query sequence in the QKey...
-		for(size_t i = refi; i < reff; i++) {
-			// Get corresponding SAKey, containing similar reference
-			// sequence & length
-			SAKey sak = qlist_.get(i);
-			// Shouldn't have identical keys in qlist_
-			assert(i == refi || qlist_.get(i) != qlist_.get(i-1));
-			// Get corresponding SANode
-			SANode *n = samap_.lookup(sak);
-			assert(n != NULL);
-			const SAVal& sav = n->payload;
-			assert(sav.repOk(*this));
-			if(sav.len > 0) {
-				nrange++;
-				satups.expand();
-				satups.back().init(sak, sav.topf, sav.topb, TSlice(salist_, sav.i, sav.len));
-				nelt += sav.len;
-#ifndef NDEBUG
-				// Shouldn't add consecutive identical entries too satups
-				if(i > refi) {
-					const SATuple b1 = satups.back();
-					const SATuple b2 = satups[satups.size()-2];
-					assert(b1.key != b2.key || b1.topf != b2.topf || b1.offs != b2.offs);
-				}
-#endif
-			}
+		if(shared_ && getLock) {
+			ThreadSafe ts(mutex_m);
+			queryQvalImpl(qv, satups, nrange, nelt);
+		} else {
+			queryQvalImpl(qv, satups, nrange, nelt);
 		}
 	}
 
@@ -522,12 +494,14 @@ class AlignmentCache {
 		bool *added,
 		bool getLock = true)
 	{
-        ThreadSafe ts(lockPtr(), shared_ && getLock);
-		assert(qk.cacheable());
-		QNode *n = qmap_.add(pool(), qk, added);
-		return (n != NULL ? &n->payload : NULL);
+		if(shared_ && getLock) {
+			ThreadSafe ts(mutex_m);
+			return addImpl(qk, added);
+		} else {
+			return addImpl(qk, added);
+		}
 	}
-
+	
 	/**
 	 * Add a new association between a read sequnce ('seq') and a
 	 * reference sequence ('')
@@ -546,8 +520,8 @@ class AlignmentCache {
 	 * ranges in this cache will become invalid and the corresponding
 	 * reads will have to be re-aligned.
 	 */
-	void clear(bool getLock = true) {
-        ThreadSafe ts(lockPtr(), shared_ && getLock);
+	void clear() {
+		ThreadSafe ts(mutex_m);
 		pool_.clear();
 		qmap_.clear();
 		qlist_.clear();
@@ -585,17 +559,9 @@ class AlignmentCache {
 	 * Return the lock object.
 	 */
 	MUTEX_T& lock() {
-	    return mutex_m;
+		return mutex_m;
 	}
 
-	/**
-	 * Return a const pointer to the lock object.  This allows us to
-	 * write const member functions that grab the lock.
-	 */
-	MUTEX_T* lockPtr() const {
-	    return const_cast<MUTEX_T*>(&mutex_m);
-	}
-
 	/**
 	 * Return true iff this cache is shared among threads.
 	 */
@@ -618,6 +584,79 @@ class AlignmentCache {
 	bool     shared_;  // true -> this cache is global
 	MUTEX_T mutex_m;    // mutex used for syncronization in case the the cache is shared.
 	uint32_t version_; // cache version
+
+private:
+
+	template <int S>
+	void queryQvalImpl(
+		const QVal& qv,
+		EList<SATuple, S>& satups,
+		size_t& nrange,
+		size_t& nelt)
+	{
+		assert(qv.repOk(*this));
+		const size_t refi = qv.offset();
+		const size_t reff = refi + qv.numRanges();
+		// For each reference sequence sufficiently similar to the
+		// query sequence in the QKey...
+		for(size_t i = refi; i < reff; i++) {
+			// Get corresponding SAKey, containing similar reference
+			// sequence & length
+			SAKey sak = qlist_.get(i);
+			// Shouldn't have identical keys in qlist_
+			assert(i == refi || qlist_.get(i) != qlist_.get(i-1));
+			// Get corresponding SANode
+			SANode *n = samap_.lookup(sak);
+			assert(n != NULL);
+			const SAVal& sav = n->payload;
+			assert(sav.repOk(*this));
+			if(sav.len > 0) {
+				nrange++;
+				satups.expand();
+				satups.back().init(sak, sav.topf, sav.topb, TSlice(salist_, sav.i, sav.len));
+				nelt += sav.len;
+#ifndef NDEBUG
+				// Shouldn't add consecutive identical entries too satups
+				if(i > refi) {
+					const SATuple b1 = satups.back();
+					const SATuple b2 = satups[satups.size()-2];
+					assert(b1.key != b2.key || b1.topf != b2.topf || b1.offs != b2.offs);
+				}
+#endif
+			}
+		}
+	}
+
+	/**
+	 * Add a new association between a read sequnce ('seq') and a
+	 * reference sequence ('')
+	 */
+	bool addOnTheFlyImpl(
+		QVal& qv,         // qval that points to the range of reference substrings
+		const SAKey& sak, // the key holding the reference substring
+		TIndexOffU topf,    // top range elt in BWT index
+		TIndexOffU botf,    // bottom range elt in BWT index
+		TIndexOffU topb,    // top range elt in BWT' index
+		TIndexOffU botb);   // bottom range elt in BWT' index
+
+	/**
+	 * Add a new query key ('qk'), usually a 2-bit encoded substring of
+	 * the read) as the key in a new Red-Black node in the qmap and
+	 * return a pointer to the node's QVal.
+	 *
+	 * The expectation is that the caller is about to set about finding
+	 * associated reference substrings, and that there will be future
+	 * calls to addOnTheFly to add associations to reference substrings
+	 * found.
+	 */
+	QVal* addImpl(
+		const QKey& qk,
+		bool *added)
+	{
+		assert(qk.cacheable());
+		QNode *n = qmap_.add(pool(), qk, added);
+		return (n != NULL ? &n->payload : NULL);
+	}
 };
 
 /**

aligner_seed.cpp

@@ -1526,7 +1526,7 @@ SeedAligner::searchSeedBi(
 		bool leaveZone = s.zones[i].first < 0;
 		//bool leaveZoneIns = zones_[i].second < 0;
 		Constraint& cons    = *zones[abs(s.zones[i].first)];
-		Constraint& insCons = *zones[abs(s.zones[i].second)];
+		//Constraint& insCons = *zones[abs(s.zones[i].second)];
 		int c = (*seq_)[off];  assert_range(0, 4, c);
 		int q = (*qual_)[off];
 		// Is it legal for us to advance on characters other than 'c'?
@@ -1598,14 +1598,14 @@ SeedAligner::searchSeedBi(
 				}
 				if(cons.canGap() && overall.canGap()) {
 					throw 1; // TODO
-					int delEx = 0;
-					if(cons.canDelete(delEx, *sc_) && overall.canDelete(delEx, *sc_)) {
-						// Try delete
-					}
-					int insEx = 0;
-					if(insCons.canInsert(insEx, *sc_) && overall.canInsert(insEx, *sc_)) {
-						// Try insert
-					}
+//					int delEx = 0;
+//					if(cons.canDelete(delEx, *sc_) && overall.canDelete(delEx, *sc_)) {
+//						// Try delete
+//					}
+//					int insEx = 0;
+//					if(insCons.canInsert(insEx, *sc_) && overall.canInsert(insEx, *sc_)) {
+//						// Try insert
+//					}
 				}
 			} // if(!bail)
 		}

aligner_seed.h

@@ -841,7 +841,7 @@ class SeedResults {
 	size_t numRepeatSeeds() const {
 		return repTot_;
 	}
-
+	
 	/**
 	 * Return fraction of seeds that align repetitively.
 	 */
@@ -1401,7 +1401,6 @@ struct SeedSearchMetrics {
 	 * SeedSearchMetrics object shread by multiple threads.
 	 */
 	void merge(const SeedSearchMetrics& m, bool getLock = false) {
-        ThreadSafe ts(&mutex_m, getLock);
 		seedsearch   += m.seedsearch;
 		nrange       += m.nrange;
 		nelt         += m.nelt;

aligner_seed2.h

@@ -2235,7 +2235,7 @@ class DescentAlignmentSelector {
 		}
 		AlnScore asc(
 			-dr.sink().bestPenalty(),  // numeric score
-			dr.query().length() - edits.size(),
+			(int)(dr.query().length() - edits.size()),
 			(int)edits.size(),         // # edits
 			ns,                        // # Ns
 			ngap);                     // # gaps

aligner_sw_common.h

@@ -209,8 +209,8 @@ struct SwMetrics {
 	 * object.  This is the only safe way to update a SwMetrics shared
 	 * by multiple threads.
 	 */
-	void merge(const SwMetrics& r, bool getLock = false) {
-        ThreadSafe ts(&mutex_m, getLock);
+	void merge(const SwMetrics& r) {
+		ThreadSafe ts(mutex_m);
 		sws        += r.sws;
 		sws10      += r.sws10;
 		sws5       += r.sws5;

aligner_sw_driver.cpp

@@ -326,11 +326,11 @@ void SwDriver::extend(
 
 		// Have to do both because whether we can get through an N depends on
 		// which direction we're coming in
-		bool fwContains = ebwtFw.contains(tmp_rdseq_);
-		tmp_rdseq_.reverse();
-		bool bwContains = ebwtBw != NULL && ebwtBw->contains(tmp_rdseq_);
-		tmp_rdseq_.reverse();
-		assert(fwContains || bwContains);
+//		bool fwContains = ebwtFw.contains(tmp_rdseq_);
+//		tmp_rdseq_.reverse();
+//		bool bwContains = ebwtBw != NULL && ebwtBw->contains(tmp_rdseq_);
+//		tmp_rdseq_.reverse();
+//		assert(fwContains || bwContains);
 	}
 #endif
 	ASSERT_ONLY(tmp_rdseq_.reverse());
@@ -471,11 +471,11 @@ void SwDriver::extend(
 
 		// Have to do both because whether we can get through an N depends on
 		// which direction we're coming in
-		bool fwContains = ebwtFw.contains(tmp_rdseq_);
-		tmp_rdseq_.reverse();
-		bool bwContains = ebwtBw != NULL && ebwtBw->contains(tmp_rdseq_);
-		tmp_rdseq_.reverse();
-		assert(fwContains || bwContains);
+//		bool fwContains = ebwtFw.contains(tmp_rdseq_);
+//		tmp_rdseq_.reverse();
+//		bool bwContains = ebwtBw != NULL && ebwtBw->contains(tmp_rdseq_);
+//		tmp_rdseq_.reverse();
+//		assert(fwContains || bwContains);
 	}
 #endif
 	assert_lt(nlex, rdlen);

aligner_swsse.h

@@ -41,8 +41,7 @@ struct SSEMetrics {
 		corerej = nrej = 0;
 	}
 
-	void merge(const SSEMetrics& o, bool getLock = false) {
-        ThreadSafe ts(&mutex_m, getLock);
+	void merge(const SSEMetrics& o) {
 		dp       += o.dp;
 		dpsat    += o.dpsat;
 		dpfail   += o.dpfail;

aligner_swsse_ee_i16.cpp

@@ -1362,7 +1362,7 @@ bool SwAligner::backtraceNucleotidesEnd2EndSseI16(
 		int readq = (*qu_)[row];
 		assert_leq(col, origCol);
 		// Get score in this cell
-		bool empty, reportedThru, canMoveThru, branch = false;
+		bool empty = false, reportedThru, canMoveThru, branch = false;
 		int cur = SSEMatrix::H;
 		if(!d.mat_.reset_[row]) {
 			d.mat_.resetRow(row);

aligner_swsse_loc_i16.cpp

@@ -1721,7 +1721,7 @@ bool SwAligner::backtraceNucleotidesLocalSseI16(
 		int readq = (*qu_)[row];
 		assert_leq(col, origCol);
 		// Get score in this cell
-		bool empty, reportedThru, canMoveThru, branch = false;
+		bool empty = false, reportedThru, canMoveThru, branch = false;
 		int cur = SSEMatrix::H;
 		if(!d.mat_.reset_[row]) {
 			d.mat_.resetRow(row);

aligner_swsse_loc_u8.cpp

@@ -1719,7 +1719,7 @@ bool SwAligner::backtraceNucleotidesLocalSseU8(
 		int readq = (*qu_)[row];
 		assert_leq(col, origCol);
 		// Get score in this cell
-		bool empty, reportedThru, canMoveThru, branch = false;
+		bool empty = false, reportedThru, canMoveThru, branch = false;
 		int cur = SSEMatrix::H;
 		if(!d.mat_.reset_[row]) {
 			d.mat_.resetRow(row);

aln_sink.h

@@ -121,8 +121,8 @@ struct ReportingMetrics {
 	 * into this object.  This is the only safe way to update a
 	 * ReportingMetrics shared by multiple threads.
 	 */
-	void merge(const ReportingMetrics& met, bool getLock = false) {
-        ThreadSafe ts(&mutex_m, getLock);
+	void merge(const ReportingMetrics& met) {
+		ThreadSafe ts(mutex_m);
 		nread         += met.nread;
 
 		npaired       += met.npaired;
@@ -847,8 +847,8 @@ class AlnSink {
 	/**
 	 * Merge given metrics in with ours by summing all individual metrics.
 	 */
-	void mergeMetrics(const ReportingMetrics& met, bool getLock = true) {
-		met_.merge(met, getLock);
+	void mergeMetrics(const ReportingMetrics& met) {
+		met_.merge(met);
 	}
 
 	/**