chimeraslayer.cpp

   1 /*
   2  *  chimeraslayer.cpp
   3  *  Mothur
   4  *
   5  *  Created by westcott on 9/25/09.
   6  *  Copyright 2009 Pschloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "chimeraslayer.h"
  11 #include "chimerarealigner.h"
  12 #include "kmerdb.hpp"
  13 #include "blastdb.hpp"
  14
  15 //***************************************************************************************************************
  16 ChimeraSlayer::ChimeraSlayer(string file, string temp, bool trim, string mode, int k, int ms, int mms, int win, float div,
  17 int minsim, int mincov, int minbs, int minsnp, int par, int it, int inc, int numw, bool r) : Chimera()  {
  18         try {
  19                 fastafile = file;
  20                 templateFileName = temp; templateSeqs = readSeqs(temp);
  21                 searchMethod = mode;
  22                 kmerSize = k;
  23                 match = ms;
  24                 misMatch = mms;
  25                 window = win;
  26                 divR = div;
  27                 minSim = minsim;
  28                 minCov = mincov;
  29                 minBS = minbs;
  30                 minSNP = minsnp;
  31                 parents = par;
  32                 iters = it;
  33                 increment = inc;
  34                 numWanted = numw;
  35                 realign = r;
  36                 trimChimera = trim;
  37
  38                 decalc = new DeCalculator();
  39
  40                 doPrep();
  41         }
  42         catch(exception& e) {
  43                 m->errorOut(e, "ChimeraSlayer", "ChimeraSlayer");
  44                 exit(1);
  45         }
  46 }
  47 //***************************************************************************************************************
  48 ChimeraSlayer::ChimeraSlayer(string file, string temp, bool trim, string name, string mode, int k, int ms, int mms, int win, float div,
  49                                                          int minsim, int mincov, int minbs, int minsnp, int par, int it, int inc, int numw, bool r) : Chimera()  {
  50         try {
  51                 fastafile = file; templateSeqs = readSeqs(fastafile);
  52                 templateFileName = temp;
  53                 searchMethod = mode;
  54                 kmerSize = k;
  55                 match = ms;
  56                 misMatch = mms;
  57                 window = win;
  58                 divR = div;
  59                 minSim = minsim;
  60                 minCov = mincov;
  61                 minBS = minbs;
  62                 minSNP = minsnp;
  63                 parents = par;
  64                 iters = it;
  65                 increment = inc;
  66                 numWanted = numw;
  67                 realign = r;
  68                 trimChimera = trim;
  69
  70                 decalc = new DeCalculator();
  71
  72                 createFilter(templateSeqs, 0.0); //just removed columns where all seqs have a gap
  73
  74                 //run filter on template
  75                 for (int i = 0; i < templateSeqs.size(); i++) { delete templateSeqs[i];  } templateSeqs.clear();
  76
  77         }
  78         catch(exception& e) {
  79                 m->errorOut(e, "ChimeraSlayer", "ChimeraSlayer");
  80                 exit(1);
  81         }
  82 }
  83 //***************************************************************************************************************
  84 int ChimeraSlayer::doPrep() {
  85         try {
  86
  87                 //read in all query seqs
  88                 vector<Sequence*> tempQuerySeqs = readSeqs(fastafile);
  89
  90                 vector<Sequence*> temp = templateSeqs;
  91                 for (int i = 0; i < tempQuerySeqs.size(); i++) {  temp.push_back(tempQuerySeqs[i]);  }
  92
  93                 createFilter(temp, 0.0); //just removed columns where all seqs have a gap
  94
  95                 for (int i = 0; i < tempQuerySeqs.size(); i++) { delete tempQuerySeqs[i];  }
  96
  97                 if (m->control_pressed) {  return 0; }
  98
  99                 //run filter on template
 100                 for (int i = 0; i < templateSeqs.size(); i++) {  if (m->control_pressed) {  return 0; }  runFilter(templateSeqs[i]);  }
 101
 102                 string  kmerDBNameLeft;
 103                 string  kmerDBNameRight;
 104
 105                 //generate the kmerdb to pass to maligner
 106                 if (searchMethod == "kmer") {
 107                         string templatePath = m->hasPath(templateFileName);
 108                         string rightTemplateFileName = templatePath + "right." + m->getRootName(m->getSimpleName(templateFileName));
 109                         databaseRight = new KmerDB(rightTemplateFileName, kmerSize);
 110
 111                         string leftTemplateFileName = templatePath + "left." + m->getRootName(m->getSimpleName(templateFileName));
 112                         databaseLeft = new KmerDB(leftTemplateFileName, kmerSize);
 113                 #ifdef USE_MPI
 114                         for (int i = 0; i < templateSeqs.size(); i++) {
 115
 116                                 if (m->control_pressed) { return 0; }
 117
 118                                 string leftFrag = templateSeqs[i]->getUnaligned();
 119                                 leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
 120
 121                                 Sequence leftTemp(templateSeqs[i]->getName(), leftFrag);
 122                                 databaseLeft->addSequence(leftTemp);
 123                         }
 124                         databaseLeft->generateDB();
 125                         databaseLeft->setNumSeqs(templateSeqs.size());
 126
 127                         for (int i = 0; i < templateSeqs.size(); i++) {
 128                                 if (m->control_pressed) { return 0; }
 129
 130                                 string rightFrag = templateSeqs[i]->getUnaligned();
 131                                 rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
 132
 133                                 Sequence rightTemp(templateSeqs[i]->getName(), rightFrag);
 134                                 databaseRight->addSequence(rightTemp);
 135                         }
 136                         databaseRight->generateDB();
 137                         databaseRight->setNumSeqs(templateSeqs.size());
 138
 139                 #else
 140                         //leftside
 141                         kmerDBNameLeft = leftTemplateFileName.substr(0,leftTemplateFileName.find_last_of(".")+1) + char('0'+ kmerSize) + "mer";
 142                         ifstream kmerFileTestLeft(kmerDBNameLeft.c_str());
 143                         bool needToGenerateLeft = true;
 144
 145                         if(kmerFileTestLeft){
 146                                 bool GoodFile = m->checkReleaseVersion(kmerFileTestLeft, m->getVersion());
 147                                 if (GoodFile) {  needToGenerateLeft = false;    }
 148                         }
 149
 150                         if(needToGenerateLeft){
 151
 152                                 for (int i = 0; i < templateSeqs.size(); i++) {
 153
 154                                         if (m->control_pressed) { return 0; }
 155
 156                                         string leftFrag = templateSeqs[i]->getUnaligned();
 157                                         leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
 158
 159                                         Sequence leftTemp(templateSeqs[i]->getName(), leftFrag);
 160                                         databaseLeft->addSequence(leftTemp);
 161                                 }
 162                                 databaseLeft->generateDB();
 163
 164                         }else {
 165                                 databaseLeft->readKmerDB(kmerFileTestLeft);
 166                         }
 167                         kmerFileTestLeft.close();
 168
 169                         databaseLeft->setNumSeqs(templateSeqs.size());
 170
 171                         //rightside
 172                         kmerDBNameRight = rightTemplateFileName.substr(0,rightTemplateFileName.find_last_of(".")+1) + char('0'+ kmerSize) + "mer";
 173                         ifstream kmerFileTestRight(kmerDBNameRight.c_str());
 174                         bool needToGenerateRight = true;
 175
 176                         if(kmerFileTestRight){
 177                                 bool GoodFile = m->checkReleaseVersion(kmerFileTestRight, m->getVersion());
 178                                 if (GoodFile) {  needToGenerateRight = false;   }
 179                         }
 180
 181                         if(needToGenerateRight){
 182
 183                                 for (int i = 0; i < templateSeqs.size(); i++) {
 184                                         if (m->control_pressed) { return 0; }
 185
 186                                         string rightFrag = templateSeqs[i]->getUnaligned();
 187                                         rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
 188
 189                                         Sequence rightTemp(templateSeqs[i]->getName(), rightFrag);
 190                                         databaseRight->addSequence(rightTemp);
 191                                 }
 192                                 databaseRight->generateDB();
 193
 194                         }else {
 195                                 databaseRight->readKmerDB(kmerFileTestRight);
 196                         }
 197                         kmerFileTestRight.close();
 198
 199                         databaseRight->setNumSeqs(templateSeqs.size());
 200                 #endif
 201                 }else if (searchMethod == "blast") {
 202
 203                         //generate blastdb
 204                         databaseLeft = new BlastDB(-1.0, -1.0, 1, -3);
 205
 206                         for (int i = 0; i < templateSeqs.size(); i++) {         databaseLeft->addSequence(*templateSeqs[i]);    }
 207                         databaseLeft->generateDB();
 208                         databaseLeft->setNumSeqs(templateSeqs.size());
 209                 }
 210
 211                 return 0;
 212
 213         }
 214         catch(exception& e) {
 215                 m->errorOut(e, "ChimeraSlayer", "doprep");
 216                 exit(1);
 217         }
 218 }
 219 //***************************************************************************************************************
 220 int ChimeraSlayer::getTemplate(Sequence* q) {
 221         try {
 222
 223                 string  kmerDBNameLeft;
 224                 string  kmerDBNameRight;
 225
 226                 //generate the kmerdb to pass to maligner
 227                 if (searchMethod == "kmer") {
 228                         string templatePath = m->hasPath(templateFileName);
 229                         string rightTemplateFileName = templatePath + "right." + m->getRootName(m->getSimpleName(templateFileName));
 230                         databaseRight = new KmerDB(rightTemplateFileName, kmerSize);
 231
 232                         string leftTemplateFileName = templatePath + "left." + m->getRootName(m->getSimpleName(templateFileName));
 233                         databaseLeft = new KmerDB(leftTemplateFileName, kmerSize);
 234 #ifdef USE_MPI
 235                         for (int i = 0; i < userTemplate.size(); i++) {
 236
 237                                 if (m->control_pressed) { return 0; }
 238
 239                                 string leftFrag = userTemplate[i]->getUnaligned();
 240                                 leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
 241
 242                                 Sequence leftTemp(userTemplate[i]->getName(), leftFrag);
 243                                 databaseLeft->addSequence(leftTemp);
 244                         }
 245                         databaseLeft->generateDB();
 246                         databaseLeft->setNumSeqs(userTemplate.size());
 247
 248                         for (int i = 0; i < userTemplate.size(); i++) {
 249                                 if (m->control_pressed) { return 0;  }
 250
 251                                 string rightFrag = userTemplate[i]->getUnaligned();
 252                                 rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
 253
 254                                 Sequence rightTemp(userTemplate[i]->getName(), rightFrag);
 255                                 databaseRight->addSequence(rightTemp);
 256                         }
 257                         databaseRight->generateDB();
 258                         databaseRight->setNumSeqs(userTemplate.size());
 259
 260 #else
 261
 262
 263                         for (int i = 0; i < userTemplate.size(); i++) {
 264
 265                                 if (m->control_pressed) { return 0; }
 266
 267                                 string leftFrag = userTemplate[i]->getUnaligned();
 268                                 leftFrag = leftFrag.substr(0, int(leftFrag.length() * 0.33));
 269
 270                                 Sequence leftTemp(userTemplate[i]->getName(), leftFrag);
 271                                 databaseLeft->addSequence(leftTemp);
 272                         }
 273                         databaseLeft->generateDB();
 274                         databaseLeft->setNumSeqs(userTemplate.size());
 275
 276                         for (int i = 0; i < userTemplate.size(); i++) {
 277                                 if (m->control_pressed) { return 0; }
 278
 279                                 string rightFrag = userTemplate[i]->getUnaligned();
 280                                 rightFrag = rightFrag.substr(int(rightFrag.length() * 0.66));
 281
 282                                 Sequence rightTemp(userTemplate[i]->getName(), rightFrag);
 283                                 databaseRight->addSequence(rightTemp);
 284                         }
 285                         databaseRight->generateDB();
 286                         databaseRight->setNumSeqs(userTemplate.size());
 287 #endif
 288                 }else if (searchMethod == "blast") {
 289
 290                         //generate blastdb
 291                         databaseLeft = new BlastDB(-1.0, -1.0, 1, -3);
 292
 293                         for (int i = 0; i < userTemplate.size(); i++) { if (m->control_pressed) { return 0; }  databaseLeft->addSequence(*userTemplate[i]);     }
 294                         databaseLeft->generateDB();
 295                         databaseLeft->setNumSeqs(userTemplate.size());
 296                 }
 297
 298                 return 0;
 299
 300         }
 301         catch(exception& e) {
 302                 m->errorOut(e, "ChimeraSlayer", "getTemplate");
 303                 exit(1);
 304         }
 305 }
 306
 307 //***************************************************************************************************************
 308 ChimeraSlayer::~ChimeraSlayer() {
 309         delete decalc;
 310         if (templateFileName != "self") {
 311                 if (searchMethod == "kmer") {  delete databaseRight;  delete databaseLeft;  }
 312                 else if (searchMethod == "blast") {  delete databaseLeft; }
 313         }else {
 314                 //delete userTemplate
 315                 for (int i = 0; i < userTemplate.size(); i++) {
 316                         delete userTemplate[i];
 317                 }
 318                 userTemplate.clear();
 319         }
 320 }
 321 //***************************************************************************************************************
 322 void ChimeraSlayer::printHeader(ostream& out) {
 323         m->mothurOutEndLine();
 324         m->mothurOut("Only reporting sequence supported by " + toString(minBS) + "% of bootstrapped results.");
 325         m->mothurOutEndLine();
 326
 327         out << "Name\tLeftParent\tRightParent\tDivQLAQRB\tPerIDQLAQRB\tBootStrapA\tDivQLBQRA\tPerIDQLBQRA\tBootStrapB\tFlag\tLeftWindow\tRightWindow\n";
 328 }
 329 //***************************************************************************************************************
 330 Sequence* ChimeraSlayer::print(ostream& out, ostream& outAcc) {
 331         try {
 332                 Sequence* trim = NULL;
 333                 if (trimChimera) { trim = new Sequence(trimQuery.getName(), trimQuery.getAligned()); }
 334
 335                 if (chimeraFlags == "yes") {
 336                         string chimeraFlag = "no";
 337                         if(  (chimeraResults[0].bsa >= minBS && chimeraResults[0].divr_qla_qrb >= divR)
 338                            ||
 339                            (chimeraResults[0].bsb >= minBS && chimeraResults[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 340
 341
 342                         if (chimeraFlag == "yes") {
 343                                 if ((chimeraResults[0].bsa >= minBS) || (chimeraResults[0].bsb >= minBS)) {
 344                                         m->mothurOut(querySeq->getName() + "\tyes"); m->mothurOutEndLine();
 345                                         outAcc << querySeq->getName() << endl;
 346
 347                                         if (trimChimera) {
 348                                                 int lengthLeft = spotMap[chimeraResults[0].winLEnd] - spotMap[chimeraResults[0].winLStart];
 349                                                 int lengthRight = spotMap[chimeraResults[0].winREnd] - spotMap[chimeraResults[0].winRStart];
 350
 351                                                 string newAligned = trim->getAligned();
 352
 353                                                 if (lengthLeft > lengthRight) { //trim right
 354                                                         for (int i = (spotMap[chimeraResults[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 355                                                 }else { //trim left
 356                                                         for (int i = 0; i < spotMap[chimeraResults[0].winLEnd]; i++) { newAligned[i] = '.'; }
 357                                                 }
 358                                                 trim->setAligned(newAligned);
 359                                         }
 360                                 }
 361                         }
 362
 363                         printBlock(chimeraResults[0], chimeraFlag, out);
 364                         out << endl;
 365                 }else {
 366                         out << querySeq->getName() << "\tno" << endl;
 367                         if (templateFileName == "self") {
 368                                 Sequence* temp = new Sequence(trimQuery.getName(), trimQuery.getAligned());
 369                                 runFilter(temp);
 370                                 userTemplate.push_back(temp);
 371                         }
 372                 }
 373
 374                 return trim;
 375
 376         }
 377         catch(exception& e) {
 378                 m->errorOut(e, "ChimeraSlayer", "print");
 379                 exit(1);
 380         }
 381 }
 382 //***************************************************************************************************************
 383 Sequence* ChimeraSlayer::print(ostream& out, ostream& outAcc, data_results leftPiece, data_results rightPiece) {
 384         try {
 385                 Sequence* trim = NULL;
 386
 387                 if (trimChimera) {
 388                         string aligned = leftPiece.trimQuery.getAligned() + rightPiece.trimQuery.getAligned();
 389                         trim = new Sequence(leftPiece.trimQuery.getName(), aligned);
 390                 }
 391
 392                 if ((leftPiece.flag == "yes") || (rightPiece.flag == "yes")) {
 393
 394                         string chimeraFlag = "no";
 395                         if (leftPiece.flag == "yes") {
 396
 397                                 if(  (leftPiece.results[0].bsa >= minBS && leftPiece.results[0].divr_qla_qrb >= divR)
 398                                         ||
 399                                         (leftPiece.results[0].bsb >= minBS && leftPiece.results[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 400                         }
 401
 402                         if (rightPiece.flag == "yes") {
 403                                 if ( (rightPiece.results[0].bsa >= minBS && rightPiece.results[0].divr_qla_qrb >= divR)
 404                                  ||
 405                                  (rightPiece.results[0].bsb >= minBS && rightPiece.results[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 406                         }
 407
 408                         bool rightChimeric = false;
 409                         bool leftChimeric = false;
 410
 411                         if (chimeraFlag == "yes") {
 412                                 //which peice is chimeric or are both
 413                                 if (rightPiece.flag == "yes") { if ((rightPiece.results[0].bsa >= minBS) || (rightPiece.results[0].bsb >= minBS)) { rightChimeric = true; } }
 414                                 if (leftPiece.flag == "yes") { if ((leftPiece.results[0].bsa >= minBS) || (leftPiece.results[0].bsb >= minBS))  { leftChimeric = true;  } }
 415
 416                                 if (rightChimeric || leftChimeric) {
 417                                         m->mothurOut(querySeq->getName() + "\tyes"); m->mothurOutEndLine();
 418                                         outAcc << querySeq->getName() << endl;
 419
 420                                         if (trimChimera) {
 421                                                 string newAligned = trim->getAligned();
 422
 423                                                 //right side is fine so keep that
 424                                                 if ((leftChimeric) && (!rightChimeric)) {
 425                                                         for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 426                                                 }else if ((!leftChimeric) && (rightChimeric)) { //leftside is fine so keep that
 427                                                         for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 428                                                 }else { //both sides are chimeric, keep longest piece
 429
 430                                                         int lengthLeftLeft = leftPiece.spotMap[leftPiece.results[0].winLEnd] - leftPiece.spotMap[leftPiece.results[0].winLStart];
 431                                                         int lengthLeftRight = leftPiece.spotMap[leftPiece.results[0].winREnd] - leftPiece.spotMap[leftPiece.results[0].winRStart];
 432
 433                                                         int longest = 1; // leftleft = 1, leftright = 2, rightleft = 3 rightright = 4
 434                                                         int length = lengthLeftLeft;
 435                                                         if (lengthLeftLeft < lengthLeftRight) { longest = 2;  length = lengthLeftRight; }
 436
 437                                                         int lengthRightLeft = rightPiece.spotMap[rightPiece.results[0].winLEnd] - rightPiece.spotMap[rightPiece.results[0].winLStart];
 438                                                         int lengthRightRight = rightPiece.spotMap[rightPiece.results[0].winREnd] - rightPiece.spotMap[rightPiece.results[0].winRStart];
 439
 440                                                         if (lengthRightLeft > length) { longest = 3; length = lengthRightLeft;  }
 441                                                         if (lengthRightRight > length) { longest = 4; }
 442
 443                                                         if (longest == 1) { //leftleft
 444                                                                 for (int i = (leftPiece.spotMap[leftPiece.results[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 445                                                         }else if (longest == 2) { //leftright
 446                                                                 //get rid of leftleft
 447                                                                 for (int i = (leftPiece.spotMap[leftPiece.results[0].winLStart]-1); i < (leftPiece.spotMap[leftPiece.results[0].winLEnd]-1); i++) { newAligned[i] = '.'; }
 448                                                                 //get rid of right
 449                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 450                                                         }else if (longest == 3) { //rightleft
 451                                                                 //get rid of left
 452                                                                 for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 453                                                                 //get rid of rightright
 454                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 455                                                         }else { //rightright
 456                                                                 //get rid of left
 457                                                                 for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 458                                                                 //get rid of rightleft
 459                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < (rightPiece.spotMap[rightPiece.results[0].winLEnd]-1); i++) { newAligned[i] = '.'; }
 460                                                         }
 461                                                 }
 462
 463                                                 trim->setAligned(newAligned);
 464                                         }
 465
 466                                 }
 467                         }
 468
 469                         printBlock(leftPiece, rightPiece, leftChimeric, rightChimeric, chimeraFlag, out);
 470                         out << endl;
 471                 }else {
 472                         out << querySeq->getName() << "\tno" << endl;
 473                         if (templateFileName == "self") {
 474                                 Sequence* temp = new Sequence(trimQuery.getName(), trimQuery.getAligned());
 475                                 runFilter(temp);
 476                                 userTemplate.push_back(temp);
 477                         }
 478                 }
 479
 480                 return trim;
 481
 482         }
 483         catch(exception& e) {
 484                 m->errorOut(e, "ChimeraSlayer", "print");
 485                 exit(1);
 486         }
 487 }
 488
 489 #ifdef USE_MPI
 490 //***************************************************************************************************************
 491 Sequence* ChimeraSlayer::print(MPI_File& out, MPI_File& outAcc, data_results leftPiece, data_results rightPiece) {
 492         try {
 493                 MPI_Status status;
 494                 bool results = false;
 495                 string outAccString = "";
 496                 string outputString = "";
 497
 498                 Sequence* trim = NULL;
 499
 500                 if (trimChimera) {
 501                         string aligned = leftPiece.trimQuery.getAligned() + rightPiece.trimQuery.getAligned();
 502                         trim = new Sequence(leftPiece.trimQuery.getName(), aligned);
 503                 }
 504
 505
 506                 if ((leftPiece.flag == "yes") || (rightPiece.flag == "yes")) {
 507
 508                         string chimeraFlag = "no";
 509                         if (leftPiece.flag == "yes") {
 510
 511                                 if(  (leftPiece.results[0].bsa >= minBS && leftPiece.results[0].divr_qla_qrb >= divR)
 512                                    ||
 513                                    (leftPiece.results[0].bsb >= minBS && leftPiece.results[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 514                         }
 515
 516                         if (rightPiece.flag == "yes") {
 517                                 if ( (rightPiece.results[0].bsa >= minBS && rightPiece.results[0].divr_qla_qrb >= divR)
 518                                         ||
 519                                         (rightPiece.results[0].bsb >= minBS && rightPiece.results[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 520                         }
 521
 522                         bool rightChimeric = false;
 523                         bool leftChimeric = false;
 524
 525                         if (chimeraFlag == "yes") {
 526                                 //which peice is chimeric or are both
 527                                 if (rightPiece.flag == "yes") { if ((rightPiece.results[0].bsa >= minBS) || (rightPiece.results[0].bsb >= minBS)) { rightChimeric = true; } }
 528                                 if (leftPiece.flag == "yes") { if ((leftPiece.results[0].bsa >= minBS) || (leftPiece.results[0].bsb >= minBS))  { leftChimeric = true;  } }
 529
 530                                 if (rightChimeric || leftChimeric) {
 531                                         cout << querySeq->getName() <<  "\tyes" << endl;
 532                                         outAccString += querySeq->getName() + "\n";
 533                                         results = true;
 534
 535                                         //write to accnos file
 536                                         int length = outAccString.length();
 537                                         char* buf2 = new char[length];
 538                                         memcpy(buf2, outAccString.c_str(), length);
 539
 540                                         MPI_File_write_shared(outAcc, buf2, length, MPI_CHAR, &status);
 541                                         delete buf2;
 542
 543                                         if (trimChimera) {
 544                                                 string newAligned = trim->getAligned();
 545
 546                                                 //right side is fine so keep that
 547                                                 if ((leftChimeric) && (!rightChimeric)) {
 548                                                         for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 549                                                 }else if ((!leftChimeric) && (rightChimeric)) { //leftside is fine so keep that
 550                                                         for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 551                                                 }else { //both sides are chimeric, keep longest piece
 552
 553                                                         int lengthLeftLeft = leftPiece.spotMap[leftPiece.results[0].winLEnd] - leftPiece.spotMap[leftPiece.results[0].winLStart];
 554                                                         int lengthLeftRight = leftPiece.spotMap[leftPiece.results[0].winREnd] - leftPiece.spotMap[leftPiece.results[0].winRStart];
 555
 556                                                         int longest = 1; // leftleft = 1, leftright = 2, rightleft = 3 rightright = 4
 557                                                         int length = lengthLeftLeft;
 558                                                         if (lengthLeftLeft < lengthLeftRight) { longest = 2;  length = lengthLeftRight; }
 559
 560                                                         int lengthRightLeft = rightPiece.spotMap[rightPiece.results[0].winLEnd] - rightPiece.spotMap[rightPiece.results[0].winLStart];
 561                                                         int lengthRightRight = rightPiece.spotMap[rightPiece.results[0].winREnd] - rightPiece.spotMap[rightPiece.results[0].winRStart];
 562
 563                                                         if (lengthRightLeft > length) { longest = 3; length = lengthRightLeft;  }
 564                                                         if (lengthRightRight > length) { longest = 4; }
 565
 566                                                         if (longest == 1) { //leftleft
 567                                                                 for (int i = (leftPiece.spotMap[leftPiece.results[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 568                                                         }else if (longest == 2) { //leftright
 569                                                                 //get rid of leftleft
 570                                                                 for (int i = (leftPiece.spotMap[leftPiece.results[0].winLStart]-1); i < (leftPiece.spotMap[leftPiece.results[0].winLEnd]-1); i++) { newAligned[i] = '.'; }
 571                                                                 //get rid of right
 572                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 573                                                         }else if (longest == 3) { //rightleft
 574                                                                 //get rid of left
 575                                                                 for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 576                                                                 //get rid of rightright
 577                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 578                                                         }else { //rightright
 579                                                                 //get rid of left
 580                                                                 for (int i = 0; i < leftPiece.spotMap[leftPiece.results[0].winREnd]; i++) { newAligned[i] = '.'; }
 581                                                                 //get rid of rightleft
 582                                                                 for (int i = (rightPiece.spotMap[rightPiece.results[0].winLStart]-1); i < (rightPiece.spotMap[rightPiece.results[0].winLEnd]-1); i++) { newAligned[i] = '.'; }
 583                                                         }
 584                                                 }
 585
 586                                                 trim->setAligned(newAligned);
 587                                         }
 588
 589                                 }
 590                         }
 591
 592                         outputString = getBlock(leftPiece, rightPiece, leftChimeric, rightChimeric, chimeraFlag);
 593                         outputString += "\n";
 594
 595                         //write to output file
 596                         int length = outputString.length();
 597                         char* buf = new char[length];
 598                         memcpy(buf, outputString.c_str(), length);
 599
 600                         MPI_File_write_shared(out, buf, length, MPI_CHAR, &status);
 601                         delete buf;
 602
 603                 }else {
 604                         outputString += querySeq->getName() + "\tno\n";
 605
 606                         //write to output file
 607                         int length = outputString.length();
 608                         char* buf = new char[length];
 609                         memcpy(buf, outputString.c_str(), length);
 610
 611                         MPI_File_write_shared(out, buf, length, MPI_CHAR, &status);
 612                         delete buf;
 613
 614                         if (template == "self") {
 615                                 Sequence temp = new Sequence(trimQuery.getName(), trimQuery.getAligned());
 616                                 runFilter(temp);
 617                                 userTemplate.push_back(temp);
 618                         }
 619                 }
 620
 621
 622                 return trim;
 623         }
 624         catch(exception& e) {
 625                 m->errorOut(e, "ChimeraSlayer", "print");
 626                 exit(1);
 627         }
 628 }
 629 //***************************************************************************************************************
 630 Sequence* ChimeraSlayer::print(MPI_File& out, MPI_File& outAcc) {
 631         try {
 632                 MPI_Status status;
 633                 bool results = false;
 634                 string outAccString = "";
 635                 string outputString = "";
 636
 637                 Sequence* trim = NULL;
 638                 if (trimChimera) { trim = new Sequence(trimQuery.getName(), trimQuery.getAligned()); }
 639
 640                 if (chimeraFlags == "yes") {
 641                         string chimeraFlag = "no";
 642                         if(  (chimeraResults[0].bsa >= minBS && chimeraResults[0].divr_qla_qrb >= divR)
 643                            ||
 644                            (chimeraResults[0].bsb >= minBS && chimeraResults[0].divr_qlb_qra >= divR) ) { chimeraFlag = "yes"; }
 645
 646
 647                         if (chimeraFlag == "yes") {
 648                                 if ((chimeraResults[0].bsa >= minBS) || (chimeraResults[0].bsb >= minBS)) {
 649                                         cout << querySeq->getName() <<  "\tyes" << endl;
 650                                         outAccString += querySeq->getName() + "\n";
 651                                         results = true;
 652
 653                                         //write to accnos file
 654                                         int length = outAccString.length();
 655                                         char* buf2 = new char[length];
 656                                         memcpy(buf2, outAccString.c_str(), length);
 657
 658                                         MPI_File_write_shared(outAcc, buf2, length, MPI_CHAR, &status);
 659                                         delete buf2;
 660
 661                                         if (trimChimera) {
 662                                                 int lengthLeft = spotMap[chimeraResults[0].winLEnd] - spotMap[chimeraResults[0].winLStart];
 663                                                 int lengthRight = spotMap[chimeraResults[0].winREnd] - spotMap[chimeraResults[0].winRStart];
 664
 665                                                 string newAligned = trim->getAligned();
 666                                                 if (lengthLeft > lengthRight) { //trim right
 667                                                         for (int i = (spotMap[chimeraResults[0].winRStart]-1); i < newAligned.length(); i++) { newAligned[i] = '.'; }
 668                                                 }else { //trim left
 669                                                         for (int i = 0; i < (spotMap[chimeraResults[0].winLEnd]-1); i++) { newAligned[i] = '.'; }
 670                                                 }
 671                                                 trim->setAligned(newAligned);
 672                                         }
 673                                 }
 674                         }
 675
 676                         outputString = getBlock(chimeraResults[0], chimeraFlag);
 677                         outputString += "\n";
 678
 679                         //write to output file
 680                         int length = outputString.length();
 681                         char* buf = new char[length];
 682                         memcpy(buf, outputString.c_str(), length);
 683
 684                         MPI_File_write_shared(out, buf, length, MPI_CHAR, &status);
 685                         delete buf;
 686
 687                 }else {
 688                         outputString += querySeq->getName() + "\tno\n";
 689
 690                         //write to output file
 691                         int length = outputString.length();
 692                         char* buf = new char[length];
 693                         memcpy(buf, outputString.c_str(), length);
 694
 695                         MPI_File_write_shared(out, buf, length, MPI_CHAR, &status);
 696                         delete buf;
 697
 698                         if (template == "self") {
 699                                 Sequence temp = new Sequence(trimQuery.getName(), trimQuery.getAligned());
 700                                 runFilter(temp);
 701                                 userTemplate.push_back(temp);
 702                         }
 703                 }
 704
 705                 return trim;
 706         }
 707         catch(exception& e) {
 708                 m->errorOut(e, "ChimeraSlayer", "print");
 709                 exit(1);
 710         }
 711 }
 712 #endif
 713
 714 //***************************************************************************************************************
 715 int ChimeraSlayer::getChimeras(Sequence* query) {
 716         try {
 717
 718                 trimQuery.setName(query->getName()); trimQuery.setAligned(query->getAligned());
 719                 printResults.trimQuery = trimQuery;
 720
 721                 chimeraFlags = "no";
 722                 printResults.flag = "no";
 723
 724                 //filter query
 725                 spotMap = runFilter(query);
 726                 printResults.spotMap = spotMap;
 727
 728                 querySeq = query;
 729
 730                 //you must create a template
 731                 vector<Sequence*> thisTemplate;
 732                 if (templateFileName != "self") { thisTemplate = templateSeqs; }
 733                 else { getTemplate(query);  thisTemplate = userTemplate; } //fills this template and creates the databases
 734
 735                 if (m->control_pressed) {  return 0;  }
 736
 737                 if (thisTemplate.size() == 0) {  return 0; } //not chimeric
 738
 739                 //referenceSeqs, numWanted, matchScore, misMatchPenalty, divR, minSimilarity
 740                 Maligner maligner(thisTemplate, numWanted, match, misMatch, divR, minSim, minCov, searchMethod, databaseLeft, databaseRight);
 741                 Slayer slayer(window, increment, minSim, divR, iters, minSNP);
 742
 743                 if (templateFileName == "self") {
 744                         if (searchMethod == "kmer") {  delete databaseRight;  delete databaseLeft;  }
 745                         else if (searchMethod == "blast") {  delete databaseLeft; }
 746                 }
 747
 748                 if (m->control_pressed) {  return 0;  }
 749
 750                 string chimeraFlag = maligner.getResults(query, decalc);
 751
 752                 if (m->control_pressed) {  return 0;  }
 753
 754                 vector<results> Results = maligner.getOutput();
 755
 756                 if (realign) {
 757                         ChimeraReAligner realigner(thisTemplate, match, misMatch);
 758                         realigner.reAlign(query, Results);
 759                 }
 760
 761                 if (chimeraFlag == "yes") {
 762
 763                         //get sequence that were given from maligner results
 764                         vector<SeqDist> seqs;
 765                         map<string, float> removeDups;
 766                         map<string, float>::iterator itDup;
 767                         map<string, string> parentNameSeq;
 768                         map<string, string>::iterator itSeq;
 769                         for (int j = 0; j < Results.size(); j++) {
 770                                 float dist = (Results[j].regionEnd - Results[j].regionStart + 1) * Results[j].queryToParentLocal;
 771                                 //only add if you are not a duplicate
 772                                 itDup = removeDups.find(Results[j].parent);
 773                                 if (itDup == removeDups.end()) { //this is not duplicate
 774                                         removeDups[Results[j].parent] = dist;
 775                                         parentNameSeq[Results[j].parent] = Results[j].parentAligned;
 776                                 }else if (dist > itDup->second) { //is this a stronger number for this parent
 777                                         removeDups[Results[j].parent] = dist;
 778                                         parentNameSeq[Results[j].parent] = Results[j].parentAligned;
 779                                 }
 780                         }
 781
 782                         for (itDup = removeDups.begin(); itDup != removeDups.end(); itDup++) {
 783                                 itSeq = parentNameSeq.find(itDup->first);
 784                                 Sequence* seq = new Sequence(itDup->first, itSeq->second);
 785
 786                                 SeqDist member;
 787                                 member.seq = seq;
 788                                 member.dist = itDup->second;
 789
 790                                 seqs.push_back(member);
 791                         }
 792
 793                         //limit number of parents to explore - default 3
 794                         if (Results.size() > parents) {
 795                                 //sort by distance
 796                                 sort(seqs.begin(), seqs.end(), compareSeqDist);
 797                                 //prioritize larger more similiar sequence fragments
 798                                 reverse(seqs.begin(), seqs.end());
 799
 800                                 for (int k = seqs.size()-1; k > (parents-1); k--)  {
 801                                         delete seqs[k].seq;
 802                                         seqs.pop_back();
 803                                 }
 804                         }
 805
 806                         //put seqs into vector to send to slayer
 807                         vector<Sequence*> seqsForSlayer;
 808
 809                         for (int k = 0; k < seqs.size(); k++) {  seqsForSlayer.push_back(seqs[k].seq);  }
 810
 811                         //mask then send to slayer...
 812                         if (seqMask != "") {
 813                                 decalc->setMask(seqMask);
 814
 815                                 //mask querys
 816                                 decalc->runMask(query);
 817
 818                                 //mask parents
 819                                 for (int k = 0; k < seqsForSlayer.size(); k++) {
 820                                         decalc->runMask(seqsForSlayer[k]);
 821                                 }
 822
 823                                 spotMap = decalc->getMaskMap();
 824                         }
 825
 826                         if (m->control_pressed) {  for (int k = 0; k < seqs.size(); k++) {  delete seqs[k].seq;   }  return 0;  }
 827
 828                         //send to slayer
 829                         chimeraFlags = slayer.getResults(query, seqsForSlayer);
 830                         if (m->control_pressed) {  return 0;  }
 831                         chimeraResults = slayer.getOutput();
 832
 833                         //free memory
 834                         for (int k = 0; k < seqs.size(); k++) {  delete seqs[k].seq;   }
 835
 836                         printResults.spotMap = spotMap;
 837                         printResults.flag = chimeraFlags;
 838                         printResults.results = chimeraResults;
 839                 }
 840
 841                 return 0;
 842         }
 843         catch(exception& e) {
 844                 m->errorOut(e, "ChimeraSlayer", "getChimeras");
 845                 exit(1);
 846         }
 847 }
 848 //***************************************************************************************************************
 849 void ChimeraSlayer::printBlock(data_struct data, string flag, ostream& out){
 850         try {
 851                 out << querySeq->getName() << '\t';
 852                 out << data.parentA.getName() << "\t" << data.parentB.getName()  << '\t';
 853
 854                 out << data.divr_qla_qrb << '\t' << data.qla_qrb << '\t' << data.bsa << '\t';
 855                 out << data.divr_qlb_qra << '\t' << data.qlb_qra << '\t' << data.bsb << '\t';
 856
 857                 out << flag << '\t' << spotMap[data.winLStart] << "-" << spotMap[data.winLEnd] << '\t' << spotMap[data.winRStart] << "-" << spotMap[data.winREnd] << '\t';
 858
 859         }
 860         catch(exception& e) {
 861                 m->errorOut(e, "ChimeraSlayer", "printBlock");
 862                 exit(1);
 863         }
 864 }
 865 //***************************************************************************************************************
 866 void ChimeraSlayer::printBlock(data_results leftdata, data_results rightdata, bool leftChimeric, bool rightChimeric, string flag, ostream& out){
 867         try {
 868
 869                 if ((leftChimeric) && (!rightChimeric)) { //print left
 870                         out << querySeq->getName() << '\t';
 871                         out << leftdata.results[0].parentA.getName() << "\t" << leftdata.results[0].parentB.getName()  << '\t';
 872
 873                         out << leftdata.results[0].divr_qla_qrb << '\t' << leftdata.results[0].qla_qrb << '\t' << leftdata.results[0].bsa << '\t';
 874                         out << leftdata.results[0].divr_qlb_qra << '\t' << leftdata.results[0].qlb_qra << '\t' << leftdata.results[0].bsb << '\t';
 875
 876                         out << flag << '\t' << leftdata.spotMap[leftdata.results[0].winLStart] << "-" << leftdata.spotMap[leftdata.results[0].winLEnd] << '\t' << leftdata.spotMap[leftdata.results[0].winRStart] << "-" << leftdata.spotMap[leftdata.results[0].winREnd] << '\t';
 877
 878                 }else if ((!leftChimeric) && (rightChimeric)) {  //print right
 879                         out << querySeq->getName() << '\t';
 880                         out << rightdata.results[0].parentA.getName() << "\t" << rightdata.results[0].parentB.getName()  << '\t';
 881
 882                         out << rightdata.results[0].divr_qla_qrb << '\t' << rightdata.results[0].qla_qrb << '\t' << rightdata.results[0].bsa << '\t';
 883                         out << rightdata.results[0].divr_qlb_qra << '\t' << rightdata.results[0].qlb_qra << '\t' << rightdata.results[0].bsb << '\t';
 884
 885                         out << flag << '\t' << rightdata.spotMap[rightdata.results[0].winLStart] << "-" << rightdata.spotMap[rightdata.results[0].winLEnd] << '\t' << rightdata.spotMap[rightdata.results[0].winRStart] << "-" << rightdata.spotMap[rightdata.results[0].winREnd] << '\t';
 886
 887                 }else  { //print both results
 888                         if (leftdata.flag == "yes") {
 889                                 out << querySeq->getName() + "_LEFT" << '\t';
 890                                 out << leftdata.results[0].parentA.getName() << "\t" << leftdata.results[0].parentB.getName()  << '\t';
 891
 892                                 out << leftdata.results[0].divr_qla_qrb << '\t' << leftdata.results[0].qla_qrb << '\t' << leftdata.results[0].bsa << '\t';
 893                                 out << leftdata.results[0].divr_qlb_qra << '\t' << leftdata.results[0].qlb_qra << '\t' << leftdata.results[0].bsb << '\t';
 894
 895                                 out << flag << '\t' << leftdata.spotMap[leftdata.results[0].winLStart] << "-" << leftdata.spotMap[leftdata.results[0].winLEnd] << '\t' << leftdata.spotMap[leftdata.results[0].winRStart] << "-" << leftdata.spotMap[leftdata.results[0].winREnd] << '\t';
 896                         }
 897
 898                         if (rightdata.flag == "yes") {
 899                                 if (leftdata.flag == "yes") { out << endl; }
 900
 901                                 out << querySeq->getName() + "_RIGHT"<< '\t';
 902                                 out << rightdata.results[0].parentA.getName() << "\t" << rightdata.results[0].parentB.getName()  << '\t';
 903
 904                                 out << rightdata.results[0].divr_qla_qrb << '\t' << rightdata.results[0].qla_qrb << '\t' << rightdata.results[0].bsa << '\t';
 905                                 out << rightdata.results[0].divr_qlb_qra << '\t' << rightdata.results[0].qlb_qra << '\t' << rightdata.results[0].bsb << '\t';
 906
 907                                 out << flag << '\t' << rightdata.spotMap[rightdata.results[0].winLStart] << "-" << rightdata.spotMap[rightdata.results[0].winLEnd] << '\t' << rightdata.spotMap[rightdata.results[0].winRStart] << "-" << rightdata.spotMap[rightdata.results[0].winREnd] << '\t';
 908
 909                         }
 910                 }
 911         }
 912         catch(exception& e) {
 913                 m->errorOut(e, "ChimeraSlayer", "printBlock");
 914                 exit(1);
 915         }
 916 }
 917 //***************************************************************************************************************
 918 string ChimeraSlayer::getBlock(data_results leftdata, data_results rightdata, bool leftChimeric, bool rightChimeric, string flag){
 919         try {
 920
 921                 string out = "";
 922
 923                 if ((leftChimeric) && (!rightChimeric)) { //get left
 924                         out += querySeq->getName() + "\t";
 925                         out += leftdata.results[0].parentA.getName() + "\t" + leftdata.results[0].parentB.getName() + "\t";
 926
 927                         out += toString(leftdata.results[0].divr_qla_qrb) + "\t" + toString(leftdata.results[0].qla_qrb) + "\t" + toString(leftdata.results[0].bsa) + "\t";
 928                         out += toString(leftdata.results[0].divr_qlb_qra) + "\t" + toString(leftdata.results[0].qlb_qra) + "\t" + toString(leftdata.results[0].bsb) + "\t";
 929
 930                         out += flag + "\t" + toString(leftdata.spotMap[leftdata.results[0].winLStart]) + "-" + toString(leftdata.spotMap[leftdata.results[0].winLEnd]) + "\t" + toString(leftdata.spotMap[leftdata.results[0].winRStart]) + "-" + toString(leftdata.spotMap[leftdata.results[0].winREnd]) + "\t";
 931
 932                 }else if ((!leftChimeric) && (rightChimeric)) {  //print right
 933                         out += querySeq->getName() + "\t";
 934                         out += rightdata.results[0].parentA.getName() + "\t" + rightdata.results[0].parentB.getName()  + "\t";
 935
 936                         out += toString(rightdata.results[0].divr_qla_qrb) + "\t" + toString(rightdata.results[0].qla_qrb) + "\t" + toString(rightdata.results[0].bsa) + "\t";
 937                         out += toString(rightdata.results[0].divr_qlb_qra) + "\t" + toString(rightdata.results[0].qlb_qra) + "\t" + toString(rightdata.results[0].bsb) + "\t";
 938
 939                         out += flag + "\t" + toString(rightdata.spotMap[rightdata.results[0].winLStart]) + "-" + toString(rightdata.spotMap[rightdata.results[0].winLEnd]) + "\t" + toString(rightdata.spotMap[rightdata.results[0].winRStart]) + "-" + toString(rightdata.spotMap[rightdata.results[0].winREnd]) + "\t";
 940
 941                 }else  { //print both results
 942
 943                         if (leftdata.flag == "yes") {
 944                                 out += querySeq->getName() + "_LEFT\t";
 945                                 out += leftdata.results[0].parentA.getName() + "\t" + leftdata.results[0].parentB.getName() + "\t";
 946
 947                                 out += toString(leftdata.results[0].divr_qla_qrb) + "\t" + toString(leftdata.results[0].qla_qrb) + "\t" + toString(leftdata.results[0].bsa) + "\t";
 948                                 out += toString(leftdata.results[0].divr_qlb_qra) + "\t" + toString(leftdata.results[0].qlb_qra) + "\t" + toString(leftdata.results[0].bsb) + "\t";
 949
 950                                 out += flag + "\t" + toString(leftdata.spotMap[leftdata.results[0].winLStart]) + "-" + toString(leftdata.spotMap[leftdata.results[0].winLEnd]) + "\t" + toString(leftdata.spotMap[leftdata.results[0].winRStart]) + "-" + toString(leftdata.spotMap[leftdata.results[0].winREnd]) + "\t";
 951                         }
 952
 953                         if (rightdata.flag == "yes") {
 954                                 if (leftdata.flag == "yes") { out += "\n"; }
 955                                 out +=  querySeq->getName() + "_RIGHT\t";
 956                                 out += rightdata.results[0].parentA.getName() + "\t" + rightdata.results[0].parentB.getName()  + "\t";
 957
 958                                 out += toString(rightdata.results[0].divr_qla_qrb) + "\t" + toString(rightdata.results[0].qla_qrb) + "\t" + toString(rightdata.results[0].bsa) + "\t";
 959                                 out += toString(rightdata.results[0].divr_qlb_qra) + "\t" + toString(rightdata.results[0].qlb_qra) + "\t" + toString(rightdata.results[0].bsb) + "\t";
 960
 961                                 out += flag + "\t" + toString(rightdata.spotMap[rightdata.results[0].winLStart]) + "-" + toString(rightdata.spotMap[rightdata.results[0].winLEnd]) + "\t" + toString(rightdata.spotMap[rightdata.results[0].winRStart]) + "-" + toString(rightdata.spotMap[rightdata.results[0].winREnd]) + "\t";
 962                         }
 963                 }
 964
 965                 return out;
 966
 967         }
 968         catch(exception& e) {
 969                 m->errorOut(e, "ChimeraSlayer", "getBlock");
 970                 exit(1);
 971         }
 972 }
 973 //***************************************************************************************************************
 974 string ChimeraSlayer::getBlock(data_struct data, string flag){
 975         try {
 976
 977                 string outputString = "";
 978
 979                 outputString += querySeq->getName() + "\t";
 980                 outputString += data.parentA.getName() + "\t" + data.parentB.getName()  + "\t";
 981
 982                 outputString += toString(data.divr_qla_qrb) + "\t" + toString(data.qla_qrb) + "\t" + toString(data.bsa) + "\t";
 983                 outputString += toString(data.divr_qlb_qra) + "\t" + toString(data.qlb_qra) + "\t" + toString(data.bsb) + "\t";
 984
 985                 outputString += flag + "\t" + toString(spotMap[data.winLStart]) + "-" + toString(spotMap[data.winLEnd]) + "\t" + toString(spotMap[data.winRStart]) + "-" + toString(spotMap[data.winREnd]) + "\t";
 986
 987                 return outputString;
 988         }
 989         catch(exception& e) {
 990                 m->errorOut(e, "ChimeraSlayer", "getBlock");
 991                 exit(1);
 992         }
 993 }
 994 //***************************************************************************************************************/
 995