VMD-L Mailing List

From: Jeff Ulrich (ulrich_at_ks.uiuc.edu)
Date: Wed Apr 30 1997 - 17:42:04 CDT

This week's VMD TECH, by Jeff Ulrich, will discuss additional
applications of the Tcl 'trace' command. A few weeks ago I hinted to
the potential for traces on vmd_pick_atom and vmd_pick_mol. This
article will be similar in spirit, though it will focus on improving
animation playback by means of the VMD globals "vmd_frame" and
"vmd_timestep." Along the way I will also discuss a few concepts
related to caching of structural information in VMD.

For the new people who just joined the list, VMDTech is a bi-weekly
article we write for vmd-l concerning some detail or topic of VMD that
we believe will be helpful to users. We plan to include much of
this material in the new version of the documentation, so please tell
us if there are any problems. Back issues of VMDTech are available
from the vmd-l mailing list archives, which can be found at
        http://www.ks.uiuc.edu/Research/vmd/mailing_list/ .

Obligatory Review
-----------------
The syntax of setting up a trace is:
        trace variable <variable-name> <operation> <procedure-name>
where <variable_name> is the variable you wish to monitor, <operation>
is 'w', 'r', or 'u' for read, write, or unset, and <procedure-name>
is the code to call when <operation> is performed on <variable-name>.

To delete a trace on a variable, use the syntax:
        trace vdelete <variable-name> <operation> <procedure-name>
and to look at existing traces on a variable, use:
        trace vinfo <variable-name>

For the remainder of this article, all traces will be explicitly
created, but I will leave trace deletions to you when they are
necessary.

Example 1 - Secondary Structure Calculations
--------------------------------------------
As many of you probably know, during animation playback secondary
structure information is calculated for the first simulation frame
only. Additional frames merely "wiggle around" the cylinders and
triangles making up the helices and beta sheets of the initial cartoon
representation. If you want VMD to actually correct for changes which
occur in the secondary structure as the simulation progresses, then it
is necessary to force additional SS calculations. You can easily do
this using the trace command. The VMD variable "vmd_frame" is a
global array which at all times contains the indices of the current
frame being displayed for each molecule. Since VMD automatically
writes to this variable upon any change in the frame number, then a
trace on this write operation will serve our purpose. To see this,
load up a molecule and its corresponding DCD or Amber file. Make a
cartoon representation of it and then cut and paste the lines below.

        trace variable vmd_frame w structure_trace
        proc structure_trace {name index op} {
                vmd_calculate_structure $index
        }

For those of you who read about traces a few weeks ago, this short
procedure is probably self-explanatory. As a quick review though, the
first line instructs VMD to call the procedure structure_trace each
time the variable vmd_frame is written to (hence the "w"). The body
of structure_trace merely calls the function that VMD uses to
calculate and display the secondary structure for the current
timestep. If you now playback your animation, the secondary structure
will be updated correctly.

These four lines will suit our purpose, but we are far from an optimal
solution. There are at least two drawbacks to our method. First of
all, the secondary structure code is slow, with significant overhead
arising from the writing of a temporary PDB file to disk. Secondly, if
you have multiple molecules loaded, structure_trace will recalculate
the SS for ALL of those molecules, whether that was your intention or
not. It would be nice if we could somehow save the secondary
structure information for every frame that we play back. Successive
simulations would then proceed much more quickly, since updating the
structure would merely be a matter of redrawing helices and triangles
on the basis of pre-computed data. It turns out that this type of
"caching" is relatively easy to do in VMD. Code to accomplish this
has been available in our script library for quite some time now. It
will be explained in more detail below.

        # start the cache for a given molecule
        proc start_sscache {{molid top}} {
            global sscache_data
            if {! [string compare $molid top]} {
                set molid [molinfo top]
            }
            global vmd_frame
            # set a trace to detect when an animation frame changes
            trace variable vmd_frame($molid) w sscache
            return
        }

        # remove the trace (need one stop for every start)
        proc stop_sscache {{molid top}} {
            if {! [string compare $molid top]} {
                set molid [molinfo top]
            }
           global vmd_frame
           trace vdelete vmd_frame($molid) w sscache
           return
        }

        # reset the whole secondary structure data cache
        proc reset_sscache {} {
          if [info exists sscache_data] {
              unset sscache_data
          }
         return
        }

        # when the frame changes, trace calls this function
        proc sscache {name index op} {
         # name == vmd_frame
         # index == molecule id of the newly changed frame
         # op == w
    
         global sscache_data

         # get the protein CA atoms
         set sel [atomselect $index "protein name CA"]

         ## get the new frame number
         # Tcl doesn't yet have it, but VMD does ...
         set frame [molinfo $index get frame]

         # see if the ss data exists in the cache
          if [info exists sscache_data($index,$frame)] {
              $sel set structure $sscache_data($index,$frame)
              return
          }

            # doesn't exist, so (re)calculate it
            vmd_calculate_structure $index
            # save the data for next time
            set sscache_data($index,$frame) [$sel get structure]

         return
        }

The procedures start_sscache and stop_sscache merely create or delete
traces on the variable vmd_frame. The syntax of this "bookwork" was
described in an earlier VMD TECH and will not be discussed here. A
particular molecule can be targeted by passing an id to these
procedures (the default is the top molecule, as usual.) The procedure
"sscache" is the real heart of the code. It checks if a secondary
structure definition for the given molecule number and frame already
exists in the Tcl array sscache_data(molecule, frame). If so, it uses
that data to redefine the "structure" keyword values. If not, then
the routine for evaluating the secondary structure is called as usual,
but the results are saved in the sscache_data array for future
reference. Once the secondary structure values are saved, the
molecule can be quickly animated and the updates can be controlled via
the animate form. Finally, the command reset_sscache resets the saved
secondary structure data for all of the molecules and all the frames.

Example 2 -- Writing Simulation-Specific Text to the Display
------------------------------------------------------------
People in our group have recently found that VMD can be used to make
informative movies for display on web sites or at conferences.
Details of the movie-making process can be found in an earlier VMD
TECH. However, unless one has access to a sophisticated movie editor
(i.e. like that which ships with IRIX 6.3), it can often seem
difficult to add text and graphics which point out important events
occuring during simulations. For instance, most movie frames are made
from screen captures of the VMD display. It would be nice to have
phrases such as "Bond rupture at 80 ps" or "Equilibration ended at 20
ps" present in the display in order to make the movie output
descriptive. A trace on vmd_frame can provide an easy way of doing
this. The code below is a simple example. It writes a simple message
to the display indicating the frame number being shown.

        trace variable vmd_frame([molinfo top]) w output_messages

        proc output_messages {args} {
                global vmd_frame
                #Deletes previous text message
                draw delete all
                draw text {0 0 0} "Timestep #$vmd_frame([molinfo top])"
        }

Of course, this output can be made more informative at the discretion
of the user. For example, suppose we are looking at a 100 frame
simulation of retinal being pulled through a window in the helices of
bacteriorhodopsin. Interesting events in this simulation might
include the role of water in the retinal rupture, the bond breakage
itself, and the fact that the window seems to attract the retinal.
The following version of output_messages will emphasize these events
with text output.
        
        proc output_messages {args} {
                global vmd_frame
                #Deletes previous text message
                draw delete all
                #Output new message corresponding to frame ID
                if {$vmd_frame(1) < 32} {
                        draw text {0 0 0} "Water encourages bond breakage"
                } elseif {$vmd_frame(1) < 60} {
                        draw text {0 0 0} "Bond breakage occurs"
                } else {
                        draw text {0 0 0} "Helix window is attractive"
                }
        }

With these text messages, the results of the simulation become much
clearer to those experiencing the playback for the first time. Frame
dependent arrows, etc. follow naturally from this example.
        

Example 3 -- Updating Representations Such as Diffusing Waters
--------------------------------------------------------------
Often one would like to alter representations in response to changes
which occur in biological systems over the course of animations. For
example, in investigating the role of diffusing waters near a flexible
side chain, it might be nice to emphasize water molecules that at any
one time are within three Angstroms of the side chain atoms. One
might do this by creating two representations in VMD. One would draw
all water molecules by lines, and a second would be a licorice
representation of those waters within three Angstroms of the target
side chain. In VMD, once these representations are defined for any
given frame, they are not updated for subsequent animation frames.
That is, if a water molecule is originally two Angstroms from the area
of interest, then it as drawn as licorice throughout the playback,
even if it ends up drifting twenty Angstroms away. In practice, it
would be nice to see waters which drift out of the three Angstrom
bubble switch from licorice to lines and vice versa. As you probably
expected, this can be done with a trace. The following example
examines waters near retinal in my hypothetical
bacteriorhodopsin/retinal system. You should change the "mol
selection" lines to suit your own purpose.

        trace variable vmd_frame([molinfo top]) w identify_close_waters

        proc identify_close_waters {args} {
                display update off
                set num [molinfo top get numreps]
                # Set up identification. I assume that the
                # water representations are the last two in
                # the list, though may be differet for your system.
                set repID_1 [expr ($num - 1)]
                set repID_2 [expr ($repID_1 - 1)]
                set molID [molinfo top]

                # Now update the representations to reflect new water
                # positions.
                mol representation Lines
                mol selection \
                  "waters and (same residue as (not within 3 of segname RET))"
                mol modrep $repID_1 $molID
                mol representation Licorice
                mol selection \
                  "waters and (same residue as (within 3 of segname RET))"
                mol modrep $repID_2 $molID
                display update on
        }

The procedure identify_close_waters updates the licorice and line
representations of water whenever a new animation frame is loaded. It
does this by recalculating distances before assigning a representation
style with the VMD selection language syntax. This method will be
slow for systems with large amounts of water, so again one might
envision caching these atom selections over time as was done with
sscache. Successive simulations would then not always have to
recalculate the positons of the water (a costly process). They would
just use indices to the relevant atoms at each time step. To do this
the code might look something like:

        trace variable vmd_frame([molinfo top]) w identify_close_waters

        proc identify_close_waters {args} {

                global far_water_data
                global near_water_data
                global vmd_frame
                display update off
                set num [molinfo top get numreps]
                # Set up identification. I assume that the
                # water representations are the last two in
                # the list, though this may be differet for your system.
                set repID_1 [expr ($num - 1)]
                set repID_2 [expr ($repID_1 - 1)]
                set molID [molinfo top]

                # Check if the waters within/beyond 3A have been
                # calculated for this frame before. If so, use the
                # results. If not, then do the calculation and store
                # the data so that next time this frame is played, we
                # need only to read in the relevant indices.

                if [info exists far_water_data($vmd_frame($molID))] {
                        set far_waters $far_water_data($vmd_frame($molID))
                        set near_waters $near_water_data($vmd_frame($molID))
                } else {
                        set far_waters [atomselect $molID \
                        "waters and (same residue as (not within 3 \
                        of segname RET))"]
                        set near_waters [atomselect $molID \
                        "waters and (same residue as (within 3 of \
                        segname RET))"]
                        set far_water_data($vmd_frame($molID)) $far_waters
                        set near_water_data($vmd_frame($molID)) $near_waters
                }

                # Form the selection string based on a given selection.
                # Yes this is a kluge but I didn't find another way.
                set sel_text "index "
                foreach i [$far_waters get index] {
                        set sel_text "$sel_text $i"
                }
                mol representation Lines
                mol selection $sel_text
                mol modrep $repID_1 $molID
                set sel_text "index "
                foreach i [$near_waters get index] {
                        set sel_text "$sel_text $i"
                }
                mol representation Licorice
                mol selection $sel_text
                mol modrep $repID_2 $molID
                display update on
        }

This got rather messy. I apologize that I am not yet the expert in
Tcl that Andrew and probably others of you are.

Example 4 -- vmd_timestep example
---------------------------------
Before signing off, I should mention at least one more VMD variable
which may be of interest to those who use VMD as a graphical front end
to molecular dynamics programs such as NAMD. This variable,
vmd_timestep, is updated whenever a new simulation timestep arrives to
VMD. As an example of its use, consider the following procedure

        trace variable vmd_timestep w timestep_temp_trace
        proc timestep_temp_trace {name id op} {
                puts "Temperature of molecule $id is $temp"
        }

This simple example prints the temperature of each incoming simulation
timestep.

Good luck,
Jeff Ulrich (ulrich_at_ks.uiuc.edu)