Difference between revisions of "BbPlane"

From PyMOLWiki
Jump to: navigation, search
(improve speed, support two_sided_lighting)
(fix H selection)
Line 69: Line 69:
             if key in coords:
             if key in coords:
                 coords[key][1] = a.coord
                 coords[key][1] = a.coord
         for a in cmd.get_model(sel + " and ((n. N extend 1 and e. H) or (r. PRO and n. CD))").atom:
         for a in cmd.get_model("(hydro or n. CD) and nbr. (" + sel + " and n. N)").atom:
             key = '/%s/%s/%s/%s' % (obj,a.segi,a.chain,a.resi)
             key = '/%s/%s/%s/%s' % (obj,a.segi,a.chain,a.resi)
             if key in coords:
             if key in coords:

Revision as of 09:42, 17 September 2010

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).


# 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
#   Modified by Thomas Holder, 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


    objSel = string: protein object or selection {default: (all)}

    color = string: color name or number {default: white}

    transp = float: transparency component (0.0--1.0) {default: 0.0}


    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)
    # format input
    transp = float(transp)
    stored.AAs = []
    coords = dict()

    # need hydrogens on peptide nitrogen
    cmd.h_add('(%s) and n. N' % objSel)
    # get the list of residue ids
    for obj in cmd.get_object_list(objSel):
        sel = obj + " and (" + objSel + ")"
        for a in cmd.get_model(sel + " and n. CA").atom:
            key = '/%s/%s/%s/%s' % (obj,a.segi,a.chain,a.resi)
            coords[key] = [a.coord,None,None]
        for a in cmd.get_model(sel + " and n. O").atom:
            key = '/%s/%s/%s/%s' % (obj,a.segi,a.chain,a.resi)
            if key in coords:
                coords[key][1] = a.coord
        for a in cmd.get_model("(hydro or n. CD) and nbr. (" + sel + " and n. N)").atom:
            key = '/%s/%s/%s/%s' % (obj,a.segi,a.chain,a.resi)
            if key in coords:
                coords[key][2] = a.coord

    # 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."

    # prepare the cgo
    obj = [

    for res in range(0, len(stored.AAs)-1):
        curIdx, nextIdx = str(stored.AAs[res]), str(stored.AAs[res+1])

        # populate the position array
        pos = [coords[curIdx][0], coords[curIdx][1], coords[nextIdx][2], coords[nextIdx][0]]

        # if the data are incomplete for any residues, ignore
        if None in pos:
            print 'peptide bond %s -> %s incomplete' % (curIdx, nextIdx)

        if cpv.distance(pos[0], pos[3]) > 4.0:
            print '%s and %s not adjacent' % (curIdx, nextIdx)

        # fill in the vertex data for the triangles; 
        for i in [0,1,2,3,2,1]:

    # finish the CGO

    # update the UI
    newName =  cmd.get_unused_name("backbonePlane")
    cmd.load_cgo(obj, newName)
    cmd.set("cgo_transparency", transp, newName)

cmd.extend("bbPlane", bbPlane)