From: John Stone (johns_at_ks.uiuc.edu)
Date: Fri Apr 30 2004 - 16:02:12 CDT

Ioana,
1) I don't understand what you were trying to accomplish with the line
'$fa global', that's clearly not valid Tcl syntax, nor is there
a 'global' operation defined for Tcl channels, at least not that
I'm aware of. The man page for the "global" keyword in Tcl is here:
  http://www.tcl.tk/man/tcl8.4/TclCmd/global.htm
You would have wanted to declare 'global fa' somewhere before
you call 'set fa [open DistBDNA.dat a]'.

If you don't have a Tcl/Tk book, I highly recommend getting one, it'll
explain all of this stuff much better than the Tcl man pages do,
There are lots of Tcl tutorials online however:
  http://resource.tcl.tk/resource/doc/start/

Regarding 2), it sounds to me like you may need to change the
"waitfor" option in the animate read command to "waitfor all"
since you're running your script in a batch mode. I suspect that
your script is terminating before it has processed all of the frames
as a result of that setup.

Hope that helps,
  John Stone
  vmd_at_ks.uiuc.edu

On Thu, Apr 29, 2004 at 09:09:49PM -0700, Ioana Cozmuta wrote:
> Hi VMD-users,
>
> I am using the script below to calculate the average distance between
> bases in a DNA strand during a MD simulation.
> I am trying to make this script work in a pbs shell via
> vmd -dispdev text -e script.vmd
>
>
> 1. At the moment I am defining the output file inside the procedure
> myBaseDist. I've tried to define the name of the output file in the main
> body of the script
> set fa [open name.file w]
> $fa global
>
> and then inside the procedure use
> global $fa
>
> but vmd does not seem to recognize the name of the file in the procedure
> this way.
> What is the correct way to make this variable global?
>
> 2. Although my test trajectory file has 10 frames, the output contains
> only 3 frames.
>
> Could anyone please help me figure out why (without getting any errors)
> the bigdcd procedure only processes 3 frames out of 10?
>
>
> Thank you,
> Ioana
>
> ************************************************
>
> #This script calculates the distance between bases of a DNA strand
>
> proc bigdcd { script args } {
> global bigdcd_frame bigdcd_proc bigdcd_firstframe vmd_frame
>
> set bigdcd_frame 0
> set bigdcd_firstframe [molinfo top get numframes]
> set bigdcd_proc $script
>
> uplevel #0 trace variable vmd_frame w bigdcd_callback
> foreach dcd $args {
> animate read dcd $dcd waitfor 0
> }
> }
>
> proc bigdcd_callback { name1 name2 op } {
> global bigdcd_frame bigdcd_proc bigdcd_firstframe vmd_frame
>
> # If we're out of frames, we're also done
> set thisframe $vmd_frame($name2)
> if { $thisframe < $bigdcd_firstframe } {
> bigdcd_done
> return
> }
>
> incr bigdcd_frame
> if { [catch {uplevel #0 $bigdcd_proc $bigdcd_frame} msg] } {
> puts stderr "bigdcd aborting at frame $bigdcd_frame\n$msg"
> bigdcd_done
> return
> }
> animate delete beg $thisframe end $thisframe
> return $msg
> }
>
> proc bigdcd_done { } {
> puts "bigdcd_done"
> after idle uplevel #0 trace vdelete vmd_frame w bigdcd_callback
> }
>
> proc dist_CM {select1 select2} {
>
> if {[$select1 num] <= 0 || [$select2 num] <= 0 } {
> error "center_of_mass: needs a selection with atoms"
> }
>
> set cm1 [veczero]
> set mass1 0
> set cm2 [veczero]
> set mass2 0
>
> foreach crd1 [$select1 get {x y z}] m1 [$select1 get mass] {
> set mass1 [expr $mass1 + $m1]
> set cm1 [vecadd $cm1 [vecscale $m1 $crd1]]
> }
>
> foreach crd2 [$select2 get {x y z}] m2 [$select2 get mass] {
> set mass2 [expr $mass2 + $m2]
> set cm2 [vecadd $cm2 [vecscale $m2 $crd2]]
> }
>
> if {$mass1 == 0 || $mass2 == 0} {
> error "center_of_mass: total mass is zero"
> }
>
> return [vecdist [vecscale [expr 1.0/$mass1] $cm1] [vecscale [expr 1.0/$mass2] $cm2]]
>
> }
>
> proc myBaseDist { frame } {
> #global fa
> set fa [open DistBDNA.dat a]
> set residues [lsort -integer -unique [[atomselect top "resid 1 to 19 or resid 21 to 39"] get resid]]
> set avg 0
> set i 1
> foreach residue $residues {
> set sel1 [atomselect top "resid $residue and not backbone"]
> set next [expr $residue + 1]
> set sel2 [atomselect top "resid $next and not backbone"]
> set dist [dist_CM $sel1 $sel2]
> set avg [expr $avg +$dist]
> set i [expr $i +1]
> }
> set N $i
> set avg [expr $avg/$N]
> puts $fa "$frame: $avg"
> close $fa
> }
>
> mol new {./m20AT/Eq20ATmod.restrt.coor} type {pdb}
> # set fa [open DistBDNA.dat a]
> # $fa global
> bigdcd myBaseDist test.dcd
>
> wait 1200
> exit

-- 
NIH Resource for Macromolecular Modeling and Bioinformatics
Beckman Institute for Advanced Science and Technology
University of Illinois, 405 N. Mathews Ave, Urbana, IL 61801
Email: johns_at_ks.uiuc.edu                 Phone: 217-244-3349              
  WWW: http://www.ks.uiuc.edu/~johns/      Fax: 217-244-6078