Commit ef8e47fc authored by Mikaël Salson's avatar Mikaël Salson

Tools: Compute C(n,k)

parent 91c41992
......@@ -217,3 +217,17 @@ int revcomp_int(int word, int size) {
string reverse(const string &text) {
return string(text.rbegin(), text.rend());
}
double nChoosek(unsigned n, unsigned k)
{
if (k > n) return 0;
if (k * 2 > n) k = n-k;
if (k == 0) return 1;
double result = n;
for( unsigned i = 2; i <= k; ++i ) {
result *= (n-i+1);
result /= i;
}
return result;
}
......@@ -149,6 +149,13 @@ string reverse(const string &text);
*/
Sequence create_sequence(string label_full, string label, string sequence, string quality);
/**
* @return the combinatorial of k among n
* @see http://stackoverflow.com/a/9331125/1192742
*/
double nChoosek(unsigned n, unsigned k);
const Sequence NULL_SEQUENCE = create_sequence("", "", "NULL", "");
bool operator==(const Sequence &s1, const Sequence &s2);
......
......@@ -275,6 +275,13 @@ void testExtractBasename() {
TEST_EXTRACT_BASENAME, extract_basename("/", true));
}
void testNChooseK() {
TAP_TEST(nChoosek(1, 10) == 0, TEST_N_CHOOSE_K, "");
TAP_TEST(nChoosek(1, 1) == 1, TEST_N_CHOOSE_K, "");
TAP_TEST(nChoosek(5, 2) == 10, TEST_N_CHOOSE_K, "");
TAP_TEST(nChoosek(8, 4) == 70, TEST_N_CHOOSE_K, "");
}
void testTools() {
testOnlineFasta1();
testFasta1();
......@@ -289,4 +296,5 @@ void testTools() {
testRevcompInt();
testExtendedNucleotides();
testExtractBasename();
testNChooseK();
}
......@@ -24,6 +24,7 @@ enum {
TEST_REVCOMP_INT,
TEST_EXTENDED_NUCL,
TEST_EXTRACT_BASENAME,
TEST_N_CHOOSE_K,
/* Storage tests */
TEST_ARRAY_KMERSTORE,
......@@ -171,6 +172,7 @@ inline void declare_tests() {
RECORD_TAP_TEST(TEST_DNA_TO_INT, "dna_to_int()");
RECORD_TAP_TEST(TEST_REVCOMP_INT, "revcomp_int()");
RECORD_TAP_TEST(TEST_EXTRACT_BASENAME, "extractBasename()");
RECORD_TAP_TEST(TEST_N_CHOOSE_K, "test nChooseK()");
RECORD_TAP_TEST(TEST_ARRAY_KMERSTORE, "Testing ArrayKmerStore");
RECORD_TAP_TEST(TEST_KMERSTORE_INSERT_ONE_SEQ, "Testing IKmerStore::insert() on one sequence");
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment