5 * Created by Pat Schloss on 12/15/08.
6 * Copyright 2008 Patrick D. Schloss. All rights reserved.
15 /**************************************************************************************************/
17 Kmer::Kmer(int size) : kmerSize(size) {
19 int power4s[9] = { 1, 4, 16, 64, 256, 1024, 4096, 16384, 65536 };
20 maxKmer = power4s[kmerSize]+1;// (int)pow(4.,k)+1;
24 /**************************************************************************************************/
26 string Kmer::getKmerString(string sequence){
27 int length = sequence.length();
28 int nKmers = length - kmerSize + 1;
29 vector<int> counts(maxKmer, 0);
31 for(int i=0;i<nKmers;i++){
32 int kmerNumber = getKmerNumber(sequence, i);
36 string kmerString = "";
37 for(int i=0;i<maxKmer;i++){
38 kmerString += getASCII(counts[i]);
44 /**************************************************************************************************/
46 int Kmer::getKmerNumber(string sequence, int index){
47 int power4s[9] = { 1, 4, 16, 64, 256, 1024, 4096, 16384, 65536 };
49 for(int i=0;i<kmerSize;i++){
50 if(toupper(sequence[i+index]) == 'A') { kmer += (0 * power4s[kmerSize-i-1]); }
51 else if(toupper(sequence[i+index]) == 'C') { kmer += (1 * power4s[kmerSize-i-1]); }
52 else if(toupper(sequence[i+index]) == 'G') { kmer += (2 * power4s[kmerSize-i-1]); }
53 else if(toupper(sequence[i+index]) == 'U') { kmer += (3 * power4s[kmerSize-i-1]); }
54 else if(toupper(sequence[i+index]) == 'T') { kmer += (3 * power4s[kmerSize-i-1]); }
55 else if(toupper(sequence[i+index]) == 'N') { return (int)power4s[kmerSize]; }
60 /**************************************************************************************************/
62 string Kmer::getKmerBases(int kmerNumber){
63 int power4s[9] = { 1, 4, 16, 64, 256, 1024, 4096, 16384, 65536 };
67 if(kmerNumber == power4s[kmerSize]){//pow(4.,7)){
68 for(int i=0;i<kmerSize;i++){
73 for(int i=0;i<kmerSize;i++){
74 int nt = (int)(kmerNumber / (float)power4s[i]) % 4;
75 if(nt == 0) { kmer = 'A' + kmer; }
76 else if(nt == 1){ kmer = 'C' + kmer; }
77 else if(nt == 2){ kmer = 'G' + kmer; }
78 else if(nt == 3){ kmer = 'T' + kmer; }
84 /**************************************************************************************************/
86 char Kmer::getASCII(int number) { return (char)(33+number); }
88 /**************************************************************************************************/
90 int Kmer::getNumber(char character) { return ((int)(character-'!')); }
92 /**************************************************************************************************/