#!/usr/bin/python
# Just splits phosphat_ester in parts and marks the connecting atom to the phosphorus with a copper (will be removed later)
from openeye.oechem import *
import os, string, sys

try:
    ifs = oemolistream()
    ifs.open(sys.argv[1])
    ofs = oemolostream()
    date = sys.argv[1][string.find(sys.argv[1], '_'):string.find(sys.argv[1], '_out')]
    ofs.open('b8.1'+date+'_ln_temp.ism')
except IndexError:
    print '\nUsage: b8.1_thio_parts_split_1.py molfile\n'
    raise SystemExit()

copper = OEGraphMol()
OEParseSmiles(copper, '[Cu]')

deprot = OEUniMolecularRxn("[P:1](=[S:2])([OH:3])>>[P:1](=[S:2])([O-:3])")

# Can't handle Phosphates in rings yet!
libgen1 = OELibraryGen("[O&!H1&!-,N&!-,S&!H1&!-:3][P&!r:2](=[S:1])([*:4])[*:5].[Cu:11].[Cu:12].[Cu:13]>>[Cu:11][A:3].[P:2]=[S:1].[Cu:12][*:4].[Cu:13][*:5]")

for mol in ifs.GetOEMols():
    deprot(mol)
    libgen1.AddStartingMaterial(mol, 0)
    libgen1.AddStartingMaterial(copper, 1)
    libgen1.AddStartingMaterial(copper, 2)
    libgen1.AddStartingMaterial(copper, 3)
    mol.Clear()
    copper.Clear()

title_list = []
for products in libgen1.GetProducts():
    if products.NumAtoms() < 76:
        num = 1
        title = products.GetTitle()+'_Z'+str(num)
        while title in title_list:
            num+=1
            title = products.GetTitle()+'_Z'+str(num)
        products.SetTitle(title)
 
        OEWriteMolecule(ofs, products)
#       print OECreateIsoSmiString(products)
ofs.close()
ifs.close()

# Remove artefacts
inf1 = open('b8.1'+date+'_ln_temp.ism', 'r')
outf1 = open('b8.1'+date+'_ln.ism', 'w')

for line in inf1.readlines():
    newline = string.replace(line, '[OH-]','[O-]')
    outf1.write(newline)

outf1.close()
inf1.close()
os.system('rm -f b8.1'+sys.argv[1][string.find(sys.argv[1], '_'):string.find(sys.argv[1], '_out')]+'_ln_temp.ism')
print '\nNOW:\nsource ~/.login\nthen invoke\nb8.2_thio_parts_split_2.py '+date+'\n'