comparison spring_package/Energy.py @ 39:172398348efd draft

"planemo upload commit 26b4018c88041ee0ca7c2976e0a012015173d7b6-dirty"
author guerler
date Fri, 22 Jan 2021 15:50:27 +0000
parents 0be0af9e695d
children
comparison
equal deleted inserted replaced
38:80a4b98121b6 39:172398348efd
1 import math 1 import math
2 from os.path import dirname, realpath
2 3
3 NTYPE = 21 4 NTYPE = 21
4 NDIST = 20 5 NDIST = 20
5 NSCALE = 2.0 6 NSCALE = 2.0
6 7
7 8
8 class Energy: 9 class Energy:
9 def __init__(self): 10 def __init__(self):
10 self.dfire = list() 11 self.dfire = list()
11 with open("spring_package/dfire/dfire.txt") as file: 12 dirPath = dirname(realpath(__file__))
13 with open("%s/Energy.data" % dirPath) as file:
12 for line in file: 14 for line in file:
13 self.dfire.append(float(line)) 15 self.dfire.append(float(line))
14 16
15 def get(self, moleculeA, moleculeB): 17 def get(self, residuesA, residuesB):
16 result = 0 18 result = 0
17 chainA = list(moleculeA.calpha.keys())[0] 19 for atomA in residuesA:
18 chainB = list(moleculeB.calpha.keys())[0] 20 indexA = self.toResCode(atomA["alignedResidue"])
19 for i in moleculeA.calpha[chainA]: 21 for atomB in residuesB:
20 atomA = moleculeA.calpha[chainA][i] 22 indexB = self.toResCode(atomB["alignedResidue"])
21 indexA = self.toResCode(atomA["residue"]) 23 dist2 = ((atomA["x"] - atomB["x"]) ** 2 +
22 for j in moleculeB.calpha[chainB]: 24 (atomA["y"] - atomB["y"]) ** 2 +
23 atomB = moleculeB.calpha[chainB][j] 25 (atomA["z"] - atomB["z"]) ** 2)
24 indexB = self.toResCode(atomB["residue"])
25 dist2 = (atomA["x"] - atomB["x"]) ** 2 + (atomA["y"] - atomB["y"]) ** 2 + (atomA["z"] - atomB["z"]) ** 2
26 dist = int((math.sqrt(dist2) * NSCALE)) 26 dist = int((math.sqrt(dist2) * NSCALE))
27 if dist < NDIST: 27 if dist < NDIST:
28 index = indexA * NTYPE * NDIST + indexB * NDIST + dist 28 index = indexA * NTYPE * NDIST + indexB * NDIST + dist
29 result = result + self.dfire[index] 29 result = result + self.dfire[index]
30 return result 30 return result
31 31
32 def getClashes(self, moleculeA, moleculeB, minDist=5.0):
33 minDist = minDist ** 2
34 clashes = 0
35 chainA = list(moleculeA.calpha.keys())[0]
36 chainB = list(moleculeB.calpha.keys())[0]
37 calphaA = moleculeA.calpha[chainA]
38 calphaB = moleculeB.calpha[chainB]
39 lenA = len(calphaA.keys())
40 lenB = len(calphaB.keys())
41 if lenA > lenB:
42 temp = calphaB
43 calphaB = calphaA
44 calphaA = temp
45 lenA = len(calphaA.keys())
46 lenB = len(calphaB.keys())
47 for i in calphaA:
48 atomA = calphaA[i]
49 for j in calphaB:
50 atomB = calphaB[j]
51 dist2 = ((atomA["x"] - atomB["x"]) ** 2 +
52 (atomA["y"] - atomB["y"]) ** 2 +
53 (atomA["z"] - atomB["z"]) ** 2)
54 if dist2 < minDist:
55 clashes = clashes + 1
56 break
57 return clashes / float(lenA)
58
32 def toResCode(self, seq): 59 def toResCode(self, seq):
33 code = dict(ALA=0, CYS=1, ASP=2, GLU=3, PHE=4, GLY=5, HIS=6, ILE=7, LYS=8, LEU=9, MET=10, 60 code = dict(A=0, C=1, D=2, E=3, F=4, G=5, H=6, I=7, K=8, L=9, M=10,
34 ASN=11, PRO=12, GLN=13, ARG=14, SER=15, THR=16, VAL=17, TRP=18, TYR=19) 61 N=11, P=12, Q=13, R=14, S=15, T=16, V=17, W=18, Y=19)
35 return code[seq] if seq in code else 20 62 return code[seq] if seq in code else 20