VMD-L Mailing List

From: Andres Morales N (andresmoralesn2_at_hotmail.com)
Date: Tue May 19 2009 - 08:51:23 CDT

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!
http://spaces.live.com/spacesapi.aspx?wx_action=create&wx_url=/friends.aspx&mkt=en-us