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