From: Andres Morales N (andresmoralesn2_at_hotmail.com)
Date: Mon May 11 2009 - 21:44:55 CDT

Hi!
 
I have a clsuter of strcutures, so I need to align them and to get the average structure. I tried with this script:
___________________________________________________

set seltext "protein"

set mol top;

set stride 1;

set frame_start 0;

set num_steps [molinfo $mol get numframes]

set frame_end $num_steps;

set frame_reference 0;

 

#Aligning Selection

 

 

set sel [atomselect $mol "$seltext"];

set reference [atomselect $mol "$seltext" frame $frame_reference]

set compare [atomselect $mol "$seltext"]

set all [atomselect $mol all]

for {set frame $frame_start} {$frame < $frame_end} {} {

$compare frame $frame

$all frame $frame

set trans_mat [measure fit $compare $reference]

$all move $trans_mat

set rmsd [measure rmsd $compare $reference]

puts "$frame $rmsd"

set frame [expr $frame + $stride];

}

display update

 

#Average Structure Calculation

set frame 0;

set indices 0;

lvarpop indices;

set newcoords 0;

lvarpop newcoords;

set num_steps2 [expr (($frame_end - $frame_start) / $stride) + 1]

puts "The number of steps is $num_steps2"

set inv_num_steps [expr 1.0 / $num_steps2]

set indices [$sel get index];

foreach indice $indices {

set suma [veczero];

puts "Index is $indice";

for {set frame $frame_start} {$frame < $frame_end} {} {

set seli [atomselect $mol "index $indice" frame $frame];

set coord [$seli get {x y z}];

foreach pieza $coord {

puts "The coordinates are $coord en frame $frame";

set suma [vecadd $suma $pieza];

}

set frame [expr $frame + $stride];

$seli delete;

}

set resul [vecscale $inv_num_steps $suma];

lappend newcoords $resul;

}

set frame 0;

$sel set {x y z} $newcoords;

$sel writepdb average_selection.pdb

display update;

____________________________________________________________________________________

 

But I did not get good results. I measured rgyr of averaged structure and it was around 8 A, so I measured rgyr of all structures from cluster, and everyone was about 11 A. So I think something is wrong with alignment. But I do not know what. I wait somebody can help me.

 

Also I was trying with the script Alignment to principal axes from http://www.ks.uiuc.edu/Research/vmd/script_library/scripts/orient/. It is only for one strcuture.

If somebody know how to use it to align many strcutures with that, please help me.

 

Thanks a lot

 

 

 

 

 

_________________________________________________________________
Discover the new Windows Vista
http://search.msn.com/results.aspx?q=windows+vista&mkt=en-US&form=QBRE