5.4 Example Program

The following code example is a simple example of how to use the oeproton library provided in QuacPac to assign partial charges to an OEGraphMol. The program opens the file foo.sdf and reads all of the molecules in the file. Each molecule has its partial charge calculated with the AM1BCC method (including AM1 optimization). Each hydrogen is then made to be implicit, and it's partial charge is added to the partial charge on it's parent heavy-atom. Finally, each molecule is written to the file foo.mol2.

 1 #!/usr/bin/env python
 2 #############################################################################
 3 #  Copyright (C) 2006, 2008,  OpenEye Scientific Software, Inc.
 4 #############################################################################
 5
 6 import os, sys
 7 from openeye.oechem import *
 8 from openeye.oequacpac import *
 9
10 def main(argv = [__name__]):
11     if len(argv) != 3:
12         OEThrow.Usage("%s <molfile> <outfile>" % argv[0])
13
14     noHydrogen = True
15     debug = False
16     mol = OEGraphMol()
17     ifs = oemolistream(argv[1])
18     ofs = oemolostream(argv[2])
19     while OEReadMolecule(ifs,mol):
20         OEAssignPartialCharges(mol,OECharges_MMFF94,noHydrogen,debug);
21         OEWriteMolecule(ofs,mol)
22
23     return 0
24
25 if __name__ == "__main__":
26     sys.exit(main(sys.argv))

Listing:5.1 Calculating Partial Charges

This second example shows how to enumerate pKa states of a molecule. Molecules are read from the input file foo.smi. The OETyperMolFunction is created so that a default output stream (std::out, SMILES format) is used, aromaticity will be called, compounds are enumerated rather than only being counted, and a maximum of 200 states per molecule will be generated. The state counter is reset for each new molecule with the OETyperMolFunction::Reset function.

 1 #!/usr/bin/env python
 2 #############################################################################
 3 #  Copyright (C) 2006, 2008,  OpenEye Scientific Software, Inc.
 4 #############################################################################
 5
 6 import os, sys
 7 from openeye.oechem import *
 8 from openeye.oequacpac import *
 9
10 def main(argv = [__name__]):
11     if len(argv) != 3:
12         OEThrow.Usage("%s <molfile> <outfile>" % argv[0])
13     ifs = oemolistream(argv[1])
14     ofs = oemolostream(argv[2])
15
16     verbose = False
17     aromatic = True
18     countOnly = False
19     maxCount = 100
20
21     tmf = OETyperMolFunction(ofs, aromatic, countOnly, maxCount)
22     mol = OEGraphMol()
23     while OEReadMolecule(ifs,mol):
24         OEEnumerateFormalCharges(mol, tmf, verbose)
25         tmf.Reset()
26
27     return 0
28
29 if __name__ == "__main__":
30     sys.exit(main(sys.argv))

Listing:5.2 Enumerating Ionization States

This third example listing expands on the previous example. Here, rather than enumerate the pKa states to std::out, they are enumerated into a stringstream. Each enumerated pKa state is then passed into tautomer enumeration. In the tautomer enumeration phase, the number of states are only counted rather than being enumerated. In both phases, the enumeration is capped at 200 states per molecule. The states of both the ionization state counter and the tautomer counter are reinitialized with their Reset functions.

 1 #!/usr/bin/env python
 2 #############################################################################
 3 #  Copyright (C) 2006, 2008,  OpenEye Scientific Software, Inc.
 4 #############################################################################
 5
 6 import os, sys
 7 from openeye.oechem import *
 8 from openeye.oequacpac import *
 9
10 def main(argv = [__name__]):
11     if len(argv) != 2:
12         OEThrow.Usage("%s <molfile>" % argv[0])
13     ifs = oemolistream(argv[1])
14
15     nulfs = oemolostream(oenul, False)
16     ofs = oemolostream()
17
18     aromatic = True
19     tyCountOnly = False
20     maxCount = 200
21     tymf = OETyperMolFunction(ofs, aromatic, tyCountOnly, maxCount)
22
23     verbose = False
24
25     taCountOnly = True
26     tamf = OETautomerMolFunction(nulfs, aromatic, taCountOnly, maxCount)
27
28
29     mol = OEGraphMol()
30     pkaState = OEGraphMol()
31
32     while OEReadMolecule(ifs,mol):
33         # enumerate pka states
34         ofs.openstring()
35         OEEnumerateFormalCharges(mol, tymf, verbose)
36         tymf.Reset()
37
38         iss = oemolistream()
39         iss.openstring(ofs.GetString())
40
41         # count tautomers of pka states
42         count = 0
43         while OEReadMolecule(iss, pkaState):
44             count += OEEnumerateTautomers(pkaState, tamf)
45             tamf.Reset()
46
47         print count,"tautomer/pka states for",mol.GetTitle()
48
49     return 0
50
51 if __name__ == "__main__":
52     sys.exit(main(sys.argv))

Listing:5.3 Counting Tautomers

				QuacPac Quality Atomic Charges, Proton Assignment and Canonicalization