| 1 | from chimera import openModels, Molecule
|
|---|
| 2 | from CGLutil.AdaptiveTree import AdaptiveTree
|
|---|
| 3 |
|
|---|
| 4 | IDEAL_BOND_DIST = 2.0
|
|---|
| 5 | MAX_BOND_DIST = IDEAL_BOND_DIST + 0.5
|
|---|
| 6 |
|
|---|
| 7 | for mol in openModels.list(modelTypes=[Molecule]):
|
|---|
| 8 | coords = []
|
|---|
| 9 | sgs = []
|
|---|
| 10 | for r in mol.residues:
|
|---|
| 11 | if r.type != "CYS":
|
|---|
| 12 | continue
|
|---|
| 13 | if "SG" not in r.atomsMap:
|
|---|
| 14 | continue
|
|---|
| 15 | sg = r.atomsMap["SG"][0]
|
|---|
| 16 | if len(sg.bonds) > 1:
|
|---|
| 17 | continue
|
|---|
| 18 | coords.append(sg.coord().data())
|
|---|
| 19 | sgs.append(sg)
|
|---|
| 20 | tree = AdaptiveTree(coords, sgs, 1.5 * MAX_BOND_DIST)
|
|---|
| 21 | for sg in sgs:
|
|---|
| 22 | if len(sg.bonds) > 1:
|
|---|
| 23 | continue
|
|---|
| 24 | sgcrd = sg.coord()
|
|---|
| 25 | others = tree.searchTree(sgcrd.data(), MAX_BOND_DIST)
|
|---|
| 26 | others = [o for o in others if o != sg and len(o.bonds) < 2]
|
|---|
| 27 | if not others:
|
|---|
| 28 | continue
|
|---|
| 29 | others.sort(lambda o1, o2: cmp(sgcrd.sqdistance(o1.coord()),
|
|---|
| 30 | sgcrd.sqdistance(o2.coord())))
|
|---|
| 31 | print "Forming disulphide bond between", sg, "and", others[0]
|
|---|
| 32 | sg.molecule.newBond(sg, others[0])
|
|---|
| 33 |
|
|---|