#!/usr/local2/bin/python2.3
from openeye.oechem import *
import os, string

# wie 2, aber =S dort vergessen!

rxn = 14
ifs = oemolistream()
ifs.open('a_kegg_filter_1.sdf')
ofs_3 = oemolostream()
ofs_3.open('b_react_ti_'+str(rxn)+'_not_filtered_ln.sdf')
ofs_4 = oemolostream()
ofs_4.open('b_react_ti_'+str(rxn)+'_not_filtered_lp.sdf')

mol_hyd = OEGraphMol()
OEParseSmiles(mol_hyd, "[OH]")
mol_wat = OEGraphMol()
OEParseSmiles(mol_wat, "O(H)H")
# Reaction 1 aromatic: OH- nucleophilic attack at a-C(=O)-X (amide etc., O exocyclic)
# Former Carbonyl-oxygen becomes OH-group, attacking OH- becomes O-
# (AROMATIC -> partly wrong annotations in KEGG)
#leaving group NOT protonated [*:4] included, because leaving group needs to be recognized only once (if there are 2!)
libgen1 = OELibraryGen("[a&!-:4][c:1]=([S:3])[!#6:2].[O:100][H:102]>>[A:4][C:1]([S:3][H:102])([O-:100])[*:2]")
#leaving group protonated - will create attacked carboxylates!(can't avoid otherwise not all pot. leaving groups recognized)
libgen2 = OELibraryGen("[c:1]=([S:3])[!#6:2].[O:100]([H:101])[H:102]>>[C:1]([S:3][H:101])([O-:100])[*+:2][H:102]")

for mol in ifs.GetOEMols():
    libgen1.AddStartingMaterial(mol, 0)
    libgen1.AddStartingMaterial(mol_hyd, 1)
    libgen2.AddStartingMaterial(mol, 0)
    libgen2.AddStartingMaterial(mol_wat, 1)
    mol.Clear()
    mol_hyd.Clear()
    mol_wat.Clear()

title_list = []
rxn_n = 1
while rxn_n < rxn:
  try:
    old_title_file = open('b_react_ti_'+str(rxn_n)+'_ln.ism','r')
    for title in old_title_file.readlines():
        title_list.append(string.split(title)[1])
    old_title_file.close()
  except IOError:
    pass
  try:
    old_title_file = open('b_react_ti_'+str(rxn_n)+'_lp.ism','r')
    for title in old_title_file.readlines():
        title_list.append(string.split(title)[1])
    old_title_file.close()
  except IOError:
    pass
  rxn_n+=1

for products in libgen1.GetProducts():
    if products.NumAtoms() < 76:
        num = 1
        title = products.GetTitle()+str(num)
        while title in title_list:
            num+=1
            title = products.GetTitle()+str(num)
        title = products.GetTitle()+str(num)
        title_list.append(title)
        products.SetTitle(title)
        OEWriteMolecule(ofs_3, products)

hydcooh = OESubSearch()
hydcooh.Init('[O-][C&X4]([SH])[O-]') # weiss nicht mehr, wozu das gut war??

for products in libgen2.GetProducts():
    if hydcooh.SingleMatch(products) != 1 and products.NumAtoms() < 76:
        num = 1
        title = products.GetTitle()+str(num)
        while title in title_list:
            num+=1
            title = products.GetTitle()+str(num)
        title_list.append(title)
        products.SetTitle(title)
        OEWriteMolecule(ofs_4, products)

ofs_3.close()
ofs_4.close()
ifs.close()