From: Jeremiah Babcock (
Date: Mon Nov 04 2013 - 21:55:13 CST

    I built a fluorophore covalently linked N-terminus to a protein and ran
simulations. I want to compare the segmental motion of the ligand to the
protein inertial frame by calculating rotation angles between principal
axes of both protein and ligand throughout my simulation timeframe. I used
the Orient plugin and designed the included script. The problem is is that
in my loop, the principal axes command is not calculating new axes for each
new trajectory frame: it calculates the first frame and afterwards it stays
constant. I know this isn't correct, as manually selecting a frame returns
different axial values for different frames. Am I doing something

# Script to find rotation angles between principal axes
package require Orient
namespace import Orient::orient

# Open a file to save inertial axes angles to
# Select protein atoms without the ligand
# Select ligand atoms
# Select protein and ligand atoms together

set file [open InertialAngles.dat w]
set protein [atomselect top "protein and not index 0 to 46"]
set ligand [atomselect top "index 0 to 46"]
set all [atomselect top "protein"]
set nfram [molinfo 0 get numframes]

# Label output file for the 4 x 1 matrix of time and axial rotation angles
put $file "Time Th_x Th_y Th_z"

# Create a loop through all DCD frames to find first the principal axes of
# protein without ligand, and then the ligand principal axes. Then find the
# angle of rotation between protein and ligand axes.

for {set j 0} {$j < $nfram } { incr j } {
    $all frame $j
    $all update

    # update frame number
    put "Frame $j"

    # convert frame number to time (ns)
    set time [expr $j*(.001)]

    # Find protein and ligand principal axes
    set Ip [Orient::calc_principalaxes $protein]
    set Il [Orient::calc_principalaxes $ligand]

    # Find rotation angles between protein and ligand principal axes
    set Ax [Orient::orient $all [lindex $Ip 0] [lindex $Il 0]]
    set Ay [Orient::orient $all [lindex $Ip 1] [lindex $Il 1]]
    set Az [Orient::orient $all [lindex $Ip 2] [lindex $Il 2]]

    # Create results list and print to file
    set type [list $time $Ax $Ay $Az]
    put $file "$type"
close $file

Thank you for your help.

Jeremiah Babcock
PhD Candidate, Molecular Biophysics
Department of Physics and Astronomy
University of Texas at San Antonio
AET 3.206, One UTSA Circle
San Antonio, TX 78249
Voice: (210) 458-6191