Working phi_psi

src/dssp.py

@@ -23,19 +23,17 @@ else:
 
     turns = get_turns(pdb_file.residues)
     get_helix(pdb_file.residues, turns)
-    # get_bends(pdb_file.residues)
+    get_bends(pdb_file.residues)
     # get_chirality(pdb_file.residues)
-    bridges = get_bridges(pdb_file.residues)
-    ladders = get_ladders(bridges, pdb_file.residues)
-    get_sheets(ladders)
+    # bridges = get_bridges(pdb_file.residues)
+    # ladders = get_ladders(bridges, pdb_file.residues)
+    # get_sheets(ladders)  
 
     # print("NBRIDGES",len(bridges))
 
-    bridges = get_bridges(pdb_file.residues)
+    # bridges = get_bridges(pdb_file.residues)
 
-    get_bonds(pdb_file.residues)
-
-    print(pdb_file.residues[1].h_bond(pdb_file.residues[2]))
+    # get_bonds(pdb_file.residues)
 
     residues = pdb_file.residues
     e_min = 0
@@ -44,4 +42,8 @@ else:
         if ene <= e_min:
             e_min = ene
             best_res = res
-    print(residues[1].resid, best_res.resid,e_min)
+    #print(residues[1].resid, best_res.resid,e_min)
+
+    #get_phi_psi(residues)
+    #print(residues[2].atoms, residues[0].resid, residues[1].resid)
+    print(get_TCO(residues[2],residues[3]))

src/structure.py

@@ -144,8 +144,10 @@ def get_helix(residues, turns):
             while(i+k in turns):
                 k+=1
                 temp_res.append(residues[i])
+                last_res_pos = residues[i].resid+k
             if(k>=1):
-                print(k,"- HELIX at", residues[i].resid)
+                helix_size=last_res_pos - residues[i].resid
+                print(turns[i].turn_type,"- HELIX at", residues[i].resid)
                 helix.append(Helix(temp_res,residues[i].resid))
             i = i+k
         else:
@@ -167,8 +169,8 @@ def get_ladders(bridges, residues):
             ladders[i] = k-1
             print("ladder", bridges[i].res_num, bridges[i+k-1].res_num)
             ladders[i] = {'start':bridges[i].res_num,
-                                           'end':bridges[i+k-1].res_num,
-                                           'bridges':temp_bridges}
+                          'end':bridges[i+k-1].res_num,
+                          'bridges':temp_bridges}
             i+=k-1
         else:
             i+=1
@@ -201,10 +203,10 @@ def res_list(ladder):
 def get_bends(residues):
     bends = {}
     for i in range(2,len(residues)-2):
-        angle = vector_angles(vectors_substr(position_vector(residues[i].atoms["CA"].coords),
-                                             position_vector(residues[i-2].atoms["CA"].coords)),
-                              vectors_substr(position_vector(residues[i+2].atoms["CA"].coords),
-                                             position_vector(residues[i].atoms["CA"].coords)))
+        angle = math.degrees(vector_angles(vectors_substr(position_vector(residues[i].atoms["CA"].coords),
+                                                          position_vector(residues[i-2].atoms["CA"].coords)),
+                                           vectors_substr(position_vector(residues[i+2].atoms["CA"].coords),
+                                                          position_vector(residues[i].atoms["CA"].coords))))
         if(angle>70):
             print("angle", residues[i].resid, angle)
             bends[residues[i].resid] = angle
@@ -240,8 +242,7 @@ def vector_angles(v1,v2):
     norm_v2 = vector_norm(v2)
     term = dot_prod/(abs(norm_v1)*abs(norm_v2))
     rad_angle = math.acos(term)
-    deg_angle = rad_angle*(180/math.pi)
-    return deg_angle
+    return rad_angle
 
 def calc_dihedral(u1, u2, u3, u4):
     """ Calculate dihedral angle method. From bioPython.PDB
@@ -299,3 +300,30 @@ def get_chirality(residues):
         
     return chirality
 
+
+def get_phi_psi(residues):
+    for i in range(len(residues)):
+        if(i==0):
+            phi = 360.0
+        else:
+            phi = calc_dihedral(residues[i-1].atoms["C"].coords,
+                                residues[i].atoms["N"].coords,
+                                residues[i].atoms["CA"].coords,
+                                residues[i].atoms["C"].coords)
+        if(i==len(residues)-1):
+            psi = 360.0
+        else:
+            psi = calc_dihedral(residues[i].atoms["N"].coords,
+                                residues[i].atoms["CA"].coords,
+                                residues[i].atoms["C"].coords,
+                                residues[i+1].atoms["N"].coords)
+
+        print("ANGLES", i, phi, psi)
+
+def get_TCO(res1, res2):
+    CO_res1 = vector_from_pos(res1.atoms["C"].coords,
+                              res1.atoms["O"].coords)
+    CO_res2 = vector_from_pos(res2.atoms["C"].coords,
+                              res2.atoms["O"].coords)
+    angle = vector_angles(CO_res1, CO_res2)
+    return(math.cos(angle))