Cart to frac

From PyMOLWiki
Revision as of 11:45, 22 November 2011 by Inchoate (talk | contribs)
Jump to navigation Jump to search

Overview

This script will convert the real space orthonormal cartesian coordinates of a selection to the fractional coordinates given the loaded cell symmetry.

Example

# load the script

run cart_to_frac.py

# fetch a protein and test

fetch 1rx1, async=0

# get the coordinates for the organic small
# molecule as fractional coordinates

print cart_to_frac("org")

The Code

import pymol
from pymol import cmd
import numpy
from numpy import cos, sin, sqrt

def cart_to_frac(objSel,quiet=0,_self=cmd):
    """
    Returns an array of fractional atomic coordinates
    for a given object or selection.

    PARAMS
      objSel -- any object or selection

      quiet -- suppress output (default, no)

      _self -- core CMD object; or none

    RETURNS
      Python list of fractional coordinates

    NOTES/EXAMPLES
      cart_to_frac org

      x = cart_to_frac("solvent", quiet=1)
    """
    a2r = numpy.pi / 180.0

    # get the model and coordinates

    m = _self.get_model(objSel)
    cart_coord = numpy.matrix(m.get_coord_list())

    # get the symmetry information
    try:
        a,b,c,alpha,beta,gamma,gp = _self.get_symmetry(objSel)
    except:
        print "Error-Failed to get symmetry. Please ensure you have a"
        print "valid object with proper crystal symmetry loaded."
        return None

    # convert to radians

    alpha = a2r * alpha
    beta  = a2r * beta
    gamma = a2r * gamma
    
    # (scaled) volume of the cell
    
    v = sqrt(1 -cos(alpha)*cos(alpha) - cos(beta)*cos(beta) - cos(gamma)*cos(gamma) + 2*cos(alpha)*cos(beta)*cos(gamma))

    tmat = numpy.matrix( [
      [ 1.0 / a, -cos(gamma)/(a*sin(gamma)), (cos(alpha)*cos(gamma)-cos(beta)) / (a*v*sin(gamma))  ],
      [ 0.0,     1.0 / b*sin(gamma),         (cos(beta) *cos(gamma)-cos(alpha))/ (b*v*sin(gamma))  ],
      [ 0.0,     0.0,                        sin(gamma) / (c*v)                                    ] ]
      )

    if not quiet:
        print (tmat * cart_coord.T).T

    # return the Nx3 results
    return (tmat * cart_coord.T).T
    
cmd.extend("cart_to_frac", cart_to_frac)