128 lines
4.8 KiB
Python
128 lines
4.8 KiB
Python
# pymol -c generate1.py
|
|
|
|
# NOTE: obsolete -- PyMOL now uses Dunbrack rotamers by default
|
|
|
|
from chempy import io
|
|
from glob import glob
|
|
|
|
sys.path.append(".")
|
|
|
|
from pymol import cmd
|
|
|
|
import pymol
|
|
from pymol import cmd
|
|
|
|
def find_torsions(res_sele):
|
|
# this routine gets a list of atom ids which correspond to the
|
|
# free torsions in the sidechain
|
|
res_sele = "("+res_sele+")"
|
|
if cmd.count_atoms("(%s and n;N,CA)"%res_sele,quiet=1)<2:
|
|
return []
|
|
result = []
|
|
done = 0
|
|
# set up the starting trunk
|
|
cmd.select("trunk","(%s and n;N,CA,C,O)"%res_sele)
|
|
first_id = cmd.identify("(%s and n;N)"%res_sele)[0]
|
|
second_id = cmd.identify("(%s and n;CA)"%res_sele)[0]
|
|
to_do_list = [(first_id,second_id)]
|
|
while len(to_do_list):
|
|
(first_id,second_id) = to_do_list.pop(0)
|
|
# find candidates for the 3rd atom in the dihedral
|
|
cmd.select("candi","(%s and (neighbor trunk))"%res_sele)
|
|
pymol.stored.list = []
|
|
cmd.iterate("candi","stored.list.append((name,ID))")
|
|
candi_list = pymol.stored.list
|
|
if len(candi_list):
|
|
# we may have a dihedral
|
|
dihe_list = []
|
|
for a in candi_list:
|
|
if cmd.select("termi","(%s and (not trunk) and (neighbor id %d))"%
|
|
(res_sele,a[1]))>0:
|
|
# apparently we do, but is it acyclic
|
|
# (i.e., can edit split the molecule?)
|
|
cmd.edit("(%s and trunk and (neighbor id %d))"%(res_sele,a[1]),
|
|
"(%s and id %d)"%(res_sele,a[1]))
|
|
if 'pkfrag2' in cmd.get_names('selections'):
|
|
# yes, so this is a valid dihedral
|
|
dihe_list.append(a)
|
|
else:
|
|
# no, so add the third atom to the trunk
|
|
cmd.select("trunk","(trunk or (%s and id %d))"%(res_sele,a[1]))
|
|
else:
|
|
# no fourth atom, so add this atom to the trunk
|
|
cmd.select("trunk","(trunk or (%s and id %d))"%(res_sele,a[1]))
|
|
if len(dihe_list):
|
|
# choose 3rd atom using alphanumeric order
|
|
dihe_list.sort()
|
|
third_id = dihe_list[0][1]
|
|
# if there is another third atom, then repeat later
|
|
if len(dihe_list)>1:
|
|
to_do_list.insert(0,(first_id,second_id))
|
|
# now choose the 4th atom, which we know exists, using a similar criterion
|
|
cmd.select("termi","(%s and (not trunk) and (neighbor id %d))"%
|
|
(res_sele,third_id))
|
|
pymol.stored.list=[]
|
|
cmd.iterate("termi","stored.list.append((name,ID))")
|
|
termi_list = pymol.stored.list
|
|
termi_list.sort()
|
|
fourth_id = termi_list[0][1]
|
|
# at this point, we should have a complete dihedral
|
|
# add the third atom into the trunk, and store the second and third for
|
|
# outward extension later on
|
|
to_do_list.append((second_id,third_id))
|
|
cmd.select("trunk","(trunk or (%s and id %d))"%(res_sele,third_id))
|
|
# cmd.show('sticks','(%s and (id %d,%d,%d,%d))'%
|
|
# (res_sele,first_id,second_id,third_id,fourth_id))
|
|
result.append((first_id,second_id,third_id,fourth_id))
|
|
return result
|
|
|
|
def measure_all_torsions(sele,prot_dict):
|
|
id_list = cmd.identify("((%s) and name ca)"%sele)
|
|
for a in id_list:
|
|
res_sele = "(alt '' and (byres id %d) and (%s))"%(a,sele)
|
|
lst = find_torsions(res_sele)
|
|
if len(lst):
|
|
cmd.iterate("(%s and n;ca)"%res_sele,"stored.resn=resn")
|
|
resn = pymol.stored.resn
|
|
print(resn)
|
|
name_list = []
|
|
for a in lst:
|
|
pymol.stored.names = []
|
|
for b in a:
|
|
cmd.iterate("(%s and id %d)"%(res_sele,b),"stored.names.append(name)")
|
|
name_list.append(tuple(pymol.stored.names))
|
|
res_dict = {}
|
|
for a in name_list:
|
|
dihe = cmd.get_dihedral("(%s and name %s)"%(res_sele,a[0]),
|
|
"(%s and name %s)"%(res_sele,a[1]),
|
|
"(%s and name %s)"%(res_sele,a[2]),
|
|
"(%s and name %s)"%(res_sele,a[3]))
|
|
res_dict[a]=dihe
|
|
if resn not in prot_dict:
|
|
prot_dict[resn]=[]
|
|
prot_dict[resn].append(res_dict)
|
|
|
|
sc_dict = {}
|
|
|
|
cmd.feedback('disable','selector','result')
|
|
cmd.feedback('disable','executive','action')
|
|
|
|
for fyle in glob("source/*"):
|
|
cmd.load(fyle)
|
|
cmd.remove("(het)")
|
|
cmd.remove("(hydro or not alt '')")
|
|
cmd.remove("(byres b>25)") # only use well-defined residues
|
|
cmd.remove("(not (byres name ca))")
|
|
cmd.remove("(not (byres name n))")
|
|
cmd.remove("(not (byres name c))")
|
|
cmd.refresh()
|
|
measure_all_torsions("all",sc_dict)
|
|
cmd.delete("all")
|
|
|
|
#for a in sc_dict.keys():
|
|
# print a,dict[a]
|
|
#print len(dict.keys())
|
|
|
|
io.pkl.toFile(sc_dict,"sc_raw.pkl")
|
|
|