VMD-L Mailing List

From: Mahad Gatti (Mahad.Gatti_at_iit.it)
Date: Wed Aug 26 2015 - 05:41:09 CDT

Hi everyone,

I'm writing a script to calculate a dihedral angle (chi torsional) in a series of dcd files, using big_dcd (working with RNA).
So far so good, I was able to write all the instructions as you will see here:

source bigdcd.tcl
proc myrmsd { frame } {

   for {set i 636} {$i<=651} {incr i} { #this range is the one I want

  set fout1 [open "dihe_$i.dat" a+]
  set molid 0

  set a1 [[atomselect top "resid $i and name O4'"] get index]
  set a2 [[atomselect top "resid $i and name C1'"] get index]
  set a3 [[atomselect top "resid $i and name N1"] get index]
  set a4 [[atomselect top "resid $i and name C2"] get index]

  label delete Dihedrals
  set molid [molinfo top]
  set curframe [molinfo top get frame]
  molinfo top set frame $frame
  label add Dihedrals $molid/$a1 $molid/$a2 $molid/$a3 $molid/$a4
# Get its value
  set curval [lindex [lindex [label list Dihedrals] 0] 4]
  molinfo top set frame $frame
  puts $fout1 [format "%f" $curval]
  close $fout1
  label delete Dihedrals
}

}
set mol [mol new mutated.prmtop type parm7 waitfor all]
# MFL DIST

mol addfile mutated.inpcrd type rst7 waitfor all
bigdcd myrmsd dcd test-1.dcd ....test-n.dcd

The problem is related to the differences between pyrimidines and purines, because the atoms involved (to be more precise the nomenclature behind) are different, thus I wanted to add a flag like an if cycle to control and choose the proper selection for both cases.

Pyrimidines = O4'-C1'-N1-C2
Purines = O4'-C1'-N9-C4

Thing is that I don't know how to use an if cycle in tcl that takes the returning value of a atomselect function (like [$a1 get resname] ) and use that information for the right selection.
My idea was to use the get resname command and check if the string is U/A/G/C, but I get errors regarding the syntax probably. I will try to write it down in pseudocode.

If {[$a1 get resname] equals (U or C)}
    then atomselect O4'-C1'-N1-C2 (pyrimidines)
else
    atomselect O4'-C1'-N9-C4 (purines)
    }

Hope that I was clear enough.

Thanks and regards.

Mahad