Rendering CG bonds & constraints with VMD
- Last Updated: Wednesday, 16 August 2017 09:24
We modified and re-wrote (check this page) a small Tcl script to read and draw Martini bonds and constraints in VMD. The bonds of the elastic networks can be drawn as well. In some cases, the script needs to read the gromacs .tpr file. For this reason, there is a latest version for use with gromacs 5 (and later) .tpr files.
The script for use with gromacs 5 can be found here: cg_bonds-v5.tcl.
The script can be found here: cg_bonds.tcl (previous version: cg_bonds.tcl.old).
The script is invoked from the VMD command line - after sourcing it (source /path/to/cg_bonds.tcl, or /path/to/cg_bonds-v5.tcl, simply type cg_bonds (options and default values listed below). You can summon an "how to" at any moment by calling the cg_bonds_usage command.
This script uses the gmxdump executable compiled with gromacs. By default, it looks for it in the /usr/bin directory; you can specify another path with a specific option (see below). The script can now parse .top and .itp files, via the -top option.
Options and default values (slightly different from the previous version!):
|-molid||"top"||VMD-defined ID of the molecule to process|
|-gmx||/usr/bin/gmxdump||absolute path to gmxdump executable (for version 5, this should point to the gmx executable, default /usr/bin/gmx)|
|-tpr||topol.tpr||absolute path to the simulation file (.tpr)|
|-top||topol.top||absolute path to the system topology files (.top linking to .itp)|
|-topoltype||"martini"||protein topology type: "martini", "elastic" or "elnedyn"|
|-net||"martini"||network to draw: "martini", "elastic" or "both"|
|-bndrcnstr||"both"||draw bonds AND/OR constraints "bonds", "constraints" or "both"|
|-cutoff||7.0||cutoff for bonds (angstroms)|
|-color||"red"||color (color name or VMD-defined ID) of elastic bonds|
|-mat||"Opaque"||material for elnedyn bonds|
|-rad||0.2||radius of elastic bonds|
|-res||6||resolution of elastic bonds|
In most of the cases, if a "classical" cut-off is used for the elastic network (0.9nm), more than 12 bonds per atom are defined and VMD refuses to draw them. BUT this script replaces bonds by cylinders, and is able to draw the entire elastic network. Note that you have to modify the default cutoff value to see all the links defined by the elastic network.
Then you can remove Martini bonds with the cg_delete_martini_bonds; and the cylinders with the cg_delete_elastic_bonds. The only option for the previous three commands is the VMD-defined ID of the molecule to process (via -molid, default is "top").
- Using tpr files from version 5 and higher is not possible.Script for use with version 5 and higher is now available!
- The cg_bonds_infos command is useless (and going to be removed in next versions probably). Removed!
- Bonds drawn for the elastic network don't follow the dynamic of the protein; the 12-bonds-limitation of VMD obliges to draw them using cylinders, not regular bonds. Please drop a mail if you have a solution to fix that...
- We need a compiled version of gromacs around and the gmxdump executable to extract bonds/constraints (which takes time); wait for the next version! It will certainly read .top and .itp files.Added! The previous method (gmxdump) is still working; and watch out, flags have slightly changed.
- The table given in output (and the script in general) is not totally reliable in some cases; hand written topologies for examples...Removed.
- According to this topic (thanks Kargar), it seems the version of VMD should be superior to 1.8.5 (excluded) to run this script properly. The current version of VMD is 1.9 (august 2011).
- Not really an issue but: cg_bonds draws ALL bonds/constraints of the system; what if we want toRender the secondary structureand draw them only for some residues?Well... Use the good representation.
- Tcl needs an empty line at the end of the files it reads (.top and .itp)!
- On the same level, Tcl seems to understand differently the simple and double quotes... Please check that before shouting.
- Script doesn't extract (yet!) harmonic bonds of "regular" elastic networks built with the martinize.py script for instance (elastic bonds defined as type 6; this comment is not valid for ElNeDyn, which defines elastic bonds as type 1).Fixed!
user@machine $ vmd protein.gro
vmd > source /home/user/scripts/cg_bonds.tcl
vmd > cg_bonds -top system.top -topoltype "elastic"
vmd > cg_bonds -gmx /home/user/bin/gmx-4.5.4/bin/gmxdump -tpr dyn.tpr -net "elastic" -cutoff 12.0 -color "orange" -mat "AOChalky" -res 12 -rad 0.1 -topoltype "elastic"
vmd > cg_delete_elastic_bonds
The few lines above will:
- load the protein.gro structure,
- source the cg_bonds.tcl script,
- draw the martini bonds (using the .top file specified in the command line) extracted from a .top containingelastic network bonds,
- draw the elastic network with the options (cutoff, color, material, resolution and radius) specified by extracting the bonds from the .tpr file,
- delete the elastic network bonds previously drawn.
All these line will be applied to the molecule with the VMD-defined ID 0 (or "top").