I want to use python script to adjust the bond angle randomly and only accept the case if the random new angle is larger than old one.
I use "mcopy" to recover the old angle if the new is smaller. But I encounter a problem: if I use "mcopy #old #new setting a" or "x" or "y", the after "adjust angle anglenew" command, the angle adjusted is not the anglenew, but a different one. I don't know why.
I hope anyone can solve my problem..
Below is the python script and init.mol2.
import os
import time
import chimera
import random
from chimera import runCommand as rc
from chimera.tkgui import saveReplyLog
import MoleculeTransform
from MoleculeTransform import transform_atom_coordinates
import numpy as np
import sys
import re
#open initial conformation init.mol2
opened = chimera.openModels.open("init.mol2")# opened as model 0 in chimera
chimera.openModels.open("init.mol2") # opened as model 1 in chimera, used to store old conformation in monte carlo simulation
rho = opened[0]
print "angleold anglenew angle_after_adjust"
for pdbi in range(0,20):
random_atom = 1
an0 = rho.atoms[random_atom-1]
an1 = rho.atoms[random_atom]
an2 = rho.atoms[random_atom+1]
angleold = chimera.angle(an0.xformCoord(),an1.xformCoord(),an2.xformCoord())
anglenew = random.random()*180
rc("adjust angle %f %s %s %s"%(anglenew,an0.oslIdent(),an1.oslIdent(),an2.oslIdent()))
angleafter = chimera.angle(an0.xformCoord(),an1.xformCoord(),an2.xformCoord())
print str(angleold)+" "+str(anglenew)+" "+str(angleafter)
if anglenew>angleold:
rc("mcopy #0 #1 settings x")
else:
rc("mcopy #1 #0 settings x")
rc("close all")
rc("stop now")
@<TRIPOS>MOLECULE
init.mol2
3 2 3 0 0
NUCLEIC_ACID
NO_CHARGES
@<TRIPOS>ATOM
1 C 1.7324 -0.5083 0.8663 C.3 1 DA 0.0000
2 C -0.6180 0.4940 -0.0060 C.3 2 DA 0.0000
3 C 1.2450 0.9460 1.8950 C.3 3 DA 0.0000
@<TRIPOS>BOND
1 1 2 1
2 2 3 1
@<TRIPOS>SUBSTRUCTURE
1 DA 1 RESIDUE 4 A DA 1 ROOT
2 DA 2 RESIDUE 4 A DA 2
3 DA 3 RESIDUE 4 A DA 1