From: jen hsin (jhsin_at_ks.uiuc.edu)
Date: Tue May 19 2009 - 09:08:41 CDT

Hi Andres,

Actually you can use VMD to do the alignment without a script.
Just use the RMSD trajectory tool plugin under the "Extensions" tab.
Or are you trying to do something specific?

-Jen

On May 19, 2009, at 8:51 AM, Andres Morales N wrote:

>
> Hi!
>
> I wrote a script to align and calculate average strcuture form
> a .dcd file.
> I wait you can look it and give me some suggestions:
>
>
> ##### Align
> ##move atoms of each strcuture from a dcd file to a position where
> alpha carbon of first residue will be located at point {0 0 0}
> set O [atomselect top "protein and resid 1 and alpha"]
> set nf [molinfo top get numframes]
> set protein [atomselect top all]
> for {set i 0 } {$i < $nf } { incr i } {
> $O frame $i
> $protein frame $i
> set ocoord [lindex [$O get {x y z}] 0]
> $protein moveby [vecscale -1 $ocoord]
> }
> ## Align vectors (between CA of resid 1 and 15) of each frame
> set P0 [atomselect top "protein and resid 15 and alpha" frame 0]
> set P0coord [lindex [$P0 get {x y z}] 0]
> set P [atomselect top "protein and resid 15 and alpha"]
> for {set i 0 } {$i < $nf } { incr i } {
> $P frame $i
> $protein frame $i
> set Pcoord [lindex [$P get {x y z}] 0]
> global M_PI
> set cosine [expr [vecdot $Pcoord $P0coord] / ( [veclength $Pcoord] *
> [veclength $P0coord])]
> set angulo1 [expr acos($cosine)*(180.0/$M_PI)]
> set A [veccross $Pcoord $P0coord]
> set M1 [transabout $A $angulo1]
> $protein move $M1
> }
>
> ## Align plane (formed between CA of resid 1 15 and 50) of each frame
> set Q0 [atomselect top "protein and resid 50 and alpha" frame 0]
> set Q0coord [lindex [$Q0 get {x y z}] 0]
> set Q [atomselect top "protein and resid 50 and alpha"]
> set N0 [veccross $P0coord $Q0coord]
> set O0 [atomselect top "protein and resid 1 and alpha" frame 0]
> set O0coord [lindex [$O0 get {x y z}] 0]
> for {set i 0 } {$i < $nf } { incr i } {
> $P frame $i
> $Q frame $i
> $protein frame $i
> set Pcoord [lindex [$P get {x y z}] 0]
> set Qcoord [lindex [$Q get {x y z}] 0]
> set N [veccross $Pcoord $Qcoord]
> global M_PI
> set cosine2 [expr [vecdot $N $N0] / ( [veclength $N] * [veclength
> $N0])]
> set angulo2 [expr acos($cosine2)*(180.0/$M_PI)]
> set B [veccross $N $N0]
> set M2 [transabout $B $angulo2]
> $protein move $M2
> }
>
> ### Average strcuture
>
> set sel [atomselect top all]
> set nf [molinfo top get numframes]
> set n [expr $nf - 1]
> set newpos [measure avpos $sel first 0 last $n]
> $sel set {x y z} $newpos
> $sel writepdb "promedio.pdb"
>
>
> Thanks for your suggestions
>
> Invite your mail contacts to join your friends list with Windows
> Live Spaces. It's easy! Try it!