On Jan 21, 2010, at 3:15 PM, Yasser Almeida Hernández wrote:
Hi... I'm performing a massive superposing of aminoacids in different conformations. The first model is a arginine "template" reference structure that is the fixed model in all matching iterations. I've loaded 236 models and i want to match the arginine guanidine group on the template guanidine. In each iteration i execute two matching, with different atoms selections. The first selection is the "right" selection: CD,NE,CZ,NH1,NH2 and the second is the "inverse" selection: CD,NE,CZ,NH2,NH1. I run the two superpositions, and i want to show that with the minimal RMSD between the guanidine groups. This is the code:
Well, I can see a variety of problems. Starting from the top...
import chimera from chimera import runCommand
right = "@CD,NE,CZ,NH1,NH2" inverse = "@CD,NE,CZ,NH2,NH1"
for i in range(2,328):
'i' will range from 2 to 327, yet you say you've loaded 236 models. I'm guessing you actually want 'i' to vary from 1 to 235 (or 236 if the "236 models" you mention don't include the template model). Therefore you should be using "range(1, 236)".
print "Superposing..." runCommand("match #"+str(i)+right+" #0"+right+"; close 1")
What is with all the "close 1" on the end of all your commands? Do you mean "wait 1" (to draw a new frame)? "Close 1" will close model 1. It makes no sense to do it over and over. Also, since you are not saving any images in this script you don't need "wait" commands either.
right_rmsd = runCommand("rmsd #"+str(i)+right+" #0"+right+"; close 1")
'runCommand' doesn't return a value. Actually it does, but it's whether or not a 'wait' is needed after the command string to update the display, not a value related to the commands' functions. It doesn't return a command-specific value pretty much for the reason exemplified by the string you used: it is two commands strung together, so which command's value should be returned? At any rate, to get the RMSD you will need to call the underlying rmsd function directly. It is in the Midas module. So this: from Midas import rmsd right_rmsd = rmsd("#"+str(i)+right, "#0"+right)
runCommand("match #"+str(i)+right+" #0"+right+"; close 1") inverse = runCommand("rmsd #"+str(i)+inverse+" #0"+right+"; close 1")
You need to assign the value (once you change to the rmsd command) to inverse_rmsd, not inverse. After the above statement, inverse no longer has the value "@CD,NE,CZ,NH2,NH1", so a lot of things will begin to go badly wrong at this point.
if right < inverse:
if right_rmsd < inverse_rmsd:
print "Minimum RMSD =",right,"\n\n"
right_rmsd, not 'right'.
runCommand("match #"+str(i)+right+" #0"+right+"; close 1")
elif right > inverse:
Again, need "_rmsd" on these. So the script does nothing if the RMSDs are exactly equal (unlikely, I know)? It seems like you could just use a simple "else:" here instead of the elif.
print "Minimum RMSD =",inverse_rmsd,"\n\n" runCommand("match #"+str(i)+right+" #0"+right+"; close 1")
In principle this script should show the matching with the minimum RMSD, but it doesn't work.
Please, help... Thanks
--Eric Eric Pettersen UCSF Computer Graphics Lab http://www.cgl.ucsf.edu