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

Source Code for Module rdkit.Chem.Pharm3D.Pharmacophore

  1  # $Id$ 
  2  # 
  3  # Copyright (C) 2004-2008 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 Geometry 
 12  import numpy 
 13  from rdkit import Chem 
 14  from rdkit.Chem import ChemicalFeatures 
 15  from rdkit.RDLogger import logger 
 16  logger = logger() 
 17   
18 -class Pharmacophore:
19 - def __init__(self, feats,initMats=True):
20 self._initializeFeats(feats) 21 nf = len(feats) 22 self._boundsMat = numpy.zeros((nf, nf), numpy.float) 23 self._boundsMat2D = numpy.zeros((nf, nf), numpy.int) 24 if initMats: 25 self._initializeMatrices()
26
27 - def _initializeFeats(self,feats):
28 self._feats = [] 29 for feat in feats: 30 if isinstance(feat,ChemicalFeatures.MolChemicalFeature): 31 pos = feat.GetPos() 32 newFeat = ChemicalFeatures.FreeChemicalFeature(feat.GetFamily(),feat.GetType(), 33 Geometry.Point3D(pos[0],pos[1],pos[2])) 34 self._feats.append(newFeat) 35 else: 36 self._feats.append(feat)
37
38 - def _initializeMatrices(self):
39 # initialize the bounds matrix with distances to start with 40 nf = len(self._feats) 41 for i in range(1, nf): 42 loci = self._feats[i].GetPos() 43 for j in range(i): 44 locj = self._feats[j].GetPos() 45 dist = loci.Distance(locj) 46 self._boundsMat[i,j] = dist 47 self._boundsMat[j,i] = dist 48 for i in range(nf): 49 for j in range(i+1,nf): 50 self._boundsMat2D[i,j] = 1000
51 52
53 - def getFeatures(self):
54 return self._feats
55
56 - def getFeature(self, i):
57 return self._feats[i]
58
59 - def getUpperBound(self, i, j):
60 if (i > j): 61 j,i = i,j 62 return self._boundsMat[i, j]
63
64 - def getLowerBound(self, i, j):
65 if (j > i): 66 j,i = i,j 67 return self._boundsMat[i,j]
68
69 - def _checkBounds(self,i,j):
70 " raises ValueError on failure " 71 nf = len(self._feats) 72 if (i < 0) or (i >= nf): 73 raise ValueError("Index out of bound") 74 if (j < 0) or (j >= nf): 75 raise ValueError("Index out of bound") 76 return True
77
78 - def setUpperBound(self, i, j, val, checkBounds=False):
79 if (checkBounds): self._checkBounds(i,j) 80 if (i > j): 81 j,i = i,j 82 self._boundsMat[i,j] = val
83
84 - def setLowerBound(self, i, j, val, checkBounds=False):
85 if (checkBounds): self._checkBounds(i,j) 86 if (j > i): 87 j,i = i,j 88 self._boundsMat[i,j] = val
89
90 - def getUpperBound2D(self, i, j):
91 if (i > j): 92 j,i = i,j 93 return self._boundsMat2D[i, j]
94
95 - def getLowerBound2D(self, i, j):
96 if (j > i): 97 j,i = i,j 98 return self._boundsMat2D[i,j]
99 100
101 - def setUpperBound2D(self, i, j, val, checkBounds=False):
102 if (checkBounds): self._checkBounds(i,j) 103 if (i > j): 104 j,i = i,j 105 self._boundsMat2D[i,j] = val
106
107 - def setLowerBound2D(self, i, j, val, checkBounds=False):
108 if (checkBounds): self._checkBounds(i,j) 109 if (j > i): 110 j,i = i,j 111 self._boundsMat2D[i,j] = val
112
113 - def __str__(self):
114 res = ' '*13 115 for i,iFeat in enumerate(self._feats): 116 res += '% 12s '%iFeat.GetFamily() 117 res += '\n' 118 for i,iFeat in enumerate(self._feats): 119 res += '% 12s '%iFeat.GetFamily() 120 for j,jFeat in enumerate(self._feats): 121 if j<i: 122 res += '% 12.3f '%self.getLowerBound(i,j) 123 elif j>i: 124 res += '% 12.3f '%self.getUpperBound(i,j) 125 else: 126 res += '% 12.3f '%0.0 127 res += '\n' 128 return res
129
130 -class ExplicitPharmacophore:
131 """ this is a pharmacophore with explicit point locations and radii 132 """
133 - def __init__(self, feats=None,radii=None):
134 if feats and radii: 135 self._initializeFeats(feats,radii)
136
137 - def _initializeFeats(self,feats,radii):
138 if len(feats)!=len(radii): 139 raise ValueError('len(feats)!=len(radii)') 140 self._feats = [] 141 self._radii = [] 142 for feat,rad in zip(feats,radii): 143 if isinstance(feat,ChemicalFeatures.MolChemicalFeature): 144 pos = feat.GetPos() 145 newFeat = ChemicalFeatures.FreeChemicalFeature(feat.GetFamily(),feat.GetType(), 146 Geometry.Point3D(pos[0],pos[1],pos[2])) 147 else: 148 newFeat = feat 149 self._feats.append(newFeat) 150 self._radii.append(rad)
151
152 - def getFeatures(self):
153 return self._feats
154
155 - def getRadii(self):
156 return self._radii
157
158 - def getFeature(self, i):
159 return self._feats[i]
160
161 - def getRadius(self, i):
162 return self._radii[i]
163
164 - def setRadius(self,i,rad):
165 self._radii[i]=rad
166
167 - def initFromString(self,text):
168 lines = text.split(r'\n') 169 self.initFromLines(lines)
170
171 - def initFromFile(self,inF):
172 self.initFromLines(inF.readlines())
173
174 - def initFromLines(self,lines):
175 from rdkit.Chem import ChemicalFeatures 176 177 import re 178 spaces = re.compile('[\ \t]+') 179 180 feats = [] 181 rads = [] 182 for lineNum,line in enumerate(lines): 183 txt = line.split('#')[0].strip() 184 if txt: 185 splitL = spaces.split(txt) 186 if len(splitL)<5: 187 logger.error('Input line %d only contains %d fields, 5 are required. Read failed.'%(lineNum,len(splitL))) 188 return 189 fName = splitL[0] 190 try: 191 xP = float(splitL[1]) 192 yP = float(splitL[2]) 193 zP = float(splitL[3]) 194 rad = float(splitL[4]) 195 except ValueError: 196 logger.error('Error parsing a number of line %d. Read failed.'%(lineNum)) 197 return 198 feats.append(ChemicalFeatures.FreeChemicalFeature(fName,fName,Geometry.Point3D(xP,yP,zP))) 199 rads.append(rad) 200 self._initializeFeats(feats,rads)
201
202 - def __str__(self):
203 res = '' 204 for feat,rad in zip(self._feats,self._radii): 205 res += '% 12s '%feat.GetFamily() 206 p = feat.GetPos() 207 res += ' % 8.4f % 8.4f % 8.4f '%(p.x,p.y,p.z) 208 res += '% 5.2f'%rad 209 res += '\n' 210 return res
211