From: Axel Kohlmeyer (
Date: Mon Aug 13 2007 - 19:08:30 CDT

On Mon, 13 Aug 2007, Jim Pfaendtner wrote:

JP> Dear VMD list,

dear jim,

JP> I am still trying to figure out why my badwater script is not
JP> working. Here is the script I am currently using:

is there any prolog that loads the psfgen plugin,
defines the topology file etc.?

check out the namd/psfgen tutorial at:

JP> set overlap [atomselect top "segid WT1 to WT12 and within 2.4 of (not
JP> segid WT1 to WT12)"]

i would make this selection:
"same residue as exwithin 2.4 of (not segid WT1 to WT12)"

i guess you want to delete whole water molecules and not just atoms...

JP> set reslist [$overlap get resid]
JP> set reslist [lsort -unique -integer $reslist]
JP> set seglist [$overlap get segid]
JP> foreach segid $seglist resid $reslist { delatom $segid $resid }
JP> writepdb tt2.pdb
JP> writepsf tt2.psf
JP> My pdb file is a solvated protein that was solvated using the VMD
JP> plugin, so the segids for the water groups are WT1 through WT12.
JP> When I run this script on my protein I dont get any errors, but the
JP> "writepsf / writepdb" commands just write empty files - it's as
JP> though my selection is grabbing all the atoms.

for delatom/writepsf/writepdb to work you'll first have to load
the psfplugin, split the segments and tell it to find the segement
coordinates etc. and then work on it.

... or try autopsf and give it the selection of
same residue as within 2.4 of (not segid WT1 to WT12)...


JP> I adapted this script from the psfgen usermanual and a previous post
JP> on the namd userlist about badwater scripts. John previously replied
JP> to this and asked if I was sure that I was getting the right
JP> selection - I am reasonably sure that I am.

JP> Thank you,
JP> Jim
JP> On Jul 16, 2007, at 3:19 PM, Jim Pfaendtner wrote:
JP> > hello,
JP> >
JP> > I am a new VMD user. I'm trying to write/adapt a TCL script to
JP> > select overlapping waters and delete them and I'm getting something
JP> > wrong.
JP> >
JP> > Here is my script:
JP> >
JP> > set badwat [atomselect top "name OH2 and within 3.0 of protein"]
JP> > foreach segid [$badwat get segid] resid [$badwat get resid]
JP> > {
JP> > delatom $segid $resid
JP> > }
JP> >
JP> > I am running it in the vmd console after loading my psf/pdb files.
JP> > The output just says:
JP> >
JP> > no segment WT1
JP> >
JP> > thanks,
JP> > Jim
JP> >

Axel Kohlmeyer
   Center for Molecular Modeling   --   University of Pennsylvania
Department of Chemistry, 231 S.34th Street, Philadelphia, PA 19104-6323
tel: 1-215-898-1582,  fax: 1-215-573-6233,  office-tel: 1-215-898-5425
If you make something idiot-proof, the universe creates a better idiot.