File: RollingHash.h

package info (click to toggle)
abyss 2.3.10-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 8,284 kB
  • sloc: cpp: 78,182; ansic: 6,512; makefile: 2,252; perl: 672; sh: 509; haskell: 412; python: 4
file content (222 lines) | stat: -rw-r--r-- 5,974 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
#ifndef ABYSS_ROLLING_HASH_H
#define ABYSS_ROLLING_HASH_H 1

#include "config.h"

#include "BloomDBG/LightweightKmer.h"
#include "BloomDBG/MaskedKmer.h"
#include "Common/Sense.h"
#include "vendor/nthash/nthash.hpp"

#include <algorithm>
#include <string>
#include <vector>
#include <cassert>
#include <boost/dynamic_bitset.hpp>
#include <cstring>

class RollingHash
{
private:

	typedef uint64_t hash_t;

	/**
	 * Determine the canonical hash value, given hash values for
	 * forward and reverse-complement of the same k-mer.
	 */
	hash_t canonicalHash(hash_t hash, hash_t rcHash) const
	{
		return (rcHash < hash) ? rcHash : hash;
	}

public:

	/**
	 * Default constructor.
	 */
	RollingHash() : m_numHashes(0), m_k(0), m_hash1(0), m_rcHash1(0) {}

	/**
	 * Constructor. Construct RollingHash object when initial k-mer
	 * is unknown.
	 * @param numHashes number of pseudo-independent hash values to compute
	 * for each k-mer
	 * @param k k-mer length
	 */
	RollingHash(unsigned numHashes, unsigned k) : m_numHashes(numHashes),
		m_k(k), m_hash1(0), m_rcHash1(0) {}

	/**
	 * Constructor. Construct RollingHash object while specifying
	 * initial k-mer to be hashed.
	 * @param kmer initial k-mer for initializing hash value(s)
	 * @param numHashes number of pseudo-independent hash values to compute
	 * for each k-mer
	 * @param k k-mer length
	 */
	RollingHash(const std::string& kmer, unsigned numHashes, unsigned k)
		: m_numHashes(numHashes), m_k(k), m_hash1(0), m_rcHash1(0)
	{
		/* init rolling hash state */
		reset(kmer);
	}

	/**
	 * Initialize hash state from sequence.
	 * @param kmer k-mer used to initialize hash state
	 */
	void reset(const std::string& kmer)
	{
		/* 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);

		if (!MaskedKmer::mask().empty())
			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
				kmer.c_str(), m_k);
	}

	/**
	 * 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)
	{
		NTC64(kmer[0], 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(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.
	 * @param prevKmer k-mer we are rolling into
	 * @param kmer current k-mer
	 */
	void rollLeft(char charIn, const char* kmer)
	{
		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);
		}
	}

	/**
	 * Get the seed hash value for the current k-mer. The seed hash
	 * value is used to calculate multiple pseudo-independant
	 * hash functions.
	 */
	size_t getHashSeed() const
	{
		return (size_t)m_hash;
	}

	/**
	 * Get hash values for current k-mer.
	 *
	 * @param hashes array for returned hash values
	 */
	void getHashes(hash_t hashes[]) const
	{
		hashes[0] = m_hash;
		for (unsigned i = 1; i < m_numHashes; ++i)
			hashes[i] = NTE64(m_hash, m_k, i);
	}

	/** Equality operator */
	bool operator==(const RollingHash& o) const
	{
		/**
		 * Note: If hash seeds are equal, then the values
		 * for all hash functions will also be equal, since
		 * the hash values are calculated from the
		 * seed in a deterministic manner. In practice seed
		 * collision is very unlikely, though!
		 */
		return m_k == o.m_k && getHashSeed() == o.getHashSeed();
	}

	/** Inequality operator */
	bool operator!=(const RollingHash& o) const
	{
		return !(*this == o);
	}

	/**
	 * 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 dir if SENSE, change last base; if ANTISENSE,
	 * change first base
	 * @param base new value for the base
	 */
	void setLastBase(char* kmer, extDirection dir, char base)
	{
		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);

		if (!MaskedKmer::mask().empty())
			m_hash = maskHash(m_hash1, m_rcHash1, MaskedKmer::mask().c_str(),
				kmer, m_k);
	}

	/**
	 * Reverse complement the hash state, so that rolling right becomes
	 * rolling left, and vice versa.  This operation is needed
	 * whenever we reverse-complement a k-mer that has an associated
	 * `RollingHash` state, so that subsequent rolling operations will
	 * produce the correct hash value.
	 */
	void reverseComplement()
	{
		std::swap(m_hash1, m_rcHash1);
	}

private:

	/** number of hash functions */
	unsigned m_numHashes;
	/** k-mer length */
	unsigned m_k;
	/** value of first hash function for current k-mer */
	hash_t m_hash1;
	/** value of first hash function for current k-mer, after
	 * reverse-complementing */
	hash_t m_rcHash1;
	/** current canonical hash value */
	hash_t m_hash;
};

#endif