Package rdkit :: Package Chem :: Module Descriptors
[hide private]
[frames] | no frames]

Source Code for Module rdkit.Chem.Descriptors

  1  # $Id$ 
  2  # 
  3  # Copyright (C) 2001-2010 greg Landrum and Rational Discovery LLC 
  4  # 
  5  #   @@ All Rights Reserved @@ 
  6  #  This file is part of the RDKit. 
  7  #  The contents are covered by the terms of the BSD license 
  8  #  which is included in the file license.txt, found at the root 
  9  #  of the RDKit source tree. 
 10  # 
 11  from rdkit import Chem 
 12  from rdkit.Chem import rdPartialCharges 
 13  import collections 
 14   
15 -def _isCallable(thing):
16 return (hasattr(collections,'Callable') and isinstance(thing,collections.Callable)) or \ 17 hasattr(thing,'__call__')
18 19 _descList=[]
20 -def _setupDescriptors(namespace):
21 global _descList,descList 22 from rdkit.Chem import GraphDescriptors,MolSurf,Lipinski,Fragments,Crippen 23 from rdkit.Chem.EState import EState_VSA 24 mods = [GraphDescriptors,MolSurf,EState_VSA,Lipinski,Crippen,Fragments] 25 26 otherMods = [Chem] 27 28 for nm,thing in namespace.items(): 29 if nm[0]!='_' and _isCallable(thing): 30 _descList.append((nm,thing)) 31 32 others = [] 33 for mod in otherMods: 34 tmp = dir(mod) 35 for name in tmp: 36 if name[0] != '_': 37 thing = getattr(mod,name) 38 if _isCallable(thing): 39 others.append(name) 40 41 for mod in mods: 42 tmp = dir(mod) 43 44 for name in tmp: 45 if name[0] != '_' and name[-1] != '_' and name not in others: 46 # filter out python reference implementations: 47 if name[:2]=='py' and name[2:] in tmp: 48 continue 49 thing = getattr(mod,name) 50 if _isCallable(thing): 51 namespace[name]=thing 52 _descList.append((name,thing)) 53 descList=_descList
54 55 from rdkit.Chem import rdMolDescriptors as _rdMolDescriptors 56 MolWt = lambda *x,**y:_rdMolDescriptors._CalcMolWt(*x,**y) 57 MolWt.version=_rdMolDescriptors._CalcMolWt_version 58 MolWt.__doc__="""The average molecular weight of the molecule 59 60 >>> MolWt(Chem.MolFromSmiles('CC')) 61 30.07 62 >>> MolWt(Chem.MolFromSmiles('[NH4+].[Cl-]')) 63 53.49... 64 65 """ 66 67 HeavyAtomMolWt=lambda x:MolWt(x,True) 68 HeavyAtomMolWt.__doc__="""The average molecular weight of the molecule ignoring hydrogens 69 70 >>> HeavyAtomMolWt(Chem.MolFromSmiles('CC')) 71 24.02... 72 >>> HeavyAtomMolWt(Chem.MolFromSmiles('[NH4+].[Cl-]')) 73 49.46 74 75 """ 76 HeavyAtomMolWt.version="1.0.0" 77 78 ExactMolWt = lambda *x,**y:_rdMolDescriptors.CalcExactMolWt(*x,**y) 79 ExactMolWt.version=_rdMolDescriptors._CalcExactMolWt_version 80 ExactMolWt.__doc__="""The exact molecular weight of the molecule 81 82 >>> ExactMolWt(Chem.MolFromSmiles('CC')) 83 30.04... 84 >>> ExactMolWt(Chem.MolFromSmiles('[13CH3]C')) 85 31.05... 86 87 """ 88 89
90 -def NumValenceElectrons(mol):
91 """ The number of valence electrons the molecule has 92 93 >>> NumValenceElectrons(Chem.MolFromSmiles('CC')) 94 14.0 95 >>> NumValenceElectrons(Chem.MolFromSmiles('C(=O)O')) 96 18.0 97 >>> NumValenceElectrons(Chem.MolFromSmiles('C(=O)[O-]')) 98 18.0 99 >>> NumValenceElectrons(Chem.MolFromSmiles('C(=O)')) 100 12.0 101 102 103 """ 104 tbl = Chem.GetPeriodicTable() 105 accum = 0.0 106 for atom in mol.GetAtoms(): 107 accum += tbl.GetNOuterElecs(atom.GetAtomicNum()) 108 accum -= atom.GetFormalCharge() 109 accum += atom.GetTotalNumHs() 110 111 return accum
112 NumValenceElectrons.version="1.0.0" 113
114 -def NumRadicalElectrons(mol):
115 """ The number of radical electrons the molecule has 116 (says nothing about spin state) 117 118 >>> NumRadicalElectrons(Chem.MolFromSmiles('CC')) 119 0.0 120 >>> NumRadicalElectrons(Chem.MolFromSmiles('C[CH3]')) 121 0.0 122 >>> NumRadicalElectrons(Chem.MolFromSmiles('C[CH2]')) 123 1.0 124 >>> NumRadicalElectrons(Chem.MolFromSmiles('C[CH]')) 125 2.0 126 >>> NumRadicalElectrons(Chem.MolFromSmiles('C[C]')) 127 3.0 128 129 130 """ 131 accum = 0.0 132 for atom in mol.GetAtoms(): 133 accum += atom.GetNumRadicalElectrons() 134 return accum
135 NumRadicalElectrons.version="1.0.0" 136
137 -def _ChargeDescriptors(mol,force=False):
138 if not force and hasattr(mol,'_chargeDescriptors'): 139 return mol._chargeDescriptors 140 chgs = rdPartialCharges.ComputeGasteigerCharges(mol) 141 minChg=500. 142 maxChg=-500. 143 for at in mol.GetAtoms(): 144 chg = float(at.GetProp('_GasteigerCharge')) 145 minChg = min(chg,minChg) 146 maxChg = max(chg,maxChg) 147 res = (minChg,maxChg) 148 mol._chargeDescriptors=res 149 return res
150 151
152 -def MaxPartialCharge(mol,force=False):
153 _,res = _ChargeDescriptors(mol,force) 154 return res
155 MaxPartialCharge.version="1.0.0" 156
157 -def MinPartialCharge(mol,force=False):
158 res,_ = _ChargeDescriptors(mol,force) 159 return res
160 MinPartialCharge.version="1.0.0" 161
162 -def MaxAbsPartialCharge(mol,force=False):
163 v1,v2 = _ChargeDescriptors(mol,force) 164 return max(abs(v1),abs(v2))
165 MaxAbsPartialCharge.version="1.0.0" 166
167 -def MinAbsPartialCharge(mol,force=False):
168 v1,v2 = _ChargeDescriptors(mol,force) 169 return min(abs(v1),abs(v2))
170 MinAbsPartialCharge.version="1.0.0" 171 172 from rdkit.Chem.EState.EState import MaxEStateIndex,MinEStateIndex,MaxAbsEStateIndex,MinAbsEStateIndex 173 174 _setupDescriptors(locals()) 175 176 177 178 #------------------------------------ 179 # 180 # doctest boilerplate 181 #
182 -def _test():
183 import doctest,sys 184 return doctest.testmod(sys.modules["__main__"],optionflags=doctest.ELLIPSIS)
185 186 if __name__ == '__main__': 187 import sys 188 failed,tried = _test() 189 sys.exit(failed) 190