From: Peter Freddolino (petefred_at_ks.uiuc.edu)
Date: Mon Mar 26 2007 - 16:00:10 CDT

Hi Alexandra,
we might need some more detail on what you mean by it not working
properly. However, the measure hbonds command requires the presence of
explicit hydrogens, so including hydrogen atoms in the selections you
use for calculations would probably be a good start.

Also, you can create your selections outside of the loop, update them
with $sel frame $i, and only delete them after the loop, rather than
making new ones for each frame, which will be slower.
Peter

alexandra.marques_at_fc.up.pt wrote:
> Hi
>
> I am trying to calculate the number of hydrogen bonds of a protein during an MD
> simulation but my script don’t work properly. Could someone please help me?
> Here is the script:
>
> set outfile [open "hbond.out" w]
> set nf [molinfo top get numframes]
> set allsel [atomselect top "all and name CA"]
> set residlist [lsort -unique -integer [$allsel get residue]]
>
> for {set i 0} {$i < $nf} {incr i} {
>
> foreach res $residlist {
> set sel1 [atomselect top "residue $res and name \"N.*\" \"O.*\" \"S.*\""]
> set sel2 [atomselect top "protein and not residue $res and name \"N.*\"
> \"O.*\" \"S.*\""]
> set hbonds [measure hbonds 3.5 60 $sel1 $sel2]
> $sel1 frame $i
> $sel2 frame $i
> $sel1 update
> $sel2 update
> $sel1 delete
> $sel2 delete
> puts $outfile "number of h-bonds in frame $i: $hbonds"
> }
> }
> close $outfile
>
> Thanks
> Alexandra
>
> -------------------------------------------------------------
> A FCUP utiliza o sistema de webmail Horde/IMP (www.horde.org)
>
> Visite: http://www.fc.up.pt/
>
>