From: Balazs Jojart (jojartb_at_gmail.com)
Date: Sun Jan 17 2010 - 13:20:20 CST

Dear Axel,
Thank you for your help!
Balazs

Axel Kohlmeyer wrote:
> On Sun, 2010-01-17 at 12:11 +0100, Balazs Jojart wrote:
>
>> Dear Users,
>> I have a problem during trajectory evaluation.
>> A would like to measure distance between two atoms in the trajectory.
>> I know, that without deleting the selection, variables I can shock my
>> computer.
>>
>
> well, currently you are only shocking people
> reading your script because of the number of
> "don't" that you have committed. there are
> quite a few good template examples that explain
> how to do this kind of thing, and you seem to
> have started from a pretty bad one.
>
>
>> Here is a part of my tcl script:
>>
>> mol new {../1.00M_GPCR_ions.pdb} type {pdb} first 0 last -1 step 1 waitfor 1
>> mol addfile {../1.00_PROD3_CASE1.dcd} type {dcd} first -1 last -1 step
>> 1 waitfor -1 0
>>
>
> what is the 0 at the end of this line for?
>
>
>> set outfile1 [open 1.00_PROD1_Asp52OD1_NA_min.txt w]
>> set nf [molinfo top get numframes]
>>
>
>
>> for {set i 1} {$i < $nf} {incr i} {
>> puts "frame $i of $nf"
>> set selatoms1 [[atomselect top "protein and resid 52 and name OD1 or
>> resname NA and resid 1"] get index]
>>
>
> ouch!!! this will produce a memory leak in _any_ case.
> the reason being, that the atomselect call returns the
> name of the new atomselect### command being created,
> but you don't store it and rather immediately use it.
> thus $selatoms contains only the indices of the atoms
> that you are interested in.
>
> what you have to do instead, is assign the return value
> of the atomselect call to a variable, and then generate
> the list of indices from that and assign it to a different
> variable.
>
>
>> set bond1 [list [lindex $selatoms1 0] [lindex $selatoms1 1]]
>>
>
> this is a pointless operation. $selatoms1 has
> the exact same content as $bond1 will have.
>
>
>> set d1 [measure bond $bond1 frame $i]
>>
>
>
>> $selatoms1 delete
>>
>
> this is where you get the error. since $selatoms contains a
> list of atom indices and not the name of an atomselect###
> command, it will try to execute those numbers, but there is
> no command of that name.
>
>
>> $bond1 delete
>>
>
> this will give you the same error. have a proper look at a
> tcl documentation or tutorial or both. variables are not
> deleted this way, but rather using the unset command. it
> is not really needed here anyways.
>
>
>> set selatoms2 [[atomselect top "protein and resid 52 and name OD1 or
>> resname NA and resid 2"] get index]
>> set bond2 [list [lindex $selatoms2 0] [lindex $selatoms2 1]]
>> set d2 [measure bond $bond2 frame $i]
>> ........
>> }
>>
>> Unfortunatelly, I obtain the following error:
>>
>
> this is not unfortunate but deserved. there are a bunch of
> issues. first of all, there is no point in creating the
> atom selections in every step. it is much better to create
> them once outside the loop and then just move them ahead
> or update them. since you want to refer to fixed atoms,
> there is not even a need to redo the selection. the order
> of their indexes will never change. so place the following
> _outside_ the loop.
>
> set sel1 [atomselect top "protein and resid 52 and name OD1 or resname
> NA and resid 1"]
> set bond1 [$sel1 list]
>
> set sel2 [atomselect top "protein and resid 52 and name OD1 or resname
> NA and resid 2"]
> set bond2 [$sel2 list]
>
> $sel1 delete
> $sel2 delete
>
> and then you can use $bond1 and $bond2 with the measure command, just
> as you did.
>
> also, you should rethink the syntax of your selection text. it may not
> do exactly what you expect. i would place parenthesis to define proper
> precedence. even if it currently works, it may fail in unexpected ways
> later on.
>
>
>>> frame 1 of 7151
>>> invalid command name "790 4926"
>>> Info) VMD for LINUXAMD64, version 1.8.7 (August 1, 2009)
>>> Info) Exiting normally.
>>>
>> Where is the problem in my script, and how can I solve it?
>>
>
> there are multiple problems, and the best solution is to study
> existing good scripts. there are examples in the VMD user's guide,
> tutorials and in the script library. i also suggest checking out
> the tcl tutorial.
>
> cheers,
> axel.
>
>> Any help will be appreciated!
>> Thank you for your help in advance,
>> Balazs
>>
>>
>
>