From: Joaquim Rui de Castro Rodrigues (
Date: Tue Nov 05 2013 - 03:46:17 CST


You must set the current frame for your "protein" and "ligand" selections inside the loop, right after $all frame $j:
$protein frame $j
$ligand frame $j

Hope this helps,
J Rui Rodrigues

De: [] Em Nome De Jeremiah Babcock []
Enviado: terša-feira, 5 de Novembro de 2013 3:55
Assunto: vmd-l: measuring angles between principal axes through multiple frames

    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 incorrect?

# 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 the
# 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