Re: tcl script for rmsd

From: Aron Broom (broomsday_at_gmail.com)
Date: Thu Oct 04 2012 - 11:18:29 CDT

As a quick note here, there is a script built-in, under the "Analysis" ->
"RMSD Trajectory Tool", that could very well do exactly what you are
interested in.

On Thu, Oct 4, 2012 at 5:43 AM, yp sun <sunyeping_at_yahoo.com.cn> wrote:

> Dear all,
>
>
> I want to analyze the NAMD simulation results, especially RMSD during
> the time. But it is difficult for me to understand the tcl script. The
> following is a tcl script from NAMD-tutorial for computing rmsd of each
> residue during the timesteps:
>
>
>
>
> proc rmsd_residue_over_time {{mol top} res} {
>
> # use frame 0 for the reference
> set reference [atomselect $mol "protein" frame 0]
> # the frame being compared
> set compare [atomselect $mol "protein"]
> #make a selection with all atoms
> set all [atomselect top all]
> #get the number of frames
> set num_steps [molinfo $mol get numframes]
> #open file for writing
> set fil [open residue_rmsd.dat w]
>
> foreach r $res {
> set rmsd($r) 0
> }
>
> #loop over all frames in the trajectory
> for {set frame 0} {$frame < $num_steps} {incr frame} {
> puts "Calculating rmsd for frame $frame ..."
> # get the correct frame
> $compare frame $frame
> $all frame $frame
> # compute the transformation
> set trans_mat [measure fit $compare $reference]
> # do the alignment
> $all move $trans_mat
>
> # compute the RMSD
> #loop through all residues
> foreach r $res {
> set ref [atomselect $mol "protein and resid $r and noh" frame 0]
> set comp [atomselect $mol "protein and resid $r and noh" frame $frame]
> set rmsd($r) [expr $rmsd($r) + [measure rmsd $comp $ref]]
> $comp delete
> $ref delete
> }
> }
>
> set ave 0
> foreach r $res {
> set rmsd($r) [expr $rmsd($r)/$num_steps]
> # print the RMSD
> puts "RMSD of residue $r is $rmsd($r)"
> puts $fil " $r \t $rmsd($r)"
> set res_b [atomselect $mol "resid $r"]
> $res_b set user $rmsd($r)
> $res_b delete
> set ave [expr $ave + $rmsd($r)]
> }
>
> set ave [expr $ave/[llength $res]]
> puts " Average rmsd per residue: $ave"
> close $filI
> }
>
>
> Could you explain in detail what are the functions of every stentences?
> For example, what role is the following loop sentence for?
>
> " foreach r $res {
> set rmsd($r) 0
> }"
>
> My simulation is done on a protein containing three chains. I feel this
> script is not suitable for calculation every residue in individual chains.
> Then how should I modify it? Thanks
>
>
> Best regards!
>
> Yeping Sun
> CAS Key Laboratory of Pathogenic Microbiology & Immunology
> INSTITUTE OF MICROBIOLOGY CHINESE ACADEMY OF SCIENCES
> NO.1 Beichen West Road,Chaoyang District,Beijing 100101,china
>

-- 
Aron Broom M.Sc
PhD Student
Department of Chemistry
University of Waterloo

This archive was generated by hypermail 2.1.6 : Tue Dec 31 2013 - 23:22:39 CST