1import sys
2import os
3import argparse
4
5import CDPL.Chem as Chem
6import CDPL.ConfGen as ConfGen
7
8
9# generates a conformer ensemble of the argument molecule using
10# the provided initialized ConfGen.ConformerGenerator instance
11def genConfEnsemble(mol: Chem.Molecule, conf_gen: ConfGen.ConformerGenerator) -> (int, int):
12 # prepare the molecule for conformer generation
13 ConfGen.prepareForConformerGeneration(mol)
14
15 # generate the conformer ensemble
16 status = conf_gen.generate(mol)
17 num_confs = conf_gen.getNumConformers()
18
19 # if sucessful, store the generated conformer ensemble as
20 # per atom 3D coordinates arrays (= the way conformers are represented in CDPKit)
21 if status == ConfGen.ReturnCode.SUCCESS or status == ConfGen.ReturnCode.TOO_MUCH_SYMMETRY:
22 conf_gen.setConformers(mol)
23 else:
24 num_confs = 0
25
26 # return status code and the number of generated conformers
27 return (status, num_confs)
28
29def main() -> None:
30 args = parseArgs()
31
32 # create reader for input molecules (format specified by file extension)
33 reader = Chem.MoleculeReader(args.in_file)
34
35 # create writer for the generated conformer ensembles (format specified by file extension)
36 writer = Chem.MolecularGraphWriter(args.out_file)
37
38 # create and initialize an instance of the class ConfGen.ConformerGenerator which
39 # will perform the actual conformer ensemble generation work
40 conf_gen = ConfGen.ConformerGenerator()
41
42 conf_gen.settings.timeout = args.max_time * 1000 # apply the -t argument
43 conf_gen.settings.minRMSD = args.min_rmsd # apply the -r argument
44 conf_gen.settings.energyWindow = args.e_window # apply the -e argument
45 conf_gen.settings.maxNumOutputConformers = args.max_confs # apply the -n argument
46
47 # dictionary mapping status codes to human readable strings
48 status_to_str = { ConfGen.ReturnCode.UNINITIALIZED : 'uninitialized',
49 ConfGen.ReturnCode.TIMEOUT : 'max. processing time exceeded',
50 ConfGen.ReturnCode.ABORTED : 'aborted',
51 ConfGen.ReturnCode.FORCEFIELD_SETUP_FAILED : 'force field setup failed',
52 ConfGen.ReturnCode.FORCEFIELD_MINIMIZATION_FAILED : 'force field structure refinement failed',
53 ConfGen.ReturnCode.FRAGMENT_LIBRARY_NOT_SET : 'fragment library not available',
54 ConfGen.ReturnCode.FRAGMENT_CONF_GEN_FAILED : 'fragment conformer generation failed',
55 ConfGen.ReturnCode.FRAGMENT_CONF_GEN_TIMEOUT : 'fragment conformer generation timeout',
56 ConfGen.ReturnCode.FRAGMENT_ALREADY_PROCESSED : 'fragment already processed',
57 ConfGen.ReturnCode.TORSION_DRIVING_FAILED : 'torsion driving failed',
58 ConfGen.ReturnCode.CONF_GEN_FAILED : 'conformer generation failed' }
59
60 # create an instance of the default implementation of the Chem.Molecule interface
61 mol = Chem.BasicMolecule()
62 i = 1
63
64 # read and process molecules one after the other until the end of input has been reached
65 try:
66 while reader.read(mol):
67 # compose a simple molecule identifier
68 mol_id = Chem.getName(mol).strip()
69
70 if mol_id == '':
71 mol_id = '#' + str(i) # fallback if name is empty
72 else:
73 mol_id = '\'%s\' (#%s)' % (mol_id, str(i))
74
75 if not args.quiet:
76 print('- Generating conformers for molecule %s...' % mol_id)
77
78 try:
79 # generate conformer ensemble for read molecule
80 status, num_confs = genConfEnsemble(mol, conf_gen)
81
82 # check for severe error reported by status code
83 if status != ConfGen.ReturnCode.SUCCESS and status != ConfGen.ReturnCode.TOO_MUCH_SYMMETRY:
84 if args.quiet:
85 print('Error: conformer ensemble generation for molecule %s failed: %s' % (mol_id, status_to_str[status]))
86 else:
87 print(' -> Conformer ensemble generation failed: %s' % status_to_str[status])
88
89 elif not args.quiet: # arrives here only if no severe error occurred
90 if status == ConfGen.ReturnCode.TOO_MUCH_SYMMETRY:
91 print(' -> Generated %s conformers (warning: too much top. symmetry - output ensemble may contain duplicates)' % str(num_confs))
92 else:
93 print(' -> Generated %s conformers' % str(num_confs))
94
95 # output generated ensemble (if available)
96 if num_confs > 0:
97 if not writer.write(mol):
98 sys.exit('Error: output of conformer ensemble for molecule %s failed' % mol_id)
99
100 except Exception as e:
101 sys.exit('Error: conformer ensemble generation or output for molecule %s failed: %s' % (mol_id, str(e)))
102
103 i += 1
104
105 except Exception as e: # handle exception raised in case of severe read errors
106 sys.exit('Error: reading molecule failed: ' + str(e))
107
108 writer.close()
109 sys.exit(0)
110
111def parseArgs() -> argparse.Namespace:
112 parser = argparse.ArgumentParser(description='Generates conformer ensembles for the given input molecules.')
113
114 parser.add_argument('-i',
115 dest='in_file',
116 required=True,
117 metavar='<file>',
118 help='Molecule input file')
119 parser.add_argument('-o',
120 dest='out_file',
121 required=True,
122 metavar='<file>',
123 help='Conformer ensemble output file')
124 parser.add_argument('-e',
125 dest='e_window',
126 required=False,
127 metavar='<float>',
128 type=float,
129 default=20.0,
130 help='Output conformer energy window (default: 20.0)')
131 parser.add_argument('-r',
132 dest='min_rmsd',
133 required=False,
134 metavar='<float>',
135 type=float,
136 default=0.5,
137 help='Output conformer RMSD threshold (default: 0.5)')
138 parser.add_argument('-t',
139 dest='max_time',
140 required=False,
141 metavar='<int>',
142 type=int,
143 default=3600,
144 help='Max. allowed molecule processing time (default: 3600 sec)')
145 parser.add_argument('-n',
146 dest='max_confs',
147 required=False,
148 metavar='<int>',
149 type=int,
150 help='Max. output ensemble size (default: 100)')
151 parser.add_argument('-q',
152 dest='quiet',
153 required=False,
154 action='store_true',
155 default=False,
156 help='Disable progress output (default: false)')
157
158 return parser.parse_args()
159
160if __name__ == '__main__':
161 main()