Difference between revisions of "Center of mass"

From PyMOLWiki
Jump to: navigation, search
(support openbabel element table for atom mass lookup)
(mass update)
Line 39: Line 39:
 
from pymol import cmd
 
from pymol import cmd
  
class AtomMass:
+
def com(selection,state=None,mass=None,object=None, quiet=1, **kwargs):
  def __init__(self):
+
  """
      from openbabel import OBElementTable
+
DESCRIPTION
      self.table = OBElementTable()
+
 
  def __getitem__(self, symbol):
+
   Places a pseudoatom at the center of mass
      return self.table.GetMass(self.table.GetAtomicNum(symbol))
 
   def __contains__(self, symbol):
 
      return self.table.GetAtomicNum(symbol) > 0
 
  
def com (selection,state=None,mass=None,name=None):
 
 
 
  """
 
 
   Author: Sean Law
 
   Author: Sean Law
 
   Michigan State University
 
   Michigan State University
 
   slaw (at) msu . edu
 
   slaw (at) msu . edu
 +
 +
SEE ALSO
 +
 +
  pseudoatom, get_com
 
   """
 
   """
   if (name == None):
+
  quiet = int(quiet)
       name=selection+"_COM"
+
   if (object == None):
   cmd.delete(name)
+
       object = cmd.get_legal_name(selection)
 
+
      object = cmd.get_unused_name(object + "_COM", 0)
 +
   cmd.delete(object)
 +
 
 
   if (state != None):
 
   if (state != None):
       x, y, z=get_com(selection,mass=mass)
+
       x, y, z=get_com(selection,mass=mass, quiet=quiet)
       print "%f %f %f" % (x, y, z)
+
       if not quiet:
       cmd.pseudoatom(name,pos=[x, y, z])
+
        print "%f %f %f" % (x, y, z)
       cmd.show("spheres",name)
+
       cmd.pseudoatom(object,pos=[x, y, z], **kwargs)
 +
       cmd.show("spheres",object)
 
   else:
 
   else:
 
       for i in range(cmd.count_states()):
 
       for i in range(cmd.count_states()):
         x, y, z=get_com(selection,mass=mass,state=i+1)
+
         x, y, z=get_com(selection,mass=mass,state=i+1, quiet=quiet)
         print "State %d:%f %f %f" % (i+1, x, y, z)
+
         if not quiet:
         cmd.pseudoatom(name,pos=[x, y, z],state=i+1)
+
            print "State %d:%f %f %f" % (i+1, x, y, z)
         cmd.show("spheres",name)
+
         cmd.pseudoatom(object,pos=[x, y, z],state=i+1, **kwargs)
 
+
         cmd.show("spheres", 'last ' + object)
  return
+
 
 
 
cmd.extend("com",com)
 
cmd.extend("com",com)
  
def get_com (selection,state=1,mass=None):
+
def get_com(selection,state=1,mass=None, quiet=1):
 +
  """
 +
DESCRIPTION
 +
 
 +
  Calculates the center of mass
  
  """
 
 
   Author: Sean Law
 
   Author: Sean Law
 
   Michigan State University
 
   Michigan State University
 
   slaw (at) msu . edu
 
   slaw (at) msu . edu
 
   """
 
   """
 +
  quiet = int(quiet)
  
  try:
 
      atmass = AtomMass()
 
  except:
 
      atmass={'H':1.00800000000000, 'C':12.0110000000000, \
 
              'N':14.0070000000000, 'O':15.9994000000000, \
 
              'P':30.9740000000000, 'S':32.0600000000000, \
 
              'F':18.9980000000000
 
          }
 
 
 
 
   totmass=0.0
 
   totmass=0.0
   if (mass!=None):
+
   if mass!=None and not quiet:
 
       print "Calculating mass-weighted COM"
 
       print "Calculating mass-weighted COM"
  
Line 101: Line 96:
 
   for a in model.atom:
 
   for a in model.atom:
 
       if (mass != None):
 
       if (mass != None):
         if (a.symbol in atmass):
+
         m = a.get_mass()
            x+= a.coord[0]*atmass[a.symbol]
+
        x+= a.coord[0]*m
            y+= a.coord[1]*atmass[a.symbol]
+
        y+= a.coord[1]*m
            z+= a.coord[2]*atmass[a.symbol]
+
        z+= a.coord[2]*m
            totmass+=atmass[a.symbol]
+
        totmass+=m
        else:
 
            print "Please add atomic mass for "+a.symbol+" to atmass"
 
 
       else:
 
       else:
 
         x+= a.coord[0]
 
         x+= a.coord[0]
Line 117: Line 110:
 
   else:
 
   else:
 
       return x/len(model.atom), y/len(model.atom), z/len(model.atom)
 
       return x/len(model.atom), y/len(model.atom), z/len(model.atom)
 +
 
cmd.extend("get_com", get_com)
 
cmd.extend("get_com", get_com)
  
 +
# vi:expandtab:sw=3
 
</source>
 
</source>
  

Revision as of 14:45, 31 March 2011

Description

This script calculates the true center-of-mass (COM) or the center-of-geometry (COG) for a given selection and returns the x, y, z values in the form of a Pseudoatom (rather than a CGO sphere). The benefit of using a Pseudoatom is that it can be selected and used in calculations. In addition, this script also iteratively goes through all states of a selection if more than one state exists and appends the corresponding COM/COG values as states into the Pseudoatom. The script itself is quite simple yet robust enough to be applied in different settings. As well, the calculation of the COM/COG is handled independently from the formation of the Pseudoatom and can be called as an independent function where applicable.

Note: In order to use the Pseudoatom in your measurements (i.e. distance), you need to invoke the distance calculation directly via the command line using the Distance function and not via the Wizard!

Usage

com selection [,state=None [,mass=None [,name=None]]]


Examples

fetch 1c3y, finish=1, multiplex=0

com 1c3y, state=1
#Create a pseudoatom representing the 1c3y COG and store it as "1c3y_COM"
#The "1c3y_COM" object will contain 1 state only

com 1c3y, state=1, name=COG
#Create a pseudoatom representing the 1c3y COG and store it as "COG"
#The "COG" object will contain 1 state only

com 1c3y, state=1, name=COM, mass=1
#Create a pseudoatom representing the 1c3y COM and store it as "COM"
#The "COM" object will contain 1 state only

com 1c3y, name=COM, mass=1
#Create a single pseudoatom containing the COM for each state found in 1c3y and store it as "COM"
#The "COM" object will contain MULTIPLE states!

PyMOL API

from pymol import cmd

def com(selection,state=None,mass=None,object=None, quiet=1, **kwargs):
   """
DESCRIPTION

   Places a pseudoatom at the center of mass

   Author: Sean Law
   Michigan State University
   slaw (at) msu . edu

SEE ALSO

   pseudoatom, get_com
   """
   quiet = int(quiet)
   if (object == None):
      object = cmd.get_legal_name(selection)
      object = cmd.get_unused_name(object + "_COM", 0)
   cmd.delete(object)

   if (state != None):
      x, y, z=get_com(selection,mass=mass, quiet=quiet)
      if not quiet:
         print "%f %f %f" % (x, y, z)
      cmd.pseudoatom(object,pos=[x, y, z], **kwargs)
      cmd.show("spheres",object)
   else:
      for i in range(cmd.count_states()):
         x, y, z=get_com(selection,mass=mass,state=i+1, quiet=quiet)
         if not quiet:
            print "State %d:%f %f %f" % (i+1, x, y, z)
         cmd.pseudoatom(object,pos=[x, y, z],state=i+1, **kwargs)
         cmd.show("spheres", 'last ' + object)

cmd.extend("com",com)

def get_com(selection,state=1,mass=None, quiet=1):
   """
DESCRIPTION

   Calculates the center of mass

   Author: Sean Law
   Michigan State University
   slaw (at) msu . edu
   """
   quiet = int(quiet)

   totmass=0.0
   if mass!=None and not quiet:
      print "Calculating mass-weighted COM"

   state=int(state)
   model = cmd.get_model(selection,state)
   x,y,z=0,0,0
   for a in model.atom:
      if (mass != None):
         m = a.get_mass()
         x+= a.coord[0]*m
         y+= a.coord[1]*m
         z+= a.coord[2]*m
         totmass+=m
      else:
         x+= a.coord[0]
         y+= a.coord[1]
         z+= a.coord[2]

   if (mass!=None):
      return x/totmass, y/totmass, z/totmass
   else:
      return x/len(model.atom), y/len(model.atom), z/len(model.atom)

cmd.extend("get_com", get_com)

# vi:expandtab:sw=3


Previous Implementation

Here is a script that calculates the center of geometry from a selection. It gets hold of the coordinates with cmd.get_model. Make sure the atoms in the selection are of equal weight.

For a sample application, see: "Convergent Evolution Examples"

## Author: Andreas Henschel 2006

from pymol import cmd
from pymol.cgo import *

def centerOfMass(selection):
   ## assumes equal weights (best called with "and name ca" suffix)
   model = cmd.get_model(selection)
   x,y,z=0,0,0
   for a in model.atom:
       x+= a.coord[0]
       y+= a.coord[1]
       z+= a.coord[2]
   return (x/len(model.atom), y/len(model.atom), z/len(model.atom))

cmd.load("/group/bioinf/Data/PDBLinks/1c7c.pdb")
cmd.select("domain", "/1c7c//A/143-283/ and name ca") ## selecting a domain

domainCenter=centerOfMass("domain")

print "Center of mass: (%.1f,%.1f,%.1f)"% domainCenter
cmd.as("cartoon", "all")
cmd.show("spheres", "domain")

## Creating a sphere CGO
com = [COLOR, 1.0, 1.0, 1.0, SPHERE]+list(domainCenter) + [3.0] ## white sphere with 3A radius
cmd.load_cgo(com, "CoM")

cmd.zoom("1c7c", 1.0)
cmd.center("domain")

#ah@bioinfws19:~/Projects/PyMOL$ pymol -qc centerOfMass4.py
#Center of mass: (-1.0,24.5,48.2)
#ah@bioinfws19:~/Projects/PyMOL$