.circleci/config.yml

@@ -2,7 +2,7 @@ version: 2
 jobs:
   build:
     docker:
-      - image: ubuntu:xenial
+      - image: ubuntu:bionic
     steps:
       - run: |
           apt-get update -qq

.github/ISSUE_TEMPLATE.rb

+# Please report
+
+- [ ] version of ABySS with `abyss-pe version`
+- [ ] distribution of Linux with `lsb_release -d`
+
+# Assembly error
+
+- [ ] complete `abyss-pe` command line
+- [ ] last 20 lines of the output of `abyss-pe`
+- [ ] number of sequenced bases
+- [ ] estimated genome size and ploidy
+- [ ] estimated sequencing depth of coverage
+
+# Build error
+
+Consider installing ABySS using [Linuxbrew](https://linuxbrew.sh) on Linux or [Homebrew](https://brew.sh) on macOS with `brew install abyss`, or using [Bioconda](https://bioconda.github.io) with `conda install abyss`.
+
+- [ ] Have you tried installing ABySS using Brew or Bioconda?
+- [ ] version of GCC or compiler with `gcc --version`
+- [ ] complete `./configure` command line
+- [ ] last 20 lines of the output of `./configure`
+- [ ] last 20 lines of the output of `make`

ABYSS/abyss.cc

@@ -73,6 +73,15 @@ static void assemble(const string& pathIn, const string& pathOut)
 	AssemblyAlgorithms::setCoverageParameters(
 			AssemblyAlgorithms::coverageHistogram(g));
 
+	if (opt::kc > 0) {
+		cout << "Minimum k-mer multiplicity kc is " << opt::kc << endl;
+		cout << "Removing low-multiplicity k-mers" << endl;
+		size_t removed = AssemblyAlgorithms::applyKmerCoverageThreshold(g, opt::kc);
+		cout << "Removed " <<  removed
+			<< " low-multiplicity k-mers, " << g.size()
+			<< " k-mers remaining" << std::endl;
+	}
+
 	cout << "Generating adjacency" << endl;
 	AssemblyAlgorithms::generateAdjacency(&g);
 

Assembly/BranchGroup.h

@@ -2,6 +2,7 @@
 #define BRANCHGROUP_H 1
 
 #include "Common/Algorithms.h"
+#include "Common/Exception.h"
 #include <algorithm> // for swap
 #include <map>
 #include <utility>
@@ -211,7 +212,7 @@ isAmbiguous(const SequenceCollectionHash& g) const
 
 namespace std {
 	template <>
-	inline void swap(BranchGroup&, BranchGroup&) { assert(false); }
+	inline void swap(BranchGroup&, BranchGroup&) NOEXCEPT { assert(false); }
 }
 
 #endif

Assembly/BranchRecordBase.h

 #ifndef ASSEMBLY_BRANCHRECORDBASE_H
 #define ASSEMBLY_BRANCHRECORDBASE_H 1
 
+#include "Common/Exception.h"
+
 #include <algorithm>
 #include <cassert>
 #include <utility>
@@ -172,7 +174,7 @@ operator Sequence() const
 
 namespace std {
 	template <>
-	inline void swap(BranchRecord& a, BranchRecord& b)
+	inline void swap(BranchRecord& a, BranchRecord& b) NOEXCEPT
 	{
 		a.swap(b);
 	}

Assembly/CoverageAlgorithm.h

@@ -112,6 +112,22 @@ void setCoverageParameters(const Histogram& h)
 	}
 }
 
+/** Remove all k-mers with multiplicity lower than the given threshold */
+static inline
+size_t applyKmerCoverageThreshold(SequenceCollectionHash& c, unsigned kc)
+{
+	if (kc == 0)
+		return 0;
+
+	for (SequenceCollectionHash::iterator it = c.begin();
+		it != c.end(); ++it) {
+		if (it->second.getMultiplicity() < kc)
+			it->second.setFlag(SF_DELETE);
+	}
+
+	return c.cleanup();
+}
+
 } // namespace AssemblyAlgorithms
 
 #endif

Assembly/Options.cc

@@ -54,6 +54,7 @@ static const char USAGE_MESSAGE[] =
 "  -t, --trim-length=N   maximum length of blunt contigs to trim [k]\n"
 "  -c, --coverage=FLOAT  remove contigs with mean k-mer coverage\n"
 "                        less than this threshold\n"
+"      --kc=N            remove all k-mers with multiplicity < N [0]\n"
 "  -b, --bubbles=N       pop bubbles shorter than N bp [3*k]\n"
 "  -b0, --no-bubbles     do not pop bubbles\n"
 "  -e, --erode=N         erode bases at the ends of blunt contigs with coverage\n"
@@ -107,6 +108,9 @@ int trimLen = -1;
 /** Coverage cutoff. */
 float coverage = -1;
 
+/** Minimum k-mer multiplicity cutoff. */
+unsigned kc = 0;
+
 /** Pop bubbles shorter than N bp. */
 int bubbleLen = -1;
 
@@ -147,7 +151,7 @@ string assemblyCmd;
 
 static const char shortopts[] = "b:c:e:E:g:k:K:mo:Q:q:s:t:v";
 
-enum { OPT_HELP = 1, OPT_VERSION, COVERAGE_HIST, OPT_DB, OPT_LIBRARY, OPT_STRAIN, OPT_SPECIES };
+enum { OPT_HELP = 1, OPT_VERSION, COVERAGE_HIST, OPT_DB, OPT_LIBRARY, OPT_STRAIN, OPT_SPECIES, OPT_KC };
 
 static const struct option longopts[] = {
 	{ "out",         required_argument, NULL, 'o' },
@@ -164,6 +168,7 @@ static const struct option longopts[] = {
 	{ "SS",          no_argument,       &opt::ss, 1 },
 	{ "no-SS",       no_argument,       &opt::ss, 0 },
 	{ "coverage",    required_argument, NULL, 'c' },
+	{ "kc",          required_argument, NULL, OPT_KC },
 	{ "coverage-hist", required_argument, NULL, COVERAGE_HIST },
 	{ "bubble-length", required_argument, NULL, 'b' },
 	{ "no-bubbles",  no_argument,       &opt::bubbleLen, 0 },
@@ -302,6 +307,9 @@ void parse(int argc, char* const* argv)
 			case OPT_SPECIES:
 				arg >> opt::metaVars[2];
 				break;
+			case OPT_KC:
+				arg >> opt::kc;
+				break;
 		}
 		if (optarg != NULL && !arg.eof()) {
 			cerr << PROGRAM ": invalid option: `-"

Assembly/Options.h

@@ -14,6 +14,7 @@ namespace opt {
 	extern unsigned erodeStrand;
 	extern unsigned trimLen;
 	extern float coverage;
+	extern unsigned kc;
 	extern unsigned bubbleLen;
 	extern unsigned ss;
 	extern bool maskCov;

Bloom/BloomFilter.h

@@ -14,6 +14,12 @@
 #include <iostream>
 #include <boost/dynamic_bitset.hpp>
 
+/*
+ * Put `BloomFilter` class in `Konnector` namespace to avoid collision with BTL
+ * `BloomFilter` class of the same name.
+ */
+namespace Konnector {
+
 /** A Bloom filter. */
 class BloomFilter
 {
@@ -165,4 +171,6 @@ class BloomFilter
 	char* m_array;
 };
 
+} // end Konnector namespace
+
 #endif

Bloom/BloomFilterWindow.h

@@ -14,12 +14,12 @@
  * A bloom filter that represents a window
  * within a larger bloom filter.
  */
-class BloomFilterWindow : public BloomFilter
+class BloomFilterWindow : public Konnector::BloomFilter
 {
 public:
 
 	/** Constructor. */
-	BloomFilterWindow() : BloomFilter() { };
+	BloomFilterWindow() : Konnector::BloomFilter() { };
 
 	/** Constructor.
 	 *
@@ -29,7 +29,7 @@ public:
 	 */
 	BloomFilterWindow(size_t fullBloomSize, size_t startBitPos,
 			size_t endBitPos, size_t hashSeed=0) :
-		BloomFilter(endBitPos - startBitPos + 1, hashSeed),
+		Konnector::BloomFilter(endBitPos - startBitPos + 1, hashSeed),
 		m_fullBloomSize(fullBloomSize),
 		m_startBitPos(startBitPos),
 		m_endBitPos(endBitPos)
@@ -63,26 +63,26 @@ public:
 	/** Return the size of the bit array. */
 	size_t size() const
 	{
-		return BloomFilter::size();
+		return Konnector::BloomFilter::size();
 	}
 
 	/** Return the number of elements with count >= max_count. */
 	size_t popcount() const
 	{
-		return BloomFilter::popcount();
+		return Konnector::BloomFilter::popcount();
 	}
 
 	/** Return the estimated false positive rate */
 	double FPR() const
 	{
-		return BloomFilter::FPR();
+		return Konnector::BloomFilter::FPR();
 	}
 
 	/** Return whether the specified bit is set. */
 	bool operator[](size_t i) const
 	{
 		if (i >= m_startBitPos && i <= m_endBitPos)
-			return BloomFilter::operator[](i - m_startBitPos);
+			return Konnector::BloomFilter::operator[](i - m_startBitPos);
 		return false;
 	}
 
@@ -96,7 +96,7 @@ public:
 	void insert(size_t i)
 	{
 		if (i >= m_startBitPos && i <= m_endBitPos)
-			BloomFilter::insert(i - m_startBitPos);
+			Konnector::BloomFilter::insert(i - m_startBitPos);
 	}
 
 	/** Add the object to this set. */
@@ -143,7 +143,7 @@ public:
 
 		if (m_size != bits) {
 			if (readOp == BITWISE_OVERWRITE) {
-				BloomFilter::resize(bits);
+				Konnector::BloomFilter::resize(bits);
 			} else {
 				std::cerr << "error: can't union/intersect bloom filters with "
 					<< "different sizes\n";

Bloom/CascadingBloomFilter.h

@@ -22,13 +22,13 @@ class CascadingBloomFilter
 	{
 		m_data.reserve(max_count);
 		for (unsigned i = 0; i < max_count; i++)
-			m_data.push_back(new BloomFilter(n, hashSeed));
+			m_data.push_back(new Konnector::BloomFilter(n, hashSeed));
 	}
 
 	/** Destructor */
 	~CascadingBloomFilter()
 	{
-		typedef std::vector<BloomFilter*>::iterator Iterator;
+		typedef std::vector<Konnector::BloomFilter*>::iterator Iterator;
 		for (Iterator i = m_data.begin(); i != m_data.end(); i++) {
 			assert(*i != NULL);
 			delete *i;
@@ -91,7 +91,7 @@ class CascadingBloomFilter
 	}
 
 	/** Get the Bloom filter for a given level */
-	BloomFilter& getBloomFilter(unsigned level)
+	Konnector::BloomFilter& getBloomFilter(unsigned level)
 	{
 		assert(m_data.at(level) != NULL);
 		return *m_data.at(level);
@@ -112,7 +112,7 @@ class CascadingBloomFilter
 
   private:
 	size_t m_hashSeed;
-	std::vector<BloomFilter*> m_data;
+	std::vector<Konnector::BloomFilter*> m_data;
 
 };
 

Bloom/bloom.cc

@@ -205,7 +205,8 @@ static const struct option longopts[] = {
 	{ NULL, 0, NULL, 0 }
 };
 
-void dieWithUsageError()
+__attribute__((noreturn))
+static void dieWithUsageError()
 {
 	cerr << "Try `" << PROGRAM
 		<< " --help' for more information.\n";
@@ -418,9 +419,9 @@ static inline void buildKonnectorBloom(size_t bits, string outputPath,
 	if (opt::windows == 0) {
 
 		if (opt::levels == 1) {
-			BloomFilter bloom(bits, opt::hashSeed);
+			Konnector::BloomFilter bloom(bits, opt::hashSeed);
 #ifdef _OPENMP
-			ConcurrentBloomFilter<BloomFilter>
+			ConcurrentBloomFilter<Konnector::BloomFilter>
 				cbf(bloom, opt::numLocks, opt::hashSeed);
 			loadFilters(cbf, argc, argv);
 #else
@@ -640,7 +641,7 @@ int combine(int argc, char** argv, BitwiseOp readOp)
 	string outputPath(argv[optind]);
 	optind++;
 
-	BloomFilter bloom;
+	Konnector::BloomFilter bloom;
 
 	for (int i = optind; i < argc; i++) {
 		string path(argv[i]);
@@ -695,7 +696,7 @@ int info(int argc, char** argv)
 		dieWithUsageError();
 	}
 
-	BloomFilter bloom;
+	Konnector::BloomFilter bloom;
 	string path = argv[optind];
 
 	if (opt::verbose)
@@ -744,9 +745,9 @@ int compare(int argc, char ** argv){
 	std::cerr << "Computing distance for 2"
 			  << " samples...\n";
 	// Get both paths and open istreams
-	BloomFilter bloomA;
+	Konnector::BloomFilter bloomA;
 	string pathA(argv[optind]);
-	BloomFilter bloomB;
+	Konnector::BloomFilter bloomB;
 	string pathB(argv[optind+1]);
 	if (opt::verbose)
 	  std::cerr << "Loading bloom filters from "
@@ -845,7 +846,7 @@ int compare(int argc, char ** argv){
 
 int memberOf(int argc, char ** argv){
 	// Initalise bloom and get globals
-	BloomFilter bloom;
+	Konnector::BloomFilter bloom;
 	parseGlobalOpts(argc, argv);
 	// Arg parser to get `m' option in case set
 	for (int c; (c = getopt_long(argc, argv,
@@ -929,11 +930,11 @@ int memberOf(int argc, char ** argv){
 /**
  * Calculate number of bases to trim from left end of sequence.
  */
-int calcLeftTrim(const Sequence& seq, unsigned k, const BloomFilter& bloom,
+int calcLeftTrim(const Sequence& seq, unsigned k, const Konnector::BloomFilter& bloom,
 	size_t minBranchLen)
 {
 	// Boost graph interface for Bloom filter
-	DBGBloom<BloomFilter> g(bloom);
+	DBGBloom<Konnector::BloomFilter> g(bloom);
 
 	// if this is the first k-mer we have found in
 	// Bloom filter, starting from the left end
@@ -1000,7 +1001,7 @@ int trim(int argc, char** argv)
 		cerr << "Loading bloom filter from `"
 			<< bloomPath << "'...\n";
 
-	BloomFilter bloom;
+	Konnector::BloomFilter bloom;
 	istream *in = openInputStream(bloomPath);
 	assert_good(*in, bloomPath);
 	bloom.read(*in);

BloomDBG/HashAgnosticCascadingBloom.h

@@ -42,12 +42,12 @@ class HashAgnosticCascadingBloom
 	{
 		m_data.reserve(levels);
 		for (unsigned i = 0; i < levels; i++)
-			m_data.push_back(new BTL::BloomFilter(size, hashes, k));
+			m_data.push_back(new BloomFilter(size, hashes, k));
 	}
 
 	/**
-	 * Constructor to load a single-level BTL::BloomFilter from
-	 * files.  This is used to make BTL::BloomFilter support the
+	 * Constructor to load a single-level BloomFilter from
+	 * files.  This is used to make BloomFilter support the
 	 * same interface as HashAgnosticCascadingBloom.
 	 */
 	HashAgnosticCascadingBloom(const string& bloomPath)
@@ -138,7 +138,7 @@ class HashAgnosticCascadingBloom
 	}
 
 	/** Get the Bloom filter for a given level */
-	BTL::BloomFilter& getBloomFilter(unsigned level)
+	BloomFilter& getBloomFilter(unsigned level)
 	{
 		assert(m_data.at(level) != NULL);
 		return *m_data.at(level);
@@ -159,7 +159,7 @@ class HashAgnosticCascadingBloom
 	void loadFilter(const string& bloomPath)
 	{
 		clear();
-		BTL::BloomFilter* bloom = new BTL::BloomFilter(bloomPath);
+		BloomFilter* bloom = new BloomFilter(bloomPath);
 		m_k = bloom->getKmerSize();
 		m_hashes = bloom->getHashNum();
 		m_data.push_back(bloom);
@@ -172,7 +172,7 @@ class HashAgnosticCascadingBloom
 	{
 		m_k = 0;
 		m_hashes = 0;
-		typedef std::vector<BTL::BloomFilter*>::iterator Iterator;
+		typedef std::vector<BloomFilter*>::iterator Iterator;
 		for (Iterator i = m_data.begin(); i != m_data.end(); i++) {
 			assert(*i != NULL);
 			delete *i;
@@ -185,7 +185,7 @@ class HashAgnosticCascadingBloom
 	/** number of hash functions */
 	unsigned m_hashes;
 	/** the array of Bloom filters */
-	std::vector<BTL::BloomFilter*> m_data;
+	std::vector<BloomFilter*> m_data;
 };
 
 #endif

BloomDBG/LightweightKmer.h

 #ifndef LIGHTWEIGHT_KMER_H
 #define LIGHTWEIGHT_KMER_H 1
 
+#include "BloomDBG/MaskedKmer.h"
+
 #include <algorithm>
 #include <cstring>
 #include <boost/shared_array.hpp>

BloomDBG/Makefile.am

@@ -26,4 +26,4 @@ abyss_bloom_dbg_SOURCES = \
 	RollingHashIterator.h \
 	SpacedSeed.h \
 	$(top_srcdir)/lib/bloomfilter/BloomFilter.hpp \
-	$(top_srcdir)/lib/rolling-hash/rolling.h
+	$(top_srcdir)/lib/nthash/nthash.hpp

BloomDBG/RollingBloomDBG.h

@@ -45,6 +45,8 @@ public:
 		: m_kmer(kmer), m_rollingHash(rollingHash) {}
 
 	const LightweightKmer& kmer() const { return m_kmer; };
+	LightweightKmer& kmer() { return m_kmer; };
+
 	const RollingHash& rollingHash() const { return m_rollingHash; }
 
 	RollingBloomDBGVertex clone() const {
@@ -63,12 +65,8 @@ public:
 
 	void setLastBase(extDirection dir, char base)
 	{
-		const unsigned k = Kmer::length();
-		if (dir == SENSE) {
-			m_rollingHash.setBase(m_kmer.c_str(), k-1, base);
-		} else {
-			m_rollingHash.setBase(m_kmer.c_str(), 0, base);
-		}
+		m_rollingHash.setLastBase(kmer().c_str(), dir, base);
+		kmer().setLastBase(dir, base);
 	}
 
 	/**

BloomDBG/RollingHash.h

@@ -2,8 +2,13 @@
 #define ABYSS_ROLLING_HASH_H 1
 
 #include "config.h"
-#include "lib/rolling-hash/rolling.h"
+
+#include "BloomDBG/LightweightKmer.h"
 #include "BloomDBG/MaskedKmer.h"
+#include "Common/Sense.h"
+#include "lib/nthash/nthash.hpp"
+
+#include <algorithm>
 #include <string>
 #include <vector>
 #include <cassert>
@@ -61,74 +66,15 @@ public:
 	 */
 	void reset(const std::string& kmer)
 	{
-		if (!MaskedKmer::mask().empty())
-			resetMasked(kmer.c_str());
-		else
-			resetUnmasked(kmer);
-	}
-
-	/**
-	 * Initialize hash values from current k-mer. When computing the hash
-	 * value, mask out "don't care" positions as per the active
-	 * k-mer mask.
-	 */
-	void resetMasked(const char* kmer)
-	{
-		const std::string& spacedSeed = MaskedKmer::mask();
-		assert(spacedSeed.length() == m_k);
-
-		/* compute first hash function for k-mer */
-		uint64_t hash1 = getFhval(m_hash1, spacedSeed.c_str(), kmer, m_k);
-
-		/* compute first hash function for reverse complement of k-mer */
-		uint64_t rcHash1 = getRhval(m_rcHash1, spacedSeed.c_str(), kmer, m_k);
-
-		m_hash = canonicalHash(hash1, rcHash1);
-	}
-
-	/**
-	 * Initialize hash values from sequence.
-	 * @param kmer k-mer used to initialize hash state
-	 */
-	void resetUnmasked(const std::string& kmer)
-	{
-		/* compute first hash function for k-mer */
-		m_hash1 = getFhval(kmer.c_str(), m_k);
-
-		/* compute first hash function for reverse complement
-		 * of k-mer */
-		m_rcHash1 = getRhval(kmer.c_str(), m_k);
+		/* compute initial hash values for forward and reverse-complement k-mer */
+		NTC64(kmer.c_str(), m_k, m_hash1, m_rcHash1);
 
+		/* get canonical hash value from forward/reverse hash values */
 		m_hash = canonicalHash(m_hash1, m_rcHash1);
-	}
 
-	/**
-	 * Compute hash values for next k-mer to the right and
-	 * update internal state.
-	 * @param kmer current k-mer
-	 * @param nextKmer k-mer we are rolling into
-	 */
-	void rollRight(const char* kmer, char charIn)
-	{
 		if (!MaskedKmer::mask().empty())
-			rollRightMasked(kmer, charIn);
-		else
-			rollRightUnmasked(kmer, charIn);
-	}
-
-	/**
-	 * Compute hash values for next k-mer to the right and
-	 * update internal state.  When computing the new hash, mask
-	 * out "don't care" positions according to the active
-	 * k-mer mask.
-	 * @param kmer current k-mer
-	 * @param nextKmer k-mer we are rolling into
-	 */
-	void rollRightMasked(const char* kmer, char charIn)
-	{
-		const std::string& spacedSeed = MaskedKmer::mask();
-		m_hash = rollHashesRight(m_hash1, m_rcHash1, spacedSeed.c_str(),
-			kmer, charIn, m_k);
+			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
+				kmer.c_str(), m_k);
 	}
 
 	/**
@@ -137,40 +83,20 @@ public:
 	 * @param kmer current k-mer
 	 * @param nextKmer k-mer we are rolling into
 	 */
-	void rollRightUnmasked(const char* kmer, char charIn)
+	void rollRight(const char* kmer, char charIn)
 	{
-		/* update first hash function */
-		rollHashesRight(m_hash1, m_rcHash1, kmer[0], charIn, m_k);
+		NTC64(kmer[0], charIn, m_k, m_hash1, m_rcHash1);
 		m_hash = canonicalHash(m_hash1, m_rcHash1);
-	}
 
-	/**
-	 * Compute hash values for next k-mer to the left and
-	 * update internal state.
-	 * @param prevKmer k-mer we are rolling into
-	 * @param kmer current k-mer
-	 */
-	void rollLeft(char charIn, const char* kmer)
-	{
-		if (!MaskedKmer::mask().empty())
-			rollLeftMasked(charIn, kmer);
-		else
-			rollLeftUnmasked(charIn, kmer);
+		if (!MaskedKmer::mask().empty()) {
+			// TODO: copying the k-mer and shifting is very inefficient;
+			// we need a specialized nthash function that rolls and masks
+			// simultaneously
+			LightweightKmer next(kmer);
+			next.shift(SENSE, charIn);
+			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
+				next.c_str(), m_k);
 		}
-
-	/**
-	 * Compute hash values for next k-mer to the left and
-	 * update internal state.  When computing the new hash, mask
-	 * out "don't care" positions according to the active
-	 * k-mer mask.
-	 * @param prevKmer k-mer we are rolling into
-	 * @param kmer current k-mer
-	 */
-	void rollLeftMasked(char charIn, const char* kmer)
-	{
-		const std::string& spacedSeed = MaskedKmer::mask();
-		m_hash = rollHashesLeft(m_hash1, m_rcHash1, spacedSeed.c_str(),
-			kmer, charIn, m_k);
 	}
 
 	/**
@@ -179,11 +105,20 @@ public:
 	 * @param prevKmer k-mer we are rolling into
 	 * @param kmer current k-mer
 	 */
-	void rollLeftUnmasked(char charIn, const char* kmer)
+	void rollLeft(char charIn, const char* kmer)
 	{
-		/* update first hash function */
-		rollHashesLeft(m_hash1, m_rcHash1, charIn, kmer[m_k-1], m_k);
+		NTC64L(kmer[m_k-1], charIn, m_k, m_hash1, m_rcHash1);
 		m_hash = canonicalHash(m_hash1, m_rcHash1);
+
+		if (!MaskedKmer::mask().empty()) {
+			// TODO: copying the k-mer and shifting is very inefficient;
+			// we need a specialized nthash function that rolls and masks
+			// simultaneously
+			LightweightKmer next(kmer);
+			next.shift(ANTISENSE, charIn);
+			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
+				next.c_str(), m_k);
+		}
 	}
 
 	/**
@@ -203,11 +138,8 @@ public:
 	 */
 	void getHashes(size_t hashes[]) const
 	{
-		uint64_t tmpHashes[MAX_HASHES];
-		multiHash(tmpHashes, m_hash, m_numHashes, m_k);
-		for (unsigned i = 0; i < m_numHashes; ++i) {
-			hashes[i] = (size_t)tmpHashes[i];
-		}
+		for (unsigned i = 0; i < m_numHashes; ++i)
+			hashes[i] = NTE64(m_hash, m_k, i);
 	}
 
 	/** Equality operator */
@@ -230,45 +162,31 @@ public:
 	}
 
 	/**
-	 * Set the base at a given position in the k-mer and update the hash
-	 * value accordingly.
-	 * @param kmer point to the k-mer char array
-	 * @param pos position of the base to be changed
-	 * @param base new value for the base
-	 */
-	void setBase(char* kmer, unsigned pos, char base)
-	{
-		if (!MaskedKmer::mask().empty())
-			setBaseMasked(kmer, pos, base);
-		else
-			setBaseUnmasked(kmer, pos, base);
-	}
-
-	/**
-	 * Set the base at a given position in the k-mer and update the hash
-	 * value accordingly.
+	 * Change the hash value to reflect a change in the first/last base of
+	 * the k-mer.
 	 * @param kmer point to the k-mer char array
-	 * @param pos position of the base to be changed
+	 * @param dir if SENSE, change last base; if ANTISENSE,
+	 * change first base
 	 * @param base new value for the base
 	 */
-	void setBaseMasked(char* kmer, unsigned pos, char base)
+	void setLastBase(char* kmer, extDirection dir, char base)
 	{
-		const std::string& spacedSeed = MaskedKmer::mask();
-		assert(spacedSeed.length() == m_k);
-		m_hash = ::setBase(m_hash1, m_rcHash1, spacedSeed.c_str(), kmer,
-			pos, base, m_k);
+		if (dir == SENSE) {
+			/* roll left to remove old last char */
+			NTC64L(kmer[m_k-1], 'A', m_k, m_hash1, m_rcHash1);
+			/* roll right to add new last char */
+			NTC64('A', base, m_k, m_hash1, m_rcHash1);
+		} else {
+			/* roll right to remove old first char */
+			NTC64(kmer[0], 'A', m_k, m_hash1, m_rcHash1);
+			/* roll left to add new first char */
+			NTC64L('A', base, m_k, m_hash1, m_rcHash1);
 		}
+		m_hash = canonicalHash(m_hash1, m_rcHash1);
 
-	/**
-	 * Set the base at a given position in the k-mer and update the hash
-	 * value accordingly.
-	 * @param kmer point to the k-mer char array
-	 * @param pos position of the base to be changed
-	 * @param base new value for the base
-	 */
-	void setBaseUnmasked(char* kmer, unsigned pos, char base)
-	{
-		m_hash = ::setBase(m_hash1, m_rcHash1, kmer, pos, base, m_k);
+		if (!MaskedKmer::mask().empty())
+			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
+				kmer, m_k);
 	}
 
 private:

BloomDBG/bloom-dbg.cc

@@ -228,7 +228,7 @@ void resumeAssemblyFromCheckpoint(int argc, char** argv,
 	HashAgnosticCascadingBloom solidKmerSet;
 
 	/* empty visited k-mers Bloom filter */
-	BTL::BloomFilter visitedKmerSet;
+	BloomFilter visitedKmerSet;
 
 	/* counters for progress messages */
 	BloomDBG::AssemblyCounters counters;

BloomDBG/bloom-dbg.h

@@ -574,7 +574,7 @@ namespace BloomDBG {
 			if (!assembledKmerSet.contains(*fwd))
 				break;
 		}
-		if (fwd.pos() > 0)
+		if (fwd != RollingHashIterator::end() && fwd.pos() > 0)
 			seq.erase(0, fwd.pos());
 
 		/* trim previously assembled k-mers from end of sequence */
@@ -584,7 +584,7 @@ namespace BloomDBG {
 			if (!assembledKmerSet.contains(*rev))
 				break;
 		}
-		if (rev.pos() > 0)
+		if (rev != RollingHashIterator::end() && rev.pos() > 0)
 			rcSeq.erase(0, rev.pos());
 
 		/* flip seq back to original orientation */
@@ -834,7 +834,7 @@ namespace BloomDBG {
 			std::ostream& out)
 	{
 		/* k-mers in previously assembled contigs */
-		BTL::BloomFilter visitedKmerSet(solidKmerSet.size(),
+		BloomFilter visitedKmerSet(solidKmerSet.size(),
 			solidKmerSet.getHashNum(), solidKmerSet.getKmerSize());
 
 		/* counters for progress messages */

ChangeLog

+2018-09-11  Ben Vandervalk <benv@bcgsc.ca>
+
+	* Release version 2.1.1
+
+	abyss-bloom-dbg:
+	* upgrade to most recent version of ntHash to reduce
+	some assembly/hashing artifacts. On a human assembly, this
+	reduced QUAST major misassemblies by 5% and increased
+	scaffold contiguity by 10%
+	* `kc` parameter now also applies to MPI assemblies (see below)
+
+	abyss-fac:
+	* change N20 and N80 to N25 and N75, respectively
+
+	ABYSS-P:
+	* add `--kc` option, with implements a hard minimum k-mer
+	multiplicity cutoff
+
+	abyss-pe:
+	* fix `zsh: no such option: pipefail` error with
+	old versions of `zsh` (fallback to `bash` instead)
+	* adding `time=1` now times *all* assembly commands
+
+	abyss-sealer:
+	* parallelize gap sealing with OpenMP (thanks to
+	@schutzekatze!)
+	* add `--gap-file` option (thanks to @schutzekatze!)
+
+	DistanceEst:
+	* add support for GFA output
+
 2018-04-13  Ben Vandervalk <benv@bcgsc.ca>
 
 	* Release version 2.1.0