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