VMD-L Mailing List

From: Axel Kohlmeyer (akohlmey_at_cmm.chem.upenn.edu)
Date: Tue Aug 21 2007 - 20:35:16 CDT

On Tue, 21 Aug 2007, J T wrote:

JT> I've also tried searching for only water oxygen contacts within a
JT> given cutoff distance to the atom selection and the problem I'm
JT> getting is that some of those oxygens are close to more than one atom
JT> from the specified atom selection. I've been playing with this for a

ouch... that is a problem that i didn't consider. nasty.
well, there are no winners in a game of losers. ;-)

JT> couple of days now. The following is a messy little hack just
JT> completed that seems to work. To note that only searching for water
JT> oxygens results in a much faster code. That will help since I plan
JT> on using this on each frame of large trajectories.
JT>
JT> After obtaining the combined list of atom indices with matching
JT> distances
JT>
JT> 1) First sort the list according to distance in decreasing order
JT> 2) Second sort the list according to residue (with -unique option)
JT> 3) Third sort the list according to distance again.

since performance is an issue, would suggest you make a
few timing tests with a few modifications suggested below.

JT>
JT> proc nwat {n sel} {
JT> set lists [measure contacts 5 [atomselect top "water and name
JT> OH2"] $sel]

here, i would check wether changing the order in
the selection makes a difference.

JT> set wlist [lindex $lists 0]
JT> set slist [lindex $lists 1]
JT> for {set i 0} {$i < [llength $wlist]} {incr i} {
JT> set satom [lindex $slist $i]
JT> set watom [lindex $wlist $i]

you can try writing the loops as:

       set wlist1 {}
       foreach satom $slist watom $wlist {

JT> set d [measure bond "$watom $satom"]
JT> lset wlist $i [list $watom $d]

... and here use:
         lappend wlist1 [list $watom $d]

JT> }
JT> set wlist2 [lsort -index 1 -decreasing $wlist]
JT> set wlist3 [lsort -unique -index 0 $wlist2]
JT> set wlist4 [lsort -index 1 $wlist3]

...and finally do this in one nested statement:

       set wlistX [lsort -index 1 [lsort -unique -index 0 [lsort -index 1 -decreasing $wlist1] ] ]

...

don't know if all of them help, but at least the change
of the loops plus the lappend should be faster since you
don't search with indexes.

cheers,
   axel.

JT> }
JT>
JT> Thanks,
JT> Jeff Tibbitt
JT> jtibbitt_at_odu.edu
JT>
JT> ------------------------------------------------------------------------
JT> ----------------------------------------------------
JT> On Aug 21, 2007, at 8:26 PM, Axel Kohlmeyer wrote:
JT>
JT> > On 8/21/07, Axel Kohlmeyer <akohlmey_at_cmm.chem.upenn.edu> wrote:
JT> >> On Tue, 21 Aug 2007, J T wrote:
JT> >>
JT> >> jeff,
JT> >>
JT> >> i suspect your problems arise from the
JT> >> fact, that you do distiguish between water
JT> >
JT> > arrggggh, typo. this should read:
JT> >
JT> > fact, that you do _not_ distiguish between water
JT> >
JT> > sorry,
JT> > axel.
JT> >
JT> >> oxygen and hydrogen and thus may have multiple
JT> >> matches with different distance.
JT> >>
JT> >> the simple workaround would be to search
JT> >> only for water oxygens. it is generally
JT> >> more consistent in many ways, too.
JT> >>
JT> >> also, i would like to be careful with the
JT> >> use of 'resid' this is not guaranteed to be
JT> >> unique (it is whatever VMD reads from the
JT> >> pdb/psf/whatever file).
JT> >>
JT> >> there are more elaborate ways to work around
JT> >> your problems, in case you absolutely need
JT> >> to include the hydrogens in your search.
JT> >>
JT> >> cheers,
JT> >> axel.
JT> >>
JT> >> JT> > > Try instead
JT> >> JT> > > set nwatlist [lsort -unique -index 1 $wlist]
JT> >> JT>
JT> >> JT> I've also tried that. Elements of the list may have the same
JT> >> residue
JT> >> JT> number, but the distance is still different, which renders the -
JT> >> JT> unique option void. I've even setting the elements of the
JT> >> list with:
JT> >> JT>
JT> >> JT> lset wlist $i [list $d $wresnum]
JT> >> JT>
JT> >> JT> instead of
JT> >> JT>
JT> >> JT> lset wlist $i "$d $wresnum"
JT> >> JT>
JT> >> JT>
JT> >> JT> Jeff Tibbitt
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> JT>
JT> >> ---------------------------------------------------------------------
JT> >> ---
JT> >> JT> ----------------------------------------------------
JT> >> JT> On Aug 21, 2007, at 5:34 PM, Nuno Loureiro Ferreira wrote:
JT> >> JT>
JT> >> JT> > J T wrote:
JT> >> JT> >>> > Instead of
JT> >> JT> >>> > set nwatlist [lsort $wlist]
JT> >> JT> >>> >
JT> >> JT> >>> > try,
JT> >> JT> >>> > set nwatlist [lsort -unique $wlist]
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >> Tried it. Since the list elements are double valued
JT> >> (distance and
JT> >> JT> >> residue number) the -unique does not work. Is a different
JT> >> way to
JT> >> JT> >> set the list up so the -unique option would filter out all
JT> >> double
JT> >> JT> >> residue numbers.
JT> >> JT> >
JT> >> JT> > Right ;-)
JT> >> JT> > Try instead
JT> >> JT> > set nwatlist [lsort -unique -index 1 $wlist]
JT> >> JT> >
JT> >> JT> > This way you will sort by uniqueness on the second element,
JT> >> in your
JT> >> JT> > case, the residue number.
JT> >> JT> >
JT> >> JT> >
JT> >> JT> >> Right now the ith list element is created with the command:
JT> >> JT> >>
JT> >> JT> >> lset wlist $i "$d $wres"
JT> >> JT> >>
JT> >> JT> >> where $d is the distance and $wres is the water residue
JT> >> number,
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >> On Aug 21, 2007, at 4:55 PM, Nuno Loureiro Ferreira wrote:
JT> >> JT> >>
JT> >> JT> >>> Hi JT
JT> >> JT> >>>
JT> >> JT> >>> Instead of
JT> >> JT> >>> set nwatlist [lsort $wlist]
JT> >> JT> >>>
JT> >> JT> >>> try,
JT> >> JT> >>> set nwatlist [lsort -unique $wlist]
JT> >> JT> >>> N.
JT> >> JT> >>>
JT> >> JT> >>> J T wrote:
JT> >> JT> >>>> Dear VMD Community,
JT> >> JT> >>>>
JT> >> JT> >>>> I've been working on a script that will select the N-closest
JT> >> JT> >>>> waters to a specified atom selection by following the
JT> >> algorithm
JT> >> JT> >>>> suggested by John Stone in an earlier post. It
JT> >> successfully
JT> >> JT> >>>> returns a sorted list of increasing distances, but there are
JT> >> JT> >>>> duplicate water residues. The script generates a list of
JT> >> all
JT> >> JT> >>>> water indices within a max cutoff distance. Then each water
JT> >> JT> >>>> index in that list is replaced with a double valued element
JT> >> JT> >>>> containing the distance and the water residue. Is there
JT> >> a way I
JT> >> JT> >>>> can use the -unique option somehow to fix this?
JT> >> JT> >>>>
JT> >> JT> >>>> Thank-you for reading,
JT> >> JT> >>>> Jeff Tibbitt
JT> >> JT> >>>> jtibbitt_at_odu.edu <mailto:jtibbitt_at_odu.edu>
JT> >> JT> >>>>
JT> >> JT> >>>>
JT> >> JT> >>>> John Stone's algorithm:
JT> >> JT> >>>> / 1) select all waters within the max cutoff distance M /
JT> >> JT> >>>> / 2) calculate the distance D of each water molecule and
JT> >> add the
JT> >> JT> >>>> index and //distance to a list or lists /
JT> >> JT> >>>> / 3) sort the list(s) by the distance (keeping the water
JT> >> index
JT> >> JT> >>>> and distance //assocation intact, if they are in separate
JT> >> JT> >>>> lists..) /
JT> >> JT> >>>> / 4) select the closest N waters from the sorted list// /
JT> >> JT> >>>> /
JT> >> JT> >>>> /
JT> >> JT> >>>> /
JT> >> JT> >>>> /
JT> >> JT> >>>> /My Script:/
JT> >> JT> >>>> proc nwat {n sel} {
JT> >> JT> >>>> set lists [measure contacts 6 [atomselect top water] $sel]
JT> >> JT> >>>> set wlist [lindex $lists 0]
JT> >> JT> >>>> set slist [lindex $lists 1]
JT> >> JT> >>>> set n [llength $wlist]
JT> >> JT> >>>> for {set i 0} {$i < $n} {incr i} {
JT> >> JT> >>>> set satom [lindex $slist $i]
JT> >> JT> >>>> set watom [lindex $wlist $i]
JT> >> JT> >>>> set wres [[atomselect top "index $watom"] get resid]
JT> >> JT> >>>> set d [measure bond "$watom $satom"]
JT> >> JT> >>>> lset wlist $i "$d $wres"
JT> >> JT> >>>> } set nwatlist [lsort $wlist]
JT> >> JT> >>>> }
JT> >> JT> >>>>
JT> >> -------------------------------------------------------------------
JT> >> JT> >>>> -----
JT> >> JT> >>>>
JT> >> JT> >>>> No virus found in this incoming message.
JT> >> JT> >>>> Checked by AVG Free Edition. Version: 7.5.484 / Virus
JT> >> Database:
JT> >> JT> >>>> 269.12.0/961 - Release Date: 19-08-2007 7:27
JT> >> JT> >>>>
JT> >> JT> >>>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >>
JT> >> JT> >
JT> >> JT>
JT> >>
JT> >> --
JT> >> =====================================================================
JT> >> ==
JT> >> Axel Kohlmeyer akohlmey_at_cmm.chem.upenn.edu http://
JT> >> www.cmm.upenn.edu
JT> >> Center for Molecular Modeling -- University of Pennsylvania
JT> >> Department of Chemistry, 231 S.34th Street, Philadelphia, PA
JT> >> 19104-6323
JT> >> tel: 1-215-898-1582, fax: 1-215-573-6233, office-tel:
JT> >> 1-215-898-5425
JT> >> =====================================================================
JT> >> ==
JT> >> If you make something idiot-proof, the universe creates a better
JT> >> idiot.
JT> >>
JT> >>
JT> >
JT> >
JT> > --
JT> > ======================================================================
JT> > =
JT> > Axel Kohlmeyer akohlmey_at_cmm.chem.upenn.edu http://
JT> > www.cmm.upenn.edu
JT> > Center for Molecular Modeling -- University of Pennsylvania
JT> > Department of Chemistry, 231 S.34th Street, Philadelphia, PA
JT> > 19104-6323
JT> > tel: 1-215-898-1582, fax: 1-215-573-6233, office-tel: 1-215-898-5425
JT> > ======================================================================
JT> > =
JT> > If you make something idiot-proof, the universe creates a better
JT> > idiot.
JT> 

-- 
=======================================================================
Axel Kohlmeyer   akohlmey_at_cmm.chem.upenn.edu   http://www.cmm.upenn.edu
   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.