Created
March 18, 2009 11:54
-
-
Save baoilleach/81077 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import math | |
import pybel | |
def sqr_dist(a, b): | |
ac = a.coords | |
bc = b.coords | |
return (ac[0]-bc[0])**2 + (ac[1]-bc[1])**2 + (ac[2]-bc[2])**2 | |
# Definitions taken from | |
# http://baoilleach.blogspot.com/2007/07/pybel-hack-that-sd-file.html | |
# I should fix them to incorporate GMC2007's comments on that page | |
HBD = pybel.Smarts("[!#6;!H0]") | |
HBA = pybel.Smarts("[$([$([#8,#16]);!$(*=N~O);" + | |
"!$(*~N=O);X1,X2]),$([#7;v3;" + | |
"!$([nH]);!$(*(-a)-a)])]") | |
HBOND_CUTOFF = 3.1 | |
SQR_HBOND_CUTOFF = HBOND_CUTOFF**2 | |
if __name__ == "__main__": | |
protein = pybel.readfile("pdb", "1HV6_macromol.pdbqt").next() | |
ligand = pybel.readfile("pdb", "1hv6_trimer_001_model_.pdb").next() | |
protein_atoms = protein.atoms | |
ligand_atoms = ligand.atoms | |
# Protein donor, ligand acceptor | |
for protein_donor in HBD.findall(protein): | |
protein_atom = protein_atoms[protein_donor[0] - 1] | |
for ligand_acceptor in HBA.findall(ligand): | |
ligand_atom = ligand_atoms[ligand_acceptor[0] - 1] | |
mysqr_dist = sqr_dist(protein_atom, ligand_atom) | |
if mysqr_dist < SQR_HBOND_CUTOFF: | |
print "H Bond between protein don (%d) and ligand acc (%d) (%.1f Ang)" % ( | |
protein_donor[0], ligand_acceptor[0], math.sqrt(mysqr_dist)) | |
# Protein acceptor, ligand donor | |
for protein_acceptor in HBA.findall(protein): | |
protein_atom = protein_atoms[protein_acceptor[0] - 1] | |
for ligand_donor in HBD.findall(ligand): | |
ligand_atom = ligand_atoms[ligand_donor[0] - 1] | |
mysqr_dist = sqr_dist(protein_atom, ligand_atom) | |
if mysqr_dist < SQR_HBOND_CUTOFF: | |
print "H Bond between protein acc (%d) and ligand don (%d) (%.1f Ang)" % ( | |
protein_acceptor[0], ligand_donor[0], math.sqrt(mysqr_dist)) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment