BbPlane
Jump to navigation
Jump to search
This script will draw a CGO plane between the backbone atoms of two neighboring residues. This is to show the planarity of the atoms. The image style this is meant to represent can be found many places, like "Introduction to Protein Structure" by Branden and Tooze (2nd ed. pp. 8).
Examples
# download the source and save as bbPlane.py
run bbPlane.py
fetch 1cll
# make planes for residues 4-9
bbPlane i. 4-10
The Source
#
# -- bbPLane.py - draws a CGO plane across the backbone atoms of
# neighboring amino acids
#
# Author: Jason Vertrees, 06/2010
# Copyright (C) Schrodinger
# Open Source License: MIT
#
from pymol.cgo import * # get constants
from pymol import cmd, stored
from chempy import cpv
def bbPlane(objSel='(all)', color='white', transp=0.0):
"""
Draws a plane across the backbone for a selection
PARAMS
objSel,
(object or selection) protein object or selection
color = string: color name or number {default: white}
transp = float: transparency component (0.0--1.0) {default: 0.0}
RETURNS
A new object representing the planes
NOTES
* You need to pass in an object or selection with at least two
amino acids. The plane spans CA_i, O_i, N-H_(i+1), and CA_(i+1)
* Avoid two-sided lighting until I rearrange the vertices
"""
# format input
transp = float(transp)
stored.AAs = []
# need hydrogens
cmd.h_add(objSel)
# get the list of residue ids
for obj in cmd.get_object_list(objSel):
cmd.iterate( obj + " and (" + objSel + ") and n. CA" , "stored.AAs.append('/"+obj+"/%s/%s/%s' % (segi,chain,resi))")
# need at least two amino acids
if len(stored.AAs) <= 1:
print "ERROR: Please provide at least two amino acids, the alpha-carbon on the 2nd is needed."
return
# prepare the cgo
obj = [
BEGIN, TRIANGLES,
COLOR,
]
obj.extend(cmd.get_color_tuple(color))
for res in range(0, len(stored.AAs)-1):
curIdx, nextIdx = str(stored.AAs[res]), str(stored.AAs[res+1])
corners = [
curIdx + "/CA",
curIdx + "/O",
nextIdx + "/N extend 1 and e. H",
nextIdx + "/CA"
]
# prolines
if cmd.count_atoms( "%s and resn PRO" % nextIdx ):
corners[2] = nextIdx + "/CD"
# if the data are incomplete for any residues, ignore
if map(cmd.count_atoms, corners) != [1,1,1,1]:
print 'peptide bond %s -> %s incomplete' % (curIdx, nextIdx)
continue
# populate the position array
pos = map(cmd.get_atom_coords, corners)
if cpv.distance(pos[0], pos[3]) > 4.0:
print '%s and %s not adjacent' % (curIdx, nextIdx)
continue
# fill in the vertex data for the triangles;
# two-sided lighting is an issue
obj.extend( [
VERTEX, pos[0][0], pos[0][1], pos[0][2],
VERTEX, pos[1][0], pos[1][1], pos[1][2],
VERTEX, pos[2][0], pos[2][1], pos[2][2],
VERTEX, pos[1][0], pos[1][1], pos[1][2],
VERTEX, pos[2][0], pos[2][1], pos[2][2],
VERTEX, pos[3][0], pos[3][1], pos[3][2],
])
# finish the CGO
obj.append(END)
# update the UI
newName = cmd.get_unused_name("backbonePlane")
cmd.load_cgo(obj, newName)
cmd.set("cgo_transparency", transp, newName)
cmd.extend("bbPlane", bbPlane)