1
#!/usr/bin/env python
2 60
"""
3
proc_tags.py
4

5
This script parses output of Psi4 calculations and writes data in SD tags.
6

7
Usage:
8
- import proc_Tags as pt
9
- pt.set_sd_tags(args)
10

11
By: Victoria T. Lim
12

13
"""
14

15 60
import openeye.oechem as oechem
16

17

18 60
def define_tag(datum, package, method, basisset):
19
    # CASE SENSITIVE
20
    # this function is case-sensitive though get_sd_list should remove sensitivity
21

22 60
    tagdict = {
23
        "QM opt energy": "QM {} Final Opt. Energy (Har) {}/{}".format(
24
            package, method, basisset),
25
        "QM opt energy scs": "QM {} Final Opt. Energy (Har) SCS-{}/{}".format(
26
            package, method, basisset),
27
        "QM opt energy initial": "QM {} Initial Opt. Energy (Har) {}/{}".format(
28
            package, method, basisset),
29
        "QM spe": "QM {} Final Single Pt. Energy (Har) {}/{}".format(
30
            package, method, basisset),
31
        "QM spe scs": "QM {} Final Single Pt. Energy (Har) SCS-{}/{}".format(
32
            package, method, basisset),
33
        "QM spe scs old": "QM {} Final SCS-Single Pt. Energy (Har) {}/{}".format( ## old version of sd tag
34
            package, method, basisset),
35
        "MM opt energy": "MM Szybki SD Energy",
36
        "original index": "Original omega conformer number",
37
        "opt runtime": "QM {} Opt. Runtime (sec) {}/{}".format(
38
            package, method, basisset),
39
        "spe runtime": "QM {} Single Pt. Runtime (sec) {}/{}".format(
40
            package, method, basisset),
41
        "opt step": "QM {} Opt. Steps {}/{}".format(package, method, basisset),
42
    }
43 60
    try:
44 60
        taglabel = tagdict[datum]
45 60
    except KeyError as exc:
46 60
        raise NameError("Error in defining string to extract SD data for "
47
            f"[ {datum} ]. Please verify that the key for define_tag is one of: ",
48
            tagdict.keys()) from exc
49

50 60
    return taglabel
51

52 60
def get_sd_list(mol, datum, package='Psi4', method=None, basisset=None, taglabel=None):
53
    """
54
    Get list of specified SD tag for all confs in mol.
55

56
    Parameters
57
    ----------
58
    mol:        OEChem molecule with all of its conformers
59
    datum:       string description of property of interest
60
        options implemented: "QM opt energy" "MM opt energy"
61
    package:    software package used for QM calculation. Psi4 or Turbomole.
62
    method:     string, for specific properties. e.g. 'mp2'
63
    basisset:   string, for specific properties. e.g. '6-31+G(d)'
64
    taglabel : string
65
        exact tag string from which to extract SD data
66

67
    Returns
68
    -------
69
    sdlist: A 1D N-length list for N conformers with property from SDTag.
70
    """
71

72 60
    if taglabel is None:
73 60
        taglabel = define_tag(datum, package, method, basisset)
74

75 60
    sd_list = []
76 60
    for j, conf in enumerate(mol.GetConfs()):
77 60
        for x in oechem.OEGetSDDataPairs(conf):
78

79
            # Case: opt did not finish --> append nan
80 60
            if "note on opt." in x.GetTag().lower(
81
            ) and "did not finish" in x.GetValue().lower():
82 0
                sd_list.append('nan')
83 0
                break
84

85
            # Case: want energy value OR want original index number
86 60
            elif taglabel.lower() in x.GetTag().lower():
87 60
                sd_list.append(x.GetValue())
88 60
                break
89

90 60
    return sd_list
91

92

93 60
def set_sd_tags(Conf, Props, calctype):
94
    """
95
    For one particular conformer, set all available SD tags based on data
96
    in Props dictionary.
97

98
    Warning
99
    -------
100
    If the exact tag already exists, and you want to add a new one then there
101
    will be duplicate tags with maybe different data. (NOT recommended).
102
    Then the function to get sd_list will only get one or the other;
103
    I think it just gets the first matching tag.
104

105
    TODO: maybe add some kind of checking to prevent duplicate tags added
106

107
    Parameters
108
    ----------
109
    Conf:       Single conformer from OEChem molecule
110
    Props:      Dictionary output from ProcessOutput function.
111
                Should contain the keys: basis, method, numSteps,
112
                initEnergy, finalEnergy, coords, time, pkg
113
    calctype: string; one of 'opt','spe','hess' for geometry optimization,
114
        single point energy calculation, or Hessian calculation
115

116
    """
117

118
    # get level of theory for setting SD tags
119 60
    method = Props['method']
120 60
    basisset = Props['basis']
121 60
    pkg = Props['package']
122

123
    # turn parameters into tag descriptions
124 60
    full_method = "{}/{}".format(method, basisset)
125 60
    cdict = {'spe': 'Single Pt.', 'opt': 'Opt.', 'hess': 'Hessian'}
126

127
    # time info can be set for all cases
128 60
    taglabel = "QM {} {} Runtime (sec) {}".format(pkg, cdict[calctype],
129
                                                  full_method)
130 60
    oechem.OEAddSDData(Conf, taglabel, str(Props['time']))
131

132
    # hessian has no other info for sd tag
133 60
    if calctype == 'hess':
134 60
        return
135

136
    # check that finalEnergy is there. if not, opt probably did not finish
137
    # make a note of that in SD tag then quit function
138 60
    if not 'finalEnergy' in Props:
139 60
        taglabel = "Note on {} {}".format(cdict[calctype], full_method)
140 60
        oechem.OEAddSDData(Conf, taglabel, "JOB DID NOT FINISH")
141 60
        return
142

143
    # Set new SD tag for conformer's final energy
144 60
    taglabel = "QM {} Final {} Energy (Har) {}".format(pkg, cdict[calctype],
145
                                                       full_method)
146 60
    oechem.OEAddSDData(Conf, taglabel, str(Props['finalEnergy']))
147

148
    # Set new SD tag for final SCS-MP2 energy if method is MP2
149 60
    if method.lower() == 'mp2':
150 60
        taglabel = "QM {} Final {} Energy (Har) SCS-{}".format(
151
            pkg, cdict[calctype], full_method)
152 60
        oechem.OEAddSDData(Conf, taglabel, str(Props['finalSCSEnergy']))
153

154
    # Add COSMO energy with outlying charge correction. Turbomole only!
155 60
    if 'ocEnergy' in Props:
156 0
        if calctype == 'spe':
157 0
            print(
158
                "Extraction of COSMO OC energy from Turbomole not yet supported for SPE calcns"
159
            )
160 0
        elif calctype == 'opt':
161 0
            taglabel = "QM {} Final {} Energy with OC correction (Har) {}".format(
162
                pkg, cdict[calctype], full_method)
163 0
            oechem.OEAddSDData(Conf, taglabel, str(Props['ocEnergy']))
164

165
    # spe has no other relevant info for sd tag
166 60
    if calctype == 'spe':
167 60
        return
168

169
    # Set new SD tag for original conformer number if not existing
170
    # !! Opt2 files should ALREADY have this !! Opt2 index is NOT orig index !!
171 60
    taglabel = "Original omega conformer number"
172 60
    if not oechem.OEHasSDData(Conf, taglabel):
173
        # if not working with confs, will have no GetIdx
174 60
        try:
175 60
            oechem.OEAddSDData(Conf, taglabel, str(Conf.GetIdx() + 1))
176 60
        except AttributeError as err:
177 60
            pass
178
    # if tag exists, append new conformer ID after the old one
179
    else:
180 0
        try:
181 0
            oldid = oechem.OEGetSDData(Conf, taglabel)
182 0
            newid = str(Conf.GetIdx() + 1)
183 0
            totid = "{}, {}".format(oldid, newid)
184 0
            oechem.OESetSDData(Conf, taglabel, totid)
185 0
        except AttributeError as err:
186 0
            pass
187

188
    # Set new SD tag for numSteps of geom. opt.
189 60
    taglabel = "QM {} {} Steps {}".format(pkg, cdict[calctype], full_method)
190 60
    oechem.OEAddSDData(Conf, taglabel, str(Props['numSteps']))
191

192
    # Set new SD tag for conformer's initial energy
193 60
    taglabel = "QM {} Initial {} Energy (Har) {}".format(
194
        pkg, cdict[calctype], full_method)
195 60
    oechem.OEAddSDData(Conf, taglabel, str(Props['initEnergy']))
196

197

198 60
def delete_tag(mol, tag):
199
    """
200
    Delete specified SD tag from all conformers of mol.
201

202
    Note: Multi-conformer molecule must be specified
203
    else will get AttributeError:
204
    'OEGraphMol' object has no attribute 'GetConfs'.
205

206
    Parameters
207
    ----------
208
    mol : multi-conformer OEChem molecule
209
    tag : string
210
        exact label of the data to delete
211

212
    """
213 60
    for j, conf in enumerate(mol.GetConfs()):
214 60
        oechem.OEDeleteSDData(conf, tag)

Read our documentation on viewing source code .

Loading