• Outreach

From: Mert Gür (gurmert_at_gmail.com)
Date: Tue Nov 27 2012 - 16:13:21 CST

Hi Mustafa,
I did not take a close look into your script. However, there is one thing I
thought is worth checking.
You get the total number of alpha carbons first and then try to select them
with the "resid" command.
Lets say you have 3 residues and the residue indices are 11 12 13 in your
pdb file. So the total number of alpha carbons you will get is 3. If you
try to select resid 3 you will end up with nothing. That might be the
problem here. You may want to check the pdb file if its residue numbering
starts with one. If it does not you have two options;
1) Renumber the residues so that they start from one
2) Use "residue" command instead of "resid". The first residue is selected
as " residue 0", so you need to subtract 1 from your numbers before you
select them.
I hope this helps.
Best,
Mert
On Nov 27, 2012 4:52 PM, "Mustafa Tekpinar" <tekpinar_at_buffalo.edu> wrote:

> Hi,
> I am trying to produce an Elastic Network Model representation by using a
> simple VMD script. However, each time I get the following error:
> "vecsub: two vectors don't have the same size"
>
> Here is the script I wrote:
> //==========================================
> #Select a CPK representation of CA atoms.
>
>
> mol modselect 0 0 name CA
> mol modstyle 0 0 CPK 1.000000 0.300000 10.000000 10.000000
>
>
>
> set R_c 10.0
> set ca_atoms [atomselect top "name CA"]
> set total_ca [\$ca_atoms num]
> puts "Total number of CA atoms are \$total_ca"
> for { set i 1 } { \$i <= \$total_ca } { incr i } {
> set coord1 [lindex [[atomselect top "resid \$i and name CA"] get {x y
> z}] 0]
> for { set j 1 } { \$j < \$i} { incr j } {
> set coord2 [lindex [[atomselect top "resid \$j and name CA"] get {x
> y z}] 0]
> set dist [veclength [vecsub \$coord1 \$coord2]]
> if {\$dist <= \$R_c} {
> if {[expr (\$i-\$j)]==1 } {
> draw color red
> draw cylinder \$coord1 \$coord2 radius 0.2
> } else {
> draw color yellow
> draw cylinder \$coord1 \$coord2 radius 0.1
> }
> }
> }
> }
>
> //==========================================
> I will really be glad if anybody can point me my mistake.
> Mustafa Tekpinar
>