From e363f37ce07e93343667d1d8a9a539ca4d535dfc Mon Sep 17 00:00:00 2001
From: Thomas Forest <thomas.forest@etu.univ-paris-diderot.fr>
Date: Sat, 14 Sep 2019 23:41:58 +0200
Subject: [PATCH] Improved bridges exhaustive search

---
 src/atom.py      |  5 +----
 src/dssp.py      |  2 +-
 src/geometry.py  |  2 +-
 src/structure.py | 46 ++++++++++++++++++++++++++++------------------
 4 files changed, 31 insertions(+), 24 deletions(-)

diff --git a/src/atom.py b/src/atom.py
index 9337f02..64fb41d 100644
--- a/src/atom.py
+++ b/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)
diff --git a/src/dssp.py b/src/dssp.py
index ab0c5e3..795d069 100644
--- a/src/dssp.py
+++ b/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))
 
diff --git a/src/geometry.py b/src/geometry.py
index b0c716f..b0b327b 100644
--- a/src/geometry.py
+++ b/src/geometry.py
@@ -3,4 +3,4 @@ class GeomObject:
     pass
 
 class Distance:
-    
+    pass
diff --git a/src/structure.py b/src/structure.py
index a6fecb5..7cf5a7c 100644
--- a/src/structure.py
+++ b/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):