- //referenceSeqs, numWanted, matchScore, misMatchPenalty, divR, minSimilarity
- maligner = new Maligner(templateSeqs, numWanted, match, misMatch, 1.01, minSim);
- slayer = new Slayer(window, increment, minSim, divR, iters);
-
- for (int i = 0; i < querySeqs.size(); i++) {
-
- string chimeraFlag = maligner->getResults(querySeqs[i]);
- //float percentIdentical = maligner->getPercentID();
- vector<results> Results = maligner->getOutput();
-
- cout << "Processing sequence: " << i+1 << endl;
-
- //for (int j = 0; j < Results.size(); j++) {
- //cout << "regionStart = " << Results[j].regionStart << "\tRegionEnd = " << Results[j].regionEnd << "\tName = " << Results[j].parent << "\tPerQP = " << Results[j].queryToParent << "\tLocalPerQP = " << Results[j].queryToParentLocal << "\tdivR = " << Results[j].divR << endl;
- //}
-
- //if (chimeraFlag == "yes") {
-
- //get sequence that were given from maligner results
- vector<SeqDist> seqs;
- map<string, float> removeDups;
- map<string, float>::iterator itDup;
- for (int j = 0; j < Results.size(); j++) {
- float dist = (Results[j].regionEnd - Results[j].regionStart + 1) * Results[j].queryToParentLocal;
- //only add if you are not a duplicate
- itDup = removeDups.find(Results[j].parent);
- if (itDup == removeDups.end()) { //this is not duplicate
- removeDups[Results[j].parent] = dist;
- }else if (dist > itDup->second) { //is this a stronger number for this parent
- removeDups[Results[j].parent] = dist;
+ for (int i = 0; i < templateSeqs.size(); i++) {
+ if (namesToAdd.count(templateSeqs[i]->getName()) != 0) {
+ thisTemplate.push_back(templateSeqs[i]);
+ }
+ }
+
+ string kmerDBNameLeft;
+ string kmerDBNameRight;
+
+ //generate the kmerdb to pass to maligner
+ if (searchMethod == "kmer") {
+ string templatePath = m->hasPath(templateFileName);
+ string rightTemplateFileName = templatePath + "right." + m->getRootName(m->getSimpleName(templateFileName));
+ databaseRight = new KmerDB(rightTemplateFileName, kmerSize);
+
+ string leftTemplateFileName = templatePath + "left." + m->getRootName(m->getSimpleName(templateFileName));
+ databaseLeft = new KmerDB(leftTemplateFileName, kmerSize);
+#ifdef USE_MPI
+ for (int i = 0; i < thisTemplate.size(); i++) {
+
+ if (m->control_pressed) { return thisTemplate; }
+
+ string leftFrag = thisTemplate[i]->getUnaligned();
+ leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
+
+ Sequence leftTemp(thisTemplate[i]->getName(), leftFrag);
+ databaseLeft->addSequence(leftTemp);
+ }
+ databaseLeft->generateDB();
+ databaseLeft->setNumSeqs(thisTemplate.size());
+
+ for (int i = 0; i < thisTemplate.size(); i++) {
+ if (m->control_pressed) { return thisTemplate; }
+
+ string rightFrag = thisTemplate[i]->getUnaligned();
+ rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
+
+ Sequence rightTemp(thisTemplate[i]->getName(), rightFrag);
+ databaseRight->addSequence(rightTemp);
+ }
+ databaseRight->generateDB();
+ databaseRight->setNumSeqs(thisTemplate.size());
+
+#else
+
+
+ for (int i = 0; i < thisTemplate.size(); i++) {
+
+ if (m->control_pressed) { return thisTemplate; }
+
+ string leftFrag = thisTemplate[i]->getUnaligned();
+ leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
+
+ Sequence leftTemp(thisTemplate[i]->getName(), leftFrag);
+ databaseLeft->addSequence(leftTemp);
+ }
+ databaseLeft->generateDB();
+ databaseLeft->setNumSeqs(thisTemplate.size());
+
+ for (int i = 0; i < thisTemplate.size(); i++) {
+ if (m->control_pressed) { return thisTemplate; }
+
+ string rightFrag = thisTemplate[i]->getUnaligned();
+ rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
+
+ Sequence rightTemp(thisTemplate[i]->getName(), rightFrag);
+ databaseRight->addSequence(rightTemp);
+ }
+ databaseRight->generateDB();
+ databaseRight->setNumSeqs(thisTemplate.size());
+#endif
+ }else if (searchMethod == "blast") {
+
+ //generate blastdb
+ databaseLeft = new BlastDB(-1.0, -1.0, 1, -3);
+
+ for (int i = 0; i < thisTemplate.size(); i++) { if (m->control_pressed) { return thisTemplate; } databaseLeft->addSequence(*thisTemplate[i]); }
+ databaseLeft->generateDB();
+ databaseLeft->setNumSeqs(thisTemplate.size());
+ }
+
+ return thisTemplate;
+
+ }
+ catch(exception& e) {
+ m->errorOut(e, "ChimeraSlayer", "getTemplate");
+ exit(1);
+ }
+}
+
+//***************************************************************************************************************
+ChimeraSlayer::~ChimeraSlayer() {
+ delete decalc;
+ if (templateFileName != "self") {
+ if (searchMethod == "kmer") { delete databaseRight; delete databaseLeft; }
+ else if (searchMethod == "blast") { delete databaseLeft; }
+ }
+}
+//***************************************************************************************************************
+void ChimeraSlayer::printHeader(ostream& out) {
+ m->mothurOutEndLine();
+ m->mothurOut("Only reporting sequence supported by " + toString(minBS) + "% of bootstrapped results.");
+ m->mothurOutEndLine();
+
+ out << "Name\tLeftParent\tRightParent\tDivQLAQRB\tPerIDQLAQRB\tBootStrapA\tDivQLBQRA\tPerIDQLBQRA\tBootStrapB\tFlag\tLeftWindow\tRightWindow\n";
+}
+//***************************************************************************************************************
+Sequence* ChimeraSlayer::print(ostream& out, ostream& outAcc) {
+ try {
+ Sequence* trim = NULL;
+ if (trimChimera) { trim = trimQuery; }
+
+ if (chimeraFlags == "yes") {
+ string chimeraFlag = "no";
+ if( (chimeraResults[0].bsa >= minBS && chimeraResults[0].divr_qla_qrb >= divR)
+ ||
+ (chimeraResults[0].bsb >= minBS && chimeraResults[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
+
+
+ if (chimeraFlag == "yes") {
+ if ((chimeraResults[0].bsa >= minBS) || (chimeraResults[0].bsb >= minBS)) {
+ m->mothurOut(querySeq->getName() + "\tyes"); m->mothurOutEndLine();
+ outAcc << querySeq->getName() << endl;
+
+ if (trimChimera) {
+ int lengthLeft = spotMap[chimeraResults[0].winLEnd] - spotMap[chimeraResults[0].winLStart];
+ int lengthRight = spotMap[chimeraResults[0].winREnd] - spotMap[chimeraResults[0].winRStart];
+
+ string newAligned = trim->getAligned();
+
+ if (lengthLeft > lengthRight) { //trim right
+ for (int i = (spotMap[chimeraResults[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
+ }else { //trim left
+ for (int i = 0; i < spotMap[chimeraResults[0].winLEnd]; i++) { newAligned[i] = '.'; }
+ }
+ trim->setAligned(newAligned);