diff --git a/src/atom.py b/src/atom.py index a4c754c..674fe44 100644 --- a/src/atom.py +++ b/src/atom.py @@ -28,11 +28,22 @@ class Residue: self.atoms[atom.atom_name] = atom self.resid = atom.res_seq_nb self.res_name = atom.res_name + self.res_letter = self.get_amino_letter(atom.res_name) self.chain_id = atom.chain_id self.insertion_code = atom.insertion_code self.indice = indice - + def get_amino_letter(self, res_name): + code3 = ['ALA', 'ARG', 'ASN', 'ASP', 'CYS', 'GLU', 'GLN', 'GLY', + 'HIS', 'ILE', 'LEU', 'LYS', 'MET', 'PHE', 'PRO', 'SER', + 'THR', 'TRP', 'TYR', 'VAL'] + + + code1 = ['A','R','N','D','C','E','Q','G','H','I','L','K','M','F','P', + 'S','T','W','Y','V'] + + return code1[code3.index(res_name)] + def h_bond(self, res2): if("H" not in res2.atoms.keys()): return(False) @@ -69,7 +80,19 @@ class Residue: return Turn(k,residues[i].resid) return False - + + def get_bends(self, residues): + i = residues.index(self) + if i >=2 and i =1 and i < len(residues)-2): + chirality = {} + angle = calc_dihedral(residues[i-1].atoms["CA"].coords, + residues[i].atoms["CA"].coords, + residues[i+1].atoms["CA"].coords, + residues[i+2].atoms["CA"].coords) + + if(angle>0 and angle<=180): + sign="+" + + if(angle<=0 and angle > -180): + sign="-" + + else: + angle = 360.0 + sign = '' + + return [angle, sign] + + def get_phi_psi(self, residues): + i = residues.index(self) + 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) + + return((phi, psi)) def get_bridges(residues): bridges = {} bridge = {} diff --git a/src/dssp.py b/src/dssp.py index 8d72452..5171b3c 100644 --- a/src/dssp.py +++ b/src/dssp.py @@ -13,22 +13,12 @@ if(len(sys.argv)<2): "call structure and parameters.") else: pdb_file = pdb.PDBFile(sys.argv[1]) - # print(pdb_file.residues[15].atoms["C"].coord_x) - #print(pdb_file.residues[2].atoms["N"].res_seq_nb, pdb_file.residues[2].atoms["N"].coord_x, pdb_file.residues[2].atoms["N"].coord_y, pdb_file.residues[2].atoms["N"].coord_z) - #print(pdb_file.residues[2].atoms["H"].coord_x) - #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])) - #print(get_bridges(pdb_file.residues)) - - + #turns = get_turns(pdb_file.residues) #get_helix(pdb_file.residues, turns) #get_bends(pdb_file.residues) # get_chirality(pdb_file.residues) - bridges = get_bridges(pdb_file.residues) - ladders = get_ladders(bridges, pdb_file.residues) - sheets = get_sheets(ladders) + # print("NBRIDGES",len(bridges)) @@ -36,32 +26,51 @@ else: # get_bonds(pdb_file.residues) - residues = pdb_file.residues - e_min = 0 - for res in residues: - ene = residues[1].h_bond(res) - if ene <= e_min: - e_min = ene - best_res = res + + # e_min = 0 + # for res in residues: + # ene = residues[1].h_bond(res) + # if ene <= e_min: + # e_min = ene + # best_res = res #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])) - for i,res in enumerate(residues): - # res.get_turns(residues, turns) - res.get_helix(residues) - #res.get_bridges(residues) - #res.get_ladders(residues, ladders) - #res.get_sheets(residues) -for i,ladder in enumerate(ladders.values()): - print(chr(65+i)) +# for i,ladder in enumerate(ladders.values()): +# print(chr(65+i)) -for ind, sheet in sheets.items(): - print(ind, sheet) - for ladder in sheet: - print(ladder['start'], ladder['end']) +# for ind, sheet in sheets.items(): +# print(ind, sheet) +# for ladder in sheet: +# print(ladder['start'], ladder['end']) +################### OUTPUT #################### +# print DSSP-style formatted header of PDB file +for elem in pdb_file.get_header() : + print(pdb_file.get_header()[elem], end="") + +# Get preliminary data for print loop +residues = pdb_file.residues +bridges = get_bridges(residues) +ladders = get_ladders(bridges, residues) +sheets = get_sheets(ladders) + +# iterating over residues +for i,res in enumerate(residues): + #res.get_turns(residues, turns) + #res.get_helix(residues) + kappa = res.get_bends(residues) + t_co = res.get_tco(residues) + alpha = res.get_chirality(residues)[0] + phi = res.get_phi_psi(residues)[0] + psi = res.get_phi_psi(residues)[1] + x_ca = res.atoms["CA"].coord_x + y_ca = res.atoms["CA"].coord_y + z_ca = res.atoms["CA"].coord_z + print(i+1, res.resid, res.chain_id, res.res_letter, round(t_co, 3), round(kappa,1), + round(alpha, 1), round(phi, 1), round(psi,1), round(x_ca,1), round(y_ca,1), round(z_ca,1)) diff --git a/src/geom.py b/src/geom.py index f69065a..e4e5b9c 100644 --- a/src/geom.py +++ b/src/geom.py @@ -1,4 +1,5 @@ -import math +import math +import numpy as np def vector_from_pos(a, b): xAB = b[0]-a[0] yAB = b[1]-a[1]