X-Git-Url: https://git.donarmstrong.com/?p=mothur.git;a=blobdiff_plain;f=kmerdb.cpp;h=9a5d23500be955b459413beee989601f6c0ec78f;hp=33761a8416524fa3a0fd10b7136ea5a097d3fd08;hb=b206f634aae1b4ce13978d203247fb64757d5482;hpb=725a3d4ff2442c79bfde0a75ed3e0904edcf03b7 diff --git a/kmerdb.cpp b/kmerdb.cpp index 33761a8..9a5d235 100644 --- a/kmerdb.cpp +++ b/kmerdb.cpp @@ -57,9 +57,12 @@ KmerDB::~KmerDB(){} vector KmerDB::findClosestSequences(Sequence* candidateSeq, int num){ try { + if (num > numSeqs) { m->mothurOut("[WARNING]: you requested " + toString(num) + " closest sequences, but the template only contains " + toString(numSeqs) + ", adjusting."); m->mothurOutEndLine(); num = numSeqs; } + vector topMatches; Kmer kmer(kmerSize); searchScore = 0; + Scores.clear(); vector matches(numSeqs, 0); // a record of the sequences with shared kmers vector timesKmerFound(kmerLocations.size()+1, 0); // a record of the kmers that we have already found @@ -75,24 +78,42 @@ vector KmerDB::findClosestSequences(Sequence* candidateSeq, int num){ } timesKmerFound[kmerNumber] = 1; // ok, we've seen the kmer now } - - vector seqMatches; - for(int i=0;i seqMatches; seqMatches.resize(numSeqs); + for(int i=0;i bestMatch) { + bestIndex = i; + bestMatch = matches[i]; + } + } + + searchScore = bestMatch; + searchScore = 100 * searchScore / (float) numKmers; // return the Sequence object corresponding to the db + topMatches.push_back(bestIndex); + Scores.push_back(searchScore); } - return topMatches; } catch(exception& e) { @@ -107,7 +128,10 @@ void KmerDB::generateDB(){ try { ofstream kmerFile; // once we have the kmerLocations folder print it out - openOutputFile(kmerDBName, kmerFile); // to a file + m->openOutputFile(kmerDBName, kmerFile); // to a file + + //output version + kmerFile << "#" << m->getVersion() << endl; for(int i=0;igetline(kmerDBFile); m->gobble(kmerDBFile); + string seqName; int seqNumber; @@ -190,6 +217,20 @@ int KmerDB::getCount(int kmer) { } } /**************************************************************************************************/ +int KmerDB::getReversed(int kmerNumber) { + try { + Kmer kmer(kmerSize); + + if (kmerNumber < 0) { return 0; } //if user gives negative number + else if (kmerNumber > maxKmer) { return 0; } //or a kmer that is bigger than maxkmer + else { return kmer.getReverseKmerNumber(kmerNumber); } // kmer is in vector range + } + catch(exception& e) { + m->errorOut(e, "KmerDB", "getReversed"); + exit(1); + } +} +/**************************************************************************************************/ vector KmerDB::getSequencesWithKmer(int kmer) { try {