From: Blake Charlebois (bdc_at_mie.utoronto.ca)
Date: Tue Jul 19 2005 - 11:28:59 CDT

Hi,

You can use Orient
(http://www.ks.uiuc.edu/Research/vmd/script_library/scripts/orient/) to get
the principal axes of your molecule, or you can use a simpler approach to
defining rotation.

I do not have what you are looking for, but I have a few procedures that may
help. I do not guarantee that they work perfectly. Also, I think vmd-l is
going to add carriage returns that will have to be removed:

proc atom_to_vector {mol_id atom_select_text {frame_num now}} {

 # e.g., set v1 [atom_to_vector top\
 # "protein and resid 1 and name CA" 0]

        set atom1 [atomselect $mol_id $atom_select_text frame $frame_num]
        set retval [lindex [$atom1 get {x y z}] 0]

        # it is a good habit to delete atom selections, though it
        # is only important in loops:
        $atom1 delete

        if {[$atom1 num]!=1} {
                puts "ERROR: expected to find 1 atom with"
                puts "ERROR: $atom_select_text,"
                puts "ERROR: found [$atom1 num]"
                return "ERROR"
                # this is not proper error trapping
        } else {
                return $retval
        }

}

proc vector_between_atoms_in_two_molecules {mol_id1 mol_id2 sel_txt1
sel_txt2\
                                                frame_num1 frame_num2} {

 # e.g., vector_between_atoms_in_two_molecules top top\
 # "protein and resid 100 and name CA"\
 # "protein and resid 190 and name CA"\
 # now now]

        set v1 [atom_to_vector $mol_id1 $sel_txt1 $frame_num1]
        set v2 [atom_to_vector $mol_id2 $sel_txt2 $frame_num2]

        return [vecsub $v1 $v2]

}

proc angle_between_vectors_in_radians {vec1 vec2} {
        set cos_of_angle [vecdot [vecnorm $vec1] [vecnorm $vec2]]

        # if cos of angle is a little bit more
        # than 1, round it to 1 to avoid an error.
        # If it is significantly more than 1, complain.
        if {$cos_of_angle>1} {
                if {$cos_of_angle<1.0000001} {
                        set cos_of_angle 1.0
                } else {
                        puts "ERROR: cannot take arccos of $cos_of_angle"
                                # tcl also outputs an error message, but
                                # this one makes it clear what went wrong
                }
        }
        if {$cos_of_angle<-1} {
                if {$cos_of_angle>-1.0000001} {
                        set cos_of_angle -1.0
                } else {
                        puts "ERROR: cannot take arccos of $cos_of_angle"
                }
        }

        return [expr {acos($cos_of_angle)}]

}

proc radians_to_degrees {angle_in_radians} {
        # pi = acos(-1) = 3.14159265359
        return [expr {$angle_in_radians*180/3.14159265359}]
}

-----Original Message-----
From: owner-vmd-l_at_ks.uiuc.edu [mailto:owner-vmd-l_at_ks.uiuc.edu] On Behalf Of
Vani Krishna
Sent: July 19, 2005 11:08 AM
To: vmd
Subject: vmd-l: measuring rotation of protein about an axis

Hello:

Is there a way to measure the rotation of a protein
about an axis for a trajectory in VMD?

Thanks
Vani

__________________________________________________
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around
http://mail.yahoo.com