Python molecule.AtomSet類代碼示例

    def createPairs(self, refnames, mobnames, slice, choice):
        # This is where the atom pairs are created
        # Get the ref atom set
        refAtms = AtomSet([])
        mobAtms = AtomSet([])
        for name in refnames:
            refnod = self.refNodes[name]
            atms = refnod.findType(Atom)
            refAtms = refAtms + atms

        for name in mobnames:
            mobnod = self.mobNodes[name]
            atms = mobnod.findType(Atom)
            mobAtms = mobAtms + atms
        # Now apply the filters to the sets
        if choice:
            refFiltAtms = refAtms.get(self.filters[choice])
            mobFiltAtms = mobAtms.get(self.filters[choice])
        if not len(refFiltAtms) == len(mobFiltAtms):
            #print 'slicing', slice
            self.warningMsg("the two sets of atoms needs to be of the same length")
            return
        ebn = self.cmdForms['pairs'].descr.entryByName
        lc = ebn['newpairs']['widget']
        for refatm, mobatm in map(None, refFiltAtms, mobFiltAtms):
            #self.pairs.append((refatm, mobatm))
            pairName = refatm.full_name() + '---' + mobatm.full_name()
            self.newPairs[pairName]=(refatm,mobatm)
            lc.add((pairName, None))

 def filterAcceptors(self, accAts):
     ntypes = ['Npl', 'Nam']
     npls = accAts.get(lambda x, ntypes=ntypes: x.babel_type=='Npl')
     nams = accAts.get(lambda x, ntypes=ntypes: x.babel_type=='Nam')
     #nAts = accAts.get(lambda x, ntypes=ntypes: x.babel_type in ntypes)
     restAts = accAts.get(lambda x, ntypes=ntypes: x.babel_type not in ntypes)
     if not restAts: restAts = AtomSet([])
     #if nAts:
     if npls:
         #for at in nAts:
         for at in npls:
             s = 0
             for b in at.bonds:
                 if b.bondOrder=='aromatic':
                     s = s + 2
                 else: s = s + b.bondOrder
             #if s<3:
             #apparently this is wrong
             if s<4:
                 restAts.append(at)
     if nams:
         #for at in nAts:
         for at in nams:
             s = 0
             for b in at.bonds:
                 if b.bondOrder=='aromatic':
                     s = s + 2
                 else: s = s + b.bondOrder
                 #s = s + b.bondOrder
             if s<3:
                 restAts.append(at)
     return restAts

def get_atoms(mol, list_of_indicies, names_to_use=['N','CA','C'], verbose=False):
    if verbose: 
        print "in get_atoms with list of indicies:" 
        print list_of_indicies
    if not len(list_of_indicies):
        raise 'invalid input: list of indicies is empty!'
    atoms = AtomSet()
    num_res = 0
    for item in list_of_indicies:
        first, second = item
        #check for valid index and for end of chain
        max_index = len(mol.chains.residues)-1
        assert first<=max_index, 'invalid start of residue range'
        assert second<=max_index, 'invalid end of residue range'
        assert second>=first, 'second index cannot be smaller than first'
        if second==max_index:
            #ie mol.chains.residues[second]==mol.chains.residues[-1]:
            these_res = mol.chains.residues[first:]
        else:
            these_res = mol.chains.residues[first:second+1]
        if verbose: print "Adding  %3d residues " %(len(these_res)),
        num_res+=len(these_res)
        if verbose: print "Now there are %d residues total" %(num_res)
        for r in these_res:
            for n in names_to_use:
                atoms.append( r.atoms.get(n)[0])
    assert len(atoms), 'invalid input: lists of indicies did not correspond to any residues!'
    if verbose: print 'returning %d atoms' %(len(atoms))
    return atoms

 def getAtoms(self,bnds):
     ats0 = AtomSet()
     for b in bnds:
         ats0.append(b.atom1)
         ats0.append(b.atom2)
     d = {}
     for a in ats0:
         d[a] = 0
     return AtomSet(d.keys())

 def doit(self, bonds):
     global var
     var=1 
     ats = AtomSet([])
     for bond in bonds:
         ats.append(bond.atom1)
         ats.append(bond.atom2)
         self.vf.removeBonds(bond.atom1, bond.atom2)
     var=0
     self.vf.GUI.VIEWER.Redraw()

 def checkForPossibleH(self, ats, blen):
     #@@FIX THIS: WHAT IS THE POINT OF THIS???
     #check that if at has all bonds, at least one is to a hydrogen
     # have to do this by element??
     probAts = AtomSet(ats.get(lambda x, blen=blen: len(x.bonds)==blen))
     #probOAts = ats.get(lambda x, blen=blen: len(x.bonds)==blen)
     #probSAts = ats.get(lambda x, blen=blen: len(x.bonds)==blen)
     if probAts:
         rAts = AtomSet([])
         for at in probAts:
             if not len(at.findHydrogens()):
                 rAts.append(at)
         if len(rAts):
             ats =  ats.subtract(rAts)
     return ats

 def setupUndoAfter(self, ats, angle,**kw):
     #no atoms, <4 atoms, 
     aSet = AtomSet(self.atomList)
     self.undoMenuString = self.name
     if len(self.atomList)==0:
         undoCmd = 'self.setTorsionGC.atomList=[]; self.setTorsionGC.update()'
     elif len(self.atomList)<4:
         #need to step back here
         undoCmd = 'self.setTorsionGC.atomList=self.setTorsionGC.atomList[:-1]; self.setTorsionGC.update()'
     elif self.origValue==self.oldValue:
         return
     else:
         restoreAngle = self.origValue
         self.undoNow = 1
         undoCmd = 'self.setTorsionGC(\''+ aSet.full_name()+ '\',' + str(restoreAngle) + ', topCommand=0)'
     self.vf.undo.addEntry((undoCmd), (self.name))

 def doit(self, ats):
     if len(ats)>2:
         if len(self.atomList):
             atSet = ats + self.atomList
         else: atSet = ats
         parent = atSet[0].parent
         parent.buildBondsByDistanceOnAtoms(atSet)
         self.atomList = AtomSet([])
         self.update(True)
     else:
         lenAts = len(self.atomList)
         last = None
         if lenAts:
             last = self.atomList[-1]
             top = self.atomList[0].top
         for at in ats:
             #check for repeats of same atom
             if lenAts and at==last:
                 continue
             #lenAts = len(self.atomList)
             #if lenAts and at==self.atomList[-1]:
             #    continue
             if lenAts and at.top!=self.atomList[-1].top:
                 msg = "intermolecular bond to %s disallowed"%(at.full_name())
                 self.warningMsg(msg)
             self.atomList.append(at)
             self.undoAtList.append(at)
             lenAts = len(self.atomList)
         self.update(True)
         #if only have one atom, there is nothing else to do
         if lenAts<2: return
         #now build bonds between pairs of atoms
         atSet = self.atomList
         if lenAts%2!=0:
             atSet = atSet[:-1]
             #all pairs of atoms will be bonded
             #so keep only the last one
             self.atomList = atSet[-1:]
             lenAts = lenAts -1
         else:
             self.vf.labelByExpression(self.atomList, negate=1, topCommand=0)
             self.atomList = AtomSet([])
         for i in range(0, lenAts, 2):
             at1 = atSet[i]
             at2 = atSet[i+1]
             self.vf.addBonds(at1, at2, origin='UserDefined', topCommand=0)
     self.update(True)

 def buildCloseContactAtoms(self, percentCutoff, ligand, comment="USER AD> "):
     pairDict = self.distanceSelector.select(ligand.allAtoms,
                     self.macro_atoms, percentCutoff=percentCutoff)
     self.pairDict = pairDict
     #reset here
     lig_close_ats = AtomSet()
     macro_close_ats = AtomSet()
     cdict = {}
     for k,v in pairDict.items():
         if len(v):
             cdict[k] = 1
         for at in v:
             if at not in macro_close_ats:
                 cdict[at] = 1
     closeAtoms = AtomSet(cdict.keys())
     lig_close_ats = closeAtoms.get(lambda x: x.top==ligand).uniq()
     #ligClAtStr = lig_close_ats.full_name()
     ligClAtStr = comment + "lig_close_ats: %d\n" %( len(lig_close_ats))
     if len(lig_close_ats):
         ligClAtStr += comment + "%s\n" %( lig_close_ats.full_name())
     macro_close_ats = closeAtoms.get(lambda x: x in self.macro_atoms).uniq()
     macroClAtStr = comment + "macro_close_ats: %d\n" %( len(macro_close_ats))
     if len(macro_close_ats):
         macroClAtStr += comment + "%s\n" %( macro_close_ats.full_name())
     #macroClAtStr = "macro_close_ats: " + len(macro_close_ats)+"\n" +macro_close_ats.full_name()
     rdict = self.results
     rdict['lig_close_atoms'] = lig_close_ats
     rdict['macro_close_atoms'] = macro_close_ats
     if self.verbose: print "macroClAtStr=", macroClAtStr
     if self.verbose: print "ligClAtStr=", ligClAtStr
     if self.verbose: print "returning "+ macroClAtStr + '==' + ligClAtStr
     return  macroClAtStr , ligClAtStr

 def set_carbon_names(self, atoms, type):
     #set carbon names explicitly
     if not atoms or not len(atoms):
         return "ERROR: set_carbon_names called with no atoms"
     assert type in ['C','A']
     if not hasattr(atoms, 'autodock_element'):
         atoms.autodock_element = atoms.element
     changed = AtomSet()
     for at in atoms:
         if at.element!='C': continue
         if at.autodock_element!=type:
             if self.rename:
                 if len(at.name)>1:
                     at.name = type + at.name[1:]
                 else:
                     at.name = type
             at.autodock_element = type
             changed.append(at)
     return changed

 def onRemoveObjectFromViewer(self, mol):
     lenAts = len(self.atomList)
     #if cmd has no atoms on its list, nothing to do
     if not lenAts:
         return
     #remove any atoms which are being deleted from viewer
     self.atomList = AtomSet(self.atomList) - mol.allAtoms
     #if some have been removed, do an update
     if lenAts!=len(self.atomList):
         self.update()
         self.extslider.set(0)

 def removeNeighbors(self, atDict):
     #filter out at-itself and at-bondedat up to 1:4
     #NB keys could be hydrogens OR donors
     for at in atDict.keys():
         closeAts = atDict[at]
         bondedAts = AtomSet([])
         for b in at.bonds:
             ###at2 = b.neighborAtom(at)
             at2 = b.atom1
             if id(at2)==id(at): at2 = b.atom2
             bondedAts.append(at2)
             #9/13 remove this:
             ##also remove 1-3
             for b2 in at2.bonds:
                 at3 = b2.atom1
                 if id(at3)==id(at2): at3 = b.atom2
                 #at3 = b2.neighborAtom(at2)
                 if id(at3)!=id(at):
                     bondedAts.append(at3)
                 #for b3 in at3.bonds:
                     #at4 = b2.neighborAtom(at3)
                     #if at4!=at and at4!=at2:
                         #bondedAts.append(at4)
         bondedAts = bondedAts.uniq()
         goodAts = []
         for i in range(len(closeAts)):
             cAt = closeAts[i]
             if cAt not in bondedAts:
                 goodAts.append(cAt)
         if len(goodAts):
             atDict[at] = goodAts
         else:
             del atDict[at]
     return atDict

 def checkAromatics(self):
     """
     checkAromatics
     """
     #this depends on userPref useProteinAromaticList
     if not len(self.cyclebonds):
         self.aromaticCs = AtomSet()
         return ""
     if self.isPeptide and self.useProteinAromaticList:
         self.aromaticCs = self.getPeptideAromatics()
         return ""
     if self.isPeptide:
         self.getPeptideBondDict()
     else:
         self.getLigandBondDict()
     counter = 0
     while counter < self.cyclecount:
         counter = counter + 1
         blist = self.bondDict[counter]
         for item in blist:
             at = item.atom1
             self._getAdjAtom(item, blist)
             #now each bond has 3 atoms specified for it: its own two and the next1 to atom1
             result = self._getNormal(item)
             item.nrmsize = result[0]
             item.nrms = result[1]
             #next find the other bond w/atom2:
             z2 = filter(lambda x,item=item, at2=item.atom2, blist=blist:x!=item and x.atom1==at2 or x.atom2==at2, blist)
             #finally, return the other atom in this bond
             item.nextbond2 = z2[0]
             if item.nextbond2==item:
                 item.nextbond2 = z2[1]
             neighbor2 = self._getnxtAtom(item.atom2,item.nextbond2)
             item.next2 = neighbor2
             #next have to check whether the normals are parallel
             #check each pair in each bond, how to keep track??
             #have to get normal at item's atom2: so have to get next2 for this bond:
         #have to loop twice to make sure neighbor has nrms
         for item in blist:
             p = item.nrms
             psize = item.nrmsize
             q = item.nextbond2.nrms
             qsize = item.nextbond2.nrmsize
             #theta is the critical test for planarity:
             #if angle between 2 nrms is 0, atoms are planar
             #NB>test is comparing theta,cos(angle), w/zero
             item.theta = Numeric.dot(p,q)/(psize*qsize)
             for p in ['next1','next2','nextbond','nextbond2']:
                 delattr(item, p)
     self.updateAromatics(self.aromaticCutOff)
     msg = '    -found '+ str(len(self.aromaticCs)) + ' aromatic carbons\n'
     return  msg

 def superimpose_cb(self):
     refAtoms = AtomSet()
     mobAtoms = AtomSet()
     for pair in self.newPairs.values():
         refAtoms.append(pair[0])
         mobAtoms.append(pair[1])
         
     apply( self.doitWrapper, (refAtoms,mobAtoms),  {} )

 def buildHydrogenBonds(self, ligand, comment="USER AD> "):
     h_pairDict = self.hydrogen_bond_builder.build(ligand.allAtoms, self.macro_atoms)
     self.h_pairDict = h_pairDict
     #keys should be from lig, values from macro 
     #sometimes are [email protected]@check [email protected]@
     h_results = {}
     for k, v in h_pairDict.items():
         h_results[k] = 1
         for at in v:
             h_results[at] = 1
     all_hb_ats = AtomSet(h_results.keys())  #all
     d = self.results
     macro_hb_ats = d['macro_hb_atoms'] = all_hb_ats.get(lambda x: x.top==self.macro)
     self.macro_hb_ats = macro_hb_ats
     # process lig
     lig_hb_ats = d['lig_hb_atoms'] = all_hb_ats.get(lambda x: x in ligand.allAtoms)
     self.lig_hb_ats = lig_hb_ats
     outS = comment + "lig_hb_atoms : %d\n"%(len(lig_hb_ats))
     for p in self.lig_hb_ats:  #intD.results['lig_hb_atoms']:
         for hb in p.hbonds:
             if hasattr(hb, 'used'): continue
             if hb.hAt is not None:
                 outS += comment + "%s,%s~%s\n"%(hb.donAt.full_name(), hb.hAt.name, hb.accAt.full_name())
             else:
                 outS += comment + "%s~%s\n"%(hb.donAt.full_name(), hb.accAt.full_name())
             hb.used = 1
             #hsg1V:B:ARG8:NH2,HH22~clean: : INI 20:N5
             #clean: : INI 20:O4,H3~hsg1V:B:ASP29:OD2
             #clean: : INI 20:O4,H3~hsg1V:B:ASP29:OD2
             #clean: : INI 20:N4,H3~hsg1V:B:GLY27:O
             #clean: : INI 20:O2,H2~hsg1V:B:ASP25:OD1
     #macroHStr = self.macro.allAtoms.get(lambda x: hasattr(x, 'hbonds') and len(x.hbonds)).full_name()
     #ligHStr = ligand.allAtoms.get(lambda x: hasattr(x, 'hbonds') and len(x.hbonds)).full_name()
     #return  macroHStr + '==' + ligHStr
     if self.verbose: 
         print  "buildHB returning:"
         print outS
     return outS