From: Sanjeev (snjv_bs_at_rediffmail.com)
Date: Sun Mar 18 2007 - 14:16:42 CDT

Dear John, I have tried to fix it with '$sel update' and the likes, but it did not help. I have now nine frames and there are valid h-bonds in all of them. The only interesting aberration I noticed is that the atomselect_n is increasing by 40 each time I run, which I do not consider to be unusual either, except that the first frame has 13 bonds/13*3=39 sized array (for acceptors, donors and hydrogens). Beyond that it beats me. I suppose the aberration is reproducible as I have run it with different PDBs. The following is the script as of now: --- hb.tcl --- set sel [atomselect top all] set n_frames [molinfo top get numframes] set PI 3.14159 puts "No. of frames: $n_frames" for {set i 0} {$i < $n_frames} {incr i} { puts "Processing frame: $i" $sel frame $i $sel update foreach {d a h} [measure hbonds 2.8 30 $sel] { puts "Sel is: $sel" set numbonds 0 foreach ctr $d { set d_ary($numbonds) $ctr set numbonds [incr numbonds] } set numbonds 0 foreach ctr $h { set h_ary($numbonds) $ctr set numbonds [incr numbonds] } set numbonds 0 foreach ctr $a { set a_ary($numbonds) $ctr set numbonds [incr numbonds] } puts "Numbonds : $numbonds" } for {set j 0} {$j < $numbonds} {incr j} { set sel [atomselect top "index $d_ary($j)" frame $i] set d_a [$sel get name] set d_r [$sel get resname] set d_d [$sel get {x y z}] set d_i [$sel get resid] set D [lindex [$sel get { x y z }] $i] set sel [atomselect top "index $a_ary($j)" frame $i] set a_a [$sel get name] set a_r [$sel get resname] set a_d [$sel get {x y z}] set a_i [$sel get resid] set A [lindex [$sel get { x y z }] $i] set sel [atomselect top "index $h_ary($j)" frame $i] set h_a [$sel get name] set h_r [$sel get resname] set h_d [$sel get {x y z}] set h_i [$sel get resid] set H [lindex [$sel get { x y z }] $i] set dist [vecdist $D $A] set angl [expr 180 - acos([vecdot [vecsub $D $H] [vecsub $H $A]]/([vecdist $D $H] * [vecdist $H $A])) * 180 / $PI] set out_1 "$i $d_a $d_r $d_i $h_a $h_r $h_i $a_a $a_r $a_i" set out_2 [format "%s %5.2f %5.1f" $out_1 $dist $angl] puts $out_2 } } unset sel unset d_ary a_ary h_ary unset n_frames --- hb.tcl --- The output is below: --- output --- vmd > source hb.tcl No. of frames: 9 Processing frame: 0 Sel is: atomselect0 Numbonds : 13 0 OG1 THR 67 HG1 THR 67 OD1 ASN 65 2.79 162.0 0 O WAT 4626 H2 WAT 4626 O1P RA3 135 2.67 176.1 0 OH TYR 98 HH TYR 98 O LYS 38 2.73 171.6 0 N ILE 86 H ILE 86 O ARG 97 2.76 169.0 0 NH1 ARG 1 1HH1 ARG 1 O PHE 5 2.73 166.4 0 NE2 GLN 14 2HE2 GLN 14 O TRP 10 2.75 162.9 0 NZ LYS 38 HZ3 LYS 38 OD1 ASN 41 2.80 163.6 0 N ALA 99 H ALA 99 O ASP 84 2.76 153.4 0 N TRP 10 H TRP 10 O THR 6 2.79 164.8 0 N CYX 23 H CYX 23 OD2 ASP 100 2.80 165.1 0 OG SER 74 HG SER 74 O PHE 76 2.70 176.8 0 NH1 ARG 77 1HH1 ARG 77 O ARG 104 2.77 157.6 0 N THR 46 H THR 46 O LEU 44 2.76 152.5 Processing frame: 1 Sel is: atomselect39 Numbonds : 0 Processing frame: 2 Sel is: atomselect39 Numbonds : 0 Processing frame: 3 Sel is: atomselect39 Numbonds : 0 Processing frame: 4 Sel is: atomselect39 Numbonds : 0 Processing frame: 5 Sel is: atomselect39 Numbonds : 0 Processing frame: 6 Sel is: atomselect39 Numbonds : 0 Processing frame: 7 Sel is: atomselect39 Numbonds : 0 Processing frame: 8 Sel is: atomselect39 Numbonds : 0 vmd > source hb.tcl No. of frames: 9 Processing frame: 0 Sel is: atomselect40 Numbonds : 13 0 OG1 THR 67 HG1 THR 67 OD1 ASN 65 2.79 162.0 0 O WAT 4626 H2 WAT 4626 O1P RA3 135 2.67 176.1 0 OH TYR 98 HH TYR 98 O LYS 38 2.73 171.6 0 N ILE 86 H ILE 86 O ARG 97 2.76 169.0 0 NH1 ARG 1 1HH1 ARG 1 O PHE 5 2.73 166.4 0 NE2 GLN 14 2HE2 GLN 14 O TRP 10 2.75 162.9 0 NZ LYS 38 HZ3 LYS 38 OD1 ASN 41 2.80 163.6 0 N ALA 99 H ALA 99 O ASP 84 2.76 153.4 0 N TRP 10 H TRP 10 O THR 6 2.79 164.8 0 N CYX 23 H CYX 23 OD2 ASP 100 2.80 165.1 0 OG SER 74 HG SER 74 O PHE 76 2.70 176.8 0 NH1 ARG 77 1HH1 ARG 77 O ARG 104 2.77 157.6 0 N THR 46 H THR 46 O LEU 44 2.76 152.5 Processing frame: 1 Sel is: atomselect79 Numbonds : 0 Processing frame: 2 Sel is: atomselect79 Numbonds : 0 Processing frame: 3 Sel is: atomselect79 Numbonds : 0 Processing frame: 4 Sel is: atomselect79 Numbonds : 0 Processing frame: 5 Sel is: atomselect79 Numbonds : 0 Processing frame: 6 Sel is: atomselect79 Numbonds : 0 Processing frame: 7 Sel is: atomselect79 Numbonds : 0 Processing frame: 8 Sel is: atomselect79 Numbonds : 0 vmd > --- It baffles me as to where the problem is. Especially now that I have use update frame as well. -Sanjeev