// molecule_smiles.cpp - Molecule's implementation of even more functions,
// notably, Structural Diagram Generation and SMILES

#include <qobject.h>
#include <qlist.h>
#include <qarray.h>
#include <qvector.h>
#include <qstack.h>
#include <qstring.h>
#include <qmessagebox.h>
#include <iostream.h>
#include <qfile.h>
#include <qtextstream.h>

#ifdef UNIX
#include <stdlib.h>
#include <time.h>
#endif

#include "render2d.h"
#include "drawable.h"
#include "molecule.h"
#include "dpoint.h"
#include "text.h"
#include "atom.h"
#include "ring.h"
#include "setofrings.h"
#include "sdg.h"
#include "defs.h"

// CleanUp(): Invoke SDG() to clean up structure
void Molecule::CleanUp() {
  SDG(true);
}

// Structure Diagram Generation - add coordinates to connectivity
// generally called after reading in a file or SMILES string which did
// not supply coordinates.  Could also be used to clean up Molecule.
// bool coord: have coordinates been set?  (true = coordinates exist, e.g.
// from file or hand drawing; false = no, strcuture supplied as connection
// table from SMILES or file)
// Method references:
// 1) Ugi I et al., Journal of Chemical research (M), 1991, 2601-2689
// 2) Christoph Steinbeck's Java implentation of above, JMDraw.
//    http://jmdraw.sourceforge.net/
void Molecule::SDG(bool coord) {
  QRect bb1;
  QPoint center1, center2;

  if (coord == true) {  // if coordinates exist, save center of bounding box
    bb1 = BoundingBoxAll();
    center1 = bb1.center();
  }

  //DPoint *t2;
  Atom *a1;

  // get unique points
  up = AllPoints();
  QVector<Atom> atoms(up.count());

  // clear "hit" flag on all atoms
  for (tmp_pt = up.first(); tmp_pt != 0; tmp_pt = up.next())
    tmp_pt->hit = false;
  // find rings (esp. find aromaticity) - do after CopyTextToDPoint()
  MakeSSSR();

  // convert "up" to JMDraw-friendly Qlist<Atom>
  // rebuild neighbors list (usually mangled by MakeSSSR)
  for (tmp_pt = up.first(); tmp_pt != 0; tmp_pt = up.next() ) {
    tmp_pt->neighbors.clear();
    for(tmp_bond = bonds.first(); tmp_bond != 0; tmp_bond = bonds.next() ) {
      if (tmp_bond->Find(tmp_pt) == true) {
	tmp_pt->neighbors.append(tmp_bond->otherPoint(tmp_pt));
	tmp_pt->bondorder.at(tmp_pt->neighbors.count() - 1) = tmp_bond->Order();
      }
    }
  }
  // first copy all DPoints
  int c1, c2, refnum;
  for (c1 = 0; c1 < up.count(); c1++) {
    tmp_pt = up.at(c1);
    a1 = new Atom(tmp_pt->element, tmp_pt->x, tmp_pt->y, tmp_pt->z);
    a1->number = tmp_pt->serial;
    a1->degree = tmp_pt->neighbors.count();
    //cout << c1 << "-degree-" << a1->degree << endl;
    atoms.insert(c1, a1);
  }
  // now build connectivity table
  for (c1 = 0; c1 < up.count(); c1++) {
    tmp_pt = up.at(c1);
    a1 = atoms.at(c1);
    for (c2 = 0; c2 < tmp_pt->neighbors.count(); c2++) {
      refnum = tmp_pt->neighbors.at(c2)->serial;
      a1->nodeTable.insert(c2, atoms.at(refnum));
      a1->bondTable[c2] = tmp_pt->bondorder[c2];
      a1->intnodeTable[c2] = refnum;
    }
  }
  /* convert this_sssr to setOfRings
  setOfRings s1;
  s1.resize(this_sssr.sssr.count());
  Ring *ring1;
  QList<DPoint> *tmp_ring;
  c1 = 0;
  for (tmp_ring = this_sssr.sssr.first(); tmp_ring != 0; 
       tmp_ring = this_sssr.sssr.next()) {
    ring1 = new Ring;
    ring1->resize(tmp_ring->count());
    c2 = 0;
    for (tmp_pt = tmp_ring->first(); tmp_pt != 0; tmp_pt = tmp_ring->next()) {
      a1 = atoms.at(tmp_pt->serial);
      ring1->insert(c2, a1);
      cout << "RA" << c2 << ":" << tmp_pt->serial << endl;
      c2++;
    }
    ring1->sort2();
    s1.insert(c1, ring1);
    c1++;
  }
  */
  // pass to SDG class.
  class SDG sdg1;
  sdg1.setAtoms(atoms);
  //sdg1.setRings(s1);
  sdg1.setBL(r->getFixedLength_bond());
  sdg1.exec();
  atoms = sdg1.getAtoms();

  cout << "SDG succeeded!" << endl;

  // convert atoms back to DPoint (essentially, just update x,y coordinates)
  for (c1 = 0; c1 < up.count(); c1++) {
    tmp_pt = up.at(c1);
    a1 = atoms.at(c1);
    tmp_pt->x = a1->x;
    tmp_pt->y = a1->y;
  }

  bb1 = BoundingBoxAll();
  int xmove = 0, ymove = 0;

  if (coord == true) { // if coordinates existed, move back into place
    center2 = bb1.center();
    xmove = center1.x() - center2.x();
    ymove = center1.y() - center2.y();
  } else { // move to top left of screen
    if (bb1.left() < 10) xmove = 10 - bb1.left();
    if (bb1.top() < 10) ymove = 10 - bb1.top();
  }
  for (tmp_pt = up.first(); tmp_pt != 0; tmp_pt = up.next()) {
    tmp_pt->x += xmove;
    tmp_pt->y += ymove;
  }

  // add hydrogens
  AddHydrogens();
}

// cheat: use Babel to make SMILES strings.
QString Molecule::ToSMILES() {
  QString smiles;

#ifdef UNIX
  QString molfile = ToMDLMolfile();
  QString babel_cmdl, babel_in, babel_out;

  srand(time(NULL));

  babel_in.setNum(rand());
  babel_in.prepend("/tmp/");
  babel_in.append(".mol");

  babel_out.setNum(rand());
  babel_out.prepend("/tmp/");
  babel_out.append(".smi");

  babel_cmdl = "babel -imdl ";
  babel_cmdl.append(babel_in);
  babel_cmdl.append(" -osmi ");
  babel_cmdl.append(babel_out);

  cout << babel_cmdl << endl;

  // create mol file
  QFile f1(babel_in);

  if (!f1.open(IO_WriteOnly)) return smiles;
  
  QTextStream t1(&f1);

  t1 << molfile;

  f1.close();

  FILE *test1;
  test1 = popen(babel_cmdl.latin1(), "r");
  int rv = pclose(test1);
  cout << "pclose():" << rv << endl;  

  QFile f2(babel_out);

  if (!f2.open(IO_ReadOnly)) return smiles;
  
  QTextStream t2(&f2);

  t2 >> smiles;

  f2.close();

  babel_cmdl = "rm " + babel_in;
  test1 = popen(babel_cmdl.latin1(), "r");
  rv = pclose(test1);
  cout << "pclose():" << rv << endl;  
  babel_cmdl = "rm " + babel_out;
  test1 = popen(babel_cmdl.latin1(), "r");
  rv = pclose(test1);
  cout << "pclose():" << rv << endl;  
#endif

  return smiles;
}

// convert this Molecule to a SMILES string
// Method references:  
// 1) Weininger D, JCICS 1988, 28, 31-36
// 2) Weininger D, Weininger A, Weininger JL, JCICS 1989, 29, 97-101
/*
QString Molecule::ToSMILES() {
  QList<DPoint> up;
  QList<Bond> broken_bonds; // copy broken bonds here, then back afterwards
  QString tmp_str, ts;
  //DPoint *tmp_pt2, *tmp_pt3;
  //int smiles_ring = 1;  // ring closure number

  // get list of unique points
  up = AllPoints();
  // clear "hit" flag on all atoms
  for (tmp_pt = up.first(); tmp_pt != 0; tmp_pt = up.next())
    tmp_pt->hit = false;
  // find rings (esp. find aromaticity) - do after CopyTextToDPoint()
  MakeSSSR();

  // break rings
  // first break bonds involved in multiple rings


  // find a starting point.  First monosubstituted carbon should do


  // unbreak rings
  if (broken_bonds.count() > 0) {
    for (tmp_bond = broken_bonds.first(); tmp_bond != 0; 
	 tmp_bond = broken_bonds.next()) {
      bonds.append(tmp_bond);
    }
  }

  return "";
}
*/

// convert SMILES string to Molecule
// (Ideally, you should call this function just after creating)
void Molecule::FromSMILES(QString sm) {
  QStringList tokens;
  QStack<DPoint> branch_tree;
  QVector<DPoint> ring_closure_array(10);
  bool ring_array_status[10]; // track which elements of ring array are used

  for (int cc = 0; cc < 10; cc++)
    ring_array_status[cc] = false;

  // tokenize
  // tokens: atoms, groups enclosed in [], (, )
  // note that numbers and symbols outside [] will break the tokenizer
  QStringList smilesTokens;
  QString prev_token, tmp_token;
  int i1;

  cout << "SMILES:" << sm << "|" << endl;

  do {
    tmp_token = "";
    if (sm[0] == '=') { // double bond
      tmp_token.append("=");
      sm.remove(0,1);
    }
    if (sm[0] == '#') { // triple bond
      tmp_token.append("#");
      sm.remove(0,1);
    }
    if (sm[0].isLetter()) {
      // extract letter token
      // look for single-letter aromatic
      if (sm[0].lower() == sm[0]) {
	tmp_token.append(sm.left(1));
	sm.remove(0,1);
      } else {
	if ( (sm[1].lower() == sm[1]) &&
	     (sm[1].isLetter()) ) { // lowercase; two-letter symbol
	  tmp_token.append(sm.left(2));
	  sm.remove(0,2);
	} else { //
	  tmp_token.append(sm.left(1));
	  sm.remove(0,1);
	}
      }
      // extract ring closure numbers
      if (sm.length() > 0) {
	do {
	  if (sm[0].isNumber()) {
	    tmp_token.append(sm.left(1));
	    sm.remove(0,1);
	  } else {
	    break;
	  }
	} while (sm.length() > 0);
      }
      smilesTokens.append(tmp_token);
    }
    if (sm[0] == '(') {
      tmp_token = "(";
      smilesTokens.append(tmp_token);
      sm.remove(0,1);
    }
    if (sm[0] == ')') {
      tmp_token = ")";
      smilesTokens.append(tmp_token);
      sm.remove(0,1);
    }
    if (sm[0] == '[') {
      i1 = sm.find("]");
      tmp_token = sm.left(i1);
      smilesTokens.append(tmp_token);
    }
    //cout << "token: " << tmp_token << endl << "left: " << sm << endl;
  } while (sm.length() > 0);

  DPoint *prev_pt = 0, *new_pt = 0;
  QString tmp_element, tmp_element_mask;
  tmp_token = "";
  prev_token = "";
  int bond_order = 0;
  bool aromatic = false, flag = false;

  for ( QStringList::Iterator it = smilesTokens.begin(); 
	it != smilesTokens.end(); 
	++it ) {
    prev_token = tmp_token;
    tmp_token = (*it).latin1();

    cout << "token: " << tmp_token << endl;

    // process tokens
    if (tmp_token == "(") { // start branch
      branch_tree.push(prev_pt);
      continue;
    }
    if (tmp_token == ")") { // end branch
      prev_pt = branch_tree.pop();
      continue;
    }
    // if not a branch, it's probably an atom
    new_pt = new DPoint;
    // calculate bond order
    aromatic = false;
    if (tmp_token[0].isLetter()) {
      if (tmp_token[0].lower() == tmp_token[0]) new_pt->aromatic = true;
    }
    if (prev_pt != 0) {
      bond_order = 1;
      if ( prev_pt->aromatic && new_pt->aromatic ) {bond_order = 4;}
    }
    if (tmp_token.left(1) == "=") {
      bond_order = 2;
      tmp_token.remove(0,1);
    }
    if (tmp_token.left(1) == "#") {
      bond_order = 3;
      tmp_token.remove(0,1);
    }
    // extract element info
    tmp_element = ""; tmp_element_mask = "";
    if (tmp_token[0].isLetter()) {
      tmp_element = tmp_token.left(1);
      tmp_element = tmp_element.upper();
    }
    if (tmp_token[1].isLetter()) tmp_element.append(tmp_token[1]);
    if (tmp_token[0] == '[') {
      int i1 = tmp_token.find("]");
      tmp_element = tmp_token.mid(1, i1 - 2);
    }
    cout << "element: " << tmp_element << endl;
    tmp_element_mask = tmp_element;
    tmp_element_mask.fill(' ');
    // add super/sub-script where appropriate

    new_pt->element = tmp_element;
    new_pt->elementmask = tmp_element_mask;
    // create Text where appropriate
    if (new_pt->element != "C") {
      Text *nt = new Text(r);
      nt->setPoint(new_pt);
      nt->setJustify(JUSTIFY_CENTER);
      nt->Highlight(false);
      nt->setText(tmp_element);
      nt->setTextMask(tmp_element_mask);
      labels.append(nt);
    }
    // create bond
    if (prev_pt != 0) {
      addBond(prev_pt, new_pt, 1, bond_order, QColor(0,0,0), true);
      prev_pt = new_pt;
    } else {
      prev_pt = new_pt;
    }
    // handle ring closure
    do {
      flag = false;
      if (tmp_token.at(tmp_token.length() - 1).isNumber()) {
	int ringnum = tmp_token.right(1).toInt();
	flag = true;
	if (ring_array_status[ringnum] == false) { // save this atom
	  ring_closure_array.insert(ringnum, new_pt);
	  ring_array_status[ringnum] = true;
	} else { // do ring closure
	  tmp_pt = ring_closure_array.at(ringnum);
	  if (tmp_pt->aromatic && new_pt->aromatic)
	    addBond(tmp_pt, new_pt, 1, 4, QColor(0,0,0), true);
	  else
	    addBond(tmp_pt, new_pt, 1, 1, QColor(0,0,0), true);
	  ring_array_status[ringnum] = false;
	}
	tmp_token.remove(tmp_token.length() - 1, 1);
      }
    } while (flag == true); 
  }

  // print atom list and connection table?
  cout << ToXML("smiles") << endl;

  SDG(false);  // generate structure coordinates
}