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

Source Code for Module rdkit.Chem.Subshape.SubshapeBuilder

  1  # $Id$ 
  2  # 
  3  # Copyright (C) 2007 by Greg Landrum  
  4  #  All rights reserved 
  5  # 
  6  from rdkit import Chem,Geometry 
  7  from rdkit.Chem import AllChem 
  8  from rdkit.Chem.Subshape import SubshapeObjects 
  9  from rdkit.Chem.Subshape import BuilderUtils 
 10  from rdkit.six.moves import cPickle 
 11  import time 
 12   
 13  #----------------------------------------------------------------------------- 
14 -class SubshapeCombineOperations(object):
15 UNION=0 16 SUM=1 17 INTERSECT=2
18 19 #-----------------------------------------------------------------------------
20 -class SubshapeBuilder(object):
21 gridDims=(20,15,10) 22 gridSpacing=0.5 23 winRad=3.0 24 nbrCount=7 25 terminalPtRadScale=0.75 26 fraction=0.25 27 stepSize=1.0 28 featFactory=None 29
30 - def SampleSubshape(self,subshape1,newSpacing):
31 ogrid=subshape1.grid 32 rgrid = Geometry.UniformGrid3D(self.gridDims[0],self.gridDims[1],self.gridDims[2], 33 newSpacing) 34 for idx in range(rgrid.GetSize()): 35 l = rgrid.GetGridPointLoc(idx) 36 v = ogrid.GetValPoint(l) 37 rgrid.SetVal(idx,v) 38 39 res = SubshapeObjects.ShapeWithSkeleton() 40 res.grid = rgrid 41 return res;
42
43 - def GenerateSubshapeShape(self,cmpd,confId=-1,addSkeleton=True,**kwargs):
44 shape = SubshapeObjects.ShapeWithSkeleton() 45 shape.grid=Geometry.UniformGrid3D(self.gridDims[0],self.gridDims[1],self.gridDims[2], 46 self.gridSpacing) 47 AllChem.EncodeShape(cmpd,shape.grid,ignoreHs=False,confId=confId) 48 if addSkeleton: 49 conf = cmpd.GetConformer(confId) 50 self.GenerateSubshapeSkeleton(shape,conf,kwargs) 51 return shape
52 - def __call__(self,cmpd,**kwargs):
53 return self.GenerateSubshapeShape(cmpd,**kwargs)
54
55 - def GenerateSubshapeSkeleton(self,shape,conf=None,terminalPtsOnly=False,skelFromConf=True):
56 if conf and skelFromConf: 57 pts = BuilderUtils.FindTerminalPtsFromConformer(conf,self.winRad,self.nbrCount) 58 else: 59 pts = BuilderUtils.FindTerminalPtsFromShape(shape,self.winRad,self.fraction) 60 pts = BuilderUtils.ClusterTerminalPts(pts,self.winRad,self.terminalPtRadScale) 61 BuilderUtils.ExpandTerminalPts(shape,pts,self.winRad) 62 if len(pts)<3: 63 raise ValueError('only found %d terminals, need at least 3'%len(pts)) 64 65 if not terminalPtsOnly: 66 pts = BuilderUtils.AppendSkeletonPoints(shape.grid,pts,self.winRad,self.stepSize) 67 for i,pt in enumerate(pts): 68 BuilderUtils.CalculateDirectionsAtPoint(pt,shape.grid,self.winRad) 69 if conf and self.featFactory: 70 BuilderUtils.AssignMolFeatsToPoints(pts,conf.GetOwningMol(),self.featFactory,self.winRad) 71 shape.skelPts=pts
72
73 - def CombineSubshapes(self,subshape1,subshape2,operation=SubshapeCombineOperations.UNION):
74 import copy 75 cs = copy.deepcopy(subshape1) 76 if operation==SubshapeCombineOperations.UNION: 77 cs.grid |= subshape2.grid 78 elif operation==SubshapeCombineOperations.SUM: 79 cs.grid += subshape2.grid 80 elif operation==SubshapeCombineOperations.INTERSECT: 81 cs.grid &= subshape2.grid 82 else: 83 raise ValueError('bad combination operation') 84 return cs
85 86 87 if __name__=='__main__': 88 from rdkit.Chem import AllChem,ChemicalFeatures 89 from rdkit.Chem.PyMol import MolViewer 90 #cmpd = Chem.MolFromSmiles('CCCc1cc(C(=O)O)ccc1') 91 #cmpd = Chem.AddHs(cmpd) 92 if 1: 93 cmpd = Chem.MolFromSmiles('C1=CC=C1C#CC1=CC=C1') 94 cmpd = Chem.AddHs(cmpd) 95 AllChem.EmbedMolecule(cmpd) 96 AllChem.UFFOptimizeMolecule(cmpd) 97 AllChem.CanonicalizeMol(cmpd) 98 print >>file('testmol.mol','w+'),Chem.MolToMolBlock(cmpd) 99 else: 100 cmpd = Chem.MolFromMolFile('testmol.mol') 101 builder=SubshapeBuilder() 102 if 1: 103 shape=builder.GenerateSubshapeShape(cmpd) 104 v = MolViewer() 105 if 1: 106 import tempfile 107 tmpFile = tempfile.mktemp('.grd') 108 v.server.deleteAll() 109 Geometry.WriteGridToFile(shape.grid,tmpFile) 110 time.sleep(1) 111 v.ShowMol(cmpd,name='testMol',showOnly=True) 112 v.server.loadSurface(tmpFile,'testGrid','',2.5) 113 v.server.resetCGO('*') 114 115 cPickle.dump(shape,file('subshape.pkl','w+')) 116 for i,pt in enumerate(shape.skelPts): 117 v.server.sphere(tuple(pt.location),.5,(1,0,1),'Pt-%d'%i) 118 if not hasattr(pt,'shapeDirs'): continue 119 momBeg = pt.location-pt.shapeDirs[0] 120 momEnd = pt.location+pt.shapeDirs[0] 121 v.server.cylinder(tuple(momBeg),tuple(momEnd),.1,(1,0,1),'v-%d'%i) 122