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