From: HongTham (hongtham0709_at_gmail.com)
Date: Thu Apr 24 2014 - 19:51:51 CDT

Hi VMD users,
I have some problem in calculating Hydrogen bond and orient that I cant
resolve them. Thanks in advance for any advices

1. Measure Hbond:
I made the Tcl script below. although in visual simulation there are some
hbonds, and by TC console as well, the file out just display {} {} {}. I
set cutoff 3.4 and angle 30.

mol new ../proion.pdb
animate delete all
mol addfile ../pro12_16.dcd waitfor all
mol addfile ../pro17_31.dcd waitfor all
mol addfile ../pro32_46.dcd waitfor all
mol addfile ../pro47_61.dcd waitfor all
mol addfile ../pro62_76.dcd waitfor all
mol addfile ../pro77_89.dcd waitfor all
mol addfile ../pro90_99.dcd waitfor all

source /home/hongtham/VMD/Nhbond.tcl

set nf [molinfo 0 get numframes]
set clist [list A C D]

foreach ch $clist {
    set fout [open "hbond$ch.dat" w]
    set sel [atomselect 0 {protein and chain $ch and resid 180 to 200 and
name "N.*" "O.*"}]

#calculation
    for { set i 0 } { $i < $nf } { incr i } {
        $sel frame $i

        puts $fout "[expr ($i*0.01)+20.1] [Nhbond $sel]"
    }
    close $fout
}
~
the source script is below:
# Count the number of H-bonds with one selection

proc Nhbond {sel {cutoff 3.4} {angle 30}} {
    return [measure hbonds 3.4 30 $sel]
}

2. measure orient
with the same selected molecule. I made the script below, i try to find
some explaination about the script but file out has nothing. Thank for any
suggestion

# analysis
    set clist [list A B C D E]
    set vec {0 0 1}
#set clist [list A ]
   # set fout [open "ori$ch.dat" w]

foreach ch $clist {
    set fout [open "ori$ch.dat" w]
    set sel [atomselect 0 "protein and chain $ch and name CA"]
    set fout [open "ori$ch.dat" w]

# distance calculaton
      for { set i 0 } { $i < $nf } { incr i } {
         $sel frame $i

 set ori [orient $sel $vec]
        puts $fout "$ori"
    }
    close $fout

and the source script is:
proc orient { sel {vec {0 0 1} }} {
    set all [atomselect [$sel molid] "all"]
    set vec [vecnorm $vec]

    set vec1 [lindex [measure inertia $sel] 1 2]
    set vec2 [vecnorm [veccross $vec1 $vec]]
    set ang [expr acos([vecdot $vec1 $vec])]
    $all move [trans center [measure center $sel] axis $vec2 $ang rad]
}

Thank you so much.
Hong Tham

Biomedical engineering department
Pukyong National University
Pusan, South Korea
email: hongtham0709_at_gmail.com