#include "alignedSequence.h"
#include "alignedStructure.h"
#include "alphabet.h"
#include "alphabetBuilder.h"
#include "coordinate3D.h"
#include "fastaReader.h"
#include "pdbReader.h"
#include "qTools.h"
#include "rmsdTools.h"
#include "sequence.h"
#include "sequenceAlignment.h"
#include "structure.h"
#include "structureAlignment.h"
#include "structureAlignmentReader.h"
#include "structureQR.h"
#include "symbol.h"
#include "symbolList.h"

#include <stdio.h>

enum {DNA, RNA, PROTEIN};

void printCoord3D(Coordinate3D* coord) {
  
  if (coord != 0) {
    printf("(%10f,%10f,%10f)\n",coord->getX(),coord->getY(),coord->getZ());
  }
  else {
    printf("0\n");
  }

  return;
}


// getStructureAlignment
//
StructureAlignment* getStructureAlignment(char* inputDir, int alphaType) {

    char* fastaFile = new char[strlen(inputDir) + 16];
    strcpy(fastaFile,inputDir);
    if (fastaFile[strlen(inputDir)] == '/') {
      strcat(fastaFile,"multiseq.fasta");
    }
    else {
      strcat(fastaFile,"/multiseq.fasta");
    }
    AlphabetBuilder* alphaBuild1 = new AlphabetBuilder();
    Alphabet* alpha1 = 0;
    switch (alphaType) {
    case DNA:
      alpha1 = alphaBuild1->getDnaAlphabet();
      break;
    case RNA:
      alpha1 = alphaBuild1->getRnaAlphabet();
      break;
    case PROTEIN:
      alpha1 = alphaBuild1->getProteinAlphabet();
      break;
    default:
      printf("Error - getStructureAlignment; invalid alphaType\n");
    }
    StructureAlignmentReader* strAlnRead1 = new StructureAlignmentReader(alpha1);
    strAlnRead1->setAlignmentFilename(fastaFile);
    StructureAlignment* structAln = strAlnRead1->getStructureAlignment();
    delete fastaFile;
    delete alphaBuild1;
    delete strAlnRead1;

    return structAln;
}


int main (int argc, char** argv) {
  
  
  // AlphabetBuilder
  
  printf("*** AlphabetBuilder Test ***\n");
  AlphabetBuilder* alphaBuild2 = new AlphabetBuilder();
  printf("Hey!\n");
  Alphabet* alpha3 = alphaBuild2->getDnaAlphabet();
  char* str2 = alpha3->toString();
  printf("alpha3: %s\n",str2);
  delete str2;
  printf("\n\n");

  // Symbol
  printf("*** Symbol Test ***\n");
  Symbol* sym0 = new Symbol('A',"ADE","Adenine");
  Symbol* sym1 = new Symbol('C',"CYT","Cytosine");
  Symbol* sym2 = new Symbol('G',"GUA","Guanine");
  Symbol* sym3 = new Symbol('U',"URA","Uracil");
  Symbol* sym4 = new Symbol('-'," - ","Gap");
  Symbol* sym5 = new Symbol('X',"XXX","Unknown");
  Symbol** symArray = new Symbol* [6];
  symArray[0] = sym0;
  symArray[1] = sym1;
  symArray[2] = sym2;
  symArray[3] = sym3;
  symArray[4] = sym4;
  symArray[5] = sym5;
  for (int i=0; i<6; i++) {
    printf("sym%d: %c, %s, %s\n",i,symArray[i]->getOne(),symArray[i]->getThree(),symArray[i]->getFull());
  }
  printf("\n\n");


  // Alphabet
  printf("*** Alphabet Test ***\n");
  Alphabet* alpha1 = new Alphabet(6,symArray);
  char* str = alpha1->toString();
  printf("alpha1: %s\n",str);
  delete str;
  delete symArray[1];
  printf("\n\n");


  // SymbolList
  printf("*** SymbolList Test ***\n");
  SymbolList* symList1 = new SymbolList(10,alpha1);
  symList1->addSymbol('C');
  symList1->addSymbol('U');
  symList1->addSymbol('A');
  symList1->addSymbol('G');
  symList1->addSymbol('G');
  symList1->addSymbol('U');
  symList1->addSymbol('A');
  symList1->addSymbol('C');
  symList1->addSymbol('G');
  symList1->addSymbol('U');
  str = symList1->toString();
  printf("symList1: %s\n",str);
  //delete str;
  printf("\n\n");
  

  // AlphabetBuilder
  
  printf("*** AlphabetBuilder Test ***\n");
  AlphabetBuilder* alphaBuild1 = new AlphabetBuilder();
  Alphabet* alpha2 = alphaBuild1->getProteinAlphabet();
  //alpha2->setBackbone("CA");
  str = alpha2->toString();
  printf("alpha2: %s\n",str);
  delete str;
  printf("\n\n");
  

  // Sequence
  printf("*** Sequence Test ***\n");
  Sequence* seq1 = new Sequence(10,alpha1);
  seq1->addSymbol('C');
  seq1->addSymbol('U');
  seq1->addSymbol('A');
  seq1->addSymbol('G');
  seq1->addSymbol('G');
  seq1->addSymbol('U');
  seq1->addSymbol('A');
  seq1->addSymbol('C');
  seq1->addSymbol('G');
  seq1->addSymbol('U');
  str = seq1->toString();
  printf("seq1: %s\n",str);
  delete str;
  printf("\n\n");

  
  // Coordinate3D
  printf("*** Coordinate3D Test ***\n");
  Coordinate3D* coord1 = new Coordinate3D(1.0,2.0,3.0);
  Coordinate3D* coord2 = new Coordinate3D(3.14,2.17,4.2);
  Coordinate3D* coord3 = new Coordinate3D(3.14,2.17,4.2);
  printf("coord1: ");
  printCoord3D(coord1);
  printf("coord2: ");
  printCoord3D(coord2);
  printf("coord3: ");
  printCoord3D(coord3);
  printf("coord1 == coord2: %d\n",coord1->equals(coord2));
  printf("coord1 == coord3: %d\n",coord1->equals(coord3));
  printf("coord2 == coord3: %d\n",coord2->equals(coord3));
  printf("\n\n");
  

  // Structure
  printf("*** Structure Test ***\n");
  
  Structure* struct1 = new Structure(10,alpha1);
  struct1->addElement('C',1.0,2.0,3.0);
  struct1->addElement('U',2.0,2.1,-3.0);
  struct1->addElement('A',3.0,2.2,3.0);
  struct1->addElement('G',4.0,2.3,-3.0);
  struct1->addElement('G',5.0,2.4,3.0);
  struct1->addElement('U',6.0,2.5,-3.0);
  struct1->addElement('A',7.0,2.6,3.0);
  struct1->addElement('C',8.0,2.7,-3.0);
  struct1->addElement('G',9.0,2.8,3.0);
  struct1->addElement('U',10.0,2.9,-3.0);
  str = struct1->toString();
  printf("struct1: %s\n",str);
  delete str;
  printf("getCoordinate(0): ");
  printCoord3D(struct1->getCoordinate(0));
  printf("getCoordinate(1): ");
  printCoord3D(struct1->getCoordinate(1));
  printf("getCoordinate(2): ");
  printCoord3D(struct1->getCoordinate(2));
  printf("getCoordinate(3): ");
  printCoord3D(struct1->getCoordinate(3));
  printf("getCoordinate(4): ");
  printCoord3D(struct1->getCoordinate(4));
  printf("getCoordinate(5): ");
  printCoord3D(struct1->getCoordinate(5));
  printf("getCoordinate(6): ");
  printCoord3D(struct1->getCoordinate(6));
  printf("getCoordinate(7): ");
  printCoord3D(struct1->getCoordinate(7));
  printf("getCoordinate(8): ");
  printCoord3D(struct1->getCoordinate(8));
  printf("getCoordinate(9): ");
  printCoord3D(struct1->getCoordinate(9));
  printf("getCoordinate(19): ");
  printCoord3D(struct1->getCoordinate(19));
  printf("\n\n");
  
  

  // AlignedSequence
  printf("*** AlignedSequence Test ***\n");
  AlignedSequence* alSeq1 = new AlignedSequence(15,alpha1);
  alSeq1->addSymbol('C');
  alSeq1->addSymbol('U');
  alSeq1->addGap();
  alSeq1->addGap();
  alSeq1->addSymbol('A');
  alSeq1->addSymbol('G');
  alSeq1->addSymbol('G');
  alSeq1->addSymbol('U');
  alSeq1->addGap();
  alSeq1->addGap();
  alSeq1->addSymbol('A');
  alSeq1->addSymbol('C');
  alSeq1->addGap();
  alSeq1->addSymbol('G');
  alSeq1->addSymbol('U');
  str = alSeq1->toString();
  printf("alSeq1: %s\n",str);
  delete str;
  printf("\n\n");
  

  // AlignedStructure
  printf("*** AlignedStructure Test ***\n");
  AlignedStructure* alStruct1 = new AlignedStructure(15,alpha1);
  alStruct1->addElement('C',1.0,2.0,3.0);
  alStruct1->addElement('U',1.1,2.1,3.1);
  alStruct1->addGap();
  alStruct1->addGap();
  alStruct1->addElement('A',1.2,2.2,3.2);
  alStruct1->addElement('G',1.3,2.3,3.3);
  alStruct1->addElement('G',1.4,2.4,3.4);
  alStruct1->addElement('U',1.5,2.5,3.5);
  alStruct1->addGap();
  alStruct1->addGap();
  alStruct1->addElement('A',1.6,2.6,3.6);
  alStruct1->addElement('C',1.7,2.7,3.7);
  alStruct1->addGap();
  alStruct1->addElement('G',1.8,2.8,3.8);
  alStruct1->addElement('U',1.9,2.9,3.9);
  str = alStruct1->toString();
  printf("alStruct1: %s\n",str);
  delete str;
  for (int i=0; i<16; i++) {
    printf("getCoordinate(%d): ", i);
    printCoord3D(alStruct1->getCoordinate(i));
  }
  for (int i=0; i<11; i++) {
    printf("getUnalignedCoordinate(%d): ", i);
    printCoord3D(alStruct1->getUnalignedCoordinate(i));
  }
  printf("\n\n");
  

  // SequenceAlignment
  printf("*** SequenceAlignment Test ***\n");

  SequenceAlignment* seqAl1 = new SequenceAlignment(10,5);
  //printf("Hey!\n");
  AlignedSequence* alSeq2 = new AlignedSequence(10,alpha2);
  alSeq2->addSymbols("MHPQEDCAT-");
  AlignedSequence* alSeq3 = new AlignedSequence(10,alpha2);
  alSeq3->addSymbols("MYP--DSV-A");
  AlignedSequence* alSeq4 = new AlignedSequence(10,alpha2);
  alSeq4->addSymbols("-HP-EDCATA");
  AlignedSequence* alSeq5 = new AlignedSequence(10,alpha2);
  alSeq5->addSymbols("MY---DC---");
  AlignedSequence* alSeq6 = new AlignedSequence(10,alpha2);
  alSeq6->addSymbols("MYP---CAT-");
  //printf("Hey!!\n");

  seqAl1->addSequence(alSeq2);
  seqAl1->addSequence(alSeq3);
  seqAl1->addSequence(alSeq4);
  seqAl1->addSequence(alSeq5);
  seqAl1->addSequence(alSeq6);
  //printf("Hey!!!\n");

  AlignedSequence* tempAlSeq;
  for (int i=0; i<5; i++) {
    tempAlSeq = seqAl1->getSequence(i);
    if (tempAlSeq != 0) {
      str = tempAlSeq->toString();
      printf("seqAl1[%d]: %s\n",i,str);
      delete str;
    }
  }
  printf("\n\n");


  // StructureAlignment
  printf("*** StructureAlignment Test ***\n");

  StructureAlignment* structAl1 = new StructureAlignment(10,3);
  AlignedStructure* alStruct2 = new AlignedStructure(10,alpha1);
  alStruct2->addElement('C',1.0,2.0,3.0);
  alStruct2->addElement('U',1.1,2.1,3.1);
  alStruct2->addGap();
  alStruct2->addGap();
  alStruct2->addElement('A',1.2,2.2,3.2);
  alStruct2->addElement('G',1.3,2.3,3.3);
  alStruct2->addElement('G',1.4,2.4,3.4);
  alStruct2->addElement('U',1.5,2.5,3.5);
  alStruct2->addGap();
  alStruct2->addGap();

  AlignedStructure* alStruct3 = new AlignedStructure(10,alpha1);
  alStruct3->addElement('C',1.0,2.0,3.0);
  alStruct3->addElement('A',1.1,2.1,3.1);
  alStruct3->addElement('A',1.1,2.1,3.1);
  alStruct3->addGap();
  alStruct3->addElement('A',1.2,2.2,3.2);
  alStruct3->addElement('G',1.3,2.3,3.3);
  alStruct3->addElement('G',1.4,2.4,3.4);
  alStruct3->addElement('U',1.5,2.5,3.5);
  alStruct3->addElement('A',1.6,2.6,3.6);
  alStruct3->addElement('U',1.7,2.7,3.7);

  AlignedStructure* alStruct4 = new AlignedStructure(10,alpha1);
  alStruct4->addElement('G',1.0,2.0,3.0);
  alStruct4->addElement('U',1.1,2.1,3.1);
  alStruct4->addElement('A',1.1,2.1,3.1);
  alStruct4->addGap();
  alStruct4->addGap();
  alStruct4->addElement('G',1.3,2.3,3.3);
  alStruct4->addElement('G',1.4,2.4,3.4);
  alStruct4->addElement('U',1.5,2.5,3.5);
  alStruct4->addGap();
  alStruct4->addElement('U',1.6,2.6,3.6);

  //printf("   About to add Structures\n");
  int check = 0;
  check = structAl1->addStructure(alStruct2);
  //printf("   check: %d\n",check);
  check = structAl1->addStructure(alStruct3);
  //printf("   check: %d\n",check);
  check = structAl1->addStructure(alStruct4);
  //printf("   check: %d\n",check);

  //printf("   Hey!!!\n");
  AlignedStructure* tempAlStruct;
  for (int i=0; i<3; i++) {
    tempAlStruct = structAl1->getStructure(i);
    if (tempAlStruct != 0) {
      str = tempAlStruct->toString();
      printf("structAl1[%d]: %s\n",i,str);
      delete str;
    }
  }
  printf("\n\n");


  // FASTAReader
  printf("*** FASTAReader Test ***\n");
  FASTAReader* fasRead1 = new FASTAReader(alpha2);
  printf("   Hey!\n");

  check = fasRead1->setFilename("data/test.fasta");
  if (check == 0) {
    printf("   File not found\n");
  }
  AlignedSequence* alSeq7 = fasRead1->getAlignedSequence(0);
  AlignedSequence* alSeq8 = fasRead1->getNextAlignedSequence();
  AlignedSequence* alSeq9 = fasRead1->getAlignedSequence(5);

  if (alSeq7 != 0) {
    str = alSeq7->toString();
    printf("   alSeq7: %s\n",str);
    delete str;
  }
  else {
    printf("   alSeq7: no AlignedSequence returned\n");
  }
  if (alSeq8 != 0) {
    str = alSeq8->toString();
    printf("   alSeq8: %s\n",str);
    delete str;
  }
  else {
    printf("   alSeq8: no AlignedSequence returned\n");
  }
  if (alSeq9 != 0) {
    str = alSeq9->toString();
    printf("   alSeq9: %s\n",str);
    delete str;
  }
  else {
    printf("   alSeq9: no AlignedSequence returned\n");
  }


  SequenceAlignment* seqAl2 = fasRead1->getSequenceAlignment();
  //printf("Hey!\n");

  if (seqAl2 == 0) {
    printf("   seqAl2: no SequenceAlignment returned\n");
  }
  else {
    printf("*** sequenceCount: %d\n",fasRead1->getSequenceCount());
    for (int i=0; i<fasRead1->getSequenceCount(); i++) {
      tempAlSeq = seqAl2->getSequence(i);
      if (tempAlSeq != 0) {
	str = tempAlSeq->toString();
	printf("   %s[%d]: %s\n",tempAlSeq->getName(),i,str);
	delete str;
      }
      else {
	printf("   seqAl2[%d]: no AlignedSequence returned\n",i);
      }
    }
  }
  printf("\n\n");


  // PDBReader
  printf("*** PDBReader Test ***\n");
  PDBReader* pdbRead1 = new PDBReader(alpha2);
  pdbRead1->setPath("data/");
  pdbRead1->setFilename("1RIQ.pdb");
  Structure* struct2 = pdbRead1->getStructure(0);
  printf("name: %s\n",struct2->getName());
  for (int i=0; i<struct2->getLength(); i++) {
    printf("%c: ",struct2->getSymbol(i)->getOne());
    printCoord3D(struct2->getCoordinate(i));
  }
  printf("\n\n");


  // StructureAlignmentReader
  printf("*** StructureAlignmentReader Test ***\n");
  StructureAlignmentReader* strAlnRead1 = new StructureAlignmentReader(alpha2);
  check = strAlnRead1->setAlignmentPath("data/");
  if (check != 1) {
    printf("   check failed\n");
  }
  check = strAlnRead1->setStructurePath("data/");
  if (check != 1) {
    printf("   check failed\n");
  }
  check = strAlnRead1->setAlignmentFilename("multiseq.fasta");
  if (check != 1) {
    printf("   check failed\n");
  }

  SequenceAlignment* seqAl3 = strAlnRead1->getSequenceAlignment();
  if (seqAl3 == 0) {
    printf("   seqAl3: no SequenceAlignment returned\n");
  }
  else {
    printf("*** sequenceCount: %d\n",seqAl3->getSequenceCount());
    for (int i=0; i<seqAl3->getSequenceCount(); i++) {
      tempAlSeq = seqAl3->getSequence(i);
      if (tempAlSeq != 0) {
	str = tempAlSeq->toString();
	printf("   %s[%d]: %s\n",tempAlSeq->getName(),i,str);
	delete str;
      }
      else {
	printf("   seqAl3[%d]: no AlignedSequence returned\n",i);
      }
    }
  }
  printf("Hey\n");
  StructureAlignment* structAl2 = strAlnRead1->getStructureAlignment();
  printf("Hey2\n");
  if (structAl2 != 0) {
  printf("Hey3\n");
    for (int i=0; i<structAl2->getStructureCount(); i++) {
      tempAlStruct = structAl2->getStructure(i);
      //printf("   Hey!\n");
      if (tempAlStruct != 0) {
	str = tempAlStruct->toString();
	printf("structAl2[%d] (%s): %s\n",i,tempAlStruct->getName(),str);
	delete str;
      }
    }
  }
  printf("Hey4\n");
  printf("\n\n");


  // QTools
  printf("*** QTools Test ***\n");
  QTools* qTools = new QTools(structAl2);
  printf("Running Q\n");
  qTools->q(1);
  FILE *out = fopen("q.txt","w");
  fprintf(out,"q:\n");
  qTools->printQ(out);
  //qTools->printQ(stdout);
  qTools->qPerResidue();
  qTools->printQPerResidue(out);
  //qTools->printQPerResidue(stdout);
  fclose(out);
  printf("\n\n");


  // StructureQR
  printf("*** StructureQR Test ***\n");
  StructureQR* structQR = new StructureQR(structAl2);
  structQR->qr();
  structQR->printColumns(stdout);
  printf("\n\n");
  
  
  // RmsdTools
  
  printf("*** RmsdTools Test **\n");
  RmsdTools* rmsdTools = new RmsdTools(structAl2);
  printf("Running RMSD\n");
  rmsdTools->rmsd();
  rmsdTools->printRmsd(stdout);
  printf("Running RMSD Per Residue\n");
  rmsdTools->rmsdPerResidue(0,1);
  rmsdTools->printRmsdPerResidue(stdout);
  printf("\n\n");
  


  /*
    int alphaType = PROTEIN;
    
    if (argc == 2) {
    char* inputDir = argv[1];
    
    StructureAlignment* structAl1 = getStructureAlignment(inputDir, alphaType);
    
    QTools* qTools = new QTools(structAl1);
    qTools->q(1);
    FILE *out = fopen("multiseq.matrix","w");
    qTools->printQ(out);
    fclose(out);
    delete structAl1;
    }
    else if (argc == 3) {
    char* params = argv[1];
    char* inputDir = argv[2];
    if (strcmp(params,"-p") == 0) {
    StructureAlignment* structAl1 = getStructureAlignment(inputDir, alphaType);
    
    QTools* qTools = new QTools(structAl1);
    qTools->qPerResidue();
    FILE *out = fopen("multiseq.qPerResidue","w");
    qTools->printQPerResidue(out);
    fclose(out);
    delete structAl1;
    }
    else if (strcmp(params,"-q") == 0) {
    StructureAlignment* structAl1 = getStructureAlignment(inputDir, alphaType);
    
    StructureQR* structQR = new StructureQR(structAl1);
    structQR->qr();
    structQR->printColumns();
    delete structAl1;
    }
    else {
    printf("Usage: qPair [-p|-q] inputdir\n");
    printf("   where inputdir is a directory containing\n");
    printf("   multiseq.fasta and pdb files for all of\n");
    printf("   the sequences in the alignment\n");
    printf("\n");
    printf(" p - calculate q per residue\n");
    printf("\n");
    printf(" q - find nonredundant set of structures\n");
    exit(1);
    }
    }
    else {
    printf("Usage: qPair [-p|-q] inputdir\n");
    printf("   where inputdir is a directory containing\n");
    printf("   multiseq.fasta and pdb files for all of\n");
    printf("   the sequences in the alignment\n");
    printf("\n");
    printf(" p - calculate q per residue\n");
    printf("\n");
    printf(" q - find nonredundant set of structures\n");
    exit(1);
    }
    
    return 0;
  */
}