Improved bridges exhaustive search

src/atom.py

@@ -47,7 +47,4 @@ class Residue:
         r_CN = self.atoms["C"].dist_atoms(res2.atoms["N"])
         # Electrostatic interaction energy, in kcal/mole
         E = q1*q2*(1/r_ON + 1/r_CH - 1/r_OH - 1/r_CN)*f
-        if(E<-0.5):
-            return(True)
-        else:
-            return(False)
+        return(E)

src/dssp.py

@@ -14,5 +14,5 @@ else:
     #print(pdb_file.residues[2].h_bond(pdb_file.residues[40]))
     #print(get_turns(pdb_file.residues))
     #print(pdb_file.residues[27].h_bond(pdb_file.residues[28]))
-    get_bridges(pdb_file.residues)
+    print(get_bridges(pdb_file.residues))
 

src/geometry.py

@@ -3,4 +3,4 @@ class GeomObject:
     pass
 
 class Distance:
-    
+    pass

src/structure.py

@@ -19,30 +19,40 @@ def get_turns(residues):
     for i,res in enumerate(residues):
         for j in range(3,6):
             if(i+j<len(residues)):
-                if(res.h_bond(residues[i+j])):
+                if(res.h_bond(residues[i+j]<-0.5)):
                     turns.append(Turn(j))
     return(turns)
 
 def get_bridges(residues):
-    # TODO : non-overlaping check to add
     bridges = []
-    for i,res in enumerate(residues):
-        for j, res in enumerate(residues):
-            if(i-1>=0 and i+1<len(residues)
-               and j-1>=0 and j+1<len(residues)):
-                if((residues[i-1].h_bond(residues[j])
-                    and residues[j].h_bond(residues[i+1]))
-                   or(residues[j-1].h_bond(residues[i])
-                      and residues[i].h_bond(residues[j+1]))):
-                    bridges.append(Bridge("para"))
+    for i in range(1,len(residues)-4):
+        E_min = 0
+        for j in range(i+2,len(residues)-1):
+            # select triplet with the minimal energy 
+
+            if(residues[i-1].h_bond(residues[j])<-0.5
+               and residues[j].h_bond(residues[i+1])<-0.5):
+                if(residues[i-1].h_bond(residues[j])+
+                   residues[j].h_bond(residues[i+1]))<E_min:
+                   E_min = residues[i-1].h_bond(residues[j])
+                   +residues[j].h_bond(residues[i+1])
+                   
+            if(residues[j-1].h_bond(residues[i])<-0.5
+               and residues[i].h_bond(residues[j+1])<-0.5):
+                if(residues[j-1].h_bond(residues[i])+
+                   residues[i].h_bond(residues[j+1]))<E_min:
+                   E_min = residues[j-1].h_bond(residues[i])
+                   +residues[i].h_bond(residues[j+1])
+                #bridges.append(Bridge("para"))
                     
-            if(i-1>=0 and i+1<len(residues)
-               and j-1>=0 and j+1<len(residues)):
-                if((residues[i].h_bond(residues[j])
-                    and residues[j].h_bond(residues[i]))
-                   or(residues[i-1].h_bond(residues[j+1])
-                      and residues[j-1].h_bond(residues[i+1]))):
-                    bridges.append(Bridge("anti"))
+                
+            if((residues[i].h_bond(residues[j])<-0.5
+                and residues[j].h_bond(residues[i])<-0.5)
+               or(residues[i-1].h_bond(residues[j+1])<-0.5
+                  and residues[j-1].h_bond(residues[i+1])<-0.5)):
+
+                #bridges.append(Bridge("anti"))
+        #bridges.append(strongest_bridge)
     return(bridges)
 
 def get_helix(turns):